schain Commit - r1548:334daf8eec15 · Jicamarca Repository

remSatsProfiles

joabAM -

r1548:334daf8eec15

parent child

Context file:

r1548:334daf8eec15

Collapse all files

The requested changes are too big and content was truncated. Show full diff

schainpy/model/data/jrodata.py +9 -4

             # Copyright (c) 2012-2020 Jicamarca Radio Observatory
             # All rights reserved.
             #
             # Distributed under the terms of the BSD 3-clause license.
             """Definition of diferent Data objects for different types of data
             Here you will find the diferent data objects for the different types
             of data, this data objects must be used as dataIn or dataOut objects in
             processing units and operations. Currently the supported data objects are:
             Voltage, Spectra, SpectraHeis, Fits, Correlation and Parameters
             """
             import copy
             import numpy
             import datetime
             import json
             import schainpy.admin
             from schainpy.utils import log
             from .jroheaderIO import SystemHeader, RadarControllerHeader
             from schainpy.model.data import _noise
             def getNumpyDtype(dataTypeCode):
                 if dataTypeCode == 0:
                     numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
                 elif dataTypeCode == 1:
                     numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
                 elif dataTypeCode == 2:
                     numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
                 elif dataTypeCode == 3:
                     numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
                 elif dataTypeCode == 4:
                     numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
                 elif dataTypeCode == 5:
                     numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
                 else:
                     raise ValueError('dataTypeCode was not defined')
                 return numpyDtype
             def getDataTypeCode(numpyDtype):
                 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
                     datatype = 0
                 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
                     datatype = 1
                 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
                     datatype = 2
                 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
                     datatype = 3
                 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
                     datatype = 4
                 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
                     datatype = 5
                 else:
                     datatype = None
                 return datatype
             def hildebrand_sekhon(data, navg):
                 """
                 This method is for the objective determination of the noise level in Doppler spectra. This
                 implementation technique is based on the fact that the standard deviation of the spectral
                 densities is equal to the mean spectral density for white Gaussian noise
                 Inputs:
                     Data    :    heights
                     navg    :    numbers of averages
                 Return:
                     mean    :    noise's level
                 """
                 sortdata = numpy.sort(data, axis=None)
                 '''
                 lenOfData = len(sortdata)
                 nums_min = lenOfData*0.5
                 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 lnoise
                 '''
                 return _noise.hildebrand_sekhon(sortdata, navg)
             class Beam:
                 def __init__(self):
                     self.codeList = []
                     self.azimuthList = []
                     self.zenithList = []
             class GenericData(object):
                 flagNoData = True
                 def copy(self, inputObj=None):
                     if inputObj == None:
                         return copy.deepcopy(self)
                     for key in list(inputObj.__dict__.keys()):
                         attribute = inputObj.__dict__[key]
                         # If this attribute is a tuple or list
                         if type(inputObj.__dict__[key]) in (tuple, list):
                             self.__dict__[key] = attribute[:]
                             continue
                         # If this attribute is another object or instance
                         if hasattr(attribute, '__dict__'):
                             self.__dict__[key] = attribute.copy()
                             continue
                         self.__dict__[key] = inputObj.__dict__[key]
                 def deepcopy(self):
                     return copy.deepcopy(self)
                 def isEmpty(self):
                     return self.flagNoData
                 def isReady(self):
                     return not self.flagNoData
             class JROData(GenericData):
                 useInputBuffer = False
                 buffer_empty = True
                 systemHeaderObj = SystemHeader()
                 radarControllerHeaderObj = RadarControllerHeader()
                 type = None
                 datatype = None  # dtype but in string
                 nProfiles = None
                 heightList = None
                 channelList = None
                 flagDiscontinuousBlock = False
                 useLocalTime = False
                 utctime = None
                 timeZone = None
                 dstFlag = None
                 errorCount = None
                 blocksize = None
                 flagDecodeData = False  # asumo q la data no esta decodificada
                 flagDeflipData = False  # asumo q la data no esta sin flip
                 flagShiftFFT = False
                 nCohInt = None
                 windowOfFilter = 1
                 C = 3e8
                 frequency = 49.92e6
                 realtime = False
                 beacon_heiIndexList = None
                 last_block = None
                 blocknow = None
                 azimuth = None
                 zenith = None
                 beam = Beam()
                 profileIndex = None
                 error = None
                 data = None
                 nmodes = None
                 metadata_list = ['heightList', 'timeZone', 'type']
                 codeList = []
                 azimuthList = []
                 elevationList = []
                 last_noise = None
                 def __str__(self):
                     return '{} - {}'.format(self.type, self.datatime())
                 def getNoise(self):
                     raise NotImplementedError
                 @property
                 def nChannels(self):
                     return len(self.channelList)
                 @property
                 def channelIndexList(self):
                     return list(range(self.nChannels))
                 @property
                 def nHeights(self):
                     return len(self.heightList)
                 def getDeltaH(self):
                     return self.heightList[1] - self.heightList[0]
                 @property
                 def ltctime(self):
                     if self.useLocalTime:
                         return self.utctime - self.timeZone * 60
                     return self.utctime
                 @property
                 def datatime(self):
                     datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
                     return datatimeValue
                 def getTimeRange(self):
                     datatime = []
                     datatime.append(self.ltctime)
                     datatime.append(self.ltctime + self.timeInterval + 1)
                     datatime = numpy.array(datatime)
                     return datatime
                 def getFmaxTimeResponse(self):
                     period = (10**-6) * self.getDeltaH() / (0.15)
                     PRF = 1. / (period * self.nCohInt)
                     fmax = PRF
                     return fmax
                 def getFmax(self):
                     PRF = 1. / (self.ippSeconds * self.nCohInt)
                     fmax = PRF
                     return fmax
                 def getVmax(self):
                     _lambda = self.C / self.frequency
                     vmax = self.getFmax() * _lambda / 2
                     return vmax
                 @property
                 def ippSeconds(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.ippSeconds
                 @ippSeconds.setter
                 def ippSeconds(self, ippSeconds):
                     '''
                     '''
                     self.radarControllerHeaderObj.ippSeconds = ippSeconds
                 @property
                 def code(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.code
                 @code.setter
                 def code(self, code):
                     '''
                     '''
                     self.radarControllerHeaderObj.code = code
                 @property
                 def nCode(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.nCode
                 @nCode.setter
                 def nCode(self, ncode):
                     '''
                     '''
                     self.radarControllerHeaderObj.nCode = ncode
                 @property
                 def nBaud(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.nBaud
                 @nBaud.setter
                 def nBaud(self, nbaud):
                     '''
                     '''
                     self.radarControllerHeaderObj.nBaud = nbaud
                 @property
                 def ipp(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.ipp
                 @ipp.setter
                 def ipp(self, ipp):
                     '''
                     '''
                     self.radarControllerHeaderObj.ipp = ipp
                 @property
                 def metadata(self):
                     '''
                     '''
                     return {attr: getattr(self, attr) for attr in self.metadata_list}
             class Voltage(JROData):
                 dataPP_POW   = None
                 dataPP_DOP   = None
                 dataPP_WIDTH = None
                 dataPP_SNR   = None
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.useLocalTime = True
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Voltage"
                     self.data = None
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
                     self.utctime = None
                     self.timeZone = 0
                     self.dstFlag = None
                     self.errorCount = None
                     self.nCohInt = None
                     self.blocksize = None
                     self.flagCohInt = False
                     self.flagDecodeData = False  # asumo q la data no esta decodificada
                     self.flagDeflipData = False  # asumo q la data no esta sin flip
                     self.flagShiftFFT = False
                     self.flagDataAsBlock = False  # Asumo que la data es leida perfil a perfil
                     self.profileIndex = 0
                     self.metadata_list = ['type', 'heightList', 'timeZone', 'nProfiles', 'channelList', 'nCohInt',
                         'code', 'nCode', 'nBaud', 'ippSeconds', 'ipp']
                 def getNoisebyHildebrand(self, channel=None, ymin_index=None, ymax_index=None):
                     """
                     Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
                     Return:
                         noiselevel
                     """
                     if channel != None:
                         data = self.data[channel,ymin_index:ymax_index]
                         nChannels = 1
                     else:
                         data = self.data[:,ymin_index:ymax_index]
                         nChannels = self.nChannels
                     noise = numpy.zeros(nChannels)
                     power = data * numpy.conjugate(data)
                     for thisChannel in range(nChannels):
                         if nChannels == 1:
                             daux = power[:].real
                         else:
                             daux = power[thisChannel, :].real
                         noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
                     return noise
                 def getNoise(self, type=1, channel=None,ymin_index=None, ymax_index=None):
                     if type == 1:
                         noise = self.getNoisebyHildebrand(channel,ymin_index, ymax_index)
                     return noise
                 def getPower(self, channel=None):
                     if channel != None:
                         data = self.data[channel]
                     else:
                         data = self.data
                     power = data * numpy.conjugate(data)
                     powerdB = 10 * numpy.log10(power.real)
                     powerdB = numpy.squeeze(powerdB)
                     return powerdB
                 @property
                 def timeInterval(self):
                     return self.ippSeconds * self.nCohInt
                 noise = property(getNoise, "I'm the 'nHeights' property.")
             class Spectra(JROData):
                 data_outlier = 0
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.data_dc = None
                     self.data_spc = None
                     self.data_cspc = None
                     self.useLocalTime = True
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Spectra"
                     self.timeZone = 0
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
                     self.pairsList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
                     self.utctime = None
                     self.nCohInt = None
                     self.nIncohInt = None
                     self.blocksize = None
                     self.nFFTPoints = None
                     self.wavelength = None
                     self.flagDecodeData = False  # asumo q la data no esta decodificada
                     self.flagDeflipData = False  # asumo q la data no esta sin flip
                     self.flagShiftFFT = False
                     self.ippFactor = 1
                     self.beacon_heiIndexList = []
                     self.noise_estimation = None
                     self.codeList = []
                     self.azimuthList = []
                     self.elevationList = []
                     self.metadata_list = ['type', 'heightList', 'timeZone', 'pairsList', 'channelList', 'nCohInt',
                         'code', 'nCode', 'nBaud', 'ippSeconds', 'ipp','nIncohInt', 'nFFTPoints', 'nProfiles']
                     self.max_nIncohInt = 1
                 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
                     """
                     Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
                     Return:
                         noiselevel
                     """
                     # if hasattr(self.nIncohInt, "__len__"): #nIncohInt is a matrix
                     #
                     #     heis =  self.data_spc.shape[2]
                     #
                     #     noise = numpy.zeros((self.nChannels, heis))
                     #     for hei in range(heis):
                     #         for channel in range(self.nChannels):
                     #             daux = self.data_spc[channel, xmin_index:xmax_index, hei]
                     #
                     #             noise[channel,hei] = hildebrand_sekhon(daux, self.nIncohInt[channel,hei])
                     #
                     # else:
                     #     noise = numpy.zeros(self.nChannels)
                     #     for channel in range(self.nChannels):
                     #         daux = self.data_spc[channel,xmin_index:xmax_index, ymin_index:ymax_index]
                     #
                     #         noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
                     noise = numpy.zeros(self.nChannels)
                     for channel in range(self.nChannels):
                         daux = self.data_spc[channel,xmin_index:xmax_index, ymin_index:ymax_index]
                         noise[channel] = hildebrand_sekhon(daux, self.max_nIncohInt)
                     return noise
                 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
                     if self.noise_estimation is not None:
                         # this was estimated by getNoise Operation defined in jroproc_spectra.py
                         return self.noise_estimation
                     else:
                         noise = self.getNoisebyHildebrand(
                             xmin_index, xmax_index, ymin_index, ymax_index)
                         return noise
                 def getFreqRangeTimeResponse(self, extrapoints=0):
                     deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
                     freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.) - deltafreq / 2
                     return freqrange
                 def getAcfRange(self, extrapoints=0):
                     deltafreq = 10. / (self.getFmax() / (self.nFFTPoints * self.ippFactor))
                     freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
                     return freqrange
                 def getFreqRange(self, extrapoints=0):
                     deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
                     freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
                     return freqrange
                 def getVelRange(self, extrapoints=0):
                     deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
                     velrange = deltav * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.)
                     if self.nmodes:
                         return velrange/self.nmodes
                     else:
                         return velrange
                 @property
                 def nPairs(self):
                     return len(self.pairsList)
                 @property
                 def pairsIndexList(self):
                     return list(range(self.nPairs))
                 @property
                 def normFactor(self):
                     pwcode = 1
                     if self.flagDecodeData:
                         pwcode = numpy.sum(self.code[0]**2)
                     #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
                     normFactor = self.nProfiles * self.nIncohInt * self.nCohInt * pwcode * self.windowOfFilter
                     return normFactor
                 @property
                 def flag_cspc(self):
                     if self.data_cspc is None:
                         return True
                     return False
                 @property
                 def flag_dc(self):
                     if self.data_dc is None:
                         return True
                     return False
                 @property
                 def timeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
                     if self.nmodes:
                         return self.nmodes*timeInterval
                     else:
                         return timeInterval
                 def getPower(self):
                     factor = self.normFactor
-                    z = numpy.divide(self.data_spc,factor)
+                    power = numpy.zeros( (self.nChannels,self.nHeights) )
+                    for ch in range(self.nChannels):
+                        if hasattr(factor,'shape'):
+                            z = numpy.divide(self.data_spc[ch],factor[ch])
+                        else:
+                            z = numpy.divide(self.data_spc[ch],factor)
                         z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
-                    avg = numpy.average(z, axis=1)
+                        avg = numpy.average(z, axis=0)
+                        power[ch,:] = 10 * numpy.log10(avg)
-                    return 10 * numpy.log10(avg)
+                    return power
                 def getCoherence(self, pairsList=None, phase=False):
                     z = []
                     if pairsList is None:
                         pairsIndexList = self.pairsIndexList
                     else:
                         pairsIndexList = []
                         for pair in pairsList:
                             if pair not in self.pairsList:
                                 raise ValueError("Pair %s is not in dataOut.pairsList" % (
                                     pair))
                             pairsIndexList.append(self.pairsList.index(pair))
                     for i in range(len(pairsIndexList)):
                         pair = self.pairsList[pairsIndexList[i]]
                         ccf = numpy.average(self.data_cspc[pairsIndexList[i], :, :], axis=0)
                         powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
                         powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
                         avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
                         if phase:
                             data = numpy.arctan2(avgcoherenceComplex.imag,
                                                  avgcoherenceComplex.real) * 180 / numpy.pi
                         else:
                             data = numpy.abs(avgcoherenceComplex)
                         z.append(data)
                     return numpy.array(z)
                 def setValue(self, value):
                     print("This property should not be initialized", value)
                     return
                 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
             class SpectraHeis(Spectra):
                 def __init__(self):
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "SpectraHeis"
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
                     self.utctime = None
                     self.blocksize = None
                     self.profileIndex = 0
                     self.nCohInt = 1
                     self.nIncohInt = 1
                 @property
                 def normFactor(self):
                     pwcode = 1
                     if self.flagDecodeData:
                         pwcode = numpy.sum(self.code[0]**2)
                     normFactor = self.nIncohInt * self.nCohInt * pwcode
                     return normFactor
                 @property
                 def timeInterval(self):
                     return self.ippSeconds * self.nCohInt * self.nIncohInt
             class Fits(JROData):
                 def __init__(self):
                     self.type = "Fits"
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
                     self.flagNoData = True
                     self.utctime = None
                     self.nCohInt = 1
                     self.nIncohInt = 1
                     self.useLocalTime = True
                     self.profileIndex = 0
                     self.timeZone = 0
                 def getTimeRange(self):
                     datatime = []
                     datatime.append(self.ltctime)
                     datatime.append(self.ltctime + self.timeInterval)
                     datatime = numpy.array(datatime)
                     return datatime
                 def getChannelIndexList(self):
                     return list(range(self.nChannels))
                 def getNoise(self, type=1):
                     if type == 1:
                         noise = self.getNoisebyHildebrand()
                     if type == 2:
                         noise = self.getNoisebySort()
                     if type == 3:
                         noise = self.getNoisebyWindow()
                     return noise
                 @property
                 def timeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
                     return timeInterval
                 @property
                 def ippSeconds(self):
                     '''
                     '''
                     return self.ipp_sec
                 noise = property(getNoise, "I'm the 'nHeights' property.")
             class Correlation(JROData):
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Correlation"
                     self.data = None
                     self.dtype = None
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
                     self.utctime = None
                     self.timeZone = 0
                     self.dstFlag = None
                     self.errorCount = None
                     self.blocksize = None
                     self.flagDecodeData = False  # asumo q la data no esta decodificada
                     self.flagDeflipData = False  # asumo q la data no esta sin flip
                     self.pairsList = None
                     self.nPoints = None
                 def getPairsList(self):
                     return self.pairsList
                 def getNoise(self, mode=2):
                     indR = numpy.where(self.lagR == 0)[0][0]
                     indT = numpy.where(self.lagT == 0)[0][0]
                     jspectra0 = self.data_corr[:, :, indR, :]
                     jspectra = copy.copy(jspectra0)
                     num_chan = jspectra.shape[0]
                     num_hei = jspectra.shape[2]
                     freq_dc = jspectra.shape[1] / 2
                     ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
                     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 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
                     noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
                     return noise
                 @property
                 def timeInterval(self):
                     return self.ippSeconds * self.nCohInt * self.nProfiles
                 def splitFunctions(self):
                     pairsList = self.pairsList
                     ccf_pairs = []
                     acf_pairs = []
                     ccf_ind = []
                     acf_ind = []
                     for l in range(len(pairsList)):
                         chan0 = pairsList[l][0]
                         chan1 = pairsList[l][1]
                         # Obteniendo pares de Autocorrelacion
                         if chan0 == chan1:
                             acf_pairs.append(chan0)
                             acf_ind.append(l)
                         else:
                             ccf_pairs.append(pairsList[l])
                             ccf_ind.append(l)
                     data_acf = self.data_cf[acf_ind]
                     data_ccf = self.data_cf[ccf_ind]
                     return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
                 @property
                 def normFactor(self):
                     acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
                     acf_pairs = numpy.array(acf_pairs)
                     normFactor = numpy.zeros((self.nPairs, self.nHeights))
                     for p in range(self.nPairs):
                         pair = self.pairsList[p]
                         ch0 = pair[0]
                         ch1 = pair[1]
                         ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
                         ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
                         normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
                     return normFactor
             class Parameters(Spectra):
                 groupList = None  # List of Pairs, Groups, etc
                 data_param = None  # Parameters obtained
                 data_pre = None  # Data Pre Parametrization
                 data_SNR = None  # Signal to Noise Ratio
                 data_outlier = None
                 abscissaList = None  # Abscissa, can be velocities, lags or time
                 utctimeInit = None  # Initial UTC time
                 paramInterval = None  # Time interval to calculate Parameters in seconds
                 useLocalTime = True
                 # Fitting
                 data_error = None  # Error of the estimation
                 constants = None
                 library = None
                 # Output signal
                 outputInterval = None  # Time interval to calculate output signal in seconds
                 data_output = None  # Out signal
                 nAvg = None
                 noise_estimation = None
                 GauSPC = None  # Fit gaussian SPC
                 max_nIncohInt = 1
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Parameters"
                     self.timeZone = 0
                 def getTimeRange1(self, interval):
                     datatime = []
                     if self.useLocalTime:
                         time1 = self.utctimeInit - self.timeZone * 60
                     else:
                         time1 = self.utctimeInit
                     datatime.append(time1)
                     datatime.append(time1 + interval)
                     datatime = numpy.array(datatime)
                     return datatime
                 @property
                 def timeInterval(self):
                     if hasattr(self, 'timeInterval1'):
                         return self.timeInterval1
                     else:
                         return self.paramInterval
                 def setValue(self, value):
                     print("This property should not be initialized")
                     return
                 def getNoise(self):
                     return self.spc_noise
                 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
             class PlotterData(object):
                 '''
                 Object to hold data to be plotted
                 '''
                 MAXNUMX = 200
                 MAXNUMY = 200
                 def __init__(self, code, exp_code, localtime=True):
                     self.key = code
                     self.exp_code = exp_code
                     self.ready = False
                     self.flagNoData = False
                     self.localtime = localtime
                     self.data = {}
                     self.meta = {}
                     self.__heights = []
                 def __str__(self):
                     dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
                     return 'Data[{}][{}]'.format(';'.join(dum), len(self.times))
                 def __len__(self):
                     return len(self.data)
                 def __getitem__(self, key):
                     if isinstance(key, int):
                         return self.data[self.times[key]]
                     elif isinstance(key, str):
                         ret = numpy.array([self.data[x][key] for x in self.times])
                         if ret.ndim > 1:
                             ret = numpy.swapaxes(ret, 0, 1)
                         return ret
                 def __contains__(self, key):
                     return key in self.data[self.min_time]
                 def setup(self):
                     '''
                     Configure object
                     '''
                     self.type = ''
                     self.ready = False
                     del self.data
                     self.data = {}
                     self.__heights = []
                     self.__all_heights = set()
                 def shape(self, key):
                     '''
                     Get the shape of the one-element data for the given key
                     '''
                     if len(self.data[self.min_time][key]):
                         return self.data[self.min_time][key].shape
                     return (0,)
                 def update(self, data, tm, meta={}):
                     '''
                     Update data object with new dataOut
                     '''
                     self.data[tm] = data
                     for key, value in meta.items():
                         setattr(self, key, value)
                 def normalize_heights(self):
                     '''
                     Ensure same-dimension of the data for different heighList
                     '''
                     H = numpy.array(list(self.__all_heights))
                     H.sort()
                     for key in self.data:
                         shape = self.shape(key)[:-1] + H.shape
                         for tm, obj in list(self.data[key].items()):
                             h = self.__heights[self.times.tolist().index(tm)]
                             if H.size == h.size:
                                 continue
                             index = numpy.where(numpy.in1d(H, h))[0]
                             dummy = numpy.zeros(shape) + numpy.nan
                             if len(shape) == 2:
                                 dummy[:, index] = obj
                             else:
                                 dummy[index] = obj
                             self.data[key][tm] = dummy
                     self.__heights = [H for tm in self.times]
                 def jsonify(self, tm, plot_name, plot_type, decimate=False):
                     '''
                     Convert data to json
                     '''
                     meta = {}
                     meta['xrange'] = []
                     dy = int(len(self.yrange)/self.MAXNUMY) + 1
                     tmp = self.data[tm][self.key]
                     shape = tmp.shape
                     if len(shape) == 2:
                         data = self.roundFloats(self.data[tm][self.key][::, ::dy].tolist())
                     elif len(shape) == 3:
                         dx = int(self.data[tm][self.key].shape[1]/self.MAXNUMX) + 1
                         data = self.roundFloats(
                             self.data[tm][self.key][::, ::dx, ::dy].tolist())
                         meta['xrange'] = self.roundFloats(self.xrange[2][::dx].tolist())
                     else:
                         data = self.roundFloats(self.data[tm][self.key].tolist())
                     ret = {
                         'plot': plot_name,
                         'code': self.exp_code,
                         'time': float(tm),
                         'data': data,
                         }
                     meta['type'] = plot_type
                     meta['interval'] = float(self.interval)
                     meta['localtime'] = self.localtime
                     meta['yrange'] = self.roundFloats(self.yrange[::dy].tolist())
                     meta.update(self.meta)
                     ret['metadata'] = meta
                     return json.dumps(ret)
                 @property
                 def times(self):
                     '''
                     Return the list of times of the current data
                     '''
                     ret = [t for t in self.data]
                     ret.sort()
                     return numpy.array(ret)
                 @property
                 def min_time(self):
                     '''
                     Return the minimun time value
                     '''
                     return self.times[0]
                 @property
                 def max_time(self):
                     '''
                     Return the maximun time value
                     '''
                     return self.times[-1]
                 # @property
                 # def heights(self):
                 #     '''
                 #     Return the list of heights of the current data
                 #     '''
                 #     return numpy.array(self.__heights[-1])
                 @staticmethod
                 def roundFloats(obj):
                     if isinstance(obj, list):
                         return list(map(PlotterData.roundFloats, obj))
                     elif isinstance(obj, float):
                         return round(obj, 2)

schainpy/model/graphics/jroplot_spectra.py +1 -1

             # Copyright (c) 2012-2020 Jicamarca Radio Observatory
             # All rights reserved.
             #
             # Distributed under the terms of the BSD 3-clause license.
             """Classes to plot Spectra data
             """
             import os
             import numpy
             from schainpy.model.graphics.jroplot_base import Plot, plt, log
             from itertools import combinations
             from matplotlib.ticker import LinearLocator
             class SpectraPlot(Plot):
                 '''
                 Plot for Spectra data
                 '''
                 CODE = 'spc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 buffering = False
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.height = 3.4 * self.nrows
                     self.cb_label = 'dB'
                     if self.showprofile:
                         self.width = 5.2 * self.ncols
                     else:
                         self.width = 4.2* self.ncols
                     self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.12})
                     self.ylabel = 'Range [km]'
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     self.update_list(dataOut)
                     data = {}
                     meta = {}
                     spc = 10*numpy.log10(dataOut.data_spc/dataOut.normFactor)
                     data['spc'] = spc
                     data['rti'] = dataOut.getPower()
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     noise = 10*numpy.log10(dataOut.getNoise()/float(norm))
                     data['noise'] = noise
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     if self.CODE == 'spc_moments':
                         data['moments'] = dataOut.moments
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     if self.CODE == 'spc_moments':
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     z = data['spc']
                     #print(z.shape, x.shape, y.shape)
                     for n, ax in enumerate(self.axes):
                         noise = self.data['noise'][n][0]
                         #print(noise)
                         if self.CODE == 'spc_moments':
                             mean = data['moments'][n, 1]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plt_profile = self.pf_axes[n].plot(
                                     data['rti'][n], y)[0]
                                 # ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
                                 #                                     color="k", linestyle="dashed", lw=1)[0]
                             if self.CODE == 'spc_moments':
                                 ax.plt_mean = ax.plot(mean, y, color='k')[0]
                         else:
                             ax.plt.set_array(z[n].T.ravel())
                             if self.showprofile:
                                 ax.plt_profile.set_data(data['rti'][n], y)
                                 #ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
                             if self.CODE == 'spc_moments':
                                 ax.plt_mean.set_data(mean, y)
                         if len(self.azimuthList) > 0 and len(self.elevationList) > 0:
                             self.titles.append('CH {}: {:2.1f}elv {:2.1f}az {:3.2f}dB'.format(self.channelList[n], noise, self.elevationList[n], self.azimuthList[n]))
                         else:
                             self.titles.append('CH {}:  {:3.2f}dB'.format(self.channelList[n], noise))
             class CrossSpectraPlot(Plot):
                 CODE = 'cspc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 zmin_coh = None
                 zmax_coh = None
                 zmin_phase = None
                 zmax_phase = None
                 realChannels = None
                 crossPairs = None
                 def setup(self):
                     self.ncols = 4
                     self.nplots = len(self.data.pairs) * 2
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.width = 3.1 * self.ncols
                     self.height = 2.6 * self.nrows
                     self.ylabel = 'Range [km]'
                     self.showprofile = False
                     self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     spc = dataOut.data_spc
                     cspc = dataOut.data_cspc
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     rawPairs = list(combinations(list(range(dataOut.nChannels)), 2))
                     meta['pairs'] = rawPairs
                     if self.crossPairs == None:
                         self.crossPairs = dataOut.pairsList
                     tmp = []
                     for n, pair in enumerate(meta['pairs']):
                         out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
                         coh = numpy.abs(out)
                         phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
                         tmp.append(coh)
                         tmp.append(phase)
                     data['cspc'] = numpy.array(tmp)
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     cspc = data['cspc']
                     for n in range(len(self.data.pairs)):
                         pair = self.crossPairs[n]
                         coh = cspc[n*2]
                         phase = cspc[n*2+1]
                         ax = self.axes[2 * n]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, coh.T,
                                                    vmin=self.zmin_coh,
                                                    vmax=self.zmax_coh,
                                                    cmap=plt.get_cmap(self.colormap_coh)
                                                    )
                         else:
                             ax.plt.set_array(coh.T.ravel())
                         self.titles.append(
                             'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
                         ax = self.axes[2 * n + 1]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, phase.T,
                                                    vmin=-180,
                                                    vmax=180,
                                                    cmap=plt.get_cmap(self.colormap_phase)
                                                    )
                         else:
                             ax.plt.set_array(phase.T.ravel())
                         self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
             class RTIPlot(Plot):
                 '''
                 Plot for RTI data
                 '''
                 CODE = 'rti'
                 colormap = 'jet'
                 plot_type = 'pcolorbuffer'
                 titles = None
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     #print("dataChannels ",self.data.channels)
                     self.nrows = len(self.data.channels)
                     self.nplots = len(self.data.channels)
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Time'
                     self.cb_label = 'dB'
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     self.titles = ['{} Channel {}'.format(
                         self.CODE.upper(), x) for x in range(self.nplots)]
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.update_list(dataOut)
                     data = {}
                     meta = {}
                     data['rti'] = dataOut.getPower()
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     noise = 10*numpy.log10(dataOut.getNoise()/float(norm))
                     data['noise'] = noise
                     return data, meta
                 def plot(self):
                     self.x = self.data.times
                     self.y = self.data.yrange
                     #print(" x, y: ",self.x, self.y)
                     self.z = self.data[self.CODE]
                     self.z = numpy.array(self.z, dtype=float)
                     self.z = numpy.ma.masked_invalid(self.z)
                     try:
                         if self.channelList != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{} Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{} Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     except:
                         if self.channelList.any() != None:
                             self.titles = ['{} Channel {}'.format(
                                 self.CODE.upper(), x) for x in self.channelList]
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     #dummy_var = self.axes #Extrañamente esto actualiza el valor axes
                     for n, ax in enumerate(self.axes):
                         self.zmin = self.zmin if self.zmin else numpy.min(self.z)
                         self.zmax = self.zmax if self.zmax else numpy.max(self.z)
                         data = self.data[-1]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data[self.CODE][n], self.y)[0]
                                 if "noise" in self.data:
                                     ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(data['noise'][n], len(self.y)), self.y,
                                                                      color="k", linestyle="dashed", lw=1)[0]
                         else:
                             ax.collections.remove(ax.collections[0])
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data[self.CODE][n], self.y)
                                 if "noise" in self.data:
                                     ax.plot_noise.set_data(numpy.repeat(data['noise'][n], len(self.y)), self.y)
             class CoherencePlot(RTIPlot):
                 '''
                 Plot for Coherence data
                 '''
                 CODE = 'coh'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = len(self.data.pairs)
                     self.nplots = len(self.data.pairs)
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Time'
                     self.plots_adjust.update({'hspace':0.6, 'left': 0.1, 'bottom': 0.1,'right':0.95})
                     if self.CODE == 'coh':
                         self.cb_label = ''
                         self.titles = [
                             'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
                     else:
                         self.cb_label = 'Degrees'
                         self.titles = [
                             'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
                 def update(self, dataOut):
                     self.update_list(dataOut)
                     data = {}
                     meta = {}
                     data['coh'] = dataOut.getCoherence()
                     meta['pairs'] = dataOut.pairsList
                     return data, meta
             class PhasePlot(CoherencePlot):
                 '''
                 Plot for Phase map data
                 '''
                 CODE = 'phase'
                 colormap = 'seismic'
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     data['phase'] = dataOut.getCoherence(phase=True)
                     meta['pairs'] = dataOut.pairsList
                     return data, meta
             class NoisePlot(Plot):
                 '''
                 Plot for noise
                 '''
                 CODE = 'noise'
                 plot_type = 'scatterbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = 1
                     self.nplots = 1
                     self.ylabel = 'Intensity [dB]'
                     self.xlabel = 'Time'
                     self.titles = ['Noise']
                     self.colorbar = False
                     self.plots_adjust.update({'right': 0.85 })
                     #if not self.titles:
                     self.titles = ['Noise Plot']
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
-                    #noise =  10*numpy.log10(dataOut.getNoise()/norm)
                     noise =  10*numpy.log10(dataOut.getNoise())
                     noise = noise.reshape(dataOut.nChannels, 1)
                     data['noise'] = noise
                     meta['yrange'] = numpy.array([])
                     return data, meta
                 def plot(self):
                     x = self.data.times
                     xmin = self.data.min_time
                     xmax = xmin + self.xrange * 60 * 60
                     Y = self.data['noise']
                     if self.axes[0].firsttime:
                         if self.ymin is None: self.ymin = numpy.nanmin(Y) - 5
                         if self.ymax is None: self.ymax = numpy.nanmax(Y) + 5
                         for ch in self.data.channels:
                             y = Y[ch]
                             self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
                         plt.legend(bbox_to_anchor=(1.18, 1.0))
                     else:
                         for ch in self.data.channels:
                             y = Y[ch]
                             self.axes[0].lines[ch].set_data(x, y)
             class PowerProfilePlot(Plot):
                 CODE = 'pow_profile'
                 plot_type = 'scatter'
                 def setup(self):
                     self.ncols = 1
                     self.nrows = 1
                     self.nplots = 1
                     self.height = 4
                     self.width = 3
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Intensity [dB]'
                     self.titles = ['Power Profile']
                     self.colorbar = False
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     data[self.CODE] = dataOut.getPower()
                     return data, meta
                 def plot(self):
                     y = self.data.yrange
                     self.y = y
                     x = self.data[-1][self.CODE]
                     if self.xmin is None: self.xmin = numpy.nanmin(x)*0.9
                     if self.xmax is None: self.xmax = numpy.nanmax(x)*1.1
                     if self.axes[0].firsttime:
                         for ch in self.data.channels:
                             self.axes[0].plot(x[ch], y, lw=1, label='Ch{}'.format(ch))
                         plt.legend()
                     else:
                         for ch in self.data.channels:
                             self.axes[0].lines[ch].set_data(x[ch], y)
             class SpectraCutPlot(Plot):
                 CODE = 'spc_cut'
                 plot_type = 'scatter'
                 buffering = False
                 heights = []
                 channelList = []
                 maintitle = "Spectra Cuts"
                 flag_setIndex = False
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.width = 4.5 * self.ncols + 2.5
                     self.height = 4.8 * self.nrows
                     self.ylabel = 'Power [dB]'
                     self.colorbar = False
                     self.plots_adjust.update({'left':0.1, 'hspace':0.3, 'right': 0.9, 'bottom':0.08})
                     if len(self.selectedHeightsList) > 0:
                         self.maintitle = "Spectra Cut"# for %d km " %(int(self.selectedHeight))
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     self.heights = dataOut.heightList
                     #print("sels: ",self.selectedHeightsList)
                     if len(self.selectedHeightsList)>0 and not self.flag_setIndex:
                         for sel_height in self.selectedHeightsList:
                             index_list = numpy.where(self.heights >= sel_height)
                             index_list = index_list[0]
                             self.height_index.append(index_list[0])
                         #print("sels i:"", self.height_index)
                         self.flag_setIndex = True
                         #print(self.height_index)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     n0 = 10*numpy.log10(dataOut.getNoise()/float(norm))
                     spc = 10*numpy.log10(dataOut.data_spc/dataOut.normFactor)
                     noise = numpy.repeat(n0,(dataOut.nFFTPoints*dataOut.nHeights)).reshape(dataOut.nChannels,dataOut.nFFTPoints,dataOut.nHeights)
                     data['spc'] = spc - noise
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0][1:]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     z = self.data[-1]['spc']
                     #print(z.shape)
                     if len(self.height_index) > 0:
                         index = self.height_index
                     else:
                         index = numpy.arange(0, len(y), int((len(y))/9))
                     #print("inde x ", index, self.axes)
                     for n, ax in enumerate(self.axes):
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.ymin = self.ymin if self.ymin else numpy.nanmin(z)
                             self.ymax = self.ymax if self.ymax else numpy.nanmax(z)
                             ax.plt = ax.plot(x, z[n, :, index].T)
                             labels = ['Range = {:2.1f}km'.format(y[i]) for i in index]
                             self.figures[0].legend(ax.plt, labels, loc='center right', prop={'size': 8})
                             ax.minorticks_on()
                             ax.grid(which='major', axis='both')
                             ax.grid(which='minor', axis='x')
                         else:
                             for i, line in enumerate(ax.plt):
                                 line.set_data(x, z[n, :, index[i]])
                         self.titles.append('CH {}'.format(self.channelList[n]))
                     plt.suptitle(self.maintitle,  fontsize=10)
             class BeaconPhase(Plot):
                 __isConfig = None
                 __nsubplots = None
                 PREFIX = 'beacon_phase'
                 def __init__(self):
                     Plot.__init__(self)
                     self.timerange = 24*60*60
                     self.isConfig = False
                     self.__nsubplots = 1
                     self.counter_imagwr = 0
                     self.WIDTH = 800
                     self.HEIGHT = 400
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.xdata = None
                     self.ydata = None
                     self.PLOT_CODE = BEACON_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                     self.filename_phase = None
                     self.figfile = None
                     self.xmin = None
                     self.xmax = None
                 def getSubplots(self):
                     ncol = 1
                     nrow = 1
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 7
                     colspan = 6
                     self.__nsubplots = 2
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH+self.WIDTHPROF,
                                       heightplot = self.HEIGHT+self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = self.getSubplots()
                     self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
                 def save_phase(self, filename_phase):
                     f = open(filename_phase,'w+')
                     f.write('\n\n')
                     f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
                     f.write('DD MM YYYY  HH MM SS   pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n' )
                     f.close()
                 def save_data(self, filename_phase, data, data_datetime):
                     f=open(filename_phase,'a')
                     timetuple_data = data_datetime.timetuple()
                     day = str(timetuple_data.tm_mday)
                     month = str(timetuple_data.tm_mon)
                     year = str(timetuple_data.tm_year)
                     hour = str(timetuple_data.tm_hour)
                     minute = str(timetuple_data.tm_min)
                     second = str(timetuple_data.tm_sec)
                     f.write(day+' '+month+' '+year+'  '+hour+' '+minute+' '+second+'   '+str(data[0])+'   '+str(data[1])+'   '+str(data[2])+'   '+str(data[3])+'\n')
                     f.close()
                 def plot(self):
                     log.warning('TODO: Not yet implemented...')
                 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
                         xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
                         timerange=None,
                         save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
                         server=None, folder=None, username=None, password=None,
                         ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
                     if dataOut.flagNoData:
                         return dataOut
                     if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
                         return
                     if pairsList == None:
                         pairsIndexList = dataOut.pairsIndexList[:10]
                     else:
                         pairsIndexList = []
                         for pair in pairsList:
                             if pair not in dataOut.pairsList:
                                 raise ValueError("Pair %s is not in dataOut.pairsList" %(pair))
                             pairsIndexList.append(dataOut.pairsList.index(pair))
                     if pairsIndexList == []:
                         return
              #         if len(pairsIndexList) > 4:
              #             pairsIndexList = pairsIndexList[0:4]
                     hmin_index = None
                     hmax_index = None
                     if hmin != None and hmax != None:
                         indexes = numpy.arange(dataOut.nHeights)
                         hmin_list = indexes[dataOut.heightList >= hmin]
                         hmax_list = indexes[dataOut.heightList <= hmax]
                         if hmin_list.any():
                             hmin_index = hmin_list[0]
                         if hmax_list.any():
                             hmax_index = hmax_list[-1]+1
                     x = dataOut.getTimeRange()
                     thisDatetime = dataOut.datatime
                     title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
                     xlabel = "Local Time"
                     ylabel = "Phase (degrees)"
                     update_figfile = False
                     nplots = len(pairsIndexList)
                     #phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
                     phase_beacon = numpy.zeros(len(pairsIndexList))
                     for i in range(nplots):
                         pair = dataOut.pairsList[pairsIndexList[i]]
                         ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
                         powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
                         powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
                         avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
                         phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
                         if dataOut.beacon_heiIndexList:
                             phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
                         else:
                             phase_beacon[i] = numpy.average(phase)
                     if not self.isConfig:
                         nplots = len(pairsIndexList)
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         if timerange != None:
                             self.timerange = timerange
                         self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
                         if ymin == None: ymin = 0
                         if ymax == None: ymax = 360
                         self.FTP_WEI = ftp_wei
                         self.EXP_CODE = exp_code
                         self.SUB_EXP_CODE = sub_exp_code
                         self.PLOT_POS = plot_pos
                         self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
                         self.isConfig = True
                         self.figfile = figfile
                         self.xdata = numpy.array([])
                         self.ydata = numpy.array([])
                         update_figfile = True
                         #open file beacon phase
                         path = '%s%03d' %(self.PREFIX, self.id)
                         beacon_file = os.path.join(path,'%s.txt'%self.name)
                         self.filename_phase = os.path.join(figpath,beacon_file)
                         #self.save_phase(self.filename_phase)
                     #store data beacon phase
                     #self.save_data(self.filename_phase, phase_beacon, thisDatetime)
                     self.setWinTitle(title)
                     title = "Phase Plot %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                     legendlabels = ["Pair (%d,%d)"%(pair[0], pair[1]) for pair in dataOut.pairsList]
                     axes = self.axesList[0]
                     self.xdata = numpy.hstack((self.xdata, x[0:1]))
                     if len(self.ydata)==0:
                         self.ydata = phase_beacon.reshape(-1,1)
                     else:
                         self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1,1)))
                     axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
                                 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
                                 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
                                 XAxisAsTime=True, grid='both'
                                 )
                     self.draw()
                     if dataOut.ltctime >= self.xmax:
                         self.counter_imagwr = wr_period
                         self.isConfig = False
                         update_figfile = True
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime,
                               update_figfile=update_figfile)
                     return dataOut
             class NoiselessSpectraPlot(Plot):
                 '''
                 Plot for Spectra data, subtracting
                 the noise in all channels, using for
                 amisr-14 data
                 '''
                 CODE = 'noiseless_spc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 buffering = False
                 channelList = []
                 last_noise = None
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.height = 3.5 * self.nrows
                     self.cb_label = 'dB'
                     if self.showprofile:
                         self.width = 5.8 * self.ncols
                     else:
                         self.width = 4.8* self.ncols
                     self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.92, 'bottom': 0.12})
                     self.ylabel = 'Range [km]'
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                 def update(self, dataOut):
                     self.update_list(dataOut)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     n0 = 10*numpy.log10(dataOut.getNoise()/float(norm))
                     spc = 10*numpy.log10(dataOut.data_spc/dataOut.normFactor)
                     noise = numpy.repeat(n0,dataOut.nHeights).reshape(dataOut.nChannels,dataOut.nHeights)
                     data['rti'] = dataOut.getPower() - noise
                     noise = numpy.repeat(n0,(dataOut.nFFTPoints*dataOut.nHeights)).reshape(dataOut.nChannels,dataOut.nFFTPoints,dataOut.nHeights)
                     data['spc'] = spc - noise
                     # data['noise'] = noise
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     z = data['spc']
                     for n, ax in enumerate(self.axes):
                         #noise = data['noise'][n]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plt_profile = self.pf_axes[n].plot(
                                     data['rti'][n], y)[0]
                         else:
                             ax.plt.set_array(z[n].T.ravel())
                             if self.showprofile:
                                 ax.plt_profile.set_data(data['rti'][n], y)
                         self.titles.append('CH {}'.format(self.channelList[n]))
             class NoiselessRTIPlot(Plot):
                 '''
                 Plot for RTI data
                 '''
                 CODE = 'noiseless_rti'
                 colormap = 'jet'
                 plot_type = 'pcolorbuffer'
                 titles = None
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 last_noise = None
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     #print("dataChannels ",self.data.channels)
                     self.nrows = len(self.data.channels)
                     self.nplots = len(self.data.channels)
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Time'
                     self.cb_label = 'dB'
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     self.titles = ['{} Channel {}'.format(
                         self.CODE.upper(), x) for x in range(self.nplots)]
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.update_list(dataOut)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     #print("Norm: ", norm)
                     #print(dataOut.nProfiles , dataOut.max_nIncohInt ,dataOut.nCohInt, dataOut.windowOfFilter)
                     n0 = 10*numpy.log10(dataOut.getNoise()/float(norm))
                     data['noise'] = n0
                     noise = numpy.repeat(n0,dataOut.nHeights).reshape(dataOut.nChannels,dataOut.nHeights)
                     data['noiseless_rti'] = dataOut.getPower() - noise
                     return data, meta
                 def plot(self):
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['noiseless_rti']
                     self.z = numpy.array(self.z, dtype=float)
                     self.z = numpy.ma.masked_invalid(self.z)
                     try:
                         if self.channelList != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{}  Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{}  Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     except:
                         if self.channelList.any() != None:
                             self.titles = ['{} Channel {}'.format(
                                 self.CODE.upper(), x) for x in self.channelList]
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     dummy_var = self.axes #Extrañamente esto actualiza el valor axes
                     #print("plot shapes ", z.shape, x.shape, y.shape)
                     for n, ax in enumerate(self.axes):
                         self.zmin = self.zmin if self.zmin else numpy.min(self.z)
                         self.zmax = self.zmax if self.zmax else numpy.max(self.z)
                         data = self.data[-1]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['noiseless_rti'][n], self.y)[0]
                         else:
                             ax.collections.remove(ax.collections[0])
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['noiseless_rti'][n], self.y)
                                 # if "noise" in self.data:
                                 #     #ax.plot_noise.set_data(numpy.repeat(data['noise'][n], len(self.y)), self.y)
                                 #     ax.plot_noise.set_data(data['noise'][n], self.y)
             class OutliersRTIPlot(Plot):
                 '''
                 Plot for data_xxxx object
                 '''
                 CODE = 'outlier_rtc' # Range Time Counts
                 colormap = 'cool'
                 plot_type = 'pcolorbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = self.data.shape('outlier_rtc')[0]
                     self.nplots = self.nrows
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     if not self.xlabel:
                         self.xlabel = 'Time'
                     self.ylabel = 'Height [km]'
                     if not self.titles:
                         self.titles = ['Outliers Ch:{}'.format(x) for x in range(self.nrows)]
                 def update(self, dataOut):
                     data = {}
                     data['outlier_rtc'] = dataOut.data_outlier
                     meta = {}
                     return data, meta
                 def plot(self):
                     # self.data.normalize_heights()
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['outlier_rtc']
                     #self.z = numpy.ma.masked_invalid(self.z)
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     for n, ax in enumerate(self.axes):
                         self.zmax = self.zmax if self.zmax is not None else numpy.max(
                             self.z[n])
                         self.zmin = self.zmin if self.zmin is not None else numpy.min(
                             self.z[n])
                         data = self.data[-1]
                         if ax.firsttime:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['outlier_rtc'][n], self.y)[0]
                                 self.pf_axes[n].set_xlabel('')
                         else:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.collections.remove(ax.collections[0])
                             ax.plt = ax.pcolormesh(x, y, z[n].T ,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['outlier_rtc'][n], self.y)
                                 self.pf_axes[n].set_xlabel('')
             class NIncohIntRTIPlot(Plot):
                 '''
                 Plot for data_xxxx object
                 '''
                 CODE = 'integrations_rtc' # Range Time Counts
                 colormap = 'BuGn'
                 plot_type = 'pcolorbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = self.data.shape('integrations_rtc')[0]
                     self.nplots = self.nrows
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     if not self.xlabel:
                         self.xlabel = 'Time'
                     self.ylabel = 'Height [km]'
                     if not self.titles:
                         self.titles = ['Integration Ch:{}'.format(x) for x in range(self.nrows)]
                 def update(self, dataOut):
                     data = {}
                     data['integrations_rtc'] = dataOut.nIncohInt
                     meta = {}
                     return data, meta
                 def plot(self):
                     # self.data.normalize_heights()
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['integrations_rtc']
                     #self.z = numpy.ma.masked_invalid(self.z)
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     for n, ax in enumerate(self.axes):
                         self.zmax = self.zmax if self.zmax is not None else numpy.max(
                             self.z[n])
                         self.zmin = self.zmin if self.zmin is not None else numpy.min(
                             self.z[n])
                         data = self.data[-1]
                         if ax.firsttime:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['integrations_rtc'][n], self.y)[0]
                                 self.pf_axes[n].set_xlabel('')
                         else:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.collections.remove(ax.collections[0])
                             ax.plt = ax.pcolormesh(x, y, z[n].T ,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['integrations_rtc'][n], self.y)
                                 self.pf_axes[n].set_xlabel('')

schainpy/model/proc/jroproc_parameters.py 0 0

	1	NO CONTENT: modified file			NO CONTENT: modified file
The requested commit or file is too big and content was truncated. Show full diff

schainpy/model/proc/jroproc_spectra.py +12 -11

             # 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.model.data import _noise
             from schainpy.utils import log
             import matplotlib.pyplot as plt
             #from scipy.optimize import curve_fit
             from schainpy.model.io.utils import getHei_index
             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):
                     # print("fft donw")
                     """
                     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, zeroPad=False):
                     #print("run spc proc")
                     try:
                         type = self.dataIn.type.decode("utf-8")
                         self.dataIn.type = type
                     except:
                         pass
                     if self.dataIn.type == "Spectra":
                         try:
                             self.dataOut.copy(self.dataIn)
                         except Exception as e:
                             print("Error dataIn ",e)
                         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
                         #print(" volts ch,prof, h: ", self.dataIn.data.shape)
                         if self.buffer is None:
                             if not zeroPad:
                                 self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nProfiles,
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                             else:
                                 self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nFFTPoints,
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                         if self.dataIn.flagDataAsBlock:
                             nVoltProfiles = self.dataIn.data.shape[1]
                             if nVoltProfiles == nProfiles or zeroPad:
                                 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 or zeroPad:
                             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
                     elif self.dataIn.type == "Parameters":
                         self.dataOut.data_spc = self.dataIn.data_spc
                         self.dataOut.data_cspc = self.dataIn.data_cspc
                         self.dataOut.data_outlier = self.dataIn.data_outlier
                         self.dataOut.nProfiles = self.dataIn.nProfiles
                         self.dataOut.nIncohInt = self.dataIn.nIncohInt
                         self.dataOut.nFFTPoints = self.dataIn.nFFTPoints
                         self.dataOut.ippFactor = self.dataIn.ippFactor
                         self.dataOut.max_nIncohInt = self.dataIn.max_nIncohInt
                         self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                         self.dataOut.ipp = self.dataIn.ipp
                         #self.dataOut.abscissaList = self.dataIn.getVelRange(1)
                         #self.dataOut.spc_noise = self.dataIn.getNoise()
                         #self.dataOut.spc_range = (self.dataIn.getFreqRange(1) , self.dataIn.getAcfRange(1) , self.dataIn.getVelRange(1))
                         # self.dataOut.normFactor = self.dataIn.normFactor
                         if hasattr(self.dataIn, 'channelList'):
                             self.dataOut.channelList = self.dataIn.channelList
                         if hasattr(self.dataIn, 'pairsList'):
                             self.dataOut.pairsList = self.dataIn.pairsList
                             self.dataOut.groupList = self.dataIn.pairsList
                         self.dataOut.flagNoData = False
                         if hasattr(self.dataIn, 'ChanDist'): #Distances of receiver channels
                             self.dataOut.ChanDist = self.dataIn.ChanDist
                         else: self.dataOut.ChanDist = None
                         #if hasattr(self.dataIn, 'VelRange'): #Velocities range
                         #    self.dataOut.VelRange = self.dataIn.VelRange
                         #else: self.dataOut.VelRange = None
                     else:
                         raise ValueError("The type of input object {} is not valid".format(
                             self.dataIn.type))
+                    #print("spc proc Done")
                 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
             class getNoiseB(Operation):
                 __slots__ =('offset','warnings', 'isConfig', 'minIndex','maxIndex','minIndexFFT','maxIndexFFT')
                 def __init__(self):
                     Operation.__init__(self)
                     self.isConfig = False
                 def setup(self, offset=None, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None, warnings=False):
                     self.warnings = warnings
                     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):
                         if self.warnings:
                             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):
                         if self.warnings:
                             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]
                     #indices relativos a los puntos de fft, puede ser de acuerdo a velocidad o frecuencia
                     minIndexFFT = 0
                     maxIndexFFT = 0
                     # validacion de velocidades
                     indminPoint = None
                     indmaxPoint = None
                     if self.dataOut.type == 'Spectra':
                         if minVel == None and maxVel == None :
                             freqrange = self.dataOut.getFreqRange(1)
                             if minFreq == None:
                                 minFreq = freqrange[0]
                             if maxFreq == None:
                                 maxFreq = freqrange[-1]
                             if (minFreq < freqrange[0]) or (minFreq > maxFreq):
                                 if self.warnings:
                                     print('minFreq: %.2f is out of the frequency range' % (minFreq))
                                     print('minFreq is setting to %.2f' % (freqrange[0]))
                                 minFreq = freqrange[0]
                             if (maxFreq > freqrange[-1]) or (maxFreq < minFreq):
                                 if self.warnings:
                                     print('maxFreq: %.2f is out of the frequency range' % (maxFreq))
                                     print('maxFreq is setting to %.2f' % (freqrange[-1]))
                                 maxFreq = freqrange[-1]
                             indminPoint = numpy.where(freqrange >= minFreq)
                             indmaxPoint = numpy.where(freqrange <= maxFreq)
                         else:
                             velrange = self.dataOut.getVelRange(1)
                             if minVel == None:
                                 minVel = velrange[0]
                             if maxVel == None:
                                 maxVel = velrange[-1]
                             if (minVel < velrange[0]) or (minVel > maxVel):
                                 if self.warnings:
                                     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):
                                 if self.warnings:
                                     print('maxVel: %.2f is out of the velocity range' % (maxVel))
                                     print('maxVel is setting to %.2f' % (velrange[-1]))
                                 maxVel = velrange[-1]
                             indminPoint = numpy.where(velrange >= minVel)
                             indmaxPoint = numpy.where(velrange <= maxVel)
                     # 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
                     #############################################################3
                     # seleccion de indices para velocidades
                     if self.dataOut.type == 'Spectra':
                         try:
                             minIndexFFT = indminPoint[0][0]
                         except:
                             minIndexFFT = 0
                         try:
                             maxIndexFFT = indmaxPoint[0][-1]
                         except:
                             maxIndexFFT = len( self.dataOut.getFreqRange(1))
                     self.minIndex, self.maxIndex, self.minIndexFFT, self.maxIndexFFT = minIndex, maxIndex, minIndexFFT, maxIndexFFT
                     self.isConfig = True
                     self.offset = 1
                     if offset!=None:
                         self.offset = 10**(offset/10)
                     #print("config getNoiseB Done")
                 def run(self, dataOut, offset=None, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None, warnings=False):
                     self.dataOut = dataOut
                     if not self.isConfig:
                         self.setup(offset, minHei, maxHei,minVel, maxVel, minFreq, maxFreq, warnings)
                     self.dataOut.noise_estimation = None
                     noise = None
                     if self.dataOut.type == 'Voltage':
                         noise = self.dataOut.getNoise(ymin_index=self.minIndex, ymax_index=self.maxIndex)
                         #print(minIndex, maxIndex,minIndexVel, maxIndexVel)
                     elif self.dataOut.type == 'Spectra':
                         #print(self.minIndex, self.maxIndex,self.minIndexFFT, self.maxIndexFFT, self.dataOut.nIncohInt)
                         noise = numpy.zeros( self.dataOut.nChannels)
                         norm = 1
                         for channel in range( self.dataOut.nChannels):
                             if not hasattr(self.dataOut.nIncohInt,'__len__'):
                                 norm = 1
                             else:
                                 norm = self.dataOut.max_nIncohInt/self.dataOut.nIncohInt[channel, self.minIndex:self.maxIndex]
                             #print("norm nIncoh: ", norm ,self.dataOut.data_spc.shape)
                             daux =  self.dataOut.data_spc[channel,self.minIndexFFT:self.maxIndexFFT, self.minIndex:self.maxIndex]
                             daux = numpy.multiply(daux, norm)
                             #print("offset: ", self.offset, 10*numpy.log10(self.offset))
                             # noise[channel] = self.getNoiseByMean(daux)/self.offset
                             #print(daux.shape, daux)
                             #noise[channel] = self.getNoiseByHS(daux, self.dataOut.max_nIncohInt)/self.offset
                             sortdata = numpy.sort(daux, axis=None)
                             noise[channel] = _noise.hildebrand_sekhon(sortdata, self.dataOut.max_nIncohInt)/self.offset
                             # data = numpy.mean(daux,axis=1)
                             # sortdata = numpy.sort(data, axis=None)
                             # noise[channel] = _noise.hildebrand_sekhon(sortdata, self.dataOut.max_nIncohInt)/self.offset
                         #noise = self.dataOut.getNoise(xmin_index=self.minIndexFFT, xmax_index=self.maxIndexFFT, ymin_index=self.minIndex, ymax_index=self.maxIndex)
                     else:
                         noise = self.dataOut.getNoise(xmin_index=self.minIndexFFT, xmax_index=self.maxIndexFFT, ymin_index=self.minIndex, ymax_index=self.maxIndex)
                     self.dataOut.noise_estimation = noise.copy() # dataOut.noise
                     #print("2: ",10*numpy.log10(self.dataOut.noise_estimation/64))
+                    #print("2: ",self.dataOut.noise_estimation)
                     #print(self.dataOut.flagNoData)
                     #print("getNoise Done", noise)
                     return self.dataOut
                 def getNoiseByMean(self,data):
                     #data debe estar ordenado
                     data = numpy.mean(data,axis=1)
                     sortdata = numpy.sort(data, axis=None)
                     #sortID=data.argsort()
                     #print(data.shape)
                     pnoise = None
                     j = 0
                     mean = numpy.mean(sortdata)
                     min = numpy.min(sortdata)
                     delta = mean - min
                     indexes = numpy.where(sortdata > (mean+delta))[0] #only array of indexes
                     #print(len(indexes))
                     if len(indexes)==0:
                         pnoise = numpy.mean(sortdata)
                     else:
                         j = indexes[0]
                         pnoise = numpy.mean(sortdata[0:j])
                     # from matplotlib import pyplot as plt
                     # plt.plot(sortdata)
                     # plt.vlines(j,(pnoise-delta),(pnoise+delta), color='r')
                     # plt.show()
                     #print("noise: ", 10*numpy.log10(pnoise))
                     return pnoise
                 def getNoiseByHS(self,data, navg):
                     #data debe estar ordenado
                     #data = numpy.mean(data,axis=1)
                     sortdata = numpy.sort(data, axis=None)
                     lenOfData = len(sortdata)
                     nums_min = lenOfData*0.2
                     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
                         #sumq -= sump**2
                         if j > nums_min:
                             rtest = float(j)/(j-1) + 1.0/navg
                             #if ((sumq*j) > (sump**2)):
                             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 lnoise
             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
                     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("runing cleanRayleigh")
                     if not self.isConfig :
                         self.setup(dataOut, min_hei,max_hei,n,timeInterval,factor_stdv)
                     tini=dataOut.utctime
                     if self.byProfiles:
                         if self.__profIndex ==  self.nIntProfiles:
                             self.__dataReady = True
                     else:
                         if (tini - self.__initime) >= self.__integrationtime:
                             self.__dataReady = True
                             self.__initime = tini
                     #if (tini.tm_min % 2) == 0 and (tini.tm_sec < 5 and self.fint==0):
                     if self.__dataReady:
                         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))
                             try:
                                 jcspc= numpy.reshape(jcspc,(int(len(jcspc)/self.nPairs),self.nPairs,self.nFFTPoints,self.nHeights))
                             except:
                                 print("no cspc")
                             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)
                             try:
                                 self.buffer2 = numpy.concatenate((self.buffer2,dataOut.data_cspc), axis=0)
                                 self.buffer3 += dataOut.data_dc
                             except:
                                 pass
                         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.max_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))
                     if raxs == 0:
                         raxs = 1
                     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_ints = None #matriz de numero de integracions (CH,HEI)
                 n = None
                 minHei_ind = None
                 maxHei_ind = None
                 navg = 1.0
                 factor = 0.0
                 dataoutliers = None # (CHANNELS, HEIGHTS)
                 def __init__(self):
                     Operation.__init__(self)
                 def setup(self, dataOut,n=None, timeInterval=None, overlapping=False, DPL=None, minHei=None, maxHei=None, avg=1,factor=0.75):
                     """
                     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.factor = factor
                     self.navg = avg
                     #self.ByLags = dataOut.ByLags   ###REDEFINIR
                     self.ByLags = False
                     self.maxProfilesInt = 1
+                    self.__nChannels = dataOut.nChannels
                     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
                     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]
                     ind_list1 = numpy.where(self.dataOut.heightList >= minHei)
                     ind_list2 = numpy.where(self.dataOut.heightList <= maxHei)
                     self.minHei_ind = ind_list1[0][0]
                     self.maxHei_ind = ind_list2[0][-1]
+                    #print("setup rem sats done")
                 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:
+                    if self.__nChannels < 2:
                         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,sortdata,navg, factor):
                     #data debe estar ordenado
                     #sortdata = numpy.sort(data, axis=None)
                     #sortID=data.argsort()
                     lenOfData = len(sortdata)
                     nums_min = lenOfData*factor
                     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
                     #print("H S done")
                     #return j,sortID
                     return j
                 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
                     #print("aes: ", self.__buffer_cspc)
                     self.__buffer_spc=numpy.array(self.__buffer_spc)
-                    if self.__buffer_cspc != None:
+                    if self.__nChannels > 1 :
                         self.__buffer_cspc=numpy.array(self.__buffer_cspc)
                     #print("FREQ_DC",self.__buffer_spc.shape,self.__buffer_cspc.shape)
                     freq_dc = int(self.__buffer_spc.shape[2] / 2)
                     #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
                     self.dataOutliers = numpy.zeros((self.nChannels,self.nHeights)) # --> almacen de outliers
                     for k in range(self.minHei_ind,self.maxHei_ind):
-                        if self.__buffer_cspc != None:
+                        if self.__nChannels > 1:
                             buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
                         outliers_IDs_cspc=[]
                         cspc_outliers_exist=False
                         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): #frecuencias en el tiempo
                                 # 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]
                                 sortdata = numpy.sort(buffer, axis=None)
                                 sortID=buffer.argsort()
                                 index = _noise.hildebrand_sekhon2(sortdata,self.navg)
                                 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1,self.factor)
                                 # fig,ax = plt.subplots()
                                 # ax.set_title(str(k)+" "+str(j))
                                 # x=range(len(sortdata))
                                 # ax.scatter(x,sortdata)
                                 # ax.axvline(index)
                                 # plt.show()
                                 indexes.append(index)
                                 #sortIDs.append(sortID)
                                 outliers_IDs=numpy.append(outliers_IDs,sortID[index:])
                             #print("Outliers: ",outliers_IDs)
                             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)
                             #print(indexmin,buffer1.shape[0], k)
                             # fig,ax = plt.subplots()
                             # ax.plot(sortdata)
                             # ax2 = ax.twinx()
                             # x=range(len(indexes))
                             # #plt.scatter(x,indexes)
                             # ax2.scatter(x,indexes)
                             # plt.show()
                             if indexmin != buffer1.shape[0]:
-                                if self.nChannels > 1:
+                                if self.__nChannels > 1:
                                     cspc_outliers_exist= True
                                 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
                                     buffer1[p,:] = numpy.NaN
                                 self.dataOutliers[i,k] = len(outliers_IDs)
                             self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
-                            if self.__buffer_cspc != None:
+                            if self.__nChannels > 1:
                                 outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
-                        if self.__buffer_cspc != None:
+                        if self.__nChannels > 1:
                             outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
                         if cspc_outliers_exist:
                             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
                                 buffer_cspc[p,:] = numpy.NaN
-                        if self.__buffer_cspc != None:
+                        if self.__nChannels > 1:
                             self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
                     nOutliers = len(outliers_IDs)
                     #print("Outliers  n: ",self.dataOutliers,nOutliers)
                     buffer=None
                     bufferH=None
                     buffer1=None
                     buffer_cspc=None
                     buffer=None
                     #data_spc = numpy.sum(self.__buffer_spc,axis=0)
                     data_spc = numpy.nansum(self.__buffer_spc,axis=0)
-                    if self.__buffer_cspc != None:
+                    if self.__nChannels > 1:
                         #data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
                         data_cspc = numpy.nansum(self.__buffer_cspc,axis=0)
                     else:
                         data_cspc = None
                     data_dc = self.__buffer_dc
                     #(CH, HEIGH)
                     self.maxProfilesInt = self.__profIndex
                     n = self.__profIndex - self.dataOutliers # n becomes a matrix
                     self.__buffer_spc = []
                     self.__buffer_cspc = []
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     #print("cleaned ",data_cspc)
                     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_ints = 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_ints = 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
                     #print("integrate", avgdata_cspc)
                     return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
                 def run(self, dataOut, n=None, DPL = None,timeInterval=None, overlapping=False, minHei=None, maxHei=None, avg=1, factor=0.75):
                     self.dataOut = dataOut
                     if n == 1:
                         return self.dataOut
                     #print("nchannels", self.dataOut.nChannels)
                     if self.dataOut.nChannels == 1:
                         self.dataOut.data_cspc = None #si es un solo canal no vale la pena acumular DATOS
                     #print("IN spc:", self.dataOut.data_spc.shape, self.dataOut.data_cspc)
                     if not self.isConfig:
                         self.setup(self.dataOut, n, timeInterval, overlapping,DPL ,minHei, maxHei, avg, factor)
                         self.isConfig = True
                     if not self.ByLags:
                         self.nProfiles=self.dataOut.nProfiles
                         self.nChannels=self.dataOut.nChannels
                         self.nHeights=self.dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(self.dataOut.utctime,
                                                                                             self.dataOut.data_spc,
                                                                                             self.dataOut.data_cspc,
                                                                                             self.dataOut.data_dc)
                     else:
                         self.nProfiles=self.dataOut.nProfiles
                         self.nChannels=self.dataOut.nChannels
                         self.nHeights=self.dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(self.dataOut.utctime,
                                                                                             self.dataOut.dataLag_spc,
                                                                                             self.dataOut.dataLag_cspc,
                                                                                             self.dataOut.dataLag_dc)
                     self.dataOut.flagNoData = True
                     if self.__dataReady:
                         if not self.ByLags:
                             if self.nChannels == 1:
                                 #print("f int", avgdata_spc.shape)
                                 self.dataOut.data_spc = avgdata_spc
                                 self.dataOut.data_cspc = None
                             else:
                                 self.dataOut.data_spc = numpy.squeeze(avgdata_spc)
                                 self.dataOut.data_cspc = numpy.squeeze(avgdata_cspc)
                             self.dataOut.data_dc = avgdata_dc
                             self.dataOut.data_outlier = self.dataOutliers
                         else:
                             self.dataOut.dataLag_spc = avgdata_spc
                             self.dataOut.dataLag_cspc = avgdata_cspc
                             self.dataOut.dataLag_dc = avgdata_dc
                             self.dataOut.data_spc=self.dataOut.dataLag_spc[:,:,:,self.dataOut.LagPlot]
                             self.dataOut.data_cspc=self.dataOut.dataLag_cspc[:,:,:,self.dataOut.LagPlot]
                             self.dataOut.data_dc=self.dataOut.dataLag_dc[:,:,self.dataOut.LagPlot]
                         self.dataOut.nIncohInt *= self.n_ints
                         self.dataOut.max_nIncohInt = self.maxProfilesInt
                         #print(self.dataOut.max_nIncohInt)
                         self.dataOut.utctime = avgdatatime
                         self.dataOut.flagNoData = False
                         #print("Faraday Integration DONE...", self.dataOut.data_cspc)
                     #print(self.dataOut.flagNoData)
                     return self.dataOut
             class removeInterference(Operation):
                 def removeInterference3(self,  min_hei = None, max_hei = None):
                     jspectra = self.dataOut.data_spc
                     #jcspectra = self.dataOut.data_cspc
                     jnoise = self.dataOut.getNoise()
                     num_incoh = self.dataOut.max_nIncohInt
                     #print(jspectra.shape)
                     num_channel, num_prof, num_hei = jspectra.shape
                     minHei = min_hei
                     maxHei = max_hei
                     ########################################################################
                     if minHei == None or (minHei < self.dataOut.heightList[0]):
                         minHei = self.dataOut.heightList[0]
                     if maxHei == None or (maxHei > self.dataOut.heightList[-1]):
                         maxHei = self.dataOut.heightList[-1]
                     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
                     ########################################################################
                     #dataOut.max_nIncohInt * dataOut.nCohInt
                     norm = self.dataOut.nIncohInt /self.dataOut.max_nIncohInt
                     #print(norm.shape)
                     # 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, :]
                         interf = jspectra[ich, :, minIndex:maxIndex]
                         #print(interf.shape)
                         inttef = interf.mean(axis=1)
                         for hei in range(num_hei):
                             temp = jspectra[ich,:, hei]
                             temp -= inttef
                             temp += jnoise[ich]*norm[ich,hei]
                             jspectra[ich,:, hei] = temp
                     # Guardar Resultados
                     self.dataOut.data_spc = jspectra
                     #self.dataOut.data_cspc = jcspectra
                     return 1
                 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_incoh = self.dataOut.max_nIncohInt
                     #print("spc: ", jspectra.shape, jcspectra)
                     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) # a half of total ranges
                         hei_interf = numpy.asmatrix(list(range(count_hei))) + num_hei - count_hei
                         hei_interf = numpy.asarray(hei_interf)[0]
                         #print(hei_interf)
                     # 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()
                         print(hei_interf[psort[list(range(offhei_interf, nhei_interf + offhei_interf))]])
                         # 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
                         print(junkspc_interf.shape)
                         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))
                         print(cinterfid)
                         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()
                         print(indAux)
                         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, minHei=None, maxHei=None):
                     self.dataOut = dataOut
                     if mode == 1:
                         self.removeInterference(interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None)
                     elif mode == 2:
                         self.removeInterference2()
                     elif mode == 3:
                         self.removeInterference3(min_hei=minHei, max_hei=maxHei)
                     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
                 incohInt = 0
                 nOutliers = 0
                 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
                     self.incohInt = 0
                     self.nOutliers = 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
                     """
                     if data_spc.all() == numpy.nan :
                         print("nan ")
                         return
                     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
                     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
                     if dataOut.flagNoData == True:
                         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)
                     self.incohInt += dataOut.nIncohInt
                     self.nOutliers += dataOut.data_outlier
                     if self.__dataReady:
                         #print("prof: ",dataOut.max_nIncohInt,self.__profIndex)
                         dataOut.data_spc = avgdata_spc
                         dataOut.data_cspc = avgdata_cspc
                         dataOut.data_dc = avgdata_dc
                         dataOut.nIncohInt = self.incohInt
                         dataOut.data_outlier = self.nOutliers
                         dataOut.utctime = avgdatatime
                         dataOut.flagNoData = False
                         dataOut.max_nIncohInt =  self.__profIndex
                         self.incohInt = 0
                         self.nOutliers = 0
                         self.__profIndex = 0
                         #print("IncohInt Done")
                     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

schainpy/model/proc/jroproc_voltage.py +394 -96

             import sys
             import numpy,math
             from scipy import interpolate
             from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
             from schainpy.model.data.jrodata import Voltage,hildebrand_sekhon
             from schainpy.utils import log
             from schainpy.model.io.utils import getHei_index
             from time import time
             import datetime
             import numpy
             #import copy
             from schainpy.model.data import _noise
             class VoltageProc(ProcessingUnit):
                 def __init__(self):
                     ProcessingUnit.__init__(self)
                     self.dataOut = Voltage()
                     self.flip = 1
                     self.setupReq = False
                 def run(self):
                     #print("running volt proc")
                     if self.dataIn.type == 'AMISR':
                         self.__updateObjFromAmisrInput()
                     if  self.dataOut.buffer_empty:
                         if self.dataIn.type == 'Voltage':
                             self.dataOut.copy(self.dataIn)
                             #print("new volts reading")
                 def __updateObjFromAmisrInput(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
                     self.dataOut.flagNoData = self.dataIn.flagNoData
                     self.dataOut.data = self.dataIn.data
                     self.dataOut.utctime = self.dataIn.utctime
                     self.dataOut.channelList = self.dataIn.channelList
                     #self.dataOut.timeInterval = self.dataIn.timeInterval
                     self.dataOut.heightList = self.dataIn.heightList
                     self.dataOut.nProfiles = self.dataIn.nProfiles
                     self.dataOut.nCohInt = self.dataIn.nCohInt
                     self.dataOut.ippSeconds = self.dataIn.ippSeconds
                     self.dataOut.frequency = self.dataIn.frequency
                     self.dataOut.azimuth = self.dataIn.azimuth
                     self.dataOut.zenith = self.dataIn.zenith
                     self.dataOut.beam.codeList = self.dataIn.beam.codeList
                     self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
                     self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
             class selectChannels(Operation):
                 def run(self, dataOut, channelList=None):
                     self.channelList = channelList
                     if self.channelList == None:
                         print("Missing channelList")
                         return dataOut
                     channelIndexList = []
                     if  type(dataOut.channelList) is not list:    #leer array desde HDF5
                         try:
                             dataOut.channelList = dataOut.channelList.tolist()
                         except Exception as e:
                             print("Select Channels: ",e)
                     for channel in self.channelList:
                         if channel not in dataOut.channelList:
                             raise ValueError("Channel %d is not in %s" %(channel, str(dataOut.channelList)))
                         index = dataOut.channelList.index(channel)
                         channelIndexList.append(index)
                     dataOut = self.selectChannelsByIndex(dataOut,channelIndexList)
                     return dataOut
                 def selectChannelsByIndex(self, dataOut, channelIndexList):
                     """
                     Selecciona un bloque de datos en base a canales segun el channelIndexList
                     Input:
                         channelIndexList    :    lista sencilla de canales a seleccionar por ej. [2,3,7]
                     Affected:
                         dataOut.data
                         dataOut.channelIndexList
                         dataOut.nChannels
                         dataOut.m_ProcessingHeader.totalSpectra
                         dataOut.systemHeaderObj.numChannels
                         dataOut.m_ProcessingHeader.blockSize
                     Return:
                         None
                     """
                     #print("selectChannelsByIndex")
                     # for channelIndex in channelIndexList:
                     #     if channelIndex not in dataOut.channelIndexList:
                     #         raise ValueError("The value %d in channelIndexList is not valid" %channelIndex)
                     if dataOut.type == 'Voltage':
                         if dataOut.flagDataAsBlock:
                             """
                             Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
                             """
                             data = dataOut.data[channelIndexList,:,:]
                         else:
                             data = dataOut.data[channelIndexList,:]
                         dataOut.data = data
                         # dataOut.channelList = [dataOut.channelList[i] for i in channelIndexList]
                         dataOut.channelList = range(len(channelIndexList))
                     elif dataOut.type == 'Spectra':
                         if hasattr(dataOut, 'data_spc'):
                             if dataOut.data_spc is None:
                                 raise ValueError("data_spc is None")
                                 return dataOut
                             else:
                                 data_spc = dataOut.data_spc[channelIndexList, :]
                                 dataOut.data_spc = data_spc
                         # if hasattr(dataOut, 'data_dc') :# and
                         #     if dataOut.data_dc is None:
                         #         raise ValueError("data_dc is None")
                         #         return dataOut
                         #     else:
                         #         data_dc = dataOut.data_dc[channelIndexList, :]
                         #         dataOut.data_dc = data_dc
                         # dataOut.channelList = [dataOut.channelList[i] for i in channelIndexList]
                         dataOut.channelList = channelIndexList
                         dataOut = self.__selectPairsByChannel(dataOut,channelIndexList)
                     return dataOut
                 def __selectPairsByChannel(self, dataOut, channelList=None):
                     #print("__selectPairsByChannel")
                     if channelList == None:
                         return
                     pairsIndexListSelected = []
                     for pairIndex in dataOut.pairsIndexList:
                         # First pair
                         if dataOut.pairsList[pairIndex][0] not in channelList:
                             continue
                         # Second pair
                         if dataOut.pairsList[pairIndex][1] not in channelList:
                             continue
                         pairsIndexListSelected.append(pairIndex)
                     if not pairsIndexListSelected:
                         dataOut.data_cspc = None
                         dataOut.pairsList = []
                         return
                     dataOut.data_cspc = dataOut.data_cspc[pairsIndexListSelected]
                     dataOut.pairsList = [dataOut.pairsList[i]
                                               for i in pairsIndexListSelected]
                     return dataOut
             class selectHeights(Operation):
                 def run(self, dataOut, minHei=None, maxHei=None, minIndex=None, maxIndex=None):
                     """
                     Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
                     minHei <= height <= maxHei
                     Input:
                         minHei    :    valor minimo de altura a considerar
                         maxHei    :    valor maximo de altura a considerar
                     Affected:
                         Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
                     Return:
 si el metodo se ejecuto con exito caso contrario devuelve 0
                     """
                     self.dataOut = dataOut
                     if minHei and maxHei:
                         if (minHei < dataOut.heightList[0]):
                             minHei = dataOut.heightList[0]
                         if (maxHei > dataOut.heightList[-1]):
                             maxHei = dataOut.heightList[-1]
                         minIndex = 0
                         maxIndex = 0
                         heights = 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)
                     self.selectHeightsByIndex(minIndex, maxIndex)
                     return dataOut
                 def selectHeightsByIndex(self, minIndex, maxIndex):
                     """
                     Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
                     minIndex <= index <= maxIndex
                     Input:
                         minIndex    :    valor de indice minimo de altura a considerar
                         maxIndex    :    valor de indice maximo de altura a considerar
                     Affected:
                         self.dataOut.data
                         self.dataOut.heightList
                     Return:
 si el metodo se ejecuto con exito caso contrario devuelve 0
                     """
                     if self.dataOut.type == 'Voltage':
                         if (minIndex < 0) or (minIndex > maxIndex):
                             raise ValueError("Height index range (%d,%d) is not valid" % (minIndex, maxIndex))
                         if (maxIndex >= self.dataOut.nHeights):
                             maxIndex = self.dataOut.nHeights
                         #voltage
                         if self.dataOut.flagDataAsBlock:
                             """
                             Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
                             """
                             data = self.dataOut.data[:,:, minIndex:maxIndex]
                         else:
                             data = self.dataOut.data[:, minIndex:maxIndex]
                         #         firstHeight = self.dataOut.heightList[minIndex]
                         self.dataOut.data = data
                         self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
                         if self.dataOut.nHeights <= 1:
                             raise ValueError("selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights))
                     elif self.dataOut.type == 'Spectra':
                         if (minIndex < 0) or (minIndex > maxIndex):
                             raise ValueError("Error selecting heights: Index range (%d,%d) is not valid" % (
                                 minIndex, maxIndex))
                         if (maxIndex >= self.dataOut.nHeights):
                             maxIndex = self.dataOut.nHeights - 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.heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
                     return 1
             class filterByHeights(Operation):
                 def run(self, dataOut, window):
                     deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                     if window == None:
                         window = (dataOut.radarControllerHeaderObj.txA/dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
                     newdelta = deltaHeight * window
                     r = dataOut.nHeights % window
                     newheights = (dataOut.nHeights-r)/window
                     if newheights <= 1:
                         raise ValueError("filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(dataOut.nHeights, window))
                     if dataOut.flagDataAsBlock:
                         """
                         Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
                         """
                         buffer = dataOut.data[:, :, 0:int(dataOut.nHeights-r)]
                         buffer = buffer.reshape(dataOut.nChannels, dataOut.nProfiles, int(dataOut.nHeights/window), window)
                         buffer = numpy.sum(buffer,3)
                     else:
                         buffer = dataOut.data[:,0:int(dataOut.nHeights-r)]
                         buffer = buffer.reshape(dataOut.nChannels,int(dataOut.nHeights/window),int(window))
                         buffer = numpy.sum(buffer,2)
                     dataOut.data = buffer
                     dataOut.heightList = dataOut.heightList[0] + numpy.arange( newheights )*newdelta
                     dataOut.windowOfFilter = window
                     return dataOut
             class setH0(Operation):
                 def run(self, dataOut, h0, deltaHeight = None):
                     if not deltaHeight:
                         deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                     nHeights = dataOut.nHeights
                     newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
                     dataOut.heightList = newHeiRange
                     return dataOut
             class deFlip(Operation):
                 def run(self, dataOut, channelList = []):
                     data = dataOut.data.copy()
                     if dataOut.flagDataAsBlock:
                         flip = self.flip
                         profileList = list(range(dataOut.nProfiles))
                         if not channelList:
                             for thisProfile in profileList:
                                 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
                                 flip *= -1.0
                         else:
                             for thisChannel in channelList:
                                 if thisChannel not in dataOut.channelList:
                                     continue
                                 for thisProfile in profileList:
                                     data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
                                     flip *= -1.0
                         self.flip = flip
                     else:
                         if not channelList:
                             data[:,:] = data[:,:]*self.flip
                         else:
                             for thisChannel in channelList:
                                 if thisChannel not in dataOut.channelList:
                                     continue
                                 data[thisChannel,:] = data[thisChannel,:]*self.flip
                         self.flip *= -1.
                     dataOut.data = data
                     return dataOut
             class setAttribute(Operation):
                 '''
                 Set an arbitrary attribute(s) to dataOut
                 '''
                 def __init__(self):
                     Operation.__init__(self)
                     self._ready = False
                 def run(self, dataOut, **kwargs):
                     for key, value in kwargs.items():
                         setattr(dataOut, key, value)
                     return dataOut
             @MPDecorator
             class printAttribute(Operation):
                 '''
                 Print an arbitrary attribute of dataOut
                 '''
                 def __init__(self):
                     Operation.__init__(self)
                 def run(self, dataOut, attributes):
                     if isinstance(attributes, str):
                         attributes = [attributes]
                     for attr in attributes:
                         if hasattr(dataOut, attr):
                             log.log(getattr(dataOut, attr), attr)
             class cleanHeightsInterf(Operation):
+                __slots__ =('heights_indx', 'repeats', 'step', 'factor', 'idate', 'idxs','config','wMask')
+                def __init__(self):
+                    self.repeats = 0
+                    self.factor=1
+                    self.wMask = None
+                    self.config = False
+                    self.idxs = None
+                    self.heights_indx = None
                 def run(self, dataOut, heightsList, repeats=0, step=0, factor=1, idate=None, startH=None, endH=None):
                     #print(dataOut.data.shape)
                     startTime = datetime.datetime.combine(idate,startH)
                     endTime =  datetime.datetime.combine(idate,endH)
                     currentTime = datetime.datetime.fromtimestamp(dataOut.utctime)
                     if currentTime < startTime or currentTime > endTime:
                         return dataOut
+                    if not self.config:
-                    wMask = numpy.asarray(factor)
-                    wMask = numpy.tile(wMask,(repeats+2))
                         #print(wMask)
                         heights = [float(hei) for hei in heightsList]
                         for r in range(repeats):
                              heights += [ (h+(step*(r+1))) for h in heights]
                         #print(heights)
                         heiList = dataOut.heightList
-                    heights_indx = [getHei_index(h,h,heiList)[0] for h in heights]
+                        self.heights_indx = [getHei_index(h,h,heiList)[0] for h in heights]
+                        self.wMask = numpy.asarray(factor)
+                        self.wMask = numpy.tile(self.wMask,(repeats+2))
+                        self.config = True
+                    """
+                    getNoisebyHildebrand(self, channel=None, ymin_index=None, ymax_index=None)
+                    """
+                    #print(self.noise =10*numpy.log10(dataOut.getNoisebyHildebrand(ymin_index=self.min_ref, ymax_index=self.max_ref)))
                     for ch in range(dataOut.data.shape[0]):
                         i = 0
-                        for hei in heights_indx:
-                            #print(dataOut.data[ch,hei])
+                        for hei in self.heights_indx:
+                            h = hei - 1
                             if dataOut.data.ndim < 3:
-                                dataOut.data[ch,hei-1] = (dataOut.data[ch,hei-1])*wMask[i]
+                                module = numpy.absolute(dataOut.data[ch,h])
+                                prev_h1 = numpy.absolute(dataOut.data[ch,h-1])
+                                dataOut.data[ch,h] = (dataOut.data[ch,h])/module * prev_h1
+                                #dataOut.data[ch,hei-1] = (dataOut.data[ch,hei-1])*self.wMask[i]
                             else:
-                                dataOut.data[ch,:,hei-1] = (dataOut.data[ch,:,hei-1])*wMask[i]
+                                module = numpy.absolute(dataOut.data[ch,:,h])
+                                prev_h1 = numpy.absolute(dataOut.data[ch,:,h-1])
+                                dataOut.data[ch,:,h] = (dataOut.data[ch,:,h])/module * prev_h1
+                                #dataOut.data[ch,:,hei-1] = (dataOut.data[ch,:,hei-1])*self.wMask[i]
                                 #print("done")
                             i += 1
                     return dataOut
             class interpolateHeights(Operation):
                 def run(self, dataOut, topLim, botLim):
                     #69 al 72 para julia
                     #82-84 para meteoros
                     if len(numpy.shape(dataOut.data))==2:
                         sampInterp = (dataOut.data[:,botLim-1] + dataOut.data[:,topLim+1])/2
                         sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
                         #dataOut.data[:,botLim:limSup+1] = sampInterp
                         dataOut.data[:,botLim:topLim+1] = sampInterp
                     else:
                         nHeights = dataOut.data.shape[2]
                         x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
                         y = dataOut.data[:,:,list(range(botLim))+list(range(topLim+1,nHeights))]
                         f = interpolate.interp1d(x, y, axis = 2)
                         xnew = numpy.arange(botLim,topLim+1)
                         ynew = f(xnew)
                         dataOut.data[:,:,botLim:topLim+1]  = ynew
                     return dataOut
             class CohInt(Operation):
                 isConfig = False
                 __profIndex = 0
                 __byTime = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer = None
                 __bufferStride = []
                 __dataReady = False
                 __profIndexStride = 0
                 __dataToPutStride = False
                 n = None
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=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 = None
                     self.__dataReady = False
                     self.byblock = byblock
                     self.stride = stride
                     if n == None and timeInterval == None:
                         raise ValueError("n or timeInterval should be specified ...")
                     if n != None:
                         self.n = n
                         self.__byTime = False
                     else:
                         self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
                         self.n = 9999
                         self.__byTime = True
                     if overlapping:
                         self.__withOverlapping = True
                         self.__buffer = None
                     else:
                         self.__withOverlapping = False
                         self.__buffer = 0
                     self.__profIndex = 0
                 def putData(self, data):
                     """
                     Add a profile to the __buffer and increase in one the __profileIndex
                     """
                     if not self.__withOverlapping:
                         self.__buffer += data.copy()
                         self.__profIndex += 1
                         return
                     #Overlapping data
                     nChannels, nHeis = data.shape
                     data = numpy.reshape(data, (1, nChannels, nHeis))
                     #If the buffer is empty then it takes the data value
                     if self.__buffer is None:
                         self.__buffer = data
                         self.__profIndex += 1
                         return
                     #If the buffer length is lower than n then stakcing the data value
                     if self.__profIndex < self.n:
                         self.__buffer = numpy.vstack((self.__buffer, data))
                         self.__profIndex += 1
                         return
                     #If the buffer length is equal to n then replacing the last buffer value with the data value
                     self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
                     self.__buffer[self.n-1] = data
                     self.__profIndex = self.n
                     return
                 def pushData(self):
                     """
                     Return the sum of the last profiles and the profiles used in the sum.
                     Affected:
                     self.__profileIndex
                     """
                     if not self.__withOverlapping:
                         data = self.__buffer
                         n = self.__profIndex
                         self.__buffer = 0
                         self.__profIndex = 0
                         return data, n
                     #Integration with Overlapping
                     data = numpy.sum(self.__buffer, axis=0)
                     # print data
                     # raise
                     n = self.__profIndex
                     return data, n
                 def byProfiles(self, data):
                     self.__dataReady = False
                     avgdata = None
                     #         n = None
                     # print data
                     # raise
                     self.putData(data)
                     if self.__profIndex == self.n:
                         avgdata, n = self.pushData()
                         self.__dataReady = True
                     return avgdata
                 def byTime(self, data, datatime):
                     self.__dataReady = False
                     avgdata = None
                     n = None
                     self.putData(data)
                     if (datatime - self.__initime) >= self.__integrationtime:
                         avgdata, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata
                 def integrateByStride(self, data, datatime):
                     # print data
                     if self.__profIndex == 0:
                         self.__buffer = [[data.copy(), datatime]]
                     else:
                         self.__buffer.append([data.copy(),datatime])
                     self.__profIndex += 1
                     self.__dataReady = False
                     if self.__profIndex == self.n * self.stride :
                         self.__dataToPutStride = True
                         self.__profIndexStride = 0
                         self.__profIndex = 0
                         self.__bufferStride = []
                         for i in range(self.stride):
                             current = self.__buffer[i::self.stride]
                             data = numpy.sum([t[0] for t in current], axis=0)
                             avgdatatime = numpy.average([t[1] for t in current])
                             # print data
                             self.__bufferStride.append((data, avgdatatime))
                     if self.__dataToPutStride:
                         self.__dataReady = True
                         self.__profIndexStride += 1
                         if self.__profIndexStride == self.stride:
                             self.__dataToPutStride = False
                         # print self.__bufferStride[self.__profIndexStride - 1]
                         # raise
                         return self.__bufferStride[self.__profIndexStride - 1]
                     return None, None
                 def integrate(self, data, datatime=None):
                     if self.__initime == None:
                         self.__initime = datatime
                     if self.__byTime:
                         avgdata = self.byTime(data, datatime)
                     else:
                         avgdata = self.byProfiles(data)
                     self.__lastdatatime = datatime
                     if avgdata is None:
                         return None, None
                     avgdatatime = self.__initime
                     deltatime = datatime - self.__lastdatatime
                     if not self.__withOverlapping:
                         self.__initime = datatime
                     else:
                         self.__initime += deltatime
                     return avgdata, avgdatatime
                 def integrateByBlock(self, dataOut):
                     times = int(dataOut.data.shape[1]/self.n)
                     avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
                     id_min = 0
                     id_max = self.n
                     for i in range(times):
                         junk = dataOut.data[:,id_min:id_max,:]
                         avgdata[:,i,:] = junk.sum(axis=1)
                         id_min += self.n
                         id_max += self.n
                     timeInterval = dataOut.ippSeconds*self.n
                     avgdatatime = (times - 1) * timeInterval + dataOut.utctime
                     self.__dataReady = True
                     return avgdata, avgdatatime
                 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
                     if not self.isConfig:
                         self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
                         self.isConfig = True
                     if dataOut.flagDataAsBlock:
                         """
                         Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
                         """
                         avgdata, avgdatatime = self.integrateByBlock(dataOut)
                         dataOut.nProfiles /= self.n
                     else:
                         if stride is None:
                             avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
                         else:
                             avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
                     #   dataOut.timeInterval *= n
                     dataOut.flagNoData = True
                     if self.__dataReady:
                         dataOut.data = avgdata
                         if not dataOut.flagCohInt:
                             dataOut.nCohInt *= self.n
                             dataOut.flagCohInt = True
                         dataOut.utctime = avgdatatime
                         # print avgdata, avgdatatime
                         # raise
                         #   dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
                         dataOut.flagNoData = False
                     return dataOut
             class Decoder(Operation):
                 isConfig = False
                 __profIndex = 0
                 code = None
                 nCode = None
                 nBaud = None
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.times = None
                     self.osamp = None
                 #         self.__setValues = False
                     self.isConfig = False
                     self.setupReq = False
                 def setup(self, code, osamp, dataOut):
                     self.__profIndex = 0
                     self.code = code
                     self.nCode = len(code)
                     self.nBaud = len(code[0])
                     if (osamp != None) and (osamp >1):
                         self.osamp = osamp
                         self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
                         self.nBaud = self.nBaud*self.osamp
                     self.__nChannels = dataOut.nChannels
                     self.__nProfiles = dataOut.nProfiles
                     self.__nHeis = dataOut.nHeights
                     if self.__nHeis < self.nBaud:
                         raise ValueError('Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud))
                     #Frequency
                     __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
                     __codeBuffer[:,0:self.nBaud] = self.code
                     self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
                     if dataOut.flagDataAsBlock:
                         self.ndatadec = self.__nHeis #- self.nBaud + 1
                         self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
                     else:
                         #Time
                         self.ndatadec = self.__nHeis #- self.nBaud + 1
                         self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
                 def __convolutionInFreq(self, data):
                     fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
                     fft_data = numpy.fft.fft(data, axis=1)
                     conv = fft_data*fft_code
                     data = numpy.fft.ifft(conv,axis=1)
                     return data
                 def __convolutionInFreqOpt(self, data):
                     raise NotImplementedError
                 def __convolutionInTime(self, data):
                     code = self.code[self.__profIndex]
                     for i in range(self.__nChannels):
                         self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
                     return self.datadecTime
                 def __convolutionByBlockInTime(self, data):
                     repetitions = int(self.__nProfiles / self.nCode)
                     junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
                     junk = junk.flatten()
                     code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
                     profilesList = range(self.__nProfiles)
                     for i in range(self.__nChannels):
                         for j in profilesList:
                             self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
                     return self.datadecTime
                 def __convolutionByBlockInFreq(self, data):
                     raise NotImplementedError("Decoder by frequency fro Blocks not implemented")
                     fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
                     fft_data = numpy.fft.fft(data, axis=2)
                     conv = fft_data*fft_code
                     data = numpy.fft.ifft(conv,axis=2)
                     return data
                 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
                     if dataOut.flagDecodeData:
                         print("This data is already decoded, recoding again ...")
                     if not self.isConfig:
                         if code is None:
                             if dataOut.code is None:
                                 raise ValueError("Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type)
                             code = dataOut.code
                         else:
                             code = numpy.array(code).reshape(nCode,nBaud)
                         self.setup(code, osamp, dataOut)
                         self.isConfig = True
                         if mode == 3:
                             sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
                         if times != None:
                             sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
                     if self.code is None:
                         print("Fail decoding: Code is not defined.")
                         return
                     self.__nProfiles = dataOut.nProfiles
                     datadec = None
                     if mode == 3:
                         mode = 0
                     if dataOut.flagDataAsBlock:
                         """
                         Decoding when data have been read as block,
                         """
                         if mode == 0:
                             datadec = self.__convolutionByBlockInTime(dataOut.data)
                         if mode == 1:
                             datadec = self.__convolutionByBlockInFreq(dataOut.data)
                     else:
                         """
                         Decoding when data have been read profile by profile
                         """
                         if mode == 0:
                             datadec = self.__convolutionInTime(dataOut.data)
                         if mode == 1:
                             datadec = self.__convolutionInFreq(dataOut.data)
                         if mode == 2:
                             datadec = self.__convolutionInFreqOpt(dataOut.data)
                     if datadec is None:
                         raise ValueError("Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode)
                     dataOut.code = self.code
                     dataOut.nCode = self.nCode
                     dataOut.nBaud = self.nBaud
                     dataOut.data = datadec
                     dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
                     dataOut.flagDecodeData = True #asumo q la data esta decodificada
                     if self.__profIndex == self.nCode-1:
                         self.__profIndex = 0
                         return dataOut
                     self.__profIndex += 1
                     return dataOut
                 #        dataOut.flagDeflipData = True #asumo q la data no esta sin flip
             class ProfileConcat(Operation):
                 isConfig = False
                 buffer = None
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.profileIndex = 0
                 def reset(self):
                     self.buffer = numpy.zeros_like(self.buffer)
                     self.start_index = 0
                     self.times = 1
                 def setup(self, data, m, n=1):
                     self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
                     self.nHeights = data.shape[1]#.nHeights
                     self.start_index = 0
                     self.times = 1
                 def concat(self, data):
                     self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
                     self.start_index = self.start_index + self.nHeights
                 def run(self, dataOut, m):
                     dataOut.flagNoData = True
                     if not self.isConfig:
                         self.setup(dataOut.data, m, 1)
                         self.isConfig = True
                     if dataOut.flagDataAsBlock:
                         raise ValueError("ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False")
                     else:
                         self.concat(dataOut.data)
                         self.times += 1
                         if self.times > m:
                             dataOut.data = self.buffer
                             self.reset()
                             dataOut.flagNoData = False
                             # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
                             deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                             xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
                             dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
                             dataOut.ippSeconds *= m
                     return dataOut
             class ProfileSelector(Operation):
                 profileIndex = None
                 # Tamanho total de los perfiles
                 nProfiles = None
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.profileIndex = 0
                 def incProfileIndex(self):
                     self.profileIndex += 1
                     if self.profileIndex >= self.nProfiles:
                         self.profileIndex = 0
                 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
                     if profileIndex < minIndex:
                         return False
                     if profileIndex > maxIndex:
                         return False
                     return True
                 def isThisProfileInList(self, profileIndex, profileList):
                     if profileIndex not in profileList:
                         return False
                     return True
                 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
                     """
                     ProfileSelector:
                     Inputs:
                         profileList        :    Index of profiles selected. Example: profileList = (0,1,2,7,8)
                         profileRangeList    :    Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
                         rangeList            :    List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
                     """
                     if rangeList is not None:
                         if type(rangeList[0]) not in (tuple, list):
                             rangeList = [rangeList]
                     dataOut.flagNoData = True
                     if dataOut.flagDataAsBlock:
                         """
                         data dimension  = [nChannels, nProfiles, nHeis]
                         """
                         if profileList != None:
                             dataOut.data = dataOut.data[:,profileList,:]
                         if profileRangeList != None:
                             minIndex = profileRangeList[0]
                             maxIndex = profileRangeList[1]
                             profileList = list(range(minIndex, maxIndex+1))
                             dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
                         if rangeList != None:
                             profileList = []
                             for thisRange in rangeList:
                                 minIndex = thisRange[0]
                                 maxIndex = thisRange[1]
                                 profileList.extend(list(range(minIndex, maxIndex+1)))
                             dataOut.data = dataOut.data[:,profileList,:]
                         dataOut.nProfiles = len(profileList)
                         dataOut.profileIndex = dataOut.nProfiles - 1
                         dataOut.flagNoData = False
                         return dataOut
                     """
                     data dimension  = [nChannels, nHeis]
                     """
                     if profileList != None:
                         if self.isThisProfileInList(dataOut.profileIndex, profileList):
                             self.nProfiles = len(profileList)
                             dataOut.nProfiles = self.nProfiles
                             dataOut.profileIndex = self.profileIndex
                             dataOut.flagNoData = False
                             self.incProfileIndex()
                         return dataOut
                     if profileRangeList != None:
                         minIndex = profileRangeList[0]
                         maxIndex = profileRangeList[1]
                         if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
                             self.nProfiles = maxIndex - minIndex + 1
                             dataOut.nProfiles = self.nProfiles
                             dataOut.profileIndex = self.profileIndex
                             dataOut.flagNoData = False
                             self.incProfileIndex()
                         return dataOut
                     if rangeList != None:
                         nProfiles = 0
                         for thisRange in rangeList:
                             minIndex = thisRange[0]
                             maxIndex = thisRange[1]
                             nProfiles += maxIndex - minIndex + 1
                         for thisRange in rangeList:
                             minIndex = thisRange[0]
                             maxIndex = thisRange[1]
                             if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
                                 self.nProfiles = nProfiles
                                 dataOut.nProfiles = self.nProfiles
                                 dataOut.profileIndex = self.profileIndex
                                 dataOut.flagNoData = False
                                 self.incProfileIndex()
                                 break
                         return dataOut
                     if beam != None: #beam is only for AMISR data
                         if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
                             dataOut.flagNoData = False
                             dataOut.profileIndex = self.profileIndex
                             self.incProfileIndex()
                         return dataOut
                     raise ValueError("ProfileSelector needs profileList, profileRangeList or rangeList parameter")
             class Reshaper(Operation):
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.__buffer = None
                     self.__nitems = 0
                 def __appendProfile(self, dataOut, nTxs):
                     if self.__buffer is None:
                         shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
                         self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
                     ini = dataOut.nHeights * self.__nitems
                     end = ini + dataOut.nHeights
                     self.__buffer[:, ini:end] = dataOut.data
                     self.__nitems += 1
                     return int(self.__nitems*nTxs)
                 def __getBuffer(self):
                     if self.__nitems == int(1./self.__nTxs):
                         self.__nitems = 0
                         return self.__buffer.copy()
                     return None
                 def __checkInputs(self, dataOut, shape, nTxs):
                     if shape is None and nTxs is None:
                         raise ValueError("Reshaper: shape of factor should be defined")
                     if nTxs:
                         if nTxs < 0:
                             raise ValueError("nTxs should be greater than 0")
                         if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
                             raise ValueError("nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs)))
                         shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
                         return shape, nTxs
                     if len(shape) != 2 and len(shape) !=  3:
                         raise ValueError("shape dimension should be equal to 2 or 3. shape = (nProfiles, nHeis) or (nChannels, nProfiles, nHeis). Actually shape = (%d, %d, %d)" %(dataOut.nChannels, dataOut.nProfiles, dataOut.nHeights))
                     if len(shape) == 2:
                         shape_tuple = [dataOut.nChannels]
                         shape_tuple.extend(shape)
                     else:
                         shape_tuple = list(shape)
                     nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
                     return shape_tuple, nTxs
                 def run(self, dataOut, shape=None, nTxs=None):
                     shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
                     dataOut.flagNoData = True
                     profileIndex = None
                     if dataOut.flagDataAsBlock:
                         dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
                         dataOut.flagNoData = False
                         profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
                     else:
                         if self.__nTxs < 1:
                             self.__appendProfile(dataOut, self.__nTxs)
                             new_data = self.__getBuffer()
                             if new_data is not None:
                                 dataOut.data = new_data
                                 dataOut.flagNoData = False
                                 profileIndex = dataOut.profileIndex*nTxs
                         else:
                             raise ValueError("nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)")
                     deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                     dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
                     dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
                     dataOut.profileIndex = profileIndex
                     dataOut.ippSeconds /= self.__nTxs
                     return dataOut
             class SplitProfiles(Operation):
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                 def run(self, dataOut, n):
                     dataOut.flagNoData = True
                     profileIndex = None
                     if dataOut.flagDataAsBlock:
                         #nchannels, nprofiles, nsamples
                         shape = dataOut.data.shape
                         if shape[2] % n != 0:
                             raise ValueError("Could not split the data, n=%d has to be multiple of %d" %(n, shape[2]))
                         new_shape = shape[0], shape[1]*n, int(shape[2]/n)
                         dataOut.data = numpy.reshape(dataOut.data, new_shape)
                         dataOut.flagNoData = False
                         profileIndex = int(dataOut.nProfiles/n) - 1
                     else:
                         raise ValueError("Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)")
                     deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                     dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
                     dataOut.nProfiles = int(dataOut.nProfiles*n)
                     dataOut.profileIndex = profileIndex
                     dataOut.ippSeconds /= n
                     return dataOut
             class CombineProfiles(Operation):
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.__remData = None
                     self.__profileIndex = 0
                 def run(self, dataOut, n):
                     dataOut.flagNoData = True
                     profileIndex = None
                     if dataOut.flagDataAsBlock:
                         #nchannels, nprofiles, nsamples
                         shape = dataOut.data.shape
                         new_shape = shape[0], shape[1]/n, shape[2]*n
                         if shape[1] % n != 0:
                             raise ValueError("Could not split the data, n=%d has to be multiple of %d" %(n, shape[1]))
                         dataOut.data = numpy.reshape(dataOut.data, new_shape)
                         dataOut.flagNoData = False
                         profileIndex = int(dataOut.nProfiles*n) - 1
                     else:
                         #nchannels, nsamples
                         if self.__remData is None:
                             newData = dataOut.data
                         else:
                             newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
                         self.__profileIndex += 1
                         if self.__profileIndex < n:
                             self.__remData = newData
                             #continue
                             return
                         self.__profileIndex = 0
                         self.__remData = None
                         dataOut.data = newData
                         dataOut.flagNoData = False
                         profileIndex = dataOut.profileIndex/n
                     deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
                     dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
                     dataOut.nProfiles = int(dataOut.nProfiles/n)
                     dataOut.profileIndex = profileIndex
                     dataOut.ippSeconds *= n
                     return dataOut
             class PulsePairVoltage(Operation):
                 '''
                 Function PulsePair(Signal Power, Velocity)
                 The real component of Lag[0] provides Intensity Information
                 The imag component of Lag[1] Phase provides Velocity Information
                 Configuration Parameters:
                 nPRF = Number of Several PRF
                 theta = Degree Azimuth angel Boundaries
                 Input:
                       self.dataOut
                       lag[N]
                 Affected:
                       self.dataOut.spc
                 '''
                 isConfig       = False
                 __profIndex    = 0
                 __initime      = None
                 __lastdatatime = None
                 __buffer       = None
                 noise          = None
                 __dataReady    = False
                 n              = None
                 __nch          = 0
                 __nHeis        = 0
                 removeDC       = False
                 ipp            = None
                 lambda_        = 0
                 def __init__(self,**kwargs):
                     Operation.__init__(self,**kwargs)
                 def setup(self, dataOut, n = None, removeDC=False):
                     '''
                     n= Numero de PRF's de entrada
                     '''
                     self.__initime        = None
                     self.__lastdatatime   = 0
                     self.__dataReady      = False
                     self.__buffer         = 0
                     self.__profIndex      = 0
                     self.noise            = None
                     self.__nch            = dataOut.nChannels
                     self.__nHeis          = dataOut.nHeights
                     self.removeDC         = removeDC
                     self.lambda_          = 3.0e8/(9345.0e6)
                     self.ippSec           = dataOut.ippSeconds
                     self.nCohInt          = dataOut.nCohInt
                     if n == None:
                         raise ValueError("n should be specified.")
                     if n != None:
                         if n<2:
                             raise ValueError("n should be greater than 2")
                     self.n       = n
                     self.__nProf = n
                     self.__buffer = numpy.zeros((dataOut.nChannels,
                                                        n,
                                                        dataOut.nHeights),
                                                       dtype='complex')
                 def putData(self,data):
                     '''
                     Add a profile to he __buffer and increase in one the __profiel Index
                     '''
                     self.__buffer[:,self.__profIndex,:]= data
                     self.__profIndex      += 1
                     return
                 def pushData(self,dataOut):
                     '''
                     Return the PULSEPAIR and the profiles used in the operation
                     Affected :  self.__profileIndex
                     '''
                     #----------------- Remove DC-----------------------------------
                     if self.removeDC==True:
                         mean    = numpy.mean(self.__buffer,1)
                         tmp     = mean.reshape(self.__nch,1,self.__nHeis)
                         dc= numpy.tile(tmp,[1,self.__nProf,1])
                         self.__buffer = self.__buffer -  dc
                     #------------------Calculo de Potencia ------------------------
                     pair0       = self.__buffer*numpy.conj(self.__buffer)
                     pair0       = pair0.real
                     lag_0       = numpy.sum(pair0,1)
                     #------------------Calculo de Ruido x canal--------------------
                     self.noise  = numpy.zeros(self.__nch)
                     for i in range(self.__nch):
                         daux         = numpy.sort(pair0[i,:,:],axis= None)
                         self.noise[i]=hildebrand_sekhon( daux ,self.nCohInt)
                     self.noise       = self.noise.reshape(self.__nch,1)
                     self.noise       = numpy.tile(self.noise,[1,self.__nHeis])
                     noise_buffer     = self.noise.reshape(self.__nch,1,self.__nHeis)
                     noise_buffer     = numpy.tile(noise_buffer,[1,self.__nProf,1])
                     #------------------ Potencia recibida= P , Potencia senal = S , Ruido= N--
                     #------------------   P= S+N  ,P=lag_0/N ---------------------------------
                     #-------------------- Power --------------------------------------------------
                     data_power       = lag_0/(self.n*self.nCohInt)
                     #------------------  Senal  ---------------------------------------------------
                     data_intensity   = pair0 - noise_buffer
                     data_intensity   = numpy.sum(data_intensity,axis=1)*(self.n*self.nCohInt)#*self.nCohInt)
                     #data_intensity   = (lag_0-self.noise*self.n)*(self.n*self.nCohInt)
                     for i in range(self.__nch):
                         for j in range(self.__nHeis):
                             if data_intensity[i][j]  < 0:
                                 data_intensity[i][j] = numpy.min(numpy.absolute(data_intensity[i][j]))
                     #----------------- Calculo de Frecuencia y Velocidad doppler--------
                     pair1            = self.__buffer[:,:-1,:]*numpy.conjugate(self.__buffer[:,1:,:])
                     lag_1            = numpy.sum(pair1,1)
                     data_freq        = (-1/(2.0*math.pi*self.ippSec*self.nCohInt))*numpy.angle(lag_1)
                     data_velocity    = (self.lambda_/2.0)*data_freq
                     #---------------- Potencia promedio estimada de la Senal-----------
                     lag_0            = lag_0/self.n
                     S                = lag_0-self.noise
                     #---------------- Frecuencia Doppler promedio ---------------------
                     lag_1            = lag_1/(self.n-1)
                     R1               = numpy.abs(lag_1)
                     #---------------- Calculo del SNR----------------------------------
                     data_snrPP       = S/self.noise
                     for i in range(self.__nch):
                         for j in range(self.__nHeis):
                             if data_snrPP[i][j]  < 1.e-20:
                                 data_snrPP[i][j] = 1.e-20
                     #----------------- Calculo del ancho espectral ----------------------
                     L                = S/R1
                     L                = numpy.where(L<0,1,L)
                     L                = numpy.log(L)
                     tmp              = numpy.sqrt(numpy.absolute(L))
                     data_specwidth   = (self.lambda_/(2*math.sqrt(2)*math.pi*self.ippSec*self.nCohInt))*tmp*numpy.sign(L)
                     n                = self.__profIndex
                     self.__buffer    = numpy.zeros((self.__nch, self.__nProf,self.__nHeis),  dtype='complex')
                     self.__profIndex = 0
                     return data_power,data_intensity,data_velocity,data_snrPP,data_specwidth,n
                 def pulsePairbyProfiles(self,dataOut):
                     self.__dataReady     =  False
                     data_power           =  None
                     data_intensity       =  None
                     data_velocity        =  None
                     data_specwidth       =  None
                     data_snrPP           =  None
                     self.putData(data=dataOut.data)
                     if self.__profIndex  == self.n:
                         data_power,data_intensity, data_velocity,data_snrPP,data_specwidth, n   = self.pushData(dataOut=dataOut)
                         self.__dataReady                   = True
                     return data_power, data_intensity, data_velocity, data_snrPP, data_specwidth
                 def pulsePairOp(self, dataOut, datatime= None):
                     if self.__initime == None:
                         self.__initime = datatime
                     data_power, data_intensity, data_velocity, data_snrPP, data_specwidth = self.pulsePairbyProfiles(dataOut)
                     self.__lastdatatime           = datatime
                     if data_power is None:
                         return None, None, None,None,None,None
                     avgdatatime    = self.__initime
                     deltatime      = datatime - self.__lastdatatime
                     self.__initime = datatime
                     return data_power, data_intensity, data_velocity, data_snrPP, data_specwidth, avgdatatime
                 def run(self, dataOut,n = None,removeDC= False, overlapping= False,**kwargs):
                     if not self.isConfig:
                         self.setup(dataOut = dataOut, n    = n , removeDC=removeDC , **kwargs)
                         self.isConfig   = True
                     data_power, data_intensity, data_velocity,data_snrPP,data_specwidth, avgdatatime = self.pulsePairOp(dataOut, dataOut.utctime)
                     dataOut.flagNoData                         = True
                     if self.__dataReady:
                         dataOut.nCohInt        *= self.n
                         dataOut.dataPP_POW      = data_intensity # S
                         dataOut.dataPP_POWER    = data_power     # P
                         dataOut.dataPP_DOP      = data_velocity
                         dataOut.dataPP_SNR      = data_snrPP
                         dataOut.dataPP_WIDTH    = data_specwidth
                         dataOut.PRFbyAngle      = self.n         #numero de PRF*cada angulo rotado que equivale a un tiempo.
                         dataOut.utctime         = avgdatatime
                         dataOut.flagNoData      = False
                     return dataOut
             # import collections
             # from scipy.stats import mode
             #
             # class Synchronize(Operation):
             #
             #     isConfig = False
             #     __profIndex = 0
             #
             #     def __init__(self, **kwargs):
             #
             #         Operation.__init__(self, **kwargs)
             # #         self.isConfig = False
             #         self.__powBuffer = None
             #         self.__startIndex = 0
             #         self.__pulseFound = False
             #
             #     def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
             #
             #         #Read data
             #
             #         powerdB = dataOut.getPower(channel = channel)
             #         noisedB = dataOut.getNoise(channel = channel)[0]
             #
             #         self.__powBuffer.extend(powerdB.flatten())
             #
             #         dataArray = numpy.array(self.__powBuffer)
             #
             #         filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
             #
             #         maxValue = numpy.nanmax(filteredPower)
             #
             #         if maxValue < noisedB + 10:
             #             #No se encuentra ningun pulso de transmision
             #             return None
             #
             #         maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
             #
             #         if len(maxValuesIndex) < 2:
             #             #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
             #             return None
             #
             #         phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
             #
             #         #Seleccionar solo valores con un espaciamiento de nSamples
             #         pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
             #
             #         if len(pulseIndex) < 2:
             #             #Solo se encontro un pulso de transmision con ancho mayor a 1
             #             return None
             #
             #         spacing = pulseIndex[1:] - pulseIndex[:-1]
             #
             #         #remover senales que se distancien menos de 10 unidades o muestras
             #         #(No deberian existir IPP menor a 10 unidades)
             #
             #         realIndex = numpy.where(spacing > 10 )[0]
             #
             #         if len(realIndex) < 2:
             #             #Solo se encontro un pulso de transmision con ancho mayor a 1
             #             return None
             #
             #         #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
             #         realPulseIndex = pulseIndex[realIndex]
             #
             #         period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
             #
             #         print "IPP = %d samples" %period
             #
             #         self.__newNSamples = dataOut.nHeights #int(period)
             #         self.__startIndex = int(realPulseIndex[0])
             #
             #         return 1
             #
             #
             #     def setup(self, nSamples, nChannels, buffer_size = 4):
             #
             #         self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
             #                                           maxlen = buffer_size*nSamples)
             #
             #         bufferList = []
             #
             #         for i in range(nChannels):
             #             bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) +  numpy.NAN,
             #                                           maxlen = buffer_size*nSamples)
             #
             #             bufferList.append(bufferByChannel)
             #
             #         self.__nSamples = nSamples
             #         self.__nChannels = nChannels
             #         self.__bufferList = bufferList
             #
             #     def run(self, dataOut, channel = 0):
             #
             #         if not self.isConfig:
             #             nSamples = dataOut.nHeights
             #             nChannels = dataOut.nChannels
             #             self.setup(nSamples, nChannels)
             #             self.isConfig = True
             #
             #         #Append new data to internal buffer
             #         for thisChannel in range(self.__nChannels):
             #             bufferByChannel = self.__bufferList[thisChannel]
             #             bufferByChannel.extend(dataOut.data[thisChannel])
             #
             #         if self.__pulseFound:
             #             self.__startIndex -= self.__nSamples
             #
             #         #Finding Tx Pulse
             #         if not self.__pulseFound:
             #             indexFound = self.__findTxPulse(dataOut, channel)
             #
             #             if indexFound == None:
             #                 dataOut.flagNoData = True
             #                 return
             #
             #             self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
             #             self.__pulseFound = True
             #             self.__startIndex = indexFound
             #
             #         #If pulse was found ...
             #         for thisChannel in range(self.__nChannels):
             #             bufferByChannel = self.__bufferList[thisChannel]
             #             #print self.__startIndex
             #             x = numpy.array(bufferByChannel)
             #             self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
             #
             #         deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
             #         dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
             # #             dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
             #
             #         dataOut.data = self.__arrayBuffer
             #
             #         self.__startIndex += self.__newNSamples
             #
             #         return
             class SSheightProfiles(Operation):
                 step          = None
                 nsamples      = None
                 bufferShape   = None
                 profileShape  = None
                 sshProfiles   = None
                 profileIndex  = None
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.isConfig = False
                 def setup(self,dataOut ,step = None , nsamples = None):
                     if step == None and nsamples == None:
                         raise ValueError("step or nheights should be specified ...")
                     self.step         = step
                     self.nsamples     = nsamples
                     self.__nChannels  = dataOut.nChannels
                     self.__nProfiles  = dataOut.nProfiles
                     self.__nHeis      = dataOut.nHeights
                     shape             = dataOut.data.shape #nchannels, nprofiles, nsamples
                     residue     =  (shape[1] - self.nsamples) % self.step
                     if residue != 0:
                         print("The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue))
                     deltaHeight      =  dataOut.heightList[1] - dataOut.heightList[0]
                     numberProfile    =  self.nsamples
                     numberSamples    =  (shape[1] - self.nsamples)/self.step
                     self.bufferShape  = int(shape[0]), int(numberSamples), int(numberProfile)  # nchannels, nsamples , nprofiles
                     self.profileShape = int(shape[0]), int(numberProfile), int(numberSamples)  # nchannels, nprofiles, nsamples
                     self.buffer       = numpy.zeros(self.bufferShape , dtype=numpy.complex)
                     self.sshProfiles  = numpy.zeros(self.profileShape, dtype=numpy.complex)
                 def run(self, dataOut, step, nsamples, code = None, repeat = None):
                     dataOut.flagNoData      = True
                     profileIndex            = None
                     #print("nProfiles, nHeights ",dataOut.nProfiles, dataOut.nHeights)
                     #print(dataOut.getFreqRange(1)/1000.)
                     #exit(1)
                     if dataOut.flagDataAsBlock:
                         dataOut.data = numpy.average(dataOut.data,axis=1)
                         #print("jee")
                     dataOut.flagDataAsBlock = False
                     if not self.isConfig:
                         self.setup(dataOut, step=step , nsamples=nsamples)
                         #print("Setup done")
                         self.isConfig = True
                     if code is not None:
                         code = numpy.array(code)
                         code_block = code
                         if repeat is not None:
                             code_block = numpy.repeat(code_block, repeats=repeat, axis=1)
                             #print(code_block.shape)
                     for i in range(self.buffer.shape[1]):
                         if code is not None:
                             self.buffer[:,i] = dataOut.data[:,i*self.step:i*self.step + self.nsamples]*code_block
                         else:
                             self.buffer[:,i]    = dataOut.data[:,i*self.step:i*self.step + self.nsamples]#*code[dataOut.profileIndex,:]
                         #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
                     for j in range(self.buffer.shape[0]):
                         self.sshProfiles[j] = numpy.transpose(self.buffer[j])
                     profileIndex  =  self.nsamples
                     deltaHeight   =  dataOut.heightList[1] - dataOut.heightList[0]
                     ippSeconds    =  (deltaHeight*1.0e-6)/(0.15)
                     #print("ippSeconds, dH: ",ippSeconds,deltaHeight)
                     try:
                         if dataOut.concat_m  is not None:
                             ippSeconds= ippSeconds/float(dataOut.concat_m)
                             #print "Profile concat %d"%dataOut.concat_m
                     except:
                         pass
                     dataOut.data            = self.sshProfiles
                     dataOut.flagNoData      = False
                     dataOut.heightList      = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
                     dataOut.nProfiles       = int(dataOut.nProfiles*self.nsamples)
                     dataOut.profileIndex    = profileIndex
                     dataOut.flagDataAsBlock = True
                     dataOut.ippSeconds      = ippSeconds
                     dataOut.step            = self.step
                     #print(numpy.shape(dataOut.data))
                     #exit(1)
                     #print("new data shape and time:", dataOut.data.shape, dataOut.utctime)
                     return dataOut
             ################################################################################3############################3
             ################################################################################3############################3
             ################################################################################3############################3
             ################################################################################3############################3
             class SSheightProfiles2(Operation):
                 '''
                 Procesa por perfiles y por bloques
                 '''
                 bufferShape   = None
                 profileShape  = None
                 sshProfiles   = None
                 profileIndex  = None
                 #nsamples      = None
                 #step          = None
                 #deltaHeight = None
                 #init_range  = None
                 __slots__ = ('step', 'nsamples', 'deltaHeight', 'init_range', 'isConfig', '__nChannels',
                             '__nProfiles', '__nHeis', 'deltaHeight', 'new_nHeights')
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.isConfig = False
                 def setup(self,dataOut ,step = None , nsamples = None):
                     if step == None and nsamples == None:
                         raise ValueError("step or nheights should be specified ...")
                     self.step         = step
                     self.nsamples     = nsamples
                     self.__nChannels  = int(dataOut.nChannels)
                     self.__nProfiles  = int(dataOut.nProfiles)
                     self.__nHeis      = int(dataOut.nHeights)
                     residue     =  (self.__nHeis - self.nsamples) % self.step
                     if residue != 0:
                         print("The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,self.__nProfiles - self.nsamples,residue))
                     self.deltaHeight      =  dataOut.heightList[1] - dataOut.heightList[0]
                     self.init_range = dataOut.heightList[0]
                     #numberProfile    =  self.nsamples
                     numberSamples    =  (self.__nHeis - self.nsamples)/self.step
                     self.new_nHeights    =  numberSamples
                     self.bufferShape  = int(self.__nChannels), int(numberSamples), int(self.nsamples)  # nchannels, nsamples , nprofiles
                     self.profileShape = int(self.__nChannels), int(self.nsamples), int(numberSamples)  # nchannels, nprofiles, nsamples
                     self.buffer       = numpy.zeros(self.bufferShape , dtype=numpy.complex)
                     self.sshProfiles  = numpy.zeros(self.profileShape, dtype=numpy.complex)
                 def getNewProfiles(self, data, code=None, repeat=None):
                     if code is not None:
                         code = numpy.array(code)
                         code_block = code
                         if repeat is not None:
                             code_block = numpy.repeat(code_block, repeats=repeat, axis=1)
                     if data.ndim < 3:
                         data = data.reshape(self.__nChannels,1,self.__nHeis )
                     #print("buff, data, :",self.buffer.shape, data.shape,self.sshProfiles.shape, code_block.shape)
                     for ch in range(self.__nChannels):
                         for i in range(int(self.new_nHeights)): #nuevas alturas
                             if code is not None:
                                 self.buffer[ch,i,:] = data[ch,:,i*self.step:i*self.step + self.nsamples]*code_block
                             else:
                                 self.buffer[ch,i,:] = data[ch,:,i*self.step:i*self.step + self.nsamples]#*code[dataOut.profileIndex,:]
                     for j in range(self.__nChannels): #en los cananles
                         self.sshProfiles[j,:,:] = numpy.transpose(self.buffer[j,:,:])
                     #print("new profs Done")
                 def run(self, dataOut, step, nsamples, code = None, repeat = None):
                     # print("running")
                     if dataOut.flagNoData == True:
                         return dataOut
                     dataOut.flagNoData = True
                     #print("init  data shape:", dataOut.data.shape)
                     #print("ch: {}  prof: {}  hs: {}".format(int(dataOut.nChannels),
                     #                int(dataOut.nProfiles),int(dataOut.nHeights)))
                     profileIndex            = None
                     # if not dataOut.flagDataAsBlock:
                     #     dataOut.nProfiles = 1
                     if not self.isConfig:
                         self.setup(dataOut, step=step , nsamples=nsamples)
                         #print("Setup done")
                         self.isConfig = True
                     dataBlock = None
                     nprof = 1
                     if dataOut.flagDataAsBlock:
                         nprof = int(dataOut.nProfiles)
                     #print("dataOut nProfiles:", dataOut.nProfiles)
                     for profile in range(nprof):
                         if dataOut.flagDataAsBlock:
                             #print("read blocks")
                             self.getNewProfiles(dataOut.data[:,profile,:], code=code, repeat=repeat)
                         else:
                             #print("read profiles")
                             self.getNewProfiles(dataOut.data, code=code, repeat=repeat) #only one channe
                         if profile == 0:
                             dataBlock = self.sshProfiles.copy()
                         else:   #by blocks
                             dataBlock = numpy.concatenate((dataBlock,self.sshProfiles), axis=1) #profile axis
                             #print("by blocks: ",dataBlock.shape, self.sshProfiles.shape)
                     profileIndex  =  self.nsamples
                     #deltaHeight   =  dataOut.heightList[1] - dataOut.heightList[0]
                     ippSeconds    =  (self.deltaHeight*1.0e-6)/(0.15)
                     dataOut.data            = dataBlock
                     #print("show me: ",self.step,self.deltaHeight, dataOut.heightList, self.new_nHeights)
                     dataOut.heightList      = numpy.arange(int(self.new_nHeights)) *self.step*self.deltaHeight + self.init_range
                     dataOut.ippSeconds      = ippSeconds
                     dataOut.step            = self.step
                     dataOut.flagNoData      = False
                     if dataOut.flagDataAsBlock:
                         dataOut.nProfiles       = int(dataOut.nProfiles*self.nsamples)
                     else:
                         dataOut.nProfiles       = int(self.nsamples)
                     dataOut.profileIndex    = dataOut.nProfiles
                     dataOut.flagDataAsBlock = True
                     dataBlock = None
                     #print("new data shape:", dataOut.data.shape, dataOut.utctime)
                     return dataOut
             class removeProfileByFaradayHS(Operation):
                 '''
                 '''
                 __buffer_data = []
                 __buffer_times = []
                 buffer = None
                 outliers_IDs_list = []
                 __slots__ = ('n','navg','profileMargin','thHistOutlier','minHei_idx','maxHei_idx','nHeights',
                     '__dh','first_utcBlock','__profIndex','init_prof','end_prof','lenProfileOut','nChannels',
                     '__count_exec','__initime','__dataReady','__ipp')
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.isConfig = False
                 def setup(self,dataOut, n=None , navg=0.8, profileMargin=50,thHistOutlier=3, minHei=None, maxHei=None):
                     if n == None and timeInterval == None:
                         raise ValueError("nprofiles or timeInterval should be specified ...")
                     if n != None:
                         self.n = n
                     self.navg = navg
                     self.profileMargin = profileMargin
                     self.thHistOutlier = thHistOutlier
                     self.__profIndex = 0
                     self.buffer = None
                     self._ipp = dataOut.ippSeconds
                     self.n_prof_released = 0
                     self.heightList = dataOut.heightList
                     self.init_prof = 0
                     self.end_prof = 0
                     self.__count_exec = 0
                     self.__profIndex = 0
                     self.first_utcBlock = None
                     self.__dh = dataOut.heightList[1] - dataOut.heightList[0]
                     minHei = minHei
                     maxHei = maxHei
                     if minHei==None :
                         minHei = dataOut.heightList[0]
                     if maxHei==None :
                         maxHei = dataOut.heightList[-1]
                     self.minHei_idx,self.maxHei_idx =  getHei_index(minHei, maxHei, dataOut.heightList)
                     self.nChannels = dataOut.nChannels
                     self.nHeights = dataOut.nHeights
                     self.test_counter = 0
                 def filterSatsProfiles(self):
                     data = self.__buffer_data
                     #print(data.shape)
                     nChannels, profiles, heights = data.shape
                     indexes=[]
                     outliers_IDs=[]
                     for c in range(nChannels):
                         for h in range(self.minHei_idx, self.maxHei_idx):
-                            power = data[c,:,h] * numpy.conjugate(data[c,:,h])
+                            power = 10* numpy.log10((data[c,:,h] * numpy.conjugate(data[c,:,h])).real)
-                            power = power.real
+                            #power = power.real
                             #power = (numpy.abs(data[c,:,h].real))
                             sortdata = numpy.sort(power, axis=None)
                             sortID=power.argsort()
-                            index = _noise.hildebrand_sekhon2(sortdata,self.navg) #0.75-> buen valor
+                            index = _noise.hildebrand_sekhon2(sortdata,self.navg)
                             indexes.append(index)
                             outliers_IDs=numpy.append(outliers_IDs,sortID[index:])
-                            # print(outliers_IDs)
+                            # print(sortdata.min(), sortdata.max(), sortdata.mean())
                             # fig,ax = plt.subplots()
                             # #ax.set_title(str(k)+" "+str(j))
                             # x=range(len(sortdata))
                             # ax.scatter(x,sortdata)
                             # ax.axvline(index)
                             # plt.grid()
                             # plt.show()
                     outliers_IDs = outliers_IDs.astype(numpy.dtype('int64'))
                     outliers_IDs = numpy.unique(outliers_IDs)
                     outs_lines = numpy.sort(outliers_IDs)
                     # #print("outliers Ids: ", outs_lines, outs_lines.shape)
                     #hist, bin_edges = numpy.histogram(outs_lines, bins=10, density=True)
                     #Agrupando el histograma de outliers,
-                    #my_bins = numpy.linspace(0,int(profiles), int(profiles/100), endpoint=False)
+                    my_bins = numpy.linspace(0,int(profiles), int(profiles/100), endpoint=False)
-                    my_bins = numpy.linspace(0,9600, 96, endpoint=False)
+                    #my_bins = numpy.linspace(0,1600, 96, endpoint=False)
                     hist, bins =  numpy.histogram(outs_lines,bins=my_bins)
                     hist_outliers_indexes = numpy.where(hist > self.thHistOutlier)       #es outlier
                     #print(hist_outliers_indexes[0])
                     bins_outliers_indexes = [int(i) for i in bins[hist_outliers_indexes]] #
                     #print(bins_outliers_indexes)
                     outlier_loc_index = []
                     # for n in range(len(bins_outliers_indexes)-1):
                     #     for k in range(bins_outliers_indexes[n]-self.profileMargin,bins_outliers_indexes[n+1]+self.profileMargin):
                     #         outlier_loc_index.append(k)
                     outlier_loc_index = [e for n in range(len(bins_outliers_indexes)-1) for e in range(bins_outliers_indexes[n]-self.profileMargin,bins_outliers_indexes[n+1]+self.profileMargin) ]
                     outlier_loc_index = numpy.asarray(outlier_loc_index)
                     #print(len(numpy.unique(outlier_loc_index)), numpy.unique(outlier_loc_index))
-                    # x, y = numpy.meshgrid(numpy.arange(profiles), self.heightList)
+                    x, y = numpy.meshgrid(numpy.arange(profiles), self.heightList)
-                    # fig, ax = plt.subplots(1,2,figsize=(8, 6))
+                    fig, ax = plt.subplots(1,2,figsize=(8, 6))
-                    # dat = data[0,:,:].real
+                    dat = data[0,:,:].real
-                    # m = numpy.nanmean(dat)
+                    dat = 10* numpy.log10((data[0,:,:] * numpy.conjugate(data[0,:,:])).real)
-                    # o = numpy.nanstd(dat)
+                    m = numpy.nanmean(dat)
-                    # #print(m, o, x.shape, y.shape)
+                    o = numpy.nanstd(dat)
-                    # c = ax[0].pcolormesh(x, y, dat.T, cmap ='YlGnBu', vmin = (m-2*o), vmax = (m+2*o))
+                    #print(m, o, x.shape, y.shape)
-                    # ax[0].vlines(outs_lines,200,600, linestyles='dashed', label = 'outs', color='w')
+                    c = ax[0].pcolormesh(x, y, dat.T, cmap ='YlGnBu', vmin = (m-2*o), vmax = (m+2*o))
-                    # fig.colorbar(c)
+                    ax[0].vlines(outs_lines,200,600, linestyles='dashed', label = 'outs', color='w')
-                    # ax[0].vlines(outlier_loc_index,650,750, linestyles='dashed', label = 'outs', color='r')
+                    fig.colorbar(c)
-                    # ax[1].hist(outs_lines,bins=my_bins)
+                    ax[0].vlines(outlier_loc_index,650,750, linestyles='dashed', label = 'outs', color='r')
+                    ax[1].hist(outs_lines,bins=my_bins)
+                    plt.show()
+                    self.outliers_IDs_list = numpy.unique(outlier_loc_index)
+                    print("outs list: ", self.outliers_IDs_list)
+                    return data
+                def filterSatsProfiles2(self):
+                    data = self.__buffer_data
+                    #print(data.shape)
+                    nChannels, profiles, heights = data.shape
+                    indexes=numpy.zeros([], dtype=int)
+                    outliers_IDs=[]
+                    for c in range(nChannels):
+                        noise_ref =10* numpy.log10((data[c,:,550:600] * numpy.conjugate(data[c,:,550:600])).real)
+                        print("Noise ",noise_ref.mean())
+                        for h in range(self.minHei_idx, self.maxHei_idx):
+                            power = 10* numpy.log10((data[c,:,h] * numpy.conjugate(data[c,:,h])).real)
+                            #power = power.real
+                            #power = (numpy.abs(data[c,:,h].real))
+                            #sortdata = numpy.sort(power, axis=None)
+                            #sortID=power.argsort()
+                            #print(sortID)
+                            th = 60 + 10
+                            index = numpy.where(power > th )
+                            if index[0].size > 10 and index[0].size < int(0.8*profiles):
+                                indexes = numpy.append(indexes, index[0])
+                            #print(index[0])
+                            #print(index[0])
+                            # fig,ax = plt.subplots()
+                            # #ax.set_title(str(k)+" "+str(j))
+                            # x=range(len(power))
+                            # ax.scatter(x,power)
+                            # #ax.axvline(index)
+                            # plt.grid()
                             # plt.show()
+                        #print(indexes)
+                    #outliers_IDs = outliers_IDs.astype(numpy.dtype('int64'))
+                    #outliers_IDs = numpy.unique(outliers_IDs)
+                    outs_lines = numpy.unique(indexes)
+                    print("outliers Ids: ", outs_lines, outs_lines.shape)
+                    #hist, bin_edges = numpy.histogram(outs_lines, bins=10, density=True)
+                    #Agrupando el histograma de outliers,
+                    my_bins = numpy.linspace(0,int(profiles), int(profiles/100), endpoint=False)
+                    #my_bins = numpy.linspace(0,1600, 96, endpoint=False)
+                    hist, bins =  numpy.histogram(outs_lines,bins=my_bins)
+                    hist_outliers_indexes = numpy.where(hist > self.thHistOutlier)       #es outlier
+                    #print(hist_outliers_indexes[0])
+                    bins_outliers_indexes = [int(i) for i in bins[hist_outliers_indexes]] #
+                    #print(bins_outliers_indexes)
+                    outlier_loc_index = []
+                    outlier_loc_index = [e for n in range(len(bins_outliers_indexes)-1) for e in range(bins_outliers_indexes[n]-self.profileMargin,bins_outliers_indexes[n+1]+self.profileMargin) ]
+                    outlier_loc_index = numpy.asarray(outlier_loc_index)
+                    outlier_loc_index = outlier_loc_index[~numpy.all(outlier_loc_index < 0)]
+                    print("outliers final: ", outlier_loc_index)
+                    x, y = numpy.meshgrid(numpy.arange(profiles), self.heightList)
+                    fig, ax = plt.subplots(1,2,figsize=(8, 6))
+                    dat = data[0,:,:].real
+                    dat = 10* numpy.log10((data[0,:,:] * numpy.conjugate(data[0,:,:])).real)
+                    m = numpy.nanmean(dat)
+                    o = numpy.nanstd(dat)
+                    #print(m, o, x.shape, y.shape)
+                    c = ax[0].pcolormesh(x, y, dat.T, cmap ='YlGnBu', vmin = (m-2*o), vmax = (m+2*o))
+                    ax[0].vlines(outs_lines,200,600, linestyles='dashed', label = 'outs', color='w')
+                    fig.colorbar(c)
+                    ax[0].vlines(outlier_loc_index,650,750, linestyles='dashed', label = 'outs', color='r')
+                    ax[1].hist(outs_lines,bins=my_bins)
+                    plt.show()
                     self.outliers_IDs_list = numpy.unique(outlier_loc_index)
+                    print("outs list: ", self.outliers_IDs_list)
                     return data
                 def cleanSpikesFFT2D(self):
                     incoh_int = 10
                     norm_img = 75
                     import matplotlib.pyplot as plt
                     import datetime
                     import cv2
                     data = self.__buffer_data
                     print("cleaning shape inpt: ",data.shape)
                     self.__buffer_data = []
                     channels , profiles, heights = data.shape
                     len_split_prof = profiles / incoh_int
                     for ch in range(channels):
                         data_10 = numpy.split(data[ch, :, self.minHei_idx:], incoh_int, axis=0) # división de los perfiles
                         print("splited data: ",len(data_10)," -> ", data_10[0].shape)
                         int_img = None
                         i_count = 0
                         n_x, n_y = data_10[0].shape
                         for s_data in data_10: #porciones de espectro
                             spectrum = numpy.fft.fft2(s_data, axes=(0,1))
                             z = numpy.abs(spectrum)
                             mg = z[2:n_y,:] #omitir dc y adjunto
                             dat = numpy.log10(mg.T)
                             i_count += 1
                             if i_count == 1:
                                 int_img = dat
                             else:
                                 int_img += dat
                                 #print(i_count)
                         min, max = int_img.min(), int_img.max()
                         int_img = ((int_img-min)*255/(max-min)).astype(numpy.uint8)
                         cv2.imshow('integrated image', int_img) #numpy.fft.fftshift(img))
                         cv2.waitKey(0)
                         #####################################################################
                         kernel_h = numpy.zeros((3,3)) #
                         kernel_h[0, :] = -2
                         kernel_h[1, :] = 3
                         kernel_h[2, :] = -2
                         kernel_5h = numpy.zeros((5,5)) #
                         kernel_5h[0, :] = -2
                         kernel_5h[1, :] = -1
                         kernel_5h[2, :] = 5
                         kernel_5h[3, :] = -1
                         kernel_5h[4, :] = -2
                         #####################################################################
                         sharp_img = cv2.filter2D(src=int_img, ddepth=-1, kernel=kernel_5h)
                         # cv2.imshow('sharp image h ', sharp_img)
                         # cv2.waitKey(0)
                         #####################################################################
                         horizontal_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (5,1)) #11
                         #####################################################################
                         detected_lines_h = cv2.morphologyEx(sharp_img, cv2.MORPH_OPEN, horizontal_kernel, iterations=1)
                         # cv2.imshow('lines horizontal', detected_lines_h) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         #####################################################################
                         ret, detected_lines_h = cv2.threshold(detected_lines_h, 200, 255, cv2.THRESH_BINARY)#
                         cv2.imshow('binary img', detected_lines_h) #numpy.fft.fftshift(detected_lines_h))
                         cv2.waitKey(0)
                         #####################################################################
                         cnts_h, h0 = cv2.findContours(detected_lines_h, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         #####################################################################
                         h_line_index = []
                         v_line_index = []
                         #cnts_h += cnts_h_s #combine large and small lines
                         # line indexes  x1, x2, y
                         for c in cnts_h:
                             #print(c)
                             if len(c) < 3: #contorno linea
                                 x1 = c[0][0][0]
                                 x2 = c[1][0][0]
                                 if x1 > 5 and x2 < (n_x-5) :
                                     start = incoh_int + (x1 * incoh_int)
                                     end = incoh_int + (x2 * incoh_int)
                                     h_line_index.append( [start, end, c[0][0][1]] )
                                 #print("x1, x2, y", c[0][0][0],c[1][0][0], c[0][0][1])
                             else:     #contorno poligono
                                 pairs = numpy.asarray([c[n][0]  for n in range(len(c))])
                                 y = numpy.unique(pairs[:,1])
                                 x = numpy.unique(pairs[:,0])
                                 #print(x)
                                 for yk in y:
                                     x0 = x[0]
                                     if x0 < 8:
                                         x0 = 10
                                     #print(x[0], x[-1], yk)
                                     h_line_index.append( [x0, x[-1], yk])
                                     #print("x1, x2, y ->p ", x[0], x[-1], yk)
                         ###################################################################
                         #print("Cleaning")
                         # #  clean Spectrum
                         spectrum = numpy.fft.fft2(data[ch,:,self.minHei_idx:], axes=(0,1))
                         z = numpy.abs(spectrum)
                         phase = numpy.angle(spectrum)
                         print("Total Horizontal", len(h_line_index))
                         if len(h_line_index) < 75 :
                             for x1, x2, y in h_line_index:
                                 print(x1, x2, y)
                                 z[x1:x2,y] =  0
                         spcCleaned = z * numpy.exp(1j*phase)
                         dat2 = numpy.log10(z[1:-1,:].T)
                         min, max =dat2.min(), dat2.max()
                         print(min, max)
                         img2 = ((dat2-min)*255/(max-min)).astype(numpy.uint8)
                         cv2.imshow('cleaned', img2) #numpy.fft.fftshift(img_cleaned))
                         cv2.waitKey(0)
                         cv2.destroyAllWindows()
                         data[ch,:,self.minHei_idx:] = numpy.fft.ifft2(spcCleaned, axes=(0,1))
                     #print("cleanOutliersByBlock Done", data.shape)
                     self.__buffer_data = data
                     return data
                 def cleanOutliersByBlock(self):
                     import matplotlib.pyplot as plt
                     import datetime
                     import cv2
                     #print(self.__buffer_data[0].shape)
                     data = self.__buffer_data#.copy()
                     print("cleaning shape inpt: ",data.shape)
                     self.__buffer_data = []
                     spectrum = numpy.fft.fft2(data[:,:,self.minHei_idx:], axes=(1,2))
                     print("spc : ",spectrum.shape)
                     (nch,nsamples, nh) = spectrum.shape
                     data2 = None
                     #print(data.shape)
                     cleanedBlock = None
                     spectrum2 = spectrum.copy()
                     for ch in range(nch):
                         dh = self.__dh
                         dt1 = (dh*1.0e-6)/(0.15)
                         dt2 = self.__buffer_times[1]-self.__buffer_times[0]
                         freqv = numpy.fft.fftfreq(nh, d=dt1)
                         freqh = numpy.fft.fftfreq(self.n, d=dt2)
                         z = numpy.abs(spectrum[ch,:,:])
                         phase = numpy.angle(spectrum[ch,:,:])
                         z1 = z[0,:]
                         #print("shape z: ", z.shape, nsamples)
                         dat = numpy.log10(z[1:nsamples,:].T)
                         pdat = numpy.log10(phase.T)
                         #print("dat mean",dat.mean())
                         mean, min, max = dat.mean(), dat.min(), dat.max()
                         img = ((dat-min)*200/(max-min)).astype(numpy.uint8)
                         # print(img.shape)
                         cv2.imshow('image', img) #numpy.fft.fftshift(img))
                         cv2.waitKey(0)
                         ''' #FUNCIONA LINEAS PEQUEÑAS
                         kernel_5h = numpy.zeros((5,3)) #
                         kernel_5h[0, :] = 2
                         kernel_5h[1, :] = 1
                         kernel_5h[2, :] = 0
                         kernel_5h[3, :] = -1
                         kernel_5h[4, :] = -2
                         sharp_imgh = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_5h)
                         cv2.imshow('sharp image  h',sharp_imgh)
                         cv2.waitKey(0)
                         horizontal_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (20,1))
                         detected_lines_h = cv2.morphologyEx(sharp_imgh, cv2.MORPH_OPEN, horizontal_kernel, iterations=1)
                         #detected_lines_h = cv2.medianBlur(detected_lines_h, 3)
                         #detected_lines_h = cv2.filter2D(src=img, ddepth=-1, kernel=kernel)
                         cv2.imshow('lines h  gray', detected_lines_h)
                         cv2.waitKey(0)
                         reth, detected_lines_h = cv2.threshold(detected_lines_h, 90, 255, cv2.THRESH_BINARY)
                         cv2.imshow('lines h ', detected_lines_h)
                         cv2.waitKey(0)
                         '''
                         '''
                         kernel_3h = numpy.zeros((3,10)) #10
                         kernel_3h[0, :] = -1
                         kernel_3h[1, :] = 2
                         kernel_3h[2, :] = -1
                         kernel_h = numpy.zeros((3,20)) #20
                         kernel_h[0, :] = -1
                         kernel_h[1, :] = 2
                         kernel_h[2, :] = -1
                         kernel_v = numpy.zeros((30,3)) #30
                         kernel_v[:, 0] = -1
                         kernel_v[:, 1] = 2
                         kernel_v[:, 2] = -1
                         kernel_4h = numpy.zeros((4,20)) #
                         kernel_4h[0, :] = 1
                         kernel_4h[1, :] = 0
                         kernel_4h[2, :] = 0
                         kernel_4h[3, :] = -1
                         kernel_5h = numpy.zeros((5,30)) #
                         kernel_5h[0, :] = 2
                         kernel_5h[1, :] = 1
                         kernel_5h[2, :] = 0
                         kernel_5h[3, :] = -1
                         kernel_5h[4, :] = -2
                         sharp_img0 = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_3h)
                         # cv2.imshow('sharp image small h',sharp_img0) # numpy.fft.fftshift(sharp_img1))
                         # cv2.waitKey(0)
                         sharp_img1 = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_h)
                         # cv2.imshow('sharp image h',sharp_img1) # numpy.fft.fftshift(sharp_img1))
                         # cv2.waitKey(0)
                         sharp_img2 = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_v)
                         # cv2.imshow('sharp image v', sharp_img2) #numpy.fft.fftshift(sharp_img2))
                         # cv2.waitKey(0)
                         sharp_imgw = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_4h)
                         # cv2.imshow('sharp image h wide', sharp_imgw) #numpy.fft.fftshift(sharp_img2))
                         # cv2.waitKey(0)
                         sharp_imgwl = cv2.filter2D(src=img, ddepth=-1, kernel=kernel_5h,	borderType = cv2.BORDER_ISOLATED)
                         cv2.imshow('sharp image h long wide', sharp_imgwl) #numpy.fft.fftshift(sharp_img2))
                         cv2.waitKey(0)
                         # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/sharp_h.jpg', sharp_img1)
                         # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/sharp_v.jpg', sharp_img2)
                         # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/input_img.jpg', img)
                         ########################small horizontal
                         horizontal_kernel_s = cv2.getStructuringElement(cv2.MORPH_RECT, (11,1)) #11
                         ######################## horizontal
                         horizontal_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (30,1)) #30
                         ######################## vertical
                         vertical_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (1,50)) #50
                         ######################## horizontal wide
                         horizontal_kernel_w = cv2.getStructuringElement(cv2.MORPH_RECT, (30,1)) # 30
                         horizontal_kernel_expand = cv2.getStructuringElement(cv2.MORPH_RECT, (3,3)) #
                         horizontal_kernel_wl = cv2.getStructuringElement(cv2.MORPH_RECT, (50,1)) #
                         detected_lines_h_s = cv2.morphologyEx(sharp_img0, cv2.MORPH_OPEN, horizontal_kernel_s, iterations=7) #7
                         detected_lines_h = cv2.morphologyEx(sharp_img1, cv2.MORPH_OPEN, horizontal_kernel, iterations=7) #7
                         detected_lines_v = cv2.morphologyEx(sharp_img2, cv2.MORPH_OPEN, vertical_kernel, iterations=7) #7
                         detected_lines_h_w = cv2.morphologyEx(sharp_imgw, cv2.MORPH_OPEN, horizontal_kernel_w, iterations=5) #5
                         detected_lines_h_wl = cv2.morphologyEx(sharp_imgwl, cv2.MORPH_OPEN, horizontal_kernel_wl, iterations=5) #
                         detected_lines_h_wl = cv2.filter2D(src=detected_lines_h_wl, ddepth=-1, kernel=horizontal_kernel_expand)
                         # cv2.imshow('lines h small gray', detected_lines_h_s) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         # cv2.imshow('lines h  gray', detected_lines_h) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         # cv2.imshow('lines v  gray', detected_lines_v) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         # cv2.imshow('lines h wide  gray', detected_lines_h_w) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         cv2.imshow('lines h long wide  gray', detected_lines_h_wl) #numpy.fft.fftshift(detected_lines_h))
                         cv2.waitKey(0)
                         reth_s, detected_lines_h_s = cv2.threshold(detected_lines_h_s, 85, 255, cv2.THRESH_BINARY)# 85
                         reth, detected_lines_h = cv2.threshold(detected_lines_h, 30, 255, cv2.THRESH_BINARY) #30
                         retv, detected_lines_v = cv2.threshold(detected_lines_v, 30, 255, cv2.THRESH_BINARY) #30
                         reth_w, detected_lines_h_w = cv2.threshold(detected_lines_h_w, 35, 255, cv2.THRESH_BINARY)#
                         reth_wl, detected_lines_h_wl = cv2.threshold(detected_lines_h_wl, 200, 255, cv2.THRESH_BINARY)#
                         cnts_h_s, h0 = cv2.findContours(detected_lines_h_s, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         cnts_h, h1 = cv2.findContours(detected_lines_h, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         cnts_v, h2 = cv2.findContours(detected_lines_v, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         cnts_h_w, h3 = cv2.findContours(detected_lines_h_w, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         cnts_h_wl, h4 = cv2.findContours(detected_lines_h_wl, cv2.RETR_LIST, cv2.CHAIN_APPROX_SIMPLE)
                         #print("horizontal ", cnts_h)
                         #print("vertical ", cnts_v)
                         # cnts, h = cv2.findContours(detected_lines_h, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
                         # print(cnts)
                         # cv2.imshow('lines h wide', detected_lines_h_w) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         cv2.imshow('lines h wide long ', detected_lines_h_wl) #numpy.fft.fftshift(detected_lines_v))
                         cv2.waitKey(0)
                         # cv2.imshow('lines h small', detected_lines_h_s) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         # cv2.imshow('lines h ', detected_lines_h) #numpy.fft.fftshift(detected_lines_h))
                         # cv2.waitKey(0)
                         # cv2.imshow('lines v ', detected_lines_v) #numpy.fft.fftshift(detected_lines_v))
                         # cv2.waitKey(0)
                         # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/lines_h.jpg', detected_lines_h)
                         # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/lines_v.jpg', detected_lines_v)
                         #cnts = cnts[0] if len(cnts) == 2 else cnts[1]
                         #y_line_index = numpy.asarray([ [c[0][0][0],c[1][0][0], c[0][0][1]]  for c in cnts_v ])
                         h_line_index = []
                         v_line_index = []
                         cnts_h += cnts_h_s #combine large and small lines
                         # line indexes  x1, x2, y
                         for c in cnts_h:
                             #print(c)
                             if len(c) < 3: #contorno linea
                                 x1 = c[0][0][0]
                                 if x1 < 8:
                                     x1 = 10
                                 h_line_index.append( [x1,c[1][0][0], c[0][0][1]] )
                                 #print("x1, x2, y", c[0][0][0],c[1][0][0], c[0][0][1])
                             else:     #contorno poligono
                                 pairs = numpy.asarray([c[n][0]  for n in range(len(c))])
                                 y = numpy.unique(pairs[:,1])
                                 x = numpy.unique(pairs[:,0])
                                 #print(x)
                                 for yk in y:
                                     x0 = x[0]
                                     if x0 < 8:
                                         x0 = 10
                                     #print(x[0], x[-1], yk)
                                     h_line_index.append( [x0, x[-1], yk])
                                     #print("x1, x2, y ->p ", x[0], x[-1], yk)
                         for c in cnts_h_w:
                             #print(c)
                             if len(c) < 3: #contorno linea
                                 x1 = c[0][0][0]
                                 if x1 < 8:
                                     x1 = 10
                                 y = c[0][0][1] - 2 # se incrementa 2 líneas x el filtro
                                 h_line_index.append( [x1,c[1][0][0],y] )
                                 #print("x1, x2, y", c[0][0][0],c[1][0][0], c[0][0][1])
                             else:     #contorno poligono
                                 pairs = numpy.asarray([c[n][0]  for n in range(len(c))])
                                 y = numpy.unique(pairs[:,1])
                                 x = numpy.unique(pairs[:,0])
                                 #print(x)
                                 for yk in y:
                                     x0 = x[0]
                                     if x0 < 8:
                                         x0 = 10
                                     h_line_index.append( [x0, x[-1], yk-2])
                         for c in cnts_h_wl:                   # # revisar
                             #print(c)
                             if len(c) < 3: #contorno linea
                                 x1 = c[0][0][0]
                                 if x1 < 8:
                                     x1 = 10
                                 y = c[0][0][1] - 2 # se incrementa 2 líneas x el filtro
                                 h_line_index.append( [x1,c[1][0][0],y] )
                                 #print("x1, x2, y", c[0][0][0],c[1][0][0], c[0][0][1])
                             else:     #contorno poligono
                                 pairs = numpy.asarray([c[n][0]  for n in range(len(c))])
                                 y = numpy.unique(pairs[:,1])
                                 x = numpy.unique(pairs[:,0])
                                 for yk in range(y[-1]-y[0]):
                                     y_k = yk +y[0]
                                     x0 = x[0]
                                     if x0 < 8:
                                         x0 = 10
                                     h_line_index.append( [x0, x[-1], y_k-2])
                         print([[c[0][0][1],c[1][0][1], c[0][0][0] ] for c in cnts_v])
                         # line indexes  y1, y2, x
                         for c in cnts_v:
                             if len(c) < 3: #contorno linea
                                 v_line_index.append( [c[0][0][1],c[1][0][1], c[0][0][0] ] )
                             else:     #contorno poligono
                                 pairs = numpy.asarray([c[n][0]  for n in range(len(c))])
                                 #print(pairs)
                                 y = numpy.unique(pairs[:,1])
                                 x = numpy.unique(pairs[:,0])
                                 for xk in x:
                                     #print(x[0], x[-1], yk)
                                     v_line_index.append( [y[0],y[-1], xk])
                         ###################################################################
                         # #  clean Horizontal
                         print("Total Horizontal", len(h_line_index))
                         if len(h_line_index) < 75 :
                             for x1, x2, y in h_line_index:
                                 #print("cleaning ",x1, x2, y)
                                 len_line = x2 - x1
                                 if y > 10 and y < (nh -10):
                                     # if y != (nh-1):
                                     #     list = [ ((z[n, y-1] + z[n,y+1])/2) for n in range(len_line)]
                                     # else:
                                     #     list = [ ((z[n, y-1] + z[n,0])/2) for n in range(len_line)]
                                     #
                                     # z[x1:x2,y] =  numpy.asarray(list)
                                     z[x1-5:x2+5,y:y+1] =  0
                         # clean vertical
                         for y1, y2, x in v_line_index:
                             len_line = y2 - y1
                             #print(x)
                             if x > 0 and x < (nsamples-2):
                                 # if x != (nsamples-1):
                                 #     list = [ ((z[x-2, n] + z[x+2,n])/2) for n in range(len_line)]
                                 # else:
                                 #     list = [ ((z[x-2, n] + z[1,n])/2) for n in range(len_line)]
                                 #
                                 # #z[x-1:x+1,y1:y2] =  numpy.asarray(list)
                                 #
                                 z[x+1,y1:y2] =  0
                         '''
                         #z[: ,[215, 217, 221, 223, 225, 340, 342, 346, 348, 350, 465, 467, 471, 473, 475]]=0
                         z[1: ,[112, 114, 118, 120, 122, 237, 239, 245, 247, 249, 362, 364, 368, 370, 372]]=0
                         # z[: ,217]=0
                         # z[: ,221]=0
                         # z[: ,223]=0
                         # z[: ,225]=0
                         dat2 = numpy.log10(z.T)
                         #print(dat2)
                         max = dat2.max()
                         #print(" min, max ", max, min)
                         img2 = ((dat2-min)*255/(max-min)).astype(numpy.uint8)
                         #img_cleaned = img2.copy()
                         #cv2.drawContours(img2, cnts_h, -1, (255,255,255), 1)
                         #cv2.drawContours(img2, cnts_v, -1, (255,255,255), 1)
                         spcCleaned = z * numpy.exp(1j*phase)
                         #print(spcCleaned)
                         # cv2.imshow('image contours', img2) #numpy.fft.fftshift(img))
                         # cv2.waitKey(0)
                         cv2.imshow('cleaned', img2) #numpy.fft.fftshift(img_cleaned))
                         cv2.waitKey(0)
                         # # cv2.imwrite('/home/soporte/Data/AMISR14/ISR/spc/spc/cleaned_{}.jpg'.format(self.test_counter), img2)
                         cv2.destroyAllWindows()
                         # self.test_counter += 1
                         #print("DC difference " ,z1 - z[0,:])
                         # m = numpy.mean(dat)
                         # o = numpy.std(dat)
                         # print("mean ", m, " std ", o)
                         # fig, ax = plt.subplots(1,2,figsize=(12, 6))
                         # #X, Y = numpy.meshgrid(numpy.sort(freqh),numpy.sort(freqv))
                         # X, Y = numpy.meshgrid(numpy.fft.fftshift(freqh),numpy.fft.fftshift(freqv))
                         #
                         # colormap = 'jet'
                         # #c = ax[0].pcolormesh(x, y, dat, cmap =colormap, vmin = (m-2*o)/2, vmax = (m+2*o))
                         # #c = ax[0].pcolormesh(X, Y, numpy.fft.fftshift(dat), cmap =colormap, vmin = 6.5, vmax = 6.8)
                         # c = ax[0].pcolormesh(X, Y, numpy.fft.fftshift(dat), cmap =colormap, vmin = (m-2*o), vmax = (m+1.5*o))
                         # fig.colorbar(c, ax=ax[0])
                         #
                         #
                         # #c = ax.pcolor( z.T , cmap ='gray', vmin = (m-2*o), vmax = (m+2*o))
                         # #date_time = datetime.datetime.fromtimestamp(self.__buffer_times[0]).strftime('%Y-%m-%d %H:%M:%S.%f')
                         # #strftime('%Y-%m-%d %H:%M:%S')
                         # #ax[0].set_title('Spectrum magnitude '+date_time)
                         # #fig.canvas.set_window_title('Spectrum magnitude {} '.format(self.n)+date_time)
                         # #print("aqui estoy2",dat2[:,:,0].shape)
                         # #c = ax[1].pcolormesh(X, Y, numpy.fft.fftshift(pdat), cmap =colormap, vmin = 4.2, vmax = 5.0)
                         # c = ax[0].pcolormesh(X, Y, numpy.fft.fftshift(dat2), cmap =colormap, vmin = (m-2*o), vmax = (m+1.5*o))
                         # #c = ax[1].pcolormesh(X, Y, numpy.fft.fftshift(pdat), cmap =colormap ) #, vmin = 0.0, vmax = 0.5)
                         # #c = ax[1].pcolormesh(x, y, dat2[:,:,0], cmap =colormap, vmin = (m-2*o)/2, vmax = (m+2*o)-1)
                         # #print("aqui estoy3")
                         # fig.colorbar(c, ax=ax[1])
                         # plt.show()
                         spectrum[ch,:,:] = spcCleaned
                     #print(data2.shape)
                     data[:,:,self.minHei_idx:] = numpy.fft.ifft2(spectrum, axes=(1,2))
                     #print("cleanOutliersByBlock Done", data.shape)
                     self.__buffer_data = data
                     return data
                 def fillBuffer(self, data, datatime):
                     if self.__profIndex == 0:
                         self.__buffer_data = data.copy()
                     else:
                         self.__buffer_data = numpy.concatenate((self.__buffer_data,data), axis=1)#en perfiles
                     self.__profIndex += 1
                     self.__buffer_times.append(datatime)
                 def getData(self, data, datatime=None):
                     if self.__profIndex  == 0:
                         self.__initime = datatime
                     self.__dataReady = False
                     self.fillBuffer(data, datatime)
                     dataBlock = None
                     if self.__profIndex == self.n:
                         #print("apnd : ",data)
                         #dataBlock = self.cleanOutliersByBlock()
                         #dataBlock = self.cleanSpikesFFT2D()
-                        dataBlock = self.filterSatsProfiles()
+                        dataBlock = self.filterSatsProfiles2()
                         self.__dataReady = True
                     return dataBlock
                     if dataBlock is None:
                         return None, None
                     return dataBlock
                 def releaseBlock(self):
                     if self.n % self.lenProfileOut != 0:
                         raise ValueError("lenProfileOut %d must be submultiple of nProfiles %d" %(self.lenProfileOut, self.n))
                         return None
                     data = self.buffer[:,self.init_prof:self.end_prof:,:]  #ch, prof, alt
                     self.init_prof = self.end_prof
                     self.end_prof += self.lenProfileOut
                     #print("data release shape: ",dataOut.data.shape, self.end_prof)
                     self.n_prof_released += 1
                     #print("f_no_data ", dataOut.flagNoData)
                     return data
                 def run(self, dataOut, n=None, navg=0.8, nProfilesOut=1, profile_margin=50,th_hist_outlier=3,minHei=None, maxHei=None):
                     #print("run op buffer 2D",dataOut.ippSeconds)
                     # self.nChannels = dataOut.nChannels
                     # self.nHeights = dataOut.nHeights
                     if not self.isConfig:
                         #print("init p idx: ", dataOut.profileIndex )
                         self.setup(dataOut,n=n, navg=navg,profileMargin=profile_margin,
                                 thHistOutlier=th_hist_outlier,minHei=minHei, maxHei=maxHei)
                         self.isConfig = True
                     dataBlock = None
                     if not dataOut.buffer_empty: #hay datos acumulados
                         if self.init_prof == 0:
                             self.n_prof_released = 0
                             self.lenProfileOut = nProfilesOut
                             dataOut.flagNoData = False
                             #print("tp 2 ",dataOut.data.shape)
                             self.init_prof = 0
                             self.end_prof = self.lenProfileOut
                             dataOut.nProfiles = self.lenProfileOut
                         if nProfilesOut == 1:
                             dataOut.flagDataAsBlock = False
                         else:
                             dataOut.flagDataAsBlock = True
                         #print("prof: ",self.init_prof)
                         dataOut.flagNoData = False
                         if numpy.isin(self.n_prof_released, self.outliers_IDs_list):
-                            #print("omitting: ", self.n_prof_released)
+                            print("omitting: ", self.n_prof_released)
                             dataOut.flagNoData = True
                         dataOut.ippSeconds = self._ipp
                         dataOut.utctime = self.first_utcBlock + self.init_prof*self._ipp
                         # print("time: ", dataOut.utctime, self.first_utcBlock, self.init_prof,self._ipp,dataOut.ippSeconds)
                         #dataOut.data = self.releaseBlock()
                         #########################################################3
                         if self.n % self.lenProfileOut != 0:
                             raise ValueError("lenProfileOut %d must be submultiple of nProfiles %d" %(self.lenProfileOut, self.n))
                             return None
                         dataOut.data = self.buffer[:,self.init_prof:self.end_prof:,:]  #ch, prof, alt
                         self.init_prof = self.end_prof
                         self.end_prof += self.lenProfileOut
-                        # #print("data release shape: ",dataOut.data.shape, self.end_prof, dataOut.flagNoData)
+                        #print("data release shape: ",dataOut.data.shape, self.end_prof, dataOut.flagNoData)
                         self.n_prof_released += 1
                         if self.end_prof >= (self.n +self.lenProfileOut):
                             self.init_prof = 0
                             self.__profIndex = 0
                             self.buffer = None
                             dataOut.buffer_empty = True
                             self.outliers_IDs_list = []
                             self.n_prof_released = 0
                             dataOut.flagNoData = False #enviar ultimo aunque sea outlier :(
                             #print("cleaning...", dataOut.buffer_empty)
                         dataOut.profileIndex    = 0 #self.lenProfileOut
                         ####################################################################
                         return dataOut
                     #print("tp 223 ",dataOut.data.shape)
                     dataOut.flagNoData = True
                     try:
                         #dataBlock = self.getData(dataOut.data.reshape(self.nChannels,1,self.nHeights), dataOut.utctime)
                         dataBlock = self.getData(numpy.reshape(dataOut.data,(self.nChannels,1,self.nHeights)), dataOut.utctime)
                         self.__count_exec +=1
                     except Exception as e:
                         print("Error getting profiles data",self.__count_exec )
                         print(e)
                         sys.exit()
                     if self.__dataReady:
                         #print("omitting: ", len(self.outliers_IDs_list))
                         self.__count_exec = 0
                         #dataOut.data =
                         #self.buffer = numpy.flip(dataBlock, axis=1)
                         self.buffer = dataBlock
                         self.first_utcBlock = self.__initime
                         dataOut.utctime = self.__initime
                         dataOut.nProfiles  = self.__profIndex
                         #dataOut.flagNoData = False
                         self.init_prof = 0
                         self.__profIndex = 0
                         self.__initime = None
                         dataBlock = None
                         self.__buffer_times = []
                         dataOut.error = False
                         dataOut.useInputBuffer = True
                         dataOut.buffer_empty = False
                         #print("1 ch: {}  prof: {}  hs: {}".format(int(dataOut.nChannels),int(dataOut.nProfiles),int(dataOut.nHeights)))
                     #print(self.__count_exec)
                     return dataOut
             class RemoveProfileSats(Operation):
                 '''
-                Omite los perfiles contaminados con señal de satelites,
+                Omite los perfiles contaminados con señal de satélites, usando una altura de referencia
                 In: minHei = min_sat_range
                     max_sat_range
                     min_hei_ref
                     max_hei_ref
                     th = diference between profiles mean, ref and sats
                 Out:
                     profile clean
                 '''
-                isConfig = False
-                min_sats = 0
-                max_sats = 999999999
-                min_ref= 0
-                max_ref= 9999999999
-                needReshape = False
-                count = 0
-                thdB = 0
-                byRanges = False
-                min_sats = None
-                max_sats = None
-                noise = 0
+                __buffer_data = []
+                __buffer_times = []
+                buffer = None
+                outliers_IDs_list = []
+                __slots__ = ('n','navg','profileMargin','thHistOutlier','minHei_idx','maxHei_idx','nHeights',
+                    'first_utcBlock','__profIndex','init_prof','end_prof','lenProfileOut','nChannels',
+                    '__count_exec','__initime','__dataReady','__ipp', 'minRef', 'maxRef', 'thdB')
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                     self.isConfig = False
+                def setup(self,dataOut, n=None , navg=0.8, profileMargin=50,thHistOutlier=3,
+                                    minHei=None, maxHei=None, minRef=None, maxRef=None, thdB=10):
+                    if n == None and timeInterval == None:
+                        raise ValueError("nprofiles or timeInterval should be specified ...")
-                def setup(self, dataOut, minHei, maxHei, minRef, maxRef, th, thdB, rangeHeiList):
+                    if n != None:
+                        self.n = n
-                    if rangeHeiList!=None:
+                    self.navg = navg
-                        self.byRanges = True
+                    self.profileMargin = profileMargin
-                    else:
+                    self.thHistOutlier = thHistOutlier
-                        if minHei==None or  maxHei==None :
+                    self.__profIndex = 0
-                            raise ValueError("Parameters heights are required")
+                    self.buffer = None
-                    if minRef==None or  maxRef==None:
+                    self._ipp = dataOut.ippSeconds
-                        raise ValueError("Parameters heights are required")
+                    self.n_prof_released = 0
+                    self.heightList = dataOut.heightList
-                    if self.byRanges:
+                    self.init_prof = 0
-                        self.min_sats = []
+                    self.end_prof = 0
-                        self.max_sats = []
+                    self.__count_exec = 0
-                        for min,max in rangeHeiList:
+                    self.__profIndex = 0
-                            a,b = getHei_index(min, max, dataOut.heightList)
+                    self.first_utcBlock = None
-                            self.min_sats.append(a)
+                    #self.__dh = dataOut.heightList[1] - dataOut.heightList[0]
-                            self.max_sats.append(b)
+                    minHei = minHei
-                    else:
+                    maxHei = maxHei
-                        self.min_sats, self.max_sats = getHei_index(minHei, maxHei, dataOut.heightList)
+                    if minHei==None :
+                        minHei = dataOut.heightList[0]
+                    if maxHei==None :
+                        maxHei = dataOut.heightList[-1]
+                    self.minHei_idx,self.maxHei_idx =  getHei_index(minHei, maxHei, dataOut.heightList)
                     self.min_ref, self.max_ref = getHei_index(minRef, maxRef, dataOut.heightList)
-                    self.th = th
+                    self.nChannels = dataOut.nChannels
+                    self.nHeights = dataOut.nHeights
+                    self.test_counter = 0
                     self.thdB = thdB
-                    self.isConfig = True
+                def filterSatsProfiles(self):
+                    data = self.__buffer_data
+                    #print(data.shape)
+                    nChannels, profiles, heights = data.shape
+                    indexes=numpy.zeros([], dtype=int)
+                    outliers_IDs=[]
+                    for c in range(nChannels):
+                        #print(self.min_ref,self.max_ref)
+                        noise_ref = 10* numpy.log10((data[c,:,self.min_ref:self.max_ref] * numpy.conjugate(data[c,:,self.min_ref:self.max_ref])).real)
+                        #print("Noise ",numpy.percentile(noise_ref,95))
+                        p95 = numpy.percentile(noise_ref,95)
+                        noise_ref = noise_ref.mean()
+                        #print("Noise ",noise_ref
-                def compareRanges(self,data, minHei,maxHei):
-                    # ref =  data[0,self.min_ref:self.max_ref] * numpy.conjugate(data[0,self.min_ref:self.max_ref])
+                        for h in range(self.minHei_idx, self.maxHei_idx):
-                    # p_ref = 10*numpy.log10(ref.real)
+                            power = 10* numpy.log10((data[c,:,h] * numpy.conjugate(data[c,:,h])).real)
-                    # m_ref =  numpy.mean(p_ref)
+                            #th = noise_ref + self.thdB
+                            th = noise_ref + 1.5*(p95-noise_ref)
+                            index = numpy.where(power > th )
+                            if index[0].size > 10 and index[0].size < int(self.navg*profiles):
+                                indexes = numpy.append(indexes, index[0])
+                            #print(index[0])
+                            #print(index[0])
-                    m_ref = self.noise
+                            # fig,ax = plt.subplots()
+                            # #ax.set_title(str(k)+" "+str(j))
+                            # x=range(len(power))
+                            # ax.scatter(x,power)
+                            # #ax.axvline(index)
+                            # plt.grid()
+                            # plt.show()
+                        #print(indexes)
-                    sats = data[0,minHei:maxHei] * numpy.conjugate(data[0,minHei:maxHei])
+                    #outliers_IDs = outliers_IDs.astype(numpy.dtype('int64'))
-                    p_sats = 10*numpy.log10(sats.real)
+                    #outliers_IDs = numpy.unique(outliers_IDs)
-                    m_sats = numpy.mean(p_sats)
-                    if m_sats > (m_ref + self.th): #and (m_sats > self.thdB):
+                    outs_lines = numpy.unique(indexes)
-                        #print("msats: ",m_sats," \tmRef: ", m_ref, "\t",(m_sats - m_ref))
-                        #print("Removing profiles...")
-                        return False
-                    return True
-                def isProfileClean(self, data):
+                    #Agrupando el histograma de outliers,
-                    '''
+                    my_bins = numpy.linspace(0,int(profiles), int(profiles/100), endpoint=True)
-                    Analiza solo 1 canal, y descarta todos...
-                    '''
-                    clean = True
-                    if self.byRanges:
+                    hist, bins =  numpy.histogram(outs_lines,bins=my_bins)
+                    hist_outliers_indexes = numpy.where(hist > self.thHistOutlier)       #es outlier
+                    hist_outliers_indexes = hist_outliers_indexes[0]
+                    if len(hist_outliers_indexes>0):
+                        hist_outliers_indexes = numpy.append(hist_outliers_indexes,hist_outliers_indexes[-1]+1)
+                    #print(hist_outliers_indexes)
+                    #print(bins, hist_outliers_indexes)
+                    bins_outliers_indexes = [int(i) for i in (bins[hist_outliers_indexes])] #
+                    outlier_loc_index = [e for n in range(len(bins_outliers_indexes)) for e in range(bins_outliers_indexes[n]-self.profileMargin,bins_outliers_indexes[n]+ profiles//100 + self.profileMargin) ]
+                    outlier_loc_index = numpy.asarray(outlier_loc_index)
+                    #print("outliers Ids: ", outlier_loc_index, outlier_loc_index.shape)
+                    outlier_loc_index = outlier_loc_index[ (outlier_loc_index >= 0) & (outlier_loc_index<profiles)]
+                    #print("outliers final: ", outlier_loc_index)
+                    # from matplotlib import pyplot as plt
+                    # x, y = numpy.meshgrid(numpy.arange(profiles), self.heightList)
+                    # fig, ax = plt.subplots(1,2,figsize=(8, 6))
+                    # dat = data[0,:,:].real
+                    # dat = 10* numpy.log10((data[0,:,:] * numpy.conjugate(data[0,:,:])).real)
+                    # m = numpy.nanmean(dat)
+                    # o = numpy.nanstd(dat)
+                    # #print(m, o, x.shape, y.shape)
+                    # c = ax[0].pcolormesh(x, y, dat.T, cmap ='YlGnBu', vmin = (m-2*o), vmax = (m+2*o))
+                    # ax[0].vlines(outs_lines,200,600, linestyles='dashed', label = 'outs', color='w')
+                    # fig.colorbar(c)
+                    # ax[0].vlines(outlier_loc_index,650,750, linestyles='dashed', label = 'outs', color='r')
+                    # ax[1].hist(outs_lines,bins=my_bins)
+                    # plt.show()
+                    self.outliers_IDs_list = outlier_loc_index
+                    #print("outs list: ", self.outliers_IDs_list)
+                    return data
+                def fillBuffer(self, data, datatime):
+                    if self.__profIndex == 0:
+                        self.__buffer_data = data.copy()
-                        for n in range(len(self.min_sats)):
-                            c = self.compareRanges(data,self.min_sats[n],self.max_sats[n])
-                            clean = clean and c
                     else:
+                        self.__buffer_data = numpy.concatenate((self.__buffer_data,data), axis=1)#en perfiles
+                    self.__profIndex += 1
+                    self.__buffer_times.append(datatime)
-                        clean =  (self.compareRanges(data, self.min_sats,self.max_sats))
+                def getData(self, data, datatime=None):
-                    return clean
+                    if self.__profIndex  == 0:
+                        self.__initime = datatime
+                    self.__dataReady = False
-                def run(self, dataOut, minHei=None, maxHei=None, minRef=None, maxRef=None, th=5, thdB=65, rangeHeiList=None):
+                    self.fillBuffer(data, datatime)
-                    dataOut.flagNoData = True
+                    dataBlock = None
+                    if self.__profIndex == self.n:
+                        #print("apnd : ",data)
+                        dataBlock = self.filterSatsProfiles()
+                        self.__dataReady = True
+                    return dataBlock
+                    if dataBlock is None:
+                        return None, None
+                    return dataBlock
+                def releaseBlock(self):
+                    if self.n % self.lenProfileOut != 0:
+                        raise ValueError("lenProfileOut %d must be submultiple of nProfiles %d" %(self.lenProfileOut, self.n))
+                        return None
+                    data = self.buffer[:,self.init_prof:self.end_prof:,:]  #ch, prof, alt
+                    self.init_prof = self.end_prof
+                    self.end_prof += self.lenProfileOut
+                    #print("data release shape: ",dataOut.data.shape, self.end_prof)
+                    self.n_prof_released += 1
+                    return data
+                def run(self, dataOut, n=None, navg=0.75, nProfilesOut=1, profile_margin=50,
+                            th_hist_outlier=3,minHei=None, maxHei=None, minRef=None, maxRef=None, thdB=10):
                     if not self.isConfig:
-                        self.setup(dataOut, minHei, maxHei, minRef, maxRef, th, thdB, rangeHeiList)
+                        #print("init p idx: ", dataOut.profileIndex )
+                        self.setup(dataOut,n=n, navg=navg,profileMargin=profile_margin,thHistOutlier=th_hist_outlier,
+                                    minHei=minHei, maxHei=maxHei, minRef=minRef, maxRef=maxRef, thdB=thdB)
                         self.isConfig = True
-                    #print(self.min_sats,self.max_sats)
-                    if dataOut.flagDataAsBlock:
-                        raise ValueError("ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False")
+                    dataBlock = None
+                    if not dataOut.buffer_empty: #hay datos acumulados
+                        if self.init_prof == 0:
+                            self.n_prof_released = 0
+                            self.lenProfileOut = nProfilesOut
+                            dataOut.flagNoData = False
+                            #print("tp 2 ",dataOut.data.shape)
+                            self.init_prof = 0
+                            self.end_prof = self.lenProfileOut
+                            dataOut.nProfiles = self.lenProfileOut
+                        if nProfilesOut == 1:
+                            dataOut.flagDataAsBlock = False
                         else:
-                        self.noise =10*numpy.log10(dataOut.getNoisebyHildebrand(ymin_index=self.min_ref, ymax_index=self.max_ref))
+                            dataOut.flagDataAsBlock = True
-                        if not self.isProfileClean(dataOut.data):
+                        #print("prof: ",self.init_prof)
+                        dataOut.flagNoData = False
+                        if numpy.isin(self.n_prof_released, self.outliers_IDs_list):
+                            #print("omitting: ", self.n_prof_released)
+                            dataOut.flagNoData = True
+                        dataOut.ippSeconds = self._ipp
+                        dataOut.utctime = self.first_utcBlock + self.init_prof*self._ipp
+                        # print("time: ", dataOut.utctime, self.first_utcBlock, self.init_prof,self._ipp,dataOut.ippSeconds)
+                        #dataOut.data = self.releaseBlock()
+                        #########################################################3
+                        if self.n % self.lenProfileOut != 0:
+                            raise ValueError("lenProfileOut %d must be submultiple of nProfiles %d" %(self.lenProfileOut, self.n))
+                            return None
+                        dataOut.data = self.buffer[:,self.init_prof:self.end_prof:,:]  #ch, prof, alt
+                        self.init_prof = self.end_prof
+                        self.end_prof += self.lenProfileOut
+                        #print("data release shape: ",dataOut.data.shape, self.end_prof, dataOut.flagNoData)
+                        self.n_prof_released += 1
+                        if self.end_prof >= (self.n +self.lenProfileOut):
+                            self.init_prof = 0
+                            self.__profIndex = 0
+                            self.buffer = None
+                            dataOut.buffer_empty = True
+                            self.outliers_IDs_list = []
+                            self.n_prof_released = 0
+                            dataOut.flagNoData = False #enviar ultimo aunque sea outlier :(
+                            #print("cleaning...", dataOut.buffer_empty)
+                        dataOut.profileIndex    = 0 #self.lenProfileOut
+                        ####################################################################
                         return dataOut
-                            #dataOut.data = numpy.full((dataOut.nChannels,dataOut.nHeights),numpy.NAN)
-                    #self.count += 1
-                    dataOut.flagNoData = False
+                    #print("tp 223 ",dataOut.data.shape)
+                    dataOut.flagNoData = True
+                    try:
+                        #dataBlock = self.getData(dataOut.data.reshape(self.nChannels,1,self.nHeights), dataOut.utctime)
+                        dataBlock = self.getData(numpy.reshape(dataOut.data,(self.nChannels,1,self.nHeights)), dataOut.utctime)
+                        self.__count_exec +=1
+                    except Exception as e:
+                        print("Error getting profiles data",self.__count_exec )
+                        print(e)
+                        sys.exit()
+                    if self.__dataReady:
+                        #print("omitting: ", len(self.outliers_IDs_list))
+                        self.__count_exec = 0
+                        #dataOut.data =
+                        #self.buffer = numpy.flip(dataBlock, axis=1)
+                        self.buffer = dataBlock
+                        self.first_utcBlock = self.__initime
+                        dataOut.utctime = self.__initime
+                        dataOut.nProfiles  = self.__profIndex
+                        #dataOut.flagNoData = False
+                        self.init_prof = 0
+                        self.__profIndex = 0
+                        self.__initime = None
+                        dataBlock = None
+                        self.__buffer_times = []
+                        dataOut.error = False
+                        dataOut.useInputBuffer = True
+                        dataOut.buffer_empty = False
+                        #print("1 ch: {}  prof: {}  hs: {}".format(int(dataOut.nChannels),int(dataOut.nProfiles),int(dataOut.nHeights)))
+                    #print(self.__count_exec)
                     return 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