schain-jc Commit - r608:3cc55c179318

data...

Julio Valdez -

r608:3cc55c179318

parent child

Context file:

r608:3cc55c179318

Collapse all files

schainpy/scripts/Meteor_JASMET30MHz.py created 644 +188 0

@@ -0,0 +1,188
	1	# DIAS 19 Y 20 FEB 2014
	2	# Comprobacion de Resultados DBS con SA
	3
	4	import os, sys
	5
	6	path = os.path.split(os.getcwd())[0]
	7	path = os.path.split(path)[0]
	8
	9	sys.path.insert(0, path)
	10
	11	from schainpy.controller import Project
	12
	13	desc = "JASMET Experiment Test"
	14	filename = "JASMETtest.xml"
	15
	16	controllerObj = Project()
	17
	18	controllerObj.setup(id = '191', name='test01', description=desc)
	19
	20	#Experimentos
	21
	22	#2014051 20 Feb 2014
	23	# path = '/home/soporte/Data/JASMET/JASMET_30/2014106'
	24	# pathFigure = '/home/soporte/workspace/Graficos/JASMET/prueba1'
	25
	26	remotefolder = "/home/wmaster/graficos"
	27	path = '/media/joscanoa/84A65E64A65E5730/soporte/Data/JASMET/JASMET_30/'
	28
	29	pathfile1 = os.path.join(os.environ['HOME'],'Pictures/JASMET30/meteor')
	30	pathfile2 = os.path.join(os.environ['HOME'],'Pictures/JASMET30/wind')
	31	pathfile3 = os.path.join(os.environ['HOME'],'Pictures/JASMET30/phase')
	32
	33	pathfig = os.path.join(os.environ['HOME'],'Pictures/JASMET30/graph')
	34
	35	startTime = '00:00:00'
	36	endTime = '23:59:59'
	37	xmin ='17.0'
	38	xmax = '34.0'
	39
	40	#------------------------------------------------------------------------------------------------
	41	readUnitConfObj = controllerObj.addReadUnit(datatype='VoltageReader',
	42	path=path,
	43	startDate='2014/04/15',
	44	endDate='2014/04/16',
	45	startTime=startTime,
	46	endTime=endTime,
	47	online=0,
	48	delay=5,
	49	walk=1)
	50
	51	opObj11 = readUnitConfObj.addOperation(name='printNumberOfBlock')
	52
	53
	54	#--------------------------------------------------------------------------------------------------
	55
	56	procUnitConfObj0 = controllerObj.addProcUnit(datatype='VoltageProc', inputId=readUnitConfObj.getId())
	57
	58	opObj00 = procUnitConfObj0.addOperation(name='selectChannels')
	59	opObj00.addParameter(name='channelList', value='0, 1, 2, 3, 4', format='intlist')
	60
	61	opObj01 = procUnitConfObj0.addOperation(name='setRadarFrequency')
	62	# opObj01.addParameter(name='frequency', value='30.e6', format='float')
	63	opObj01.addParameter(name='frequency', value='50.e6', format='float')
	64
	65	#opObj11 = procUnitConfObj0.addOperation(name='Decoder', optype='other')
	66	#--------------------------------------------------------------------------------------------------
	67
	68	procUnitConfObj1 = controllerObj.addProcUnit(datatype='ParametersProc', inputId=procUnitConfObj0.getId())
	69	procUnitConfObj1.addParameter(name='nSeconds', value='100', format='int')
	70
	71	opObj10 = procUnitConfObj1.addOperation(name='MeteorDetection')
	72	# opObj10.addParameter(name='predefinedPhaseShifts', value='-89.5, 41.5, 0.0, -138.0, -85.5', format='floatlist')
	73	opObj10.addParameter(name='predefinedPhaseShifts', value='0, 0, 0, 0, 0', format='floatlist')
	74	opObj10.addParameter(name='cohDetection', value='0', format='bool')
	75	opObj10.addParameter(name='noise_multiple', value='4', format='int')
	76	opObj10.addParameter(name='SNRThresh', value='5', format='float')
	77	opObj10.addParameter(name='phaseThresh', value='20', format='float')
	78	opObj10.addParameter(name='azimuth', value='45', format='float')
	79	opObj10.addParameter(name='hmin', value='68', format='float')
	80	opObj10.addParameter(name='hmax', value='112', format='float')
	81	opObj10.addParameter(name='saveAll', value='1', format='bool')
	82
	83	opObj12 = procUnitConfObj1.addOperation(name='HDF5Writer', optype='other')
	84	opObj12.addParameter(name='path', value=pathfile1)
	85	opObj12.addParameter(name='blocksPerFile', value='1000', format='int')
	86	opObj12.addParameter(name='metadataList',value='type,inputUnit,heightList,paramInterval',format='list')
	87	opObj12.addParameter(name='dataList',value='data_param',format='list')
	88	opObj12.addParameter(name='mode',value='0',format='int')
	89	#Tiene que ser de 3 dimensiones, append en lugar de aumentar una dimension
	90
	91	opObj13 = procUnitConfObj1.addOperation(name='SkyMapPlot', optype='other')
	92	opObj13.addParameter(name='id', value='1', format='int')
	93	opObj13.addParameter(name='wintitle', value='Sky Map', format='str')
	94	opObj13.addParameter(name='save', value='1', format='bool')
	95	opObj13.addParameter(name='figpath', value=pathfig, format='str')
	96	opObj13.addParameter(name='ftp', value='1', format='int')
	97	opObj13.addParameter(name='exp_code', value='15', format='int')
	98	opObj13.addParameter(name='sub_exp_code', value='1', format='int')
	99
	100
	101	#--------------------------------------------------------------------------------------------------
	102	procUnitConfObj2 = controllerObj.addProcUnit(datatype='ParametersProc', inputId=procUnitConfObj1.getId())
	103
	104	opObj22 = procUnitConfObj2.addOperation(name='WindProfiler', optype='other')
	105	opObj22.addParameter(name='technique', value='Meteors', format='str')
	106	opObj22.addParameter(name='nHours', value='0.25', format='float')
	107	opObj22.addParameter(name='SNRThresh', value='12.0', format='float')
	108
	109	opObj23 = procUnitConfObj2.addOperation(name='WindProfilerPlot', optype='other')
	110	opObj23.addParameter(name='id', value='2', format='int')
	111	opObj23.addParameter(name='wintitle', value='Wind Profiler', format='str')
	112	opObj23.addParameter(name='save', value='1', format='bool')
	113	opObj23.addParameter(name='figpath', value = pathfig, format='str')
	114	opObj23.addParameter(name='zmin', value='-120', format='int')
	115	opObj23.addParameter(name='zmax', value='120', format='int')
	116	# opObj12.addParameter(name='zmin_ver', value='-0.8', format='float')
	117	# opObj12.addParameter(name='zmax_ver', value='0.8', format='float')
	118	# opObj23.addParameter(name='SNRmin', value='-10', format='int')
	119	# opObj23.addParameter(name='SNRmax', value='60', format='int')
	120	# opObj23.addParameter(name='SNRthresh', value='0', format='float')
	121	opObj23.addParameter(name='xmin', value=xmin, format='float')
	122	opObj23.addParameter(name='xmax', value=xmax, format='float')
	123	opObj23.addParameter(name='ftp', value='1', format='int')
	124	opObj23.addParameter(name='exp_code', value='15', format='int')
	125	opObj23.addParameter(name='sub_exp_code', value='1', format='int')
	126
	127	opObj24 = procUnitConfObj2.addOperation(name='HDF5Writer', optype='other')
	128	opObj24.addParameter(name='path', value=pathfile2)
	129	opObj24.addParameter(name='blocksPerFile', value='1000', format='int')
	130	opObj24.addParameter(name='metadataList',value='type,inputUnit,outputInterval',format='list')
	131	opObj24.addParameter(name='dataList',value='data_output,utctime',format='list')
	132
	133	#--------------------------------------------------------------------------------------------------
	134	procUnitConfObj3 = controllerObj.addProcUnit(datatype='ParametersProc', inputId=procUnitConfObj1.getId())
	135
	136
	137	opObj31 = procUnitConfObj3.addOperation(name='PhaseCalibration', optype='other')
	138	opObj31.addParameter(name='nHours', value='0.25', format='float')
	139	opObj31.addParameter(name='distances', value='-15, -15, 12, 12', format='intlist')
	140	# opObj31.addParameter(name='distances', value='-25, -25, 20, 20', format='intlist')
	141	opObj31.addParameter(name='pairs', value='(0,3),(1,2)', format='pairslist')
	142	opObj31.addParameter(name='hmin', value='68', format='float')
	143	opObj31.addParameter(name='hmax', value='112', format='float')
	144
	145	opObj32 = procUnitConfObj3.addOperation(name='PhasePlot', optype='other')
	146	opObj32.addParameter(name='id', value='201', format='int')
	147	opObj32.addParameter(name='wintitle', value='PhaseCalibration', format='str')
	148	# opObj32.addParameter(name='timerange', value='300', format='int')
	149	opObj32.addParameter(name='xmin', value=xmin, format='float')
	150	opObj32.addParameter(name='xmax', value=xmax, format='float')
	151	# opObj32.addParameter(name='xmin', value='19', format='float')
	152	# opObj32.addParameter(name='xmax', value='20', format='float')
	153	opObj32.addParameter(name='ymin', value='-180', format='float')
	154	opObj32.addParameter(name='ymax', value='180', format='float')
	155	opObj32.addParameter(name='figpath', value=pathfig, format='str')
	156	opObj32.addParameter(name='ftp', value='1', format='int')
	157	opObj32.addParameter(name='exp_code', value='15', format='int')
	158	opObj32.addParameter(name='sub_exp_code', value='1', format='int')
	159
	160	opObj33 = procUnitConfObj3.addOperation(name='HDF5Writer', optype='other')
	161	opObj33.addParameter(name='path', value=pathfile3)
	162	opObj33.addParameter(name='blocksPerFile', value='1000', format='int')
	163	opObj33.addParameter(name='metadataList',value='type,inputUnit,outputInterval',format='list')
	164	opObj33.addParameter(name='dataList',value='data_output,utctime',format='list')
	165	# opObj25.addParameter(name='mode',value='1,0,0',format='intlist')
	166
	167	#--------------------------------------------------------------------------------------------------
	168
	169	procUnitConfObj4 = controllerObj.addProcUnit(name='SendToServer')
	170	procUnitConfObj4.addParameter(name='server', value='jro-app.igp.gob.pe', format='str')
	171	procUnitConfObj4.addParameter(name='username', value='wmaster', format='str')
	172	procUnitConfObj4.addParameter(name='password', value='mst2010vhf', format='str')
	173	procUnitConfObj4.addParameter(name='localfolder', value=pathfig, format='str')
	174	procUnitConfObj4.addParameter(name='remotefolder', value=remotefolder, format='str')
	175	procUnitConfObj4.addParameter(name='ext', value='.png', format='str')
	176	procUnitConfObj4.addParameter(name='period', value=20, format='int')
	177	procUnitConfObj4.addParameter(name='protocol', value='ftp', format='str')
	178
	179	#--------------------------------------------------------------------------------------------------
	180
	181	print "Escribiendo el archivo XML"
	182	controllerObj.writeXml(filename)
	183	print "Leyendo el archivo XML"
	184	controllerObj.readXml(filename)
	185
	186	controllerObj.createObjects()
	187	controllerObj.connectObjects()
	188	controllerObj.run() No newline at end of file

schainpy/scripts/Meteor_JASMET_PreProc.py created 644 +56 0

@@ -0,0 +1,56
	1	"""
	2	Se debe verficar que el disco de datos se encuentra montado en el sistema
	3	"""
	4	import os, sys
	5
	6	path = os.path.split(os.getcwd())[0]
	7	path = os.path.split(path)[0]
	8
	9	sys.path.insert(0, path)
	10
	11	from schainpy.controller import Project
	12
	13	desc = "Meteor Experiment Test"
	14	filename = "meteor20130812.xml"
	15
	16	controllerObj = Project()
	17	controllerObj.setup(id = '191', name='meteor_test01', description=desc)
	18
	19	path = '/home/dsuarez/.gvfs/data on 10.10.20.13/Jasmet50'
	20
	21
	22	readUnitConfObj = controllerObj.addReadUnit(datatype='Voltage',
	23	path=path,
	24	startDate='2014/04/16',
	25	endDate='2014/04/16',
	26	startTime='00:00:00',
	27	endTime='23:59:59',
	28	online=0,
	29	walk=1)
	30
	31	opObj11 = readUnitConfObj.addOperation(name='printNumberOfBlock')
	32
	33	procUnitConfObj0 = controllerObj.addProcUnit(datatype='Voltage', inputId=readUnitConfObj.getId())
	34
	35
	36	opObj11 = procUnitConfObj0.addOperation(name='Decoder', optype='other')
	37
	38
	39	opObj11 = procUnitConfObj0.addOperation(name='CohInt', optype='other')
	40	opObj11.addParameter(name='n', value='2', format='int')
	41
	42	opObj11 = procUnitConfObj0.addOperation(name='VoltageWriter', optype='other')
	43	opObj11.addParameter(name='path', value='/home/jasmet/jasmet30_abril')
	44	opObj11.addParameter(name='blocksPerFile', value='100', format='int')
	45	opObj11.addParameter(name='profilesPerBlock', value='200', format='int')
	46
	47
	48
	49	print "Escribiendo el archivo XML"
	50	controllerObj.writeXml(filename)
	51	print "Leyendo el archivo XML"
	52	controllerObj.readXml(filename)
	53
	54	controllerObj.createObjects()
	55	controllerObj.connectObjects()
	56	controllerObj.run()

schainpy/model/data/jrodata.py +14 -5

             '''
             $Author: murco $
             $Id: JROData.py 173 2012-11-20 15:06:21Z murco $
             '''
             import copy
             import numpy
             import datetime
             from jroheaderIO import SystemHeader, RadarControllerHeader
             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:
                     -1        :    any error
                     anoise    :    noise's level
                 """
                 sortdata = numpy.sort(data,axis=None)
                 lenOfData = len(sortdata)
                 nums_min = lenOfData/10
                 if (lenOfData/10) > 2:
                     nums_min = lenOfData/10
                 else:
                     nums_min = 2
                 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
                 stdv = numpy.sqrt((sumq - lnoise**2)/(j - 1))
                 return lnoise
             class Beam:
                 def __init__(self):
                     self.codeList = []
                     self.azimuthList = []
                     self.zenithList = []
             class GenericData(object):
                 flagNoData = True
                 def __init__(self):
                     raise ValueError, "This class has not been implemented"
                 def copy(self, inputObj=None):
                     if inputObj == None:
                         return copy.deepcopy(self)
                     for key in inputObj.__dict__.keys():
                         self.__dict__[key] = inputObj.__dict__[key]
                 def deepcopy(self):
                     return copy.deepcopy(self)
                 def isEmpty(self):
                     return self.flagNoData
             class JROData(GenericData):
             #    m_BasicHeader = BasicHeader()
             #    m_ProcessingHeader = ProcessingHeader()
                 systemHeaderObj = SystemHeader()
                 radarControllerHeaderObj = RadarControllerHeader()
             #    data = None
                 type = None
                 datatype = None #dtype but in string
             #     dtype = None
             #    nChannels = None
             #    nHeights = None
                 nProfiles = None
                 heightList = None
                 channelList = None
                 flagDiscontinuousBlock = False
                 useLocalTime = False
                 utctime = None
                 timeZone = None
                 dstFlag = None
                 errorCount = None
                 blocksize = None
             #     nCode = None
             #
             #     nBaud = None
             #
             #     code = None
                 flagDecodeData = False #asumo q la data no esta decodificada
                 flagDeflipData = False #asumo q la data no esta sin flip
                 flagShiftFFT = False
             #     ippSeconds = None
             #     timeInterval = None
                 nCohInt = None
             #     noise = None
                 windowOfFilter = 1
                 #Speed of ligth
                 C = 3e8
                 frequency = 49.92e6
                 realtime = False
                 beacon_heiIndexList = None
                 last_block = None
                 blocknow = None
                 azimuth = None
                 zenith = None
                 beam = Beam()
                 profileIndex = None
                 def __init__(self):
                     raise ValueError, "This class has not been implemented"
                 def getNoise(self):
                     raise ValueError, "Not implemented"
                 def getNChannels(self):
                     return len(self.channelList)
                 def getChannelIndexList(self):
                     return range(self.nChannels)
                 def getNHeights(self):
                     return len(self.heightList)
                 def getHeiRange(self, extrapoints=0):
                     heis = self.heightList
             #        deltah = self.heightList[1] - self.heightList[0]
             #
             #        heis.append(self.heightList[-1])
                     return heis
                 def getltctime(self):
                     if self.useLocalTime:
                         return self.utctime - self.timeZone*60
                     return self.utctime
                 def getDatatime(self):
                     datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
                     return datatimeValue
                 def getTimeRange(self):
                     datatime = []
                     datatime.append(self.ltctime)
                     datatime.append(self.ltctime + self.timeInterval+60)
                     datatime = numpy.array(datatime)
                     return datatime
                 def getFmax(self):
                     PRF = 1./(self.ippSeconds * self.nCohInt)
                     fmax = PRF/2.
                     return fmax
                 def getVmax(self):
                     _lambda = self.C/self.frequency
                     vmax = self.getFmax() * _lambda
                     return vmax
                 def get_ippSeconds(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.ippSeconds
                 def set_ippSeconds(self, ippSeconds):
                     '''
                     '''
                     self.radarControllerHeaderObj.ippSeconds = ippSeconds
                     return
                 def get_dtype(self):
                     '''
                     '''
                     return getNumpyDtype(self.datatype)
                 def set_dtype(self, numpyDtype):
                     '''
                     '''
                     self.datatype = getDataTypeCode(numpyDtype)
                 def get_code(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.code
                 def set_code(self, code):
                     '''
                     '''
                     self.radarControllerHeaderObj.code = code
                     return
                 def get_ncode(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.nCode
                 def set_ncode(self, nCode):
                     '''
                     '''
                     self.radarControllerHeaderObj.nCode = nCode
                     return
                 def get_nbaud(self):
                     '''
                     '''
                     return self.radarControllerHeaderObj.nBaud
                 def set_nbaud(self, nBaud):
                     '''
                     '''
                     self.radarControllerHeaderObj.nBaud = nBaud
                     return
             #     def getTimeInterval(self):
             #
             #         raise IOError, "This method should be implemented inside each Class"
                 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
                 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
                 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
                 #noise = property(getNoise, "I'm the 'nHeights' property.")
                 datatime = property(getDatatime, "I'm the 'datatime' property")
                 ltctime = property(getltctime, "I'm the 'ltctime' property")
                 ippSeconds = property(get_ippSeconds, set_ippSeconds)
                 dtype = property(get_dtype, set_dtype)
             #     timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
                 code = property(get_code, set_code)
                 nCode = property(get_ncode, set_ncode)
                 nBaud = property(get_nbaud, set_nbaud)
             class Voltage(JROData):
                 #data es un numpy array de 2 dmensiones (canales, alturas)
                 data = None
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.useLocalTime = True
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Voltage"
                     self.data = None
             #         self.dtype = None
             #        self.nChannels = 0
             #        self.nHeights = 0
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
             #        self.channelIndexList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
                     self.utctime = None
                     self.timeZone = None
                     self.dstFlag = None
                     self.errorCount = None
                     self.nCohInt = 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.flagShiftFFT = False
                     self.flagDataAsBlock = False    #Asumo que la data es leida perfil a perfil
                     self.profileIndex = 0
                 def getNoisebyHildebrand(self, channel = None):
                     """
                     Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
                     Return:
                         noiselevel
                     """
                     if channel != None:
                         data = self.data[channel]
                         nChannels = 1
                     else:
                         data = self.data
                         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):
                     if type == 1:
                         noise = self.getNoisebyHildebrand(channel)
                     return 10*numpy.log10(noise)
                 def getPower(self, channel = None):
                     if channel != None:
                         data = self.data[channel]
                     else:
                         data = self.data
                     power = data * numpy.conjugate(data)
                     return 10*numpy.log10(power.real)
                 def getTimeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt
                     return timeInterval
                 noise = property(getNoise, "I'm the 'nHeights' property.")
                 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
             class Spectra(JROData):
                 #data es un numpy array de 2 dmensiones (canales, perfiles, alturas)
                 data_spc = None
                 #data es un numpy array de 2 dmensiones (canales, pares, alturas)
                 data_cspc = None
                 #data es un numpy array de 2 dmensiones (canales, alturas)
                 data_dc = None
                 nFFTPoints = None
             #     nPairs = None
                 pairsList = None
                 nIncohInt = None
                 wavelength = None #Necesario para cacular el rango de velocidad desde la frecuencia
                 nCohInt = None #se requiere para determinar el valor de timeInterval
                 ippFactor = None
                 profileIndex = 0
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.useLocalTime = True
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Spectra"
             #        self.data = None
             #         self.dtype = None
             #        self.nChannels = 0
             #        self.nHeights = 0
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
             #        self.channelIndexList = 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.noise = None
                     self.beacon_heiIndexList = []
                     self.noise_estimation = None
                 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
                     """
                     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)
                     return noise
                 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
                     if self.noise_estimation != None:
                         return self.noise_estimation #this was estimated by getNoise Operation defined in jroproc_spectra.py
                     else:
                         noise = self.getNoisebyHildebrand(xmin_index, xmax_index, ymin_index, ymax_index)
                         return noise
                 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.) - deltav/2
                     return velrange
                 def getNPairs(self):
                     return len(self.pairsList)
                 def getPairsIndexList(self):
                     return range(self.nPairs)
                 def getNormFactor(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
                 def getFlagCspc(self):
                     if self.data_cspc == None:
                         return True
                     return False
                 def getFlagDc(self):
                     if self.data_dc == None:
                         return True
                     return False
                 def getTimeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles
                     return timeInterval
                 nPairs = property(getNPairs, "I'm the 'nPairs' property.")
                 pairsIndexList = property(getPairsIndexList, "I'm the 'pairsIndexList' property.")
                 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
                 flag_cspc = property(getFlagCspc)
                 flag_dc = property(getFlagDc)
                 noise = property(getNoise, "I'm the 'nHeights' property.")
                 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
             class SpectraHeis(Spectra):
                 data_spc = None
                 data_cspc = None
                 data_dc = None
                 nFFTPoints = None
             #     nPairs = None
                 pairsList = None
                 nCohInt = None
                 nIncohInt = None
                 def __init__(self):
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "SpectraHeis"
             #         self.dtype = None
             #        self.nChannels = 0
             #        self.nHeights = 0
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
             #        self.channelIndexList = None
                     self.flagNoData = True
                     self.flagDiscontinuousBlock = False
             #         self.nPairs = 0
                     self.utctime = None
                     self.blocksize = None
                     self.profileIndex = 0
                     self.nCohInt = 1
                     self.nIncohInt = 1
                 def getNormFactor(self):
                     pwcode = 1
                     if self.flagDecodeData:
                         pwcode = numpy.sum(self.code[0]**2)
                     normFactor = self.nIncohInt*self.nCohInt*pwcode
                     return normFactor
                 def getTimeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
                     return timeInterval
                 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
                 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
             class Fits(JROData):
                 heightList = None
                 channelList = None
                 flagNoData = True
                 flagDiscontinuousBlock = False
                 useLocalTime = False
                 utctime = None
                 timeZone = None
             #     ippSeconds = None
             #     timeInterval = None
                 nCohInt = None
                 nIncohInt = None
                 noise = None
                 windowOfFilter = 1
                 #Speed of ligth
                 C = 3e8
                 frequency = 49.92e6
                 realtime = False
                 def __init__(self):
                     self.type = "Fits"
                     self.nProfiles = None
                     self.heightList = None
                     self.channelList = None
             #         self.channelIndexList = None
                     self.flagNoData = True
                     self.utctime = None
                     self.nCohInt = 1
                     self.nIncohInt = 1
                     self.useLocalTime = True
                     self.profileIndex = 0
             #         self.utctime = None
             #         self.timeZone = None
             #         self.ltctime = None
             #         self.timeInterval = None
             #         self.header = None
             #         self.data_header = None
             #         self.data = None
             #         self.datatime = None
             #         self.flagNoData = False
             #         self.expName = ''
             #         self.nChannels = None
             #         self.nSamples = None
             #         self.dataBlocksPerFile = None
             #         self.comments = ''
             #
                 def getltctime(self):
                     if self.useLocalTime:
                         return self.utctime - self.timeZone*60
                     return self.utctime
                 def getDatatime(self):
                     datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
                     return datatime
                 def getTimeRange(self):
                     datatime = []
                     datatime.append(self.ltctime)
                     datatime.append(self.ltctime + self.timeInterval)
                     datatime = numpy.array(datatime)
                     return datatime
                 def getHeiRange(self):
                     heis = self.heightList
                     return heis
                 def isEmpty(self):
                     return self.flagNoData
                 def getNHeights(self):
                     return len(self.heightList)
                 def getNChannels(self):
                     return len(self.channelList)
                 def getChannelIndexList(self):
                     return range(self.nChannels)
                 def getNoise(self, type = 1):
                     #noise = numpy.zeros(self.nChannels)
                     if type == 1:
                         noise = self.getNoisebyHildebrand()
                     if type == 2:
                         noise = self.getNoisebySort()
                     if type == 3:
                         noise = self.getNoisebyWindow()
                     return noise
                 def getTimeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
                     return timeInterval
                 datatime = property(getDatatime, "I'm the 'datatime' property")
                 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
                 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
                 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
                 noise = property(getNoise, "I'm the 'nHeights' property.")
                 datatime = property(getDatatime, "I'm the 'datatime' property")
                 ltctime = property(getltctime, "I'm the 'ltctime' property")
                 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
             class Correlation(JROData):
                 noise = None
                 SNR = None
                 pairsAutoCorr = None    #Pairs of Autocorrelation
                 #--------------------------------------------------
                 data_corr = None
                 data_volt = None
                 lagT = None # each element value is a profileIndex
                 lagR = None # each element value is in km
                 pairsList = None
                 calculateVelocity = None
                 nPoints = None
                 nAvg = None
                 bufferSize = None
                 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 = None
                     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 getLagTRange(self, extrapoints=0):
                     lagTRange = self.lagT
                     diff = lagTRange[1] - lagTRange[0]
                     extra = numpy.arange(1,extrapoints + 1)*diff + lagTRange[-1]
                     lagTRange = numpy.hstack((lagTRange, extra))
                     return lagTRange
                 def getLagRRange(self, extrapoints=0):
                     return self.lagR
                 def getPairsList(self):
                     return self.pairsList
                 def getCalculateVelocity(self):
                     return self.calculateVelocity
                 def getNPoints(self):
                     return self.nPoints
                 def getNAvg(self):
                     return self.nAvg
                 def getBufferSize(self):
                     return self.bufferSize
                 def getPairsAutoCorr(self):
                     pairsList = self.pairsList
                     pairsAutoCorr = numpy.zeros(self.nChannels, dtype = 'int')*numpy.nan
                     for l in range(len(pairsList)):
                         firstChannel = pairsList[l][0]
                         secondChannel = pairsList[l][1]
                         #Obteniendo pares de Autocorrelacion
                         if firstChannel == secondChannel:
                             pairsAutoCorr[firstChannel] = int(l)
                     pairsAutoCorr = pairsAutoCorr.astype(int)
                     return pairsAutoCorr
                 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[range(0,1)] = ind_vel[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(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
                 def getTimeInterval(self):
                     timeInterval = self.ippSeconds * self.nCohInt * self.nPoints
                     return timeInterval
                 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
             #     pairsList = property(getPairsList, "I'm the 'pairsList' property.")
             #     nPoints = property(getNPoints, "I'm the 'nPoints' property.")
                 calculateVelocity = property(getCalculateVelocity, "I'm the 'calculateVelocity' property.")
                 nAvg = property(getNAvg, "I'm the 'nAvg' property.")
                 bufferSize = property(getBufferSize, "I'm the 'bufferSize' property.")
             class Parameters(JROData):
                 #Information from previous data
                 inputUnit = None        #Type of data to be processed
                 operation = None        #Type of operation to parametrize
                 normFactor = None       #Normalization Factor
                 groupList = None        #List of Pairs, Groups, etc
                 #Parameters
                 data_param = None       #Parameters obtained
                 data_pre = None         #Data Pre Parametrization
                 data_SNR = None         #Signal to Noise Ratio
-                heightRange = None      #Heights
+            #     heightRange = None      #Heights
-                abscissaRange = None    #Abscissa, can be velocities, lags or time
+                abscissaList = None    #Abscissa, can be velocities, lags or time
                 noise = None            #Noise Potency
-            #     initUtcTime = None      #Initial UTC time
+                utctimeInit = None      #Initial UTC time
                 paramInterval = None    #Time interval to calculate Parameters in seconds
                 #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
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     self.radarControllerHeaderObj = RadarControllerHeader()
                     self.systemHeaderObj = SystemHeader()
                     self.type = "Parameters"
                 def getTimeRange1(self):
                     datatime = []
-                    datatime.append(self.ltctime)
+                    if self.useLocalTime:
-                    datatime.append(self.ltctime + self.outputInterval - 1)
+                        time1 = self.utctimeInit - self.timeZone*60
+                    else:
+                        time1 = utctimeInit
+            #         datatime.append(self.utctimeInit)
+            #         datatime.append(self.utctimeInit + self.outputInterval)
+                    datatime.append(time1)
+                    datatime.append(time1 + self.outputInterval)
                     datatime = numpy.array(datatime)
                     return datatime

schainpy/model/graphics/figure.py +2 -2

             import os
             import numpy
             import time, datetime
             import mpldriver
             from schainpy.model.proc.jroproc_base import Operation
             def isRealtime(utcdatatime):
                 utcnow = time.mktime(time.localtime())
                 delta = abs(utcnow - utcdatatime) # abs
                 if delta >= 30.:
                     return False
                 return True
             class Figure(Operation):
                 __driver = mpldriver
                 __isConfigThread = False
                 fig = None
                 id = None
                 wintitle = None
                 width = None
                 height = None
                 nplots = None
                 timerange = None
                 axesObjList = []
                 WIDTH = None
                 HEIGHT = None
                 PREFIX = 'fig'
                 xmin = None
                 xmax = None
                 counter_imagwr = 0
                 figfile = None
                 def __init__(self):
                     raise ValueError, "This method is not implemented"
                 def __del__(self):
                     self.__driver.closeFigure()
                 def getFilename(self, name, ext='.png'):
                     path = '%s%03d' %(self.PREFIX, self.id)
                     filename = '%s_%s%s' %(self.PREFIX, name, ext)
                     return os.path.join(path, filename)
                 def getAxesObjList(self):
                     return self.axesObjList
                 def getSubplots(self):
                     raise ValueError, "Abstract method: This method should be defined"
                 def getScreenDim(self, widthplot, heightplot):
                     nrow, ncol = self.getSubplots()
                     widthscreen = widthplot*ncol
                     heightscreen = heightplot*nrow
                     return widthscreen, heightscreen
                 def getTimeLim(self, x, xmin=None, xmax=None, timerange=None):
                     if self.xmin != None and self.xmax != None:
                         if timerange == None:
                             timerange = self.xmax - self.xmin
                         xmin = self.xmin + timerange
                         xmax = self.xmax + timerange
                         return xmin, xmax
                     if timerange == None and (xmin==None or xmax==None):
                         timerange = 14400   #seconds
                         #raise ValueError, "(timerange) or (xmin & xmax) should be defined"
                     if timerange != None:
-                        txmin = x[0] - x[0] % min(timerange/10, 10*60)
+                        txmin = x[0] #- x[0] % min(timerange/10, 10*60)
                     else:
-                        txmin = x[0] - x[0] % 10*60
+                        txmin = x[0] #- x[0] % 10*60
                     thisdatetime = datetime.datetime.utcfromtimestamp(txmin)
                     thisdate = datetime.datetime.combine(thisdatetime.date(), datetime.time(0,0,0))
                     if timerange != None:
                         xmin = (thisdatetime - thisdate).seconds/(60*60.)
                         xmax = xmin + timerange/(60*60.)
                     mindt = thisdate + datetime.timedelta(hours=xmin) - datetime.timedelta(seconds=time.timezone)
                     xmin_sec = time.mktime(mindt.timetuple())
                     maxdt = thisdate + datetime.timedelta(hours=xmax) - datetime.timedelta(seconds=time.timezone)
                     xmax_sec = time.mktime(maxdt.timetuple())
                     return xmin_sec, xmax_sec
                 def init(self, id, nplots, wintitle):
                     raise ValueError, "This method has been replaced with createFigure"
                 def createFigure(self, id, wintitle, widthplot=None, heightplot=None, show=True):
                     """
                     Crea la figura de acuerdo al driver y parametros seleccionados seleccionados.
                     Las dimensiones de la pantalla es calculada a partir de los atributos self.WIDTH
                     y self.HEIGHT y el numero de subplots (nrow, ncol)
                     Input:
                         id    :    Los parametros necesarios son
                         wintitle    :
                     """
                     if widthplot == None:
                         widthplot = self.WIDTH
                     if heightplot == None:
                         heightplot = self.HEIGHT
                     self.id = id
                     self.wintitle = wintitle
                     self.widthscreen, self.heightscreen = self.getScreenDim(widthplot, heightplot)
                     self.fig = self.__driver.createFigure(id=self.id,
                                                           wintitle=self.wintitle,
                                                           width=self.widthscreen,
                                                           height=self.heightscreen,
                                                           show=show)
                     self.axesObjList = []
                     self.counter_imagwr = 0
                 def setDriver(self, driver=mpldriver):
                     self.__driver = driver
                 def setTitle(self, title):
                     self.__driver.setTitle(self.fig, title)
                 def setWinTitle(self, title):
                     self.__driver.setWinTitle(self.fig, title=title)
                 def setTextFromAxes(self, text):
                     raise ValueError, "Este metodo ha sido reemplazaado con el metodo setText de la clase Axes"
                 def makeAxes(self, nrow, ncol, xpos, ypos, colspan, rowspan):
                     raise ValueError, "Este metodo ha sido reemplazaado con el metodo addAxes"
                 def addAxes(self, *args):
                     """
                     Input:
                         *args    :    Los parametros necesarios son
                                       nrow, ncol, xpos, ypos, colspan, rowspan
                     """
                     axesObj = Axes(self.fig, *args)
                     self.axesObjList.append(axesObj)
                 def saveFigure(self, figpath, figfile, *args):
                     filename = os.path.join(figpath, figfile)
                     fullpath = os.path.split(filename)[0]
                     if not os.path.exists(fullpath):
                         subpath = os.path.split(fullpath)[0]
                         if not os.path.exists(subpath):
                             os.mkdir(subpath)
                         os.mkdir(fullpath)
                     self.__driver.saveFigure(self.fig, filename, *args)
                 def save(self, figpath, figfile=None, save=True, ftp=False, wr_period=1, thisDatetime=None, update_figfile=True):
                     self.counter_imagwr += 1
                     if self.counter_imagwr < wr_period:
                         return
                     self.counter_imagwr = 0
                     if save:
                         if not figfile:
                             if not thisDatetime:
                                 raise ValueError, "Saving figure: figfile or thisDatetime should be defined"
                                 return
                             str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
                             figfile = self.getFilename(name = str_datetime)
                         if self.figfile == None:
                             self.figfile = figfile
                         if update_figfile:
                             self.figfile = figfile
                         # store png plot to local folder
                         self.saveFigure(figpath, self.figfile)
                     if not ftp:
                         return
                     if not thisDatetime:
                         return
                     # store png plot to FTP server according to RT-Web format
                     ftp_filename = self.getNameToFtp(thisDatetime, self.FTP_WEI, self.EXP_CODE, self.SUB_EXP_CODE, self.PLOT_CODE, self.PLOT_POS)
             #         ftp_filename = os.path.join(figpath, name)
                     self.saveFigure(figpath, ftp_filename)
                 def getNameToFtp(self, thisDatetime, FTP_WEI, EXP_CODE, SUB_EXP_CODE, PLOT_CODE, PLOT_POS):
                     YEAR_STR = '%4.4d'%thisDatetime.timetuple().tm_year
                     DOY_STR = '%3.3d'%thisDatetime.timetuple().tm_yday
                     FTP_WEI = '%2.2d'%FTP_WEI
                     EXP_CODE = '%3.3d'%EXP_CODE
                     SUB_EXP_CODE = '%2.2d'%SUB_EXP_CODE
                     PLOT_CODE = '%2.2d'%PLOT_CODE
                     PLOT_POS = '%2.2d'%PLOT_POS
                     name = YEAR_STR + DOY_STR + FTP_WEI + EXP_CODE + SUB_EXP_CODE + PLOT_CODE + PLOT_POS
                     return name
                 def draw(self):
                     self.__driver.draw(self.fig)
                 def run(self):
                     raise ValueError, "This method is not implemented"
                 def close(self, show=False):
                     self.__driver.closeFigure(show=show, fig=self.fig)
                 axesList = property(getAxesObjList)
             class Axes:
                 __driver = mpldriver
                 fig = None
                 ax = None
                 plot = None
                 __missing = 1E30
                 __firsttime = None
                 __showprofile = False
                 xmin = None
                 xmax = None
                 ymin = None
                 ymax = None
                 zmin = None
                 zmax = None
                 x_buffer = None
                 z_buffer = None
                 decimationx = None
                 decimationy = None
                 __MAXNUMX = 300
                 __MAXNUMY = 150
                 def __init__(self, *args):
                     """
                     Input:
                         *args    :    Los parametros necesarios son
                                       fig, nrow, ncol, xpos, ypos, colspan, rowspan
                     """
                     ax = self.__driver.createAxes(*args)
                     self.fig = args[0]
                     self.ax = ax
                     self.plot = None
                     self.__firsttime = True
                     self.idlineList = []
                     self.x_buffer = numpy.array([])
                     self.z_buffer = numpy.array([])
                 def setText(self, text):
                     self.__driver.setAxesText(self.ax, text)
                 def setXAxisAsTime(self):
                     pass
                 def pline(self, x, y,
                            xmin=None, xmax=None,
                            ymin=None, ymax=None,
                            xlabel='', ylabel='',
                            title='',
                            **kwargs):
                     """
                     Input:
                         x        :
                         y        :
                         xmin    :
                         xmax    :
                         ymin    :
                         ymax    :
                         xlabel    :
                         ylabel    :
                         title    :
                         **kwargs :    Los parametros aceptados son
                                       ticksize
                                       ytick_visible
                     """
                     if self.__firsttime:
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         self.plot = self.__driver.createPline(self.ax, x, y,
                                                               xmin, xmax,
                                                               ymin, ymax,
                                                               xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title,
                                                               **kwargs)
                         self.idlineList.append(0)
                         self.__firsttime = False
                         return
                     self.__driver.pline(self.plot, x, y, xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title)
                 def addpline(self, x, y, idline, **kwargs):
                     lines = self.ax.lines
                     if idline in self.idlineList:
                         self.__driver.set_linedata(self.ax, x, y, idline)
                     if  idline not in(self.idlineList):
                         self.__driver.addpline(self.ax, x, y, **kwargs)
                         self.idlineList.append(idline)
                     return
                 def pmultiline(self, x, y,
                                xmin=None, xmax=None,
                                ymin=None, ymax=None,
                                xlabel='', ylabel='',
                                title='',
                                **kwargs):
                     if self.__firsttime:
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         self.plot = self.__driver.createPmultiline(self.ax, x, y,
                                                               xmin, xmax,
                                                               ymin, ymax,
                                                               xlabel=xlabel,
                                                               ylabel=ylabel,
                                                               title=title,
                                                               **kwargs)
                         self.__firsttime = False
                         return
                     self.__driver.pmultiline(self.plot, x, y, xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title)
                 def pmultilineyaxis(self, x, y,
                                xmin=None, xmax=None,
                                ymin=None, ymax=None,
                                xlabel='', ylabel='',
                                title='',
                                **kwargs):
                     if self.__firsttime:
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         self.plot = self.__driver.createPmultilineYAxis(self.ax, x, y,
                                                               xmin, xmax,
                                                               ymin, ymax,
                                                               xlabel=xlabel,
                                                               ylabel=ylabel,
                                                               title=title,
                                                               **kwargs)
                         if self.xmin == None: self.xmin = xmin
                         if self.xmax == None: self.xmax = xmax
                         if self.ymin == None: self.ymin = ymin
                         if self.ymax == None: self.ymax = ymax
                         self.__firsttime = False
                         return
                     self.__driver.pmultilineyaxis(self.plot, x, y, xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title)
                 def pcolor(self, x, y, z,
                            xmin=None, xmax=None,
                            ymin=None, ymax=None,
                            zmin=None, zmax=None,
                            xlabel='', ylabel='',
                            title='', rti = False, colormap='jet',
                            **kwargs):
                     """
                     Input:
                         x        :
                         y        :
                         x        :
                         xmin    :
                         xmax    :
                         ymin    :
                         ymax    :
                         zmin    :
                         zmax    :
                         xlabel    :
                         ylabel    :
                         title    :
                         **kwargs :    Los parametros aceptados son
                                       ticksize=9,
                                       cblabel=''
                                       rti = True or False
                     """
                     if self.__firsttime:
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmin == None: zmin = numpy.nanmin(z)
                         if zmax == None: zmax = numpy.nanmax(z)
                         self.plot = self.__driver.createPcolor(self.ax, x, y, z,
                                                                xmin, xmax,
                                                                ymin, ymax,
                                                                zmin, zmax,
                                                                xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title,
                                                                 colormap=colormap,
                                                                 **kwargs)
                         if self.xmin == None: self.xmin = xmin
                         if self.xmax == None: self.xmax = xmax
                         if self.ymin == None: self.ymin = ymin
                         if self.ymax == None: self.ymax = ymax
                         if self.zmin == None: self.zmin = zmin
                         if self.zmax == None: self.zmax = zmax
                         self.__firsttime = False
                         return
                     if rti:
                         self.__driver.addpcolor(self.ax, x, y, z, self.zmin, self.zmax,
                                                 xlabel=xlabel,
                                                 ylabel=ylabel,
                                                 title=title,
                                                 colormap=colormap)
                         return
                     self.__driver.pcolor(self.plot, z,
                                          xlabel=xlabel,
                                          ylabel=ylabel,
                                          title=title)
                 def pcolorbuffer(self, x, y, z,
                            xmin=None, xmax=None,
                            ymin=None, ymax=None,
                            zmin=None, zmax=None,
                            xlabel='', ylabel='',
                            title='', rti = True, colormap='jet',
                            maxNumX = None, maxNumY = None,
                            **kwargs):
                     if maxNumX == None:
                         maxNumX = self.__MAXNUMX
                     if maxNumY == None:
                         maxNumY = self.__MAXNUMY
                     if self.__firsttime:
                         self.z_buffer = z
                         self.x_buffer = numpy.hstack((self.x_buffer, x))
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmin == None: zmin = numpy.nanmin(z)
                         if zmax == None: zmax = numpy.nanmax(z)
                         self.plot = self.__driver.createPcolor(self.ax, self.x_buffer, y, z,
                                                                xmin, xmax,
                                                                ymin, ymax,
                                                                zmin, zmax,
                                                                xlabel=xlabel,
                                                                 ylabel=ylabel,
                                                                 title=title,
                                                                 colormap=colormap,
                                                                 **kwargs)
                         if self.xmin == None: self.xmin = xmin
                         if self.xmax == None: self.xmax = xmax
                         if self.ymin == None: self.ymin = ymin
                         if self.ymax == None: self.ymax = ymax
                         if self.zmin == None: self.zmin = zmin
                         if self.zmax == None: self.zmax = zmax
                         self.__firsttime = False
                         return
                     self.x_buffer = numpy.hstack((self.x_buffer, x[-1]))
                     self.z_buffer = numpy.hstack((self.z_buffer, z))
                     if self.decimationx == None:
                         deltax = float(self.xmax - self.xmin)/maxNumX
                         deltay = float(self.ymax - self.ymin)/maxNumY
                         resolutionx = self.x_buffer[2]-self.x_buffer[0]
                         resolutiony = y[1]-y[0]
                         self.decimationx = numpy.ceil(deltax / resolutionx)
                         self.decimationy = numpy.ceil(deltay / resolutiony)
                     z_buffer = self.z_buffer.reshape(-1,len(y))
                     x_buffer = self.x_buffer[::self.decimationx]
                     y_buffer = y[::self.decimationy]
                     z_buffer = z_buffer[::self.decimationx, ::self.decimationy]
                     #===================================================
                     x_buffer, y_buffer, z_buffer = self.__fillGaps(x_buffer, y_buffer, z_buffer)
                     self.__driver.addpcolorbuffer(self.ax, x_buffer, y_buffer, z_buffer, self.zmin, self.zmax,
                                             xlabel=xlabel,
                                             ylabel=ylabel,
                                             title=title,
                                             colormap=colormap)
                 def polar(self, x, y,
                            title='', xlabel='',ylabel='',**kwargs):
                     if self.__firsttime:
                         self.plot = self.__driver.createPolar(self.ax, x, y, title = title, xlabel = xlabel, ylabel = ylabel)
                         self.__firsttime = False
                         self.x_buffer = x
                         self.y_buffer = y
                         return
                     self.x_buffer = numpy.hstack((self.x_buffer,x))
                     self.y_buffer = numpy.hstack((self.y_buffer,y))
                     self.__driver.polar(self.plot, self.x_buffer, self.y_buffer, xlabel=xlabel,
                                                         ylabel=ylabel,
                                                         title=title)
                 def __fillGaps(self, x_buffer, y_buffer, z_buffer):
                     deltas = x_buffer[1:] - x_buffer[0:-1]
                     x_median = numpy.median(deltas)
                     index = numpy.where(deltas >= 2*x_median)
                     if len(index[0]) != 0:
                         z_buffer[index[0],::] = self.__missing
                         z_buffer = numpy.ma.masked_inside(z_buffer,0.99*self.__missing,1.01*self.__missing)
                     return x_buffer, y_buffer, z_buffer
   No newline at end of file

schainpy/model/graphics/jroplot_parameters.py +211 -13

             import os
             import datetime
             import numpy
             from figure import Figure, isRealtime
             from plotting_codes import *
             class MomentsPlot(Figure):
                 isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
                 PREFIX = 'prm'
                 def __init__(self):
                     self.isConfig = False
                     self.__nsubplots = 1
                     self.WIDTH = 280
                     self.HEIGHT = 250
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.counter_imagwr = 0
                     self.PLOT_CODE = MOMENTS_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                 def getSubplots(self):
                     ncol = int(numpy.sqrt(self.nplots)+0.9)
                     nrow = int(self.nplots*1./ncol + 0.9)
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 1
                     colspan = 1
                     if showprofile:
                         ncolspan = 3
                         colspan = 2
                         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()
                     counter = 0
                     for y in range(nrow):
                         for x in range(ncol):
                             if counter >= self.nplots:
                                 break
                             self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
                             if showprofile:
                                 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
                             counter += 1
                 def run(self, dataOut, id, wintitle="", channelList=None, showprofile=True,
                         xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=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, realtime=False):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         ymin            :    None,
                         ymax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
                     if dataOut.flagNoData:
                         return None
                     if realtime:
                         if not(isRealtime(utcdatatime = dataOut.utctime)):
                             print 'Skipping this plot function'
                             return
                     if channelList == None:
                         channelIndexList = dataOut.channelIndexList
                     else:
                         channelIndexList = []
                         for channel in channelList:
                             if channel not in dataOut.channelList:
                                 raise ValueError, "Channel %d is not in dataOut.channelList"
                             channelIndexList.append(dataOut.channelList.index(channel))
                     factor = dataOut.normFactor
                     x = dataOut.abscissaList
                     y = dataOut.heightList
                     z = dataOut.data_pre[channelIndexList,:,:]/factor
                     z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                     avg = numpy.average(z, axis=1)
                     noise = dataOut.noise/factor
                     zdB = 10*numpy.log10(z)
                     avgdB = 10*numpy.log10(avg)
                     noisedB = 10*numpy.log10(noise)
                     #thisDatetime = dataOut.datatime
                     thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
                     title = wintitle + " Parameters"
                     xlabel = "Velocity (m/s)"
                     ylabel = "Range (Km)"
                     if not self.isConfig:
                         nplots = len(channelIndexList)
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmin == None: zmin = numpy.nanmin(avgdB)*0.9
                         if zmax == None: zmax = numpy.nanmax(avgdB)*0.9
                         self.FTP_WEI = ftp_wei
                         self.EXP_CODE = exp_code
                         self.SUB_EXP_CODE = sub_exp_code
                         self.PLOT_POS = plot_pos
                         self.isConfig = True
                     self.setWinTitle(title)
                     for i in range(self.nplots):
                         str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
                         title = "Channel %d: %4.2fdB: %s" %(dataOut.channelList[i]+1, noisedB[i], str_datetime)
                         axes = self.axesList[i*self.__nsubplots]
                         axes.pcolor(x, y, zdB[i,:,:],
                                     xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
                                     xlabel=xlabel, ylabel=ylabel, title=title,
                                     ticksize=9, cblabel='')
                                     #Mean Line
                         mean = dataOut.data_param[i, 1, :]
                         axes.addpline(mean, y, idline=0, color="black", linestyle="solid", lw=1)
                         if self.__showprofile:
                             axes = self.axesList[i*self.__nsubplots +1]
                             axes.pline(avgdB[i], y,
                                     xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
                                     xlabel='dB', ylabel='', title='',
                                     ytick_visible=False,
                                     grid='x')
                             noiseline = numpy.repeat(noisedB[i], len(y))
                             axes.addpline(noiseline, y, idline=1, color="black", linestyle="dashed", lw=2)
                     self.draw()
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime)
             class SkyMapPlot(Figure):
                 __isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
-                PREFIX = 'prm'
+                PREFIX = 'mmap'
                 def __init__(self):
                     self.__isConfig = False
                     self.__nsubplots = 1
             #         self.WIDTH = 280
             #         self.HEIGHT = 250
                     self.WIDTH = 600
                     self.HEIGHT = 600
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.counter_imagwr = 0
-                    self.PLOT_CODE = SKYMAP_CODE
+                    self.PLOT_CODE = MSKYMAP_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                 def getSubplots(self):
                     ncol = int(numpy.sqrt(self.nplots)+0.9)
                     nrow = int(self.nplots*1./ncol + 0.9)
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=False, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 1
                     colspan = 1
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH, #+ self.WIDTHPROF,
                                       heightplot = self.HEIGHT,# + self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = 1,1
                     counter = 0
                     x = 0
                     y = 0
                     self.addAxes(1, 1, 0, 0, 1, 1, True)
                 def run(self, dataOut, id, wintitle="", channelList=None, showprofile=False,
                         xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=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, realtime=False):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         ymin            :    None,
                         ymax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
-                    arrayParameters = dataOut.data_param
+                    arrayParameters = dataOut.data_param[0,:]
                     error = arrayParameters[:,-1]
                     indValid = numpy.where(error == 0)[0]
                     finalMeteor = arrayParameters[indValid,:]
                     finalAzimuth = finalMeteor[:,4]
                     finalZenith = finalMeteor[:,5]
                     x = finalAzimuth*numpy.pi/180
                     y = finalZenith
                     #thisDatetime = dataOut.datatime
                     thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
                     title = wintitle + " Parameters"
                     xlabel = "Zonal Zenith Angle (deg) "
                     ylabel = "Meridional Zenith Angle (deg)"
                     if not self.__isConfig:
                         nplots = 1
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         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.firstdate = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
                         self.__isConfig = True
                     self.setWinTitle(title)
                     i = 0
                     str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
                     axes = self.axesList[i*self.__nsubplots]
                     nevents = axes.x_buffer.shape[0] + x.shape[0]
                     title = "Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n" %(self.firstdate,str_datetime,nevents)
                     axes.polar(x, y,
                                 title=title, xlabel=xlabel, ylabel=ylabel,
                                 ticksize=9, cblabel='')
                     self.draw()
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime)
             class WindProfilerPlot(Figure):
                 __isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
                 PREFIX = 'wind'
                 def __init__(self):
-                    self.timerange = 2*60*60
+                    self.timerange = None
                     self.__isConfig = False
                     self.__nsubplots = 1
                     self.WIDTH = 800
                     self.HEIGHT = 150
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.counter_imagwr = 0
                     self.PLOT_CODE = WIND_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                     self.tmin = None
                     self.tmax = None
                     self.xmin = None
                     self.xmax = None
                     self.figfile = None
                 def getSubplots(self):
                     ncol = 1
                     nrow = self.nplots
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 1
                     colspan = 1
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH + self.WIDTHPROF,
                                       heightplot = self.HEIGHT + self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = self.getSubplots()
                     counter = 0
                     for y in range(nrow):
                         if counter >= self.nplots:
                             break
                         self.addAxes(nrow, ncol*ncolspan, y, 0, colspan, 1)
                         counter += 1
-                def run(self, dataOut, id, wintitle="", channelList=None,
+                def run(self, dataOut, id, wintitle="", channelList=None, showprofile='False',
                         xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
                         zmax_ver = None, zmin_ver = None, SNRmin = None, SNRmax = None,
                         timerange=None, SNRthresh = None,
                         save=False, figpath='./', lastone=0,figfile=None, ftp=False, wr_period=1, show=True,
                         server=None, folder=None, username=None, password=None,
                         ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         ymin            :    None,
                         ymax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
                     if channelList == None:
                         channelIndexList = dataOut.channelIndexList
                     else:
                         channelIndexList = []
                         for channel in channelList:
                             if channel not in dataOut.channelList:
                                 raise ValueError, "Channel %d is not in dataOut.channelList"
                             channelIndexList.append(dataOut.channelList.index(channel))
-                    if timerange != None:
+            #         if timerange != None:
-                        self.timerange = timerange
+            #             self.timerange = timerange
+            #
-                    tmin = None
+            #         tmin = None
-                    tmax = None
+            #         tmax = None
                     x = dataOut.getTimeRange1()
             #         y = dataOut.heightList
                     y = dataOut.heightList
                     z = dataOut.data_output.copy()
                     nplots = z.shape[0]    #Number of wind dimensions estimated
                     nplotsw = nplots
                     #If there is a SNR function defined
                     if dataOut.data_SNR != None:
                         nplots += 1
                         SNR = dataOut.data_SNR
                         SNRavg = numpy.average(SNR, axis=0)
                         SNRdB = 10*numpy.log10(SNR)
                         SNRavgdB = 10*numpy.log10(SNRavg)
                         if SNRthresh == None: SNRthresh = -5.0
                         ind = numpy.where(SNRavg < 10**(SNRthresh/10))[0]
                         for i in range(nplotsw):
                             z[i,ind] = numpy.nan
-                    showprofile = False
+            #         showprofile = False
             #        thisDatetime = dataOut.datatime
                     thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
                     title = wintitle + "Wind"
                     xlabel = ""
                     ylabel = "Range (Km)"
                     if not self.__isConfig:
                         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 = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmax == None: zmax = numpy.nanmax(abs(z[range(2),:]))
                         #if numpy.isnan(zmax): zmax = 50
                         if zmin == None: zmin = -zmax
                         if nplotsw == 3:
                             if zmax_ver == None: zmax_ver = numpy.nanmax(abs(z[2,:]))
                             if zmin_ver == None: zmin_ver = -zmax_ver
                         if dataOut.data_SNR != None:
                             if SNRmin == None:  SNRmin = numpy.nanmin(SNRavgdB)
                             if SNRmax == None:  SNRmax = numpy.nanmax(SNRavgdB)
                         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.setWinTitle(title)
                     if ((self.xmax - x[1]) < (x[1]-x[0])):
                         x[1] = self.xmax
                     strWind = ['Zonal', 'Meridional', 'Vertical']
                     strCb = ['Velocity (m/s)','Velocity (m/s)','Velocity (cm/s)']
                     zmaxVector = [zmax, zmax, zmax_ver]
                     zminVector = [zmin, zmin, zmin_ver]
                     windFactor = [1,1,100]
                     for i in range(nplotsw):
                         title = "%s Wind: %s" %(strWind[i], thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                         axes = self.axesList[i*self.__nsubplots]
                         z1 = z[i,:].reshape((1,-1))*windFactor[i]
                         axes.pcolorbuffer(x, y, z1,
                                     xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zminVector[i], zmax=zmaxVector[i],
                                     xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
                                     ticksize=9, cblabel=strCb[i], cbsize="1%", colormap="RdBu_r" )
                     if dataOut.data_SNR != None:
                         i += 1
                         title = "Signal Noise Ratio (SNR): %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                         axes = self.axesList[i*self.__nsubplots]
                         SNRavgdB = SNRavgdB.reshape((1,-1))
                         axes.pcolorbuffer(x, y, SNRavgdB,
                                     xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=SNRmin, zmax=SNRmax,
                                     xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
                                     ticksize=9, cblabel='', cbsize="1%", colormap="jet")
                     self.draw()
                     if x[1] >= self.axesList[0].xmax:
                         self.counter_imagwr = wr_period
                         self.__isConfig = False
                         self.figfile = None
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime,
                               update_figfile=False)
             class ParametersPlot(Figure):
                 __isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
                 PREFIX = 'prm'
                 def __init__(self):
                     self.timerange = 2*60*60
                     self.__isConfig = False
                     self.__nsubplots = 1
                     self.WIDTH = 800
                     self.HEIGHT = 150
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.counter_imagwr = 0
                     self.PLOT_CODE = PARMS_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                     self.tmin = None
                     self.tmax = None
                     self.xmin = None
                     self.xmax = None
                     self.figfile = None
                 def getSubplots(self):
                     ncol = 1
                     nrow = self.nplots
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 1
                     colspan = 1
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH + self.WIDTHPROF,
                                       heightplot = self.HEIGHT + self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = self.getSubplots()
                     counter = 0
                     for y in range(nrow):
                         for x in range(ncol):
                             if counter >= self.nplots:
                                 break
                             self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
                             if showprofile:
                                 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
                             counter += 1
                 def run(self, dataOut, id, wintitle="", channelList=None, showprofile=False,
                         xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,timerange=None,
                         parameterIndex = None, onlyPositive = False,
-                        SNRthresh = -numpy.inf, SNR = True, SNRmin = None, SNRmax = None,
+                        SNRthresh = -numpy.inf, SNR = True, SNRmin = None, SNRmax = None, onlySNR = False,
                         DOP = True,
                         zlabel = "", parameterName = "", parameterObject = "data_param",
                         save=False, figpath='./', lastone=0,figfile=None, ftp=False, wr_period=1, show=True,
                         server=None, folder=None, username=None, password=None,
                         ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         ymin            :    None,
                         ymax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
                     data_param = getattr(dataOut, parameterObject)
                     if channelList == None:
                         channelIndexList = numpy.arange(data_param.shape[0])
                     else:
                         channelIndexList = numpy.array(channelList)
                     nchan = len(channelIndexList) #Number of channels being plotted
                     if nchan < 1:
                         return
                     nGraphsByChannel = 0
                     if SNR:
                         nGraphsByChannel += 1
                     if DOP:
                         nGraphsByChannel += 1
                     if nGraphsByChannel < 1:
                         return
                     nplots = nGraphsByChannel*nchan
                     if timerange != None:
                         self.timerange = timerange
                     #tmin = None
                     #tmax = None
                     if parameterIndex == None:
                         parameterIndex = 1
                     x = dataOut.getTimeRange1()
                     y = dataOut.heightList
                     z = data_param[channelIndexList,parameterIndex,:].copy()
                     zRange = dataOut.abscissaList
             #         nChannels = z.shape[0]    #Number of wind dimensions estimated
             #        thisDatetime = dataOut.datatime
                     if dataOut.data_SNR != None:
                         SNRarray = dataOut.data_SNR[channelIndexList,:]
                         SNRdB = 10*numpy.log10(SNRarray)
             #             SNRavgdB = 10*numpy.log10(SNRavg)
                         ind = numpy.where(SNRdB < 10**(SNRthresh/10))
                         z[ind] = numpy.nan
                     thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
                     title = wintitle + " Parameters Plot" #: %s" %(thisDatetime.strftime("%d-%b-%Y"))
                     xlabel = ""
                     ylabel = "Range (Km)"
+                    if (SNR and not onlySNR): nplots = 2*nplots
                     if onlyPositive:
                         colormap = "jet"
                         zmin = 0
                     else: colormap = "RdBu_r"
                     if not self.__isConfig:
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmin == None: zmin = numpy.nanmin(zRange)
                         if zmax == None: zmax = numpy.nanmax(zRange)
                         if SNR:
                             if SNRmin == None:  SNRmin = numpy.nanmin(SNRdB)
                             if SNRmax == None:  SNRmax = numpy.nanmax(SNRdB)
                         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.setWinTitle(title)
                     if ((self.xmax - x[1]) < (x[1]-x[0])):
                         x[1] = self.xmax
                     for i in range(nchan):
+                        if (SNR and not onlySNR): j = 2*i
+                        else: j = i
                         j = nGraphsByChannel*i
+                        if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
+                            title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
+                        if not onlySNR:
+                            axes = self.axesList[j*self.__nsubplots]
+                            z1 = z[i,:].reshape((1,-1))
+                            axes.pcolorbuffer(x, y, z1,
+                                    xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
+                                    xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,colormap=colormap,
+                                    ticksize=9, cblabel=zlabel, cbsize="1%")
                         if DOP:
                             title = "%s Channel %d: %s" %(parameterName, channelIndexList[i]+1, thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                             if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
                                 title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
                             axes = self.axesList[j]
                             z1 = z[i,:].reshape((1,-1))
                             axes.pcolorbuffer(x, y, z1,
                                         xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
                                         xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,colormap=colormap,
                                         ticksize=9, cblabel=zlabel, cbsize="1%")
                         if SNR:
                             title = "Channel %d Signal Noise Ratio (SNR): %s" %(channelIndexList[i]+1, thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
+                            axes = self.axesList[(j)*self.__nsubplots]
+                            if not onlySNR:
+                                axes = self.axesList[(j + 1)*self.__nsubplots]
                             axes = self.axesList[(j + nGraphsByChannel-1)]
                             z1 = SNRdB[i,:].reshape((1,-1))
                             axes.pcolorbuffer(x, y, z1,
                                     xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=SNRmin, zmax=SNRmax,
                                     xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,colormap="jet",
                                     ticksize=9, cblabel=zlabel, cbsize="1%")
                     self.draw()
                     if x[1] >= self.axesList[0].xmax:
                         self.counter_imagwr = wr_period
                         self.__isConfig = False
                         self.figfile = None
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime,
                               update_figfile=False)
             class SpectralFittingPlot(Figure):
                 __isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
                 PREFIX = 'prm'
                 N = None
                 ippSeconds = None
                 def __init__(self):
                     self.__isConfig = False
                     self.__nsubplots = 1
                     self.PLOT_CODE = SPECFIT_CODE
                     self.WIDTH = 450
                     self.HEIGHT = 250
                     self.WIDTHPROF = 0
                     self.HEIGHTPROF = 0
                 def getSubplots(self):
                     ncol = int(numpy.sqrt(self.nplots)+0.9)
                     nrow = int(self.nplots*1./ncol + 0.9)
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=False, show=True):
                     showprofile = False
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 5
                     colspan = 4
                     if showprofile:
                         ncolspan = 5
                         colspan = 4
                         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()
                     counter = 0
                     for y in range(nrow):
                         for x in range(ncol):
                             if counter >= self.nplots:
                                 break
                             self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
                             if showprofile:
                                 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
                             counter += 1
                 def run(self, dataOut, id, cutHeight=None, fit=False, wintitle="", channelList=None, showprofile=True,
                         xmin=None, xmax=None, ymin=None, ymax=None,
                         save=False, figpath='./', figfile=None, show=True):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
                     if cutHeight==None:
                         h=270
                     heightindex = numpy.abs(cutHeight - dataOut.heightList).argmin()
                     cutHeight = dataOut.heightList[heightindex]
                     factor = dataOut.normFactor
                     x = dataOut.abscissaList[:-1]
                     #y = dataOut.getHeiRange()
                     z = dataOut.data_pre[:,:,heightindex]/factor
                     z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                     avg = numpy.average(z, axis=1)
                     listChannels = z.shape[0]
                     #Reconstruct Function
                     if fit==True:
                         groupArray = dataOut.groupList
                         listChannels = groupArray.reshape((groupArray.size))
                         listChannels.sort()
                         spcFitLine = numpy.zeros(z.shape)
                         constants = dataOut.constants
                         nGroups = groupArray.shape[0]
                         nChannels = groupArray.shape[1]
                         nProfiles = z.shape[1]
                         for f in range(nGroups):
                             groupChann = groupArray[f,:]
                             p = dataOut.data_param[f,:,heightindex]
             #                 p = numpy.array([ 89.343967,0.14036615,0.17086219,18.89835291,1.58388365,1.55099167])
                             fitLineAux = dataOut.library.modelFunction(p, constants)*nProfiles
                             fitLineAux = fitLineAux.reshape((nChannels,nProfiles))
                             spcFitLine[groupChann,:] = fitLineAux
             #             spcFitLine = spcFitLine/factor
                         z = z[listChannels,:]
                         spcFitLine = spcFitLine[listChannels,:]
                         spcFitLinedB = 10*numpy.log10(spcFitLine)
                     zdB = 10*numpy.log10(z)
                     #thisDatetime = dataOut.datatime
                     thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
                     title = wintitle + " Doppler Spectra: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
                     xlabel = "Velocity (m/s)"
                     ylabel = "Spectrum"
                     if not self.__isConfig:
                         nplots = listChannels.size
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         if xmin == None: xmin = numpy.nanmin(x)
                         if xmax == None: xmax = numpy.nanmax(x)
                         if ymin == None: ymin = numpy.nanmin(zdB)
                         if ymax == None: ymax = numpy.nanmax(zdB)+2
                         self.__isConfig = True
                     self.setWinTitle(title)
                     for i in range(self.nplots):
             #             title = "Channel %d: %4.2fdB" %(dataOut.channelList[i]+1, noisedB[i])
                         title = "Height %4.1f km\nChannel %d:" %(cutHeight, listChannels[i]+1)
                         axes = self.axesList[i*self.__nsubplots]
                         if fit == False:
                             axes.pline(x, zdB[i,:],
                                         xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
                                         xlabel=xlabel, ylabel=ylabel, title=title
                                         )
                         if fit == True:
                             fitline=spcFitLinedB[i,:]
                             y=numpy.vstack([zdB[i,:],fitline] )
                             legendlabels=['Data','Fitting']
                             axes.pmultilineyaxis(x, y,
                                         xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
                                         xlabel=xlabel, ylabel=ylabel, title=title,
                                         legendlabels=legendlabels, marker=None,
                                         linestyle='solid', grid='both')
                     self.draw()
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime)
             class EWDriftsPlot(Figure):
                 __isConfig = None
                 __nsubplots = None
                 WIDTHPROF = None
                 HEIGHTPROF = None
                 PREFIX = 'drift'
                 def __init__(self):
                     self.timerange = 2*60*60
                     self.isConfig = False
                     self.__nsubplots = 1
                     self.WIDTH = 800
                     self.HEIGHT = 150
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.counter_imagwr = 0
                     self.PLOT_CODE = EWDRIFT_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                     self.tmin = None
                     self.tmax = None
                     self.xmin = None
                     self.xmax = None
                     self.figfile = None
                 def getSubplots(self):
                     ncol = 1
                     nrow = self.nplots
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 1
                     colspan = 1
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH + self.WIDTHPROF,
                                       heightplot = self.HEIGHT + self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = self.getSubplots()
                     counter = 0
                     for y in range(nrow):
                         if counter >= self.nplots:
                             break
                         self.addAxes(nrow, ncol*ncolspan, y, 0, colspan, 1)
                         counter += 1
                 def run(self, dataOut, id, wintitle="", channelList=None,
                         xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
                         zmaxVertical = None, zminVertical = None, zmaxZonal = None, zminZonal = None,
                         timerange=None, SNRthresh = -numpy.inf, SNRmin = None, SNRmax = None, SNR_1 = False,
                         save=False, figpath='./', lastone=0,figfile=None, ftp=False, wr_period=1, show=True,
                         server=None, folder=None, username=None, password=None,
                         ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
                     """
                     Input:
                         dataOut         :
                         id        :
                         wintitle        :
                         channelList     :
                         showProfile     :
                         xmin            :    None,
                         xmax            :    None,
                         ymin            :    None,
                         ymax            :    None,
                         zmin            :    None,
                         zmax            :    None
                     """
                     if timerange != None:
                         self.timerange = timerange
                     tmin = None
                     tmax = None
                     x = dataOut.getTimeRange1()
             #         y = dataOut.heightList
                     y = dataOut.heightList
                     z = dataOut.data_output
                     nplots = z.shape[0]    #Number of wind dimensions estimated
                     nplotsw = nplots
                     #If there is a SNR function defined
                     if dataOut.data_SNR != None:
                         nplots += 1
                         SNR = dataOut.data_SNR
                         if SNR_1:
                             SNR += 1
                         SNRavg = numpy.average(SNR, axis=0)
                         SNRdB = 10*numpy.log10(SNR)
                         SNRavgdB = 10*numpy.log10(SNRavg)
                         ind = numpy.where(SNRavg < 10**(SNRthresh/10))[0]
                         for i in range(nplotsw):
                             z[i,ind] = numpy.nan
                     showprofile = False
             #        thisDatetime = dataOut.datatime
                     thisDatetime = datetime.datetime.utcfromtimestamp(x[1])
                     title = wintitle + " EW Drifts"
                     xlabel = ""
                     ylabel = "Height (Km)"
                     if not self.__isConfig:
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
                         if ymin == None: ymin = numpy.nanmin(y)
                         if ymax == None: ymax = numpy.nanmax(y)
                         if zmaxZonal == None: zmaxZonal = numpy.nanmax(abs(z[0,:]))
                         if zminZonal == None: zminZonal = -zmaxZonal
                         if zmaxVertical == None: zmaxVertical = numpy.nanmax(abs(z[1,:]))
                         if zminVertical == None: zminVertical = -zmaxVertical
                         if dataOut.data_SNR != None:
                             if SNRmin == None:  SNRmin = numpy.nanmin(SNRavgdB)
                             if SNRmax == None:  SNRmax = numpy.nanmax(SNRavgdB)
                         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.setWinTitle(title)
                     if ((self.xmax - x[1]) < (x[1]-x[0])):
                         x[1] = self.xmax
                     strWind = ['Zonal','Vertical']
                     strCb = 'Velocity (m/s)'
                     zmaxVector = [zmaxZonal, zmaxVertical]
                     zminVector = [zminZonal, zminVertical]
                     for i in range(nplotsw):
                         title = "%s Drifts: %s" %(strWind[i], thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                         axes = self.axesList[i*self.__nsubplots]
                         z1 = z[i,:].reshape((1,-1))
                         axes.pcolorbuffer(x, y, z1,
                                     xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zminVector[i], zmax=zmaxVector[i],
                                     xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
                                     ticksize=9, cblabel=strCb, cbsize="1%", colormap="RdBu_r")
                     if dataOut.data_SNR != None:
                         i += 1
                         if SNR_1:
                             title = "Signal Noise Ratio + 1 (SNR+1): %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                         else:
                             title = "Signal Noise Ratio (SNR): %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                         axes = self.axesList[i*self.__nsubplots]
                         SNRavgdB = SNRavgdB.reshape((1,-1))
                         axes.pcolorbuffer(x, y, SNRavgdB,
                                     xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=SNRmin, zmax=SNRmax,
                                     xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
                                     ticksize=9, cblabel='', cbsize="1%", colormap="jet")
                     self.draw()
                     if x[1] >= self.axesList[0].xmax:
                         self.counter_imagwr = wr_period
                         self.__isConfig = False
                         self.figfile = None
+            class PhasePlot(Figure):
+                __isConfig = None
+                __nsubplots = None
+                PREFIX = 'mphase'
+                def __init__(self):
+                    self.timerange = 24*60*60
+                    self.__isConfig = False
+                    self.__nsubplots = 1
+                    self.counter_imagwr = 0
+                    self.WIDTH = 600
+                    self.HEIGHT = 300
+                    self.WIDTHPROF = 120
+                    self.HEIGHTPROF = 0
+                    self.xdata = None
+                    self.ydata = None
+                    self.PLOT_CODE = MPHASE_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
+                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 run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
+                        xmin=None, xmax=None, ymin=None, ymax=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 timerange != None:
+                        self.timerange = timerange
+                    tmin = None
+                    tmax = None
+                    x = dataOut.getTimeRange1()
+                    y = dataOut.getHeiRange()
+                    #thisDatetime = dataOut.datatime
+                    thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[1])
+                    title = wintitle + " Phase of Beacon Signal" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
+                    xlabel = "Local Time"
+                    ylabel = "Phase"
+                    #phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
+                    phase_beacon = dataOut.data_output
+                    if not self.__isConfig:
+                        self.nplots = phase_beacon.size
+                        self.setup(id=id,
+                                   nplots=self.nplots,
+                                   wintitle=wintitle,
+                                   showprofile=showprofile,
+                                   show=show)
+                        tmin, tmax = self.getTimeLim(x, xmin, xmax)
+                        if ymin == None: ymin = numpy.nanmin(phase_beacon) - 10.0
+                        if ymax == None: ymax = numpy.nanmax(phase_beacon) + 10.0
+                        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([])
+                        #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 Offset %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
+                    legendlabels = ["phase %d"%(chan) for chan in numpy.arange(self.nplots)]
+                    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=tmin, xmax=tmax, 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 x[1] >= self.axesList[0].xmax:
+                        self.counter_imagwr = wr_period
+                        del self.xdata
+                        del self.ydata
+                        self.__isConfig = False
+                    if self.figfile == None:
+                        str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
+                        self.figfile = self.getFilename(name = str_datetime)
+                    if figpath != '':
+                        self.counter_imagwr += 1
+                        if (self.counter_imagwr>=wr_period):
+                            # store png plot to local folder
+                            self.saveFigure(figpath, self.figfile)
+                            # store png plot to FTP server according to RT-Web format
+                            name = self.getNameToFtp(thisDatetime, self.FTP_WEI, self.EXP_CODE, self.SUB_EXP_CODE, self.PLOT_CODE, self.PLOT_POS)
+                            ftp_filename = os.path.join(figpath, name)
+                            self.saveFigure(figpath, ftp_filename)
+                            self.counter_imagwr = 0
+                        self.figfile = None
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime,
                               update_figfile=False)

schainpy/model/graphics/mpldriver.py +5 -3

             import numpy
             import datetime
             import sys
             import matplotlib
             if 'linux' in sys.platform:
                 matplotlib.use("TKAgg")
             if 'darwin' in sys.platform:
                 matplotlib.use("WXAgg")
             #Qt4Agg', 'GTK', 'GTKAgg', 'ps', 'agg', 'cairo', 'MacOSX', 'GTKCairo', 'WXAgg', 'template', 'TkAgg', 'GTK3Cairo', 'GTK3Agg', 'svg', 'WebAgg', 'CocoaAgg', 'emf', 'gdk', 'WX'
             import matplotlib.pyplot
             from mpl_toolkits.axes_grid1 import make_axes_locatable
             from matplotlib.ticker import *
             ###########################################
             #Actualizacion de las funciones del driver
             ###########################################
             def createFigure(id, wintitle, width, height, facecolor="w", show=True):
                 matplotlib.pyplot.ioff()
                 fig = matplotlib.pyplot.figure(num=id, facecolor=facecolor)
                 fig.canvas.manager.set_window_title(wintitle)
                 fig.canvas.manager.resize(width, height)
                 matplotlib.pyplot.ion()
                 if show:
                     matplotlib.pyplot.show()
                 return fig
             def closeFigure(show=False, fig=None):
                 matplotlib.pyplot.ioff()
                 matplotlib.pyplot.pause(0.1)
                 if show:
                     matplotlib.pyplot.show()
                 if fig != None:
                     matplotlib.pyplot.close(fig)
                     matplotlib.pyplot.pause(0.1)
                     matplotlib.pyplot.ion()
                     return
                 matplotlib.pyplot.close("all")
                 matplotlib.pyplot.pause(0.1)
                 matplotlib.pyplot.ion()
                 return
             def saveFigure(fig, filename):
                 matplotlib.pyplot.ioff()
                 fig.savefig(filename)
                 matplotlib.pyplot.ion()
             def setWinTitle(fig, title):
                 fig.canvas.manager.set_window_title(title)
             def setTitle(fig, title):
                 fig.suptitle(title)
             def createAxes(fig, nrow, ncol, xpos, ypos, colspan, rowspan, polar=False):
                 matplotlib.pyplot.ioff()
                 matplotlib.pyplot.figure(fig.number)
                 axes = matplotlib.pyplot.subplot2grid((nrow, ncol),
                                                     (xpos, ypos),
                                                     colspan=colspan,
                                                     rowspan=rowspan,
                                                     polar=polar)
                 matplotlib.pyplot.ion()
                 return axes
             def setAxesText(ax, text):
                 ax.annotate(text,
                             xy = (.1, .99),
                             xycoords = 'figure fraction',
                             horizontalalignment = 'left',
                             verticalalignment = 'top',
                             fontsize = 10)
             def printLabels(ax, xlabel, ylabel, title):
                 ax.set_xlabel(xlabel, size=11)
                 ax.set_ylabel(ylabel, size=11)
                 ax.set_title(title, size=8)
             def createPline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='',
                             ticksize=9, xtick_visible=True, ytick_visible=True,
                             nxticks=4, nyticks=10,
                             grid=None,color='blue'):
                 """
                 Input:
                     grid    :    None, 'both', 'x', 'y'
                 """
                 matplotlib.pyplot.ioff()
                 ax.set_xlim([xmin,xmax])
                 ax.set_ylim([ymin,ymax])
                 printLabels(ax, xlabel, ylabel, title)
                 ######################################################
                 if (xmax-xmin)<=1:
                     xtickspos = numpy.linspace(xmin,xmax,nxticks)
                     xtickspos = numpy.array([float("%.1f"%i) for i in xtickspos])
                     ax.set_xticks(xtickspos)
                 else:
                     xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
             #         xtickspos = numpy.arange(nxticks)*float(xmax-xmin)/float(nxticks) + int(xmin)
                     ax.set_xticks(xtickspos)
                 for tick in ax.get_xticklabels():
                     tick.set_visible(xtick_visible)
                 for tick in ax.xaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 ######################################################
                 for tick in ax.get_yticklabels():
                     tick.set_visible(ytick_visible)
                 for tick in ax.yaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 ax.plot(x, y, color=color)
                 iplot = ax.lines[-1]
                 ######################################################
                 if '0.' in matplotlib.__version__[0:2]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if '1.0.' in matplotlib.__version__[0:4]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if grid != None:
                     ax.grid(b=True, which='major', axis=grid)
                 matplotlib.pyplot.tight_layout()
                 matplotlib.pyplot.ion()
                 return iplot
             def set_linedata(ax, x, y, idline):
                 ax.lines[idline].set_data(x,y)
             def pline(iplot, x, y, xlabel='', ylabel='', title=''):
                 ax = iplot.get_axes()
                 printLabels(ax, xlabel, ylabel, title)
                 set_linedata(ax, x, y, idline=0)
             def addpline(ax, x, y, color, linestyle, lw):
                 ax.plot(x,y,color=color,linestyle=linestyle,lw=lw)
             def createPcolor(ax, x, y, z, xmin, xmax, ymin, ymax, zmin, zmax,
                              xlabel='', ylabel='', title='', ticksize = 9,
                              colormap='jet',cblabel='', cbsize="5%",
                              XAxisAsTime=False):
                 matplotlib.pyplot.ioff()
                 divider = make_axes_locatable(ax)
                 ax_cb = divider.new_horizontal(size=cbsize, pad=0.05)
                 fig = ax.get_figure()
                 fig.add_axes(ax_cb)
                 ax.set_xlim([xmin,xmax])
                 ax.set_ylim([ymin,ymax])
                 printLabels(ax, xlabel, ylabel, title)
                 imesh = ax.pcolormesh(x,y,z.T, vmin=zmin, vmax=zmax, cmap=matplotlib.pyplot.get_cmap(colormap))
                 cb =  matplotlib.pyplot.colorbar(imesh, cax=ax_cb)
                 cb.set_label(cblabel)
             #    for tl in ax_cb.get_yticklabels():
             #        tl.set_visible(True)
                 for tick in ax.yaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 for tick in ax.xaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 for tick in cb.ax.get_yticklabels():
                     tick.set_fontsize(ticksize)
                 ax_cb.yaxis.tick_right()
                 if '0.' in matplotlib.__version__[0:2]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return imesh
                 if '1.0.' in matplotlib.__version__[0:4]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return imesh
                 matplotlib.pyplot.tight_layout()
                 if XAxisAsTime:
                     func = lambda x, pos: ('%s') %(datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
                     ax.xaxis.set_major_formatter(FuncFormatter(func))
                     ax.xaxis.set_major_locator(LinearLocator(7))
                 matplotlib.pyplot.ion()
                 return imesh
             def pcolor(imesh, z, xlabel='', ylabel='', title=''):
                 z = z.T
                 ax = imesh.get_axes()
                 printLabels(ax, xlabel, ylabel, title)
                 imesh.set_array(z.ravel())
             def addpcolor(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', colormap='jet'):
                 printLabels(ax, xlabel, ylabel, title)
                 ax.pcolormesh(x,y,z.T,vmin=zmin,vmax=zmax, cmap=matplotlib.pyplot.get_cmap(colormap))
             def addpcolorbuffer(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', colormap='jet'):
                 printLabels(ax, xlabel, ylabel, title)
                 ax.collections.remove(ax.collections[0])
                 ax.pcolormesh(x,y,z.T,vmin=zmin,vmax=zmax, cmap=matplotlib.pyplot.get_cmap(colormap))
             def createPmultiline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
                             ticksize=9, xtick_visible=True, ytick_visible=True,
                             nxticks=4, nyticks=10,
                             grid=None):
                 """
                 Input:
                     grid    :    None, 'both', 'x', 'y'
                 """
                 matplotlib.pyplot.ioff()
                 lines = ax.plot(x.T, y)
                 leg = ax.legend(lines, legendlabels, loc='upper right')
                 leg.get_frame().set_alpha(0.5)
                 ax.set_xlim([xmin,xmax])
                 ax.set_ylim([ymin,ymax])
                 printLabels(ax, xlabel, ylabel, title)
                 xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
                 ax.set_xticks(xtickspos)
                 for tick in ax.get_xticklabels():
                     tick.set_visible(xtick_visible)
                 for tick in ax.xaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 for tick in ax.get_yticklabels():
                     tick.set_visible(ytick_visible)
                 for tick in ax.yaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 iplot = ax.lines[-1]
                 if '0.' in matplotlib.__version__[0:2]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if '1.0.' in matplotlib.__version__[0:4]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if grid != None:
                     ax.grid(b=True, which='major', axis=grid)
                 matplotlib.pyplot.tight_layout()
                 matplotlib.pyplot.ion()
                 return iplot
             def pmultiline(iplot, x, y, xlabel='', ylabel='', title=''):
                 ax = iplot.get_axes()
                 printLabels(ax, xlabel, ylabel, title)
                 for i in range(len(ax.lines)):
                     line = ax.lines[i]
                     line.set_data(x[i,:],y)
             def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
                             ticksize=9, xtick_visible=True, ytick_visible=True,
                             nxticks=4, nyticks=10, marker='.', markersize=10, linestyle="None",
                             grid=None, XAxisAsTime=False):
                 """
                 Input:
                     grid    :    None, 'both', 'x', 'y'
                 """
                 matplotlib.pyplot.ioff()
             #    lines = ax.plot(x, y.T, marker=marker,markersize=markersize,linestyle=linestyle)
                 lines = ax.plot(x, y.T, linestyle=linestyle, marker=marker, markersize=markersize)
                 leg = ax.legend(lines, legendlabels, loc='upper left', bbox_to_anchor=(1.01, 1.00), numpoints=1, handlelength=1.5, \
                                 handletextpad=0.5, borderpad=0.5, labelspacing=0.5, borderaxespad=0.)
                 for label in leg.get_texts(): label.set_fontsize(9)
                 ax.set_xlim([xmin,xmax])
                 ax.set_ylim([ymin,ymax])
                 printLabels(ax, xlabel, ylabel, title)
             #    xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
             #    ax.set_xticks(xtickspos)
                 for tick in ax.get_xticklabels():
                     tick.set_visible(xtick_visible)
                 for tick in ax.xaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 for tick in ax.get_yticklabels():
                     tick.set_visible(ytick_visible)
                 for tick in ax.yaxis.get_major_ticks():
                     tick.label.set_fontsize(ticksize)
                 iplot = ax.lines[-1]
                 if '0.' in matplotlib.__version__[0:2]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if '1.0.' in matplotlib.__version__[0:4]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if grid != None:
                     ax.grid(b=True, which='major', axis=grid)
                 matplotlib.pyplot.tight_layout()
                 if XAxisAsTime:
                     func = lambda x, pos: ('%s') %(datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
                     ax.xaxis.set_major_formatter(FuncFormatter(func))
                     ax.xaxis.set_major_locator(LinearLocator(7))
                 matplotlib.pyplot.ion()
                 return iplot
             def pmultilineyaxis(iplot, x, y, xlabel='', ylabel='', title=''):
                 ax = iplot.get_axes()
                 printLabels(ax, xlabel, ylabel, title)
                 for i in range(len(ax.lines)):
                     line = ax.lines[i]
                     line.set_data(x,y[i,:])
             def createPolar(ax, x, y,
                             xlabel='', ylabel='', title='', ticksize = 9,
                              colormap='jet',cblabel='', cbsize="5%",
                              XAxisAsTime=False):
                 matplotlib.pyplot.ioff()
                 ax.plot(x,y,'bo', markersize=5)
             #     ax.set_rmax(90)
                 ax.set_ylim(0,90)
                 ax.set_yticks(numpy.arange(0,90,20))
-                ax.text(0, -110, ylabel, rotation='vertical', va ='center', ha = 'center' ,size='11')
+            #     ax.text(0, -110, ylabel, rotation='vertical', va ='center', ha = 'center' ,size='11')
+            #     ax.text(0, 50, ylabel, rotation='vertical', va ='center', ha = 'left' ,size='11')
             #     ax.text(100, 100, 'example', ha='left', va='center', rotation='vertical')
-                printLabels(ax, xlabel, '', title)
+                ax.yaxis.labelpad = 230
+                printLabels(ax, xlabel, ylabel, title)
                 iplot = ax.lines[-1]
                 if '0.' in matplotlib.__version__[0:2]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
                 if '1.0.' in matplotlib.__version__[0:4]:
                     print "The matplotlib version has to be updated to 1.1 or newer"
                     return iplot
             #     if grid != None:
             #         ax.grid(b=True, which='major', axis=grid)
                 matplotlib.pyplot.tight_layout()
                 matplotlib.pyplot.ion()
                 return iplot
             def polar(iplot, x, y, xlabel='', ylabel='', title=''):
                 ax = iplot.get_axes()
             #     ax.text(0, -110, ylabel, rotation='vertical', va ='center', ha = 'center',size='11')
-                printLabels(ax, xlabel, '', title)
+                printLabels(ax, xlabel, ylabel, title)
                 set_linedata(ax, x, y, idline=0)
             def draw(fig):
                 if type(fig) == 'int':
                     raise ValueError, "Error drawing: Fig parameter should be a matplotlib figure object figure"
                 fig.canvas.draw()

schainpy/model/graphics/plotting_codes.py +7 -6

             '''
             @author: roj-idl71
             '''
             #USED IN jroplot_spectra.py
             RTI_CODE = 0            #Range time intensity (RTI).
             SPEC_CODE = 1           #Spectra (and Cross-spectra) information.
             CROSS_CODE = 2          #Cross-Correlation information.
             COH_CODE = 3            #Coherence map.
             BASE_CODE = 4           #Base lines graphic.
             ROW_CODE = 5            #Row Spectra.
             TOTAL_CODE = 6          #Total Power.
             DRIFT_CODE = 7          #Drifts graphics.
             HEIGHT_CODE = 8         #Height profile.
             PHASE_CODE = 9          #Signal Phase.
             POWER_CODE = 16
             NOISE_CODE = 17
             BEACON_CODE = 18
             #USED IN jroplot_parameters.py
-            MOMENTS_CODE = 20
-            SKYMAP_CODE = 21
             WIND_CODE = 22
-            PARMS_CODE = 23
+            MSKYMAP_CODE = 23
-            SPECFIT_CODE = 24
+            MPHASE_CODE = 24
-            EWDRIFT_CODE = 25
+            MOMENTS_CODE = 25
+            PARMS_CODE = 26
+            SPECFIT_CODE = 27
+            EWDRIFT_CODE = 28

schainpy/model/io/jroIO_HDF5.py +134 -32

             import numpy
             import time
             import os
             import h5py
             import re
+            import tables
             from schainpy.model.data.jrodata import *
             from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
             from schainpy.model.io.jroIO_base import *
             class HDF5Reader(ProcessingUnit):
                 ext = ".hdf5"
                 optchar = "D"
                 timezone = None
                 secStart = None
                 secEnd = None
                 fileIndex = None
                 blockIndex = None
                 blocksPerFile = None
                 path = None
                 #List of Files
                 filenameList = None
                 datetimeList = None
                 #Hdf5 File
                 fpMetadata = None
                 pathMeta = None
                 listMetaname = None
                 listMeta = None
                 listDataname = None
                 listData = None
                 listShapes = None
                 fp = None
                 #dataOut reconstruction
                 dataOut = None
                 nRecords = None
                 def __init__(self):
                     self.dataOut = self.__createObjByDefault()
                     return
                 def __createObjByDefault(self):
                     dataObj = Parameters()
                     return dataObj
                 def setup(self,path=None,
                                 startDate=None,
                                 endDate=None,
                                 startTime=datetime.time(0,0,0),
                                 endTime=datetime.time(23,59,59),
                                 walk=True,
                                 timezone='ut',
                                 all=0,
                                 online=False,
                                 ext=None):
                     if ext==None:
                         ext = self.ext
                     self.timezone = timezone
             #         self.all = all
             #         self.online = online
                     self.path = path
                     startDateTime = datetime.datetime.combine(startDate,startTime)
                     endDateTime = datetime.datetime.combine(endDate,endTime)
                     secStart = (startDateTime-datetime.datetime(1970,1,1)).total_seconds()
                     secEnd = (endDateTime-datetime.datetime(1970,1,1)).total_seconds()
                     self.secStart = secStart
                     self.secEnd = secEnd
                     if not(online):
                         #Busqueda de archivos offline
                         self.__searchFilesOffline(path, startDate, endDate, ext, startTime, endTime, secStart, secEnd, walk)
                     else:
                         self.__searchFilesOnline(path, walk)
                     if not(self.filenameList):
                         print "There is no files into the folder: %s"%(path)
                         sys.exit(-1)
             #         self.__getExpParameters()
                     self.fileIndex = -1
                     self.__setNextFileOffline()
                     self.__readMetadata()
                     self.blockIndex = 0
                     return
                 def __searchFilesOffline(self,
                                         path,
                                         startDate,
                                         endDate,
                                         ext,
                                         startTime=datetime.time(0,0,0),
                                         endTime=datetime.time(23,59,59),
                                         secStart = 0,
                                         secEnd = numpy.inf,
                                         walk=True):
             #         self.__setParameters(path, startDate, endDate, startTime, endTime, walk)
             #
             #         self.__checkPath()
             #
             #         self.__findDataForDates()
             #
             #         self.__selectDataForTimes()
             #
             #         for i in range(len(self.filenameList)):
             #             print "%s" %(self.filenameList[i])
                     pathList = []
                     if not walk:
                         #pathList.append(path)
                         multi_path = path.split(',')
                         for single_path in multi_path:
                             pathList.append(single_path)
                     else:
                         #dirList = []
                         multi_path = path.split(',')
                         for single_path in multi_path:
                             dirList = []
                             for thisPath in os.listdir(single_path):
                                 if not os.path.isdir(os.path.join(single_path,thisPath)):
                                     continue
-                                if not isRadarFolder(thisPath):
+                                if not isDoyFolder(thisPath):
                                     continue
                                 dirList.append(thisPath)
                             if not(dirList):
                                 return None, None
                             thisDate = startDate
                             while(thisDate <= endDate):
                                 year = thisDate.timetuple().tm_year
                                 doy = thisDate.timetuple().tm_yday
                                 matchlist = fnmatch.filter(dirList, '?' + '%4.4d%3.3d' % (year,doy) + '*')
                                 if len(matchlist) == 0:
                                     thisDate += datetime.timedelta(1)
                                     continue
                                 for match in matchlist:
                                     pathList.append(os.path.join(single_path,match))
                                 thisDate += datetime.timedelta(1)
                     if pathList == []:
                         print "Any folder was found for the date range: %s-%s" %(startDate, endDate)
                         return None, None
                     print "%d folder(s) was(were) found for the date range: %s - %s" %(len(pathList), startDate, endDate)
                     filenameList = []
                     datetimeList = []
                     pathDict = {}
                     filenameList_to_sort = []
                     for i in range(len(pathList)):
                         thisPath = pathList[i]
                         fileList = glob.glob1(thisPath, "*%s" %ext)
                         fileList.sort()
                         pathDict.setdefault(fileList[0])
                         pathDict[fileList[0]] = i
                         filenameList_to_sort.append(fileList[0])
                     filenameList_to_sort.sort()
                     for file in filenameList_to_sort:
                         thisPath = pathList[pathDict[file]]
                         fileList = glob.glob1(thisPath, "*%s" %ext)
                         fileList.sort()
                         for file in fileList:
                             filename = os.path.join(thisPath,file)
                             thisDatetime = self.__isFileinThisTime(filename, secStart, secEnd)
                             if not(thisDatetime):
                                 continue
                             filenameList.append(filename)
                             datetimeList.append(thisDatetime)
                     if not(filenameList):
                         print "Any file was found for the time range %s - %s" %(startTime, endTime)
                         return None, None
                     print "%d file(s) was(were) found for the time range: %s - %s" %(len(filenameList), startTime, endTime)
                     print
                     for i in range(len(filenameList)):
                         print "%s -> [%s]" %(filenameList[i], datetimeList[i].ctime())
                     self.filenameList = filenameList
                     self.datetimeList = datetimeList
                     return pathList, filenameList
                 def __isFileinThisTime(self, filename, startSeconds, endSeconds):
                     """
                     Retorna 1 si el archivo de datos se encuentra dentro del rango de horas especificado.
                     Inputs:
                         filename            :    nombre completo del archivo de datos en formato Jicamarca (.r)
                         startTime          :    tiempo inicial del rango seleccionado en formato datetime.time
                         endTime            :    tiempo final del rango seleccionado en formato datetime.time
                     Return:
                         Boolean    :    Retorna True si el archivo de datos contiene datos en el rango de
                                         fecha especificado, de lo contrario retorna False.
                     Excepciones:
                         Si el archivo no existe o no puede ser abierto
                         Si la cabecera no puede ser leida.
                     """
                     try:
                         fp = fp = h5py.File(filename,'r')
                     except IOError:
                         traceback.print_exc()
                         raise IOError, "The file %s can't be opened" %(filename)
                     grp = fp['Data']
                     timeAux = grp['time']
                     time0 = timeAux[:][0].astype(numpy.float)   #Time Vector
                     fp.close()
                     if self.timezone == 'lt':
                         time0 -= 5*3600
                     boolTimer = numpy.logical_and(time0 >= startSeconds,time0 < endSeconds)
                     if not (numpy.any(boolTimer)):
                         return None
                     thisDatetime = datetime.datetime.utcfromtimestamp(time0[0])
                     return thisDatetime
                 def __checkPath(self):
                     if os.path.exists(self.path):
                         self.status = 1
                     else:
                         self.status = 0
                         print 'Path:%s does not exists'%self.path
                     return
                 def __setNextFileOffline(self):
                     idFile = self.fileIndex
                     idFile += 1
                     if not(idFile < len(self.filenameList)):
                         print "No more Files"
                         return 0
                     filename = self.filenameList[idFile]
                     filePointer = h5py.File(filename,'r')
                     self.flagIsNewFile = 1
                     self.fileIndex = idFile
                     self.filename = filename
                     self.fp = filePointer
                     print "Setting the file: %s"%self.filename
                     self.__readMetadata()
                     self.__setBlockList()
             #         self.nRecords = self.fp['Data'].attrs['blocksPerFile']
                     self.nRecords = self.fp['Data'].attrs['nRecords']
                     self.blockIndex = 0
                     return 1
                 def __setBlockList(self):
                     '''
                     self.fp
                     self.startDateTime
                     self.endDateTime
                     self.blockList
                     self.blocksPerFile
                     '''
                     filePointer = self.fp
                     secStart = self.secStart
                     secEnd = self.secEnd
                     grp = filePointer['Data']
                     timeVector = grp['time'].value.astype(numpy.float)[0]
                     if self.timezone == 'lt':
                         timeVector -= 5*3600
                     ind = numpy.where(numpy.logical_and(timeVector >= secStart , timeVector < secEnd))[0]
                     self.blockList = ind
                     self.blocksPerFile = len(ind)
                     return
                 def __readMetadata(self):
                     '''
                     self.pathMeta
                     self.listShapes
                     self.listMetaname
                     self.listMeta
                     '''
                     grp = self.fp['Data']
                     pathMeta = os.path.join(self.path, grp.attrs['metadata'])
                     if pathMeta == self.pathMeta:
                         return
                     else:
                         self.pathMeta = pathMeta
                     filePointer = h5py.File(self.pathMeta,'r')
                     groupPointer = filePointer['Metadata']
                     listMetaname = []
                     listMetadata = []
                     for item in groupPointer.items():
                         name = item[0]
                         if name=='array dimensions':
                             table = groupPointer[name][:]
                             listShapes = {}
                             for shapes in table:
                                 listShapes[shapes[0]] = numpy.array([shapes[1],shapes[2],shapes[3],shapes[4]])
                         else:
                             data = groupPointer[name].value
                             listMetaname.append(name)
                             listMetadata.append(data)
                             if name=='type':
                                 self.__initDataOut(data)
                     filePointer.close()
                     self.listShapes = listShapes
                     self.listMetaname = listMetaname
                     self.listMeta = listMetadata
                     return
                 def __readData(self):
                     grp = self.fp['Data']
                     listdataname = []
                     listdata = []
                     for item in grp.items():
                         name = item[0]
                         if name == 'time':
                             listdataname.append('utctime')
                             timeAux = grp[name].value.astype(numpy.float)[0]
                             listdata.append(timeAux)
                             continue
                         listdataname.append(name)
                         array = self.__setDataArray(self.nRecords, grp[name],self.listShapes[name])
                         listdata.append(array)
                     self.listDataname = listdataname
                     self.listData = listdata
                     return
                 def __setDataArray(self, nRecords, dataset, shapes):
                     nChannels = shapes[0]    #Dimension 0
                     nPoints = shapes[1]      #Dimension 1, number of Points or Parameters
                     nSamples = shapes[2]     #Dimension 2, number of samples or ranges
                     mode = shapes[3]
             #         if nPoints>1:
             #             arrayData = numpy.zeros((nRecords,nChannels,nPoints,nSamples))
             #         else:
             #             arrayData = numpy.zeros((nRecords,nChannels,nSamples))
             #
             #         chn = 'channel'
             #
             #         for i in range(nChannels):
             #
             #             data = dataset[chn + str(i)].value
             #
             #             if nPoints>1:
             #                 data = numpy.rollaxis(data,2)
             #
             #             arrayData[:,i,:] = data
                     arrayData = numpy.zeros((nRecords,nChannels,nPoints,nSamples))
                     doSqueeze = False
                     if mode == 0:
                         strds = 'channel'
                         nDatas = nChannels
                         newShapes = (nRecords,nPoints,nSamples)
                         if nPoints == 1:
                             doSqueeze = True
                             axisSqueeze = 2
                     else:
                         strds = 'param'
                         nDatas = nPoints
                         newShapes = (nRecords,nChannels,nSamples)
                         if nChannels == 1:
                             doSqueeze = True
                             axisSqueeze = 1
                     for i in range(nDatas):
                         data = dataset[strds + str(i)].value
                         data = data.reshape(newShapes)
                         if mode == 0:
                             arrayData[:,i,:,:] = data
                         else:
                             arrayData[:,:,i,:] = data
                     if doSqueeze:
                             arrayData = numpy.squeeze(arrayData, axis=axisSqueeze)
                     return arrayData
                 def __initDataOut(self, type):
             #         if type =='Parameters':
             #             self.dataOut = Parameters()
             #         elif type =='Spectra':
             #             self.dataOut = Spectra()
             #         elif type =='Voltage':
             #             self.dataOut = Voltage()
             #         elif type =='Correlation':
             #             self.dataOut = Correlation()
                     return
                 def __setDataOut(self):
                     listMeta = self.listMeta
                     listMetaname = self.listMetaname
                     listDataname = self.listDataname
                     listData = self.listData
                     blockIndex = self.blockIndex
                     blockList = self.blockList
                     for i in range(len(listMeta)):
                         setattr(self.dataOut,listMetaname[i],listMeta[i])
                     for j in range(len(listData)):
                         if listDataname[j]=='utctime':
             #                 setattr(self.dataOut,listDataname[j],listData[j][blockList[blockIndex]])
                             setattr(self.dataOut,'utctimeInit',listData[j][blockList[blockIndex]])
                             continue
                         setattr(self.dataOut,listDataname[j],listData[j][blockList[blockIndex],:])
                     return self.dataOut.data_param
                 def getData(self):
             #         if self.flagNoMoreFiles:
             #             self.dataOut.flagNoData = True
             #             print 'Process finished'
             #             return 0
             #
                     if self.blockIndex==self.blocksPerFile:
                          if not( self.__setNextFileOffline() ):
                             self.dataOut.flagNoData = True
                             return 0
             #
             #         if self.datablock == None: # setear esta condicion cuando no hayan datos por leers
             #             self.dataOut.flagNoData = True
             #             return 0
                     self.__readData()
                     self.__setDataOut()
                     self.dataOut.flagNoData = False
                     self.blockIndex += 1
                     return
                 def run(self, **kwargs):
                     if not(self.isConfig):
                         self.setup(**kwargs)
             #             self.setObjProperties()
                         self.isConfig = True
                     self.getData()
                     return
             class HDF5Writer(Operation):
                 ext = ".hdf5"
                 optchar = "D"
                 metaoptchar = "M"
                 metaFile = None
+                filename = None
                 path = None
                 setFile = None
                 fp = None
                 grp = None
                 ds = None
                 firsttime = True
                 #Configurations
                 blocksPerFile = None
                 blockIndex = None
                 dataOut = None
                 #Data Arrays
                 dataList = None
                 metadataList = None
                 arrayDim = None
                 tableDim = None
             #     dtype = [('arrayName', 'S20'),('nChannels', 'i'), ('nPoints', 'i'), ('nSamples', 'i'),('mode', 'b')]
                 dtype = [('arrayName', 'S20'),('nDimensions', 'i'), ('dim2', 'i'), ('dim1', 'i'),('dim0', 'i'),('mode', 'b')]
                 mode = None
                 nDatas = None    #Number of datasets to be stored per array
                 nDims = None  #Number Dimensions in each dataset
+                nDimsForDs = None
                 def __init__(self):
                     Operation.__init__(self)
                     self.isConfig = False
                     return
                 def setup(self, dataOut, **kwargs):
                     self.path = kwargs['path']
                     if kwargs.has_key('ext'):
                         self.ext = kwargs['ext']
                     if kwargs.has_key('blocksPerFile'):
                         self.blocksPerFile = kwargs['blocksPerFile']
                     else:
                         self.blocksPerFile = 10
                     self.metadataList = kwargs['metadataList']
                     self.dataList = kwargs['dataList']
                     self.dataOut = dataOut
                     if kwargs.has_key('mode'):
                         mode = kwargs['mode']
                         if type(mode) == int:
                             mode = numpy.zeros(len(self.dataList)) + mode
                     else:
                         mode = numpy.zeros(len(self.dataList))
                     self.mode = mode
                     arrayDim = numpy.zeros((len(self.dataList),5))
                     #Table dimensions
                     dtype0 = self.dtype
                     tableList = []
                     for i in range(len(self.dataList)):
                         dataAux = getattr(self.dataOut, self.dataList[i])
                         if type(dataAux)==float or type(dataAux)==int:
                             arrayDim[i,0] = 1
                         else:
                             arrayDim0 = dataAux.shape
                             arrayDim[i,0] = len(arrayDim0)
                             arrayDim[i,4] = mode[i]
                             if len(arrayDim0) == 3:
                                 arrayDim[i,1:-1] = numpy.array(arrayDim0)
                             elif len(arrayDim0) == 2:
                                 arrayDim[i,2:-1] = numpy.array(arrayDim0) #nHeights
                             elif len(arrayDim0) == 1:
                                 arrayDim[i,3] = arrayDim0
                             elif len(arrayDim0) == 0:
                                 arrayDim[i,0] = 1
                                 arrayDim[i,3] = 1
                         table = numpy.array((self.dataList[i],) + tuple(arrayDim[i,:]),dtype = dtype0)
                         tableList.append(table)
                     self.arrayDim = arrayDim
                     self.tableDim = numpy.array(tableList, dtype = dtype0)
                     self.blockIndex = 0
                     return
                 def putMetadata(self):
                     fp = self.createMetadataFile()
                     self.writeMetadata(fp)
                     fp.close()
                     return
                 def createMetadataFile(self):
                     ext = self.ext
                     path = self.path
                     setFile = self.setFile
                     timeTuple = time.localtime(self.dataOut.utctime)
                     subfolder = ''
+                    fullpath = os.path.join( path, subfolder )
+                    if not( os.path.exists(fullpath) ):
+                        os.mkdir(fullpath)
+                        setFile = -1 #inicializo mi contador de seteo
+                    subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
                     fullpath = os.path.join( path, subfolder )
                     if not( os.path.exists(fullpath) ):
                         os.mkdir(fullpath)
                         setFile = -1 #inicializo mi contador de seteo
                     else:
                         filesList = os.listdir( fullpath )
-                        if len( filesList ) > 0:
                         filesList = sorted( filesList, key=str.lower )
+                        if len( filesList ) > 0:
+                            filesList = [k for k in filesList if 'M' in k]
                             filen = filesList[-1]
                             # el filename debera tener el siguiente formato
                             # 0 1234 567 89A BCDE (hex)
                             # x YYYY DDD SSS .ext
                             if isNumber( filen[8:11] ):
                                 setFile = int( filen[8:11] ) #inicializo mi contador de seteo al seteo del ultimo file
                             else:
                                 setFile = -1
                         else:
                             setFile = -1 #inicializo mi contador de seteo
                     setFile += 1
                     file = '%s%4.4d%3.3d%3.3d%s' % (self.metaoptchar,
                                                     timeTuple.tm_year,
                                                     timeTuple.tm_yday,
                                                     setFile,
                                                     ext )
                     filename = os.path.join( path, subfolder, file )
                     self.metaFile = file
                     #Setting HDF5 File
                     fp = h5py.File(filename,'w')
                     return fp
                 def writeMetadata(self, fp):
                     grp = fp.create_group("Metadata")
                     grp.create_dataset('array dimensions', data = self.tableDim, dtype = self.dtype)
                     for i in range(len(self.metadataList)):
                         grp.create_dataset(self.metadataList[i], data=getattr(self.dataOut, self.metadataList[i]))
                     return
                 def setNextFile(self):
                     ext = self.ext
                     path = self.path
                     setFile = self.setFile
                     mode = self.mode
-                    if self.fp != None:
-                        self.fp.close()
                     timeTuple = time.localtime(self.dataOut.utctime)
                     subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
                     fullpath = os.path.join( path, subfolder )
-                    if not( os.path.exists(fullpath) ):
-                        os.mkdir(fullpath)
+                    if os.path.exists(fullpath):
-                        setFile = -1 #inicializo mi contador de seteo
-                    else:
                         filesList = os.listdir( fullpath )
+                        filesList = [k for k in filesList if 'D' in k]
                         if len( filesList ) > 0:
                             filesList = sorted( filesList, key=str.lower )
                             filen = filesList[-1]
                             # el filename debera tener el siguiente formato
                             # 0 1234 567 89A BCDE (hex)
                             # x YYYY DDD SSS .ext
                             if isNumber( filen[8:11] ):
                                 setFile = int( filen[8:11] ) #inicializo mi contador de seteo al seteo del ultimo file
                             else:
                                 setFile = -1
                         else:
                             setFile = -1 #inicializo mi contador de seteo
                     setFile += 1
                     file = '%s%4.4d%3.3d%3.3d%s' % (self.optchar,
                                                     timeTuple.tm_year,
                                                     timeTuple.tm_yday,
                                                     setFile,
                                                     ext )
                     filename = os.path.join( path, subfolder, file )
                     #Setting HDF5 File
                     fp = h5py.File(filename,'w')
                     grp = fp.create_group("Data")
                     grp.attrs['metadata'] = self.metaFile
             #         grp.attrs['blocksPerFile'] = 0
                     ds = []
                     data = []
+                    nDimsForDs = []
                     nDatas = numpy.zeros(len(self.dataList))
                     nDims = self.arrayDim[:,0]
+                    nDim1 = self.arrayDim[:,2]
+                    nDim0 = self.arrayDim[:,3]
                     for i in range(len(self.dataList)):
                         if nDims[i]==1:
-                            ds0 = grp.create_dataset(self.dataList[i], (1,1), maxshape=(1,None) , chunks = True, dtype='S20')
+            #                 ds0 = grp.create_dataset(self.dataList[i], (1,1), maxshape=(1,self.blocksPerFile) , chunks = True, dtype='S20')
+                            ds0 = grp.create_dataset(self.dataList[i], (1,1), maxshape=(1,self.blocksPerFile) , chunks = True, dtype=numpy.float64)
                             ds.append(ds0)
                             data.append([])
+                            nDimsForDs.append(nDims[i])
                         else:
                             if mode[i]==0:
                                 strMode = "channel"
                                 nDatas[i] = self.arrayDim[i,1]
                             else:
                                 strMode = "param"
                                 nDatas[i] = self.arrayDim[i,2]
                             if nDims[i]==2:
                                 nDatas[i] = self.arrayDim[i,2]
                             grp0 = grp.create_group(self.dataList[i])
                             for j in range(int(nDatas[i])):
                                 tableName = strMode + str(j)
                                 if nDims[i] == 3:
-                                    ds0 = grp0.create_dataset(tableName, (1,1,1) , maxshape=(None,None,None), chunks=True)
+                                    ds0 = grp0.create_dataset(tableName, (nDim1[i],nDim0[i],1) , data = numpy.zeros((nDim1[i],nDim0[i],1)) ,maxshape=(None,nDim0[i],None), chunks=True)
                                 else:
-                                    ds0 = grp0.create_dataset(tableName, (1,1) , maxshape=(None,None), chunks=True)
+                                    ds0 = grp0.create_dataset(tableName, (1,nDim0[i]), data = numpy.zeros((1,nDim0[i])) , maxshape=(None,nDim0[i]), chunks=True)
                                 ds.append(ds0)
                                 data.append([])
+                                nDimsForDs.append(nDims[i])
                     self.nDatas = nDatas
                     self.nDims = nDims
+                    self.nDimsForDs = nDimsForDs
                     #Saving variables
                     print 'Writing the file: %s'%filename
+                    self.filename = filename
                     self.fp = fp
                     self.grp = grp
+                    self.grp.attrs.modify('nRecords', 1)
                     self.ds = ds
                     self.data = data
                     self.setFile = setFile
                     self.firsttime = True
                     self.blockIndex = 0
                     return
                 def putData(self):
-                    self.setBlock()
-                    self.writeBlock()
+                    if not self.firsttime:
+                        self.fp.flush()
+                        self.fp.close()
+                        self.readBlock()
                     if self.blockIndex == self.blocksPerFile:
                         self.setNextFile()
+                    self.setBlock()
+                    self.writeBlock()
                     return
+                def readBlock(self):
+                    '''
+                    data Array configured
+                    self.data
+                    '''
+                    ds = self.ds
+                            #Setting HDF5 File
+                    fp = h5py.File(self.filename,'r+')
+                    grp = fp["Data"]
+                    ind = 0
+            #         grp.attrs['blocksPerFile'] = 0
+                    for i in range(len(self.dataList)):
+                        if self.nDims[i]==1:
+                            ds0 = grp[self.dataList[i]]
+                            ds[ind] = ds0
+                            ind += 1
+                        else:
+                            if self.mode[i]==0:
+                                strMode = "channel"
+                            else:
+                                strMode = "param"
+                            grp0 = grp[self.dataList[i]]
+                            for j in range(int(self.nDatas[i])):
+                                tableName = strMode + str(j)
+                                ds0 = grp0[tableName]
+                                ds[ind] = ds0
+                                ind += 1
+                    self.fp = fp
+                    self.grp = grp
+                    self.ds = ds
+                    return
                 def setBlock(self):
                     '''
                     data Array configured
                     self.data
                     '''
                     #Creating Arrays
                     data = self.data
                     nDatas = self.nDatas
                     nDims = self.nDims
                     mode = self.mode
                     ind = 0
                     for i in range(len(self.dataList)):
                         dataAux = getattr(self.dataOut,self.dataList[i])
                         if nDims[i] == 1:
-                            data[ind] = numpy.array([str(dataAux)]).reshape((1,1))
+            #                 data[ind] = numpy.array([str(dataAux)]).reshape((1,1))
-                            if not self.firsttime:
+                            data[ind] = dataAux
-                                data[ind] = numpy.hstack((self.ds[ind][:], self.data[ind]))
+            #                 if not self.firsttime:
+            #                     data[ind] = numpy.hstack((self.ds[ind][:], self.data[ind]))
                             ind += 1
                         else:
                             for j in range(int(nDatas[i])):
                                 if (mode[i] == 0) or (nDims[i] == 2):       #In case division per channel or Dimensions is only 1
                                     data[ind] = dataAux[j,:]
                                 else:
                                     data[ind] = dataAux[:,j,:]
-                                if nDims[i] == 3:
+            #                     if nDims[i] == 3:
-                                    data[ind] = data[ind].reshape((data[ind].shape[0],data[ind].shape[1],1))
+            #                         data[ind] = data[ind].reshape((data[ind].shape[0],data[ind].shape[1],1))
-                                    if not self.firsttime:
+            #                         if not self.firsttime:
-                                        data[ind] = numpy.dstack((self.ds[ind][:], data[ind]))
+            #                             data[ind] = numpy.dstack((self.ds[ind][:], data[ind]))
-                                else:
+            #                     else:
-                                    data[ind] = data[ind].reshape((1,data[ind].shape[0]))
+            #                         data[ind] = data[ind].reshape((1,data[ind].shape[0]))
-                                    if not self.firsttime:
+            #                         if not self.firsttime:
-                                        data[ind] = numpy.vstack((self.ds[ind][:], data[ind]))
+            #                             data[ind] = numpy.vstack((self.ds[ind][:], data[ind]))
                                 ind += 1
                     self.data = data
                     return
                 def writeBlock(self):
                     '''
                     Saves the block in the HDF5 file
                     '''
                     for i in range(len(self.ds)):
+                        if self.firsttime:
+            #                 self.ds[i].resize(self.data[i].shape)
+            #                 self.ds[i][self.blockIndex,:] = self.data[i]
+                            if type(self.data[i]) == numpy.ndarray:
+                                nDims1 = len(self.ds[i].shape)
+                                if nDims1 == 3:
+                                    self.data[i] = self.data[i].reshape((self.data[i].shape[0],self.data[i].shape[1],1))
                                     self.ds[i].resize(self.data[i].shape)
                             self.ds[i][:] = self.data[i]
+                        else:
+                            if self.nDimsForDs[i] == 1:
+                                self.ds[i].resize((self.ds[i].shape[0], self.ds[i].shape[1] + 1))
+                                self.ds[i][0,-1] = self.data[i]
+                            elif self.nDimsForDs[i] == 2:
+                                self.ds[i].resize((self.ds[i].shape[0] + 1,self.ds[i].shape[1]))
+                                self.ds[i][self.blockIndex,:] = self.data[i]
+                            elif self.nDimsForDs[i] == 3:
+                                dataShape = self.data[i].shape
+                                dsShape = self.ds[i].shape
+                                if dataShape[0]==dsShape[0]:
+                                    self.ds[i].resize((self.ds[i].shape[0],self.ds[i].shape[1],self.ds[i].shape[2]+1))
+                                    self.ds[i][:,:,-1] = self.data[i]
+                                else:
+                                    self.ds[i].resize((self.ds[i].shape[0] + dataShape[0],self.ds[i].shape[1],self.ds[i].shape[2]))
+                                    self.ds[i][dsShape[0]:,:,0] = self.data[i]
+            #                 self.ds[i].append(self.data[i])
+            #                 self.fp.flush()
+            #         if not self.firsttime:
+            #             self.fp.root.Data._v_attrs.nRecords = self.blockIndex
-                    self.blockIndex += 1
+            #         if self.firsttime:
+            #             self.fp.close()
-                    self.grp.attrs.modify('nRecords', self.blockIndex)
+            #             self.readBlock2()
+                    self.blockIndex += 1
                     self.firsttime = False
                     return
                 def run(self, dataOut, **kwargs):
                     if not(self.isConfig):
                         self.setup(dataOut, **kwargs)
                         self.isConfig = True
                         self.putMetadata()
                         self.setNextFile()
                     self.putData()
                     return

schainpy/model/proc/jroproc_parameters.py +407 -49

             import numpy
             import math
             from scipy import optimize
             from scipy import interpolate
             from scipy import signal
             from scipy import stats
             import re
             import datetime
             import copy
             import sys
             import importlib
             import itertools
             from jroproc_base import ProcessingUnit, Operation
             from schainpy.model.data.jrodata import Parameters
             class ParametersProc(ProcessingUnit):
                 nSeconds = None
                 def __init__(self):
                     ProcessingUnit.__init__(self)
             #         self.objectDict = {}
                     self.buffer = None
                     self.firstdatatime = None
                     self.profIndex = 0
                     self.dataOut = Parameters()
                 def __updateObjFromInput(self):
                     self.dataOut.inputUnit = self.dataIn.type
                     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.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.nHeights = self.dataIn.nHeights
             #         self.dataOut.nChannels = self.dataIn.nChannels
                     self.dataOut.nBaud = self.dataIn.nBaud
                     self.dataOut.nCode = self.dataIn.nCode
                     self.dataOut.code = self.dataIn.code
             #        self.dataOut.nProfiles = self.dataOut.nFFTPoints
                     self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
                     self.dataOut.utctime = self.firstdatatime
                     self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
                     self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
             #        self.dataOut.nCohInt = self.dataIn.nCohInt
             #        self.dataOut.nIncohInt = 1
                     self.dataOut.ippSeconds = self.dataIn.ippSeconds
             #        self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
                     self.dataOut.timeInterval = self.dataIn.timeInterval
                     self.dataOut.heightList = self.dataIn.getHeiRange()
                     self.dataOut.frequency = self.dataIn.frequency
                 def run(self, nSeconds = None, nProfiles = None):
                     if self.firstdatatime == None:
                         self.firstdatatime = self.dataIn.utctime
                     #----------------------    Voltage Data    ---------------------------
                     if self.dataIn.type == "Voltage":
                         self.dataOut.flagNoData = True
                         if nSeconds != None:
                             self.nSeconds = nSeconds
                             self.nProfiles= int(numpy.floor(nSeconds/(self.dataIn.ippSeconds*self.dataIn.nCohInt)))
                         if self.buffer == None:
                                 self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                            self.nProfiles,
                                                            self.dataIn.nHeights),
                                                            dtype='complex')
                         self.buffer[:,self.profIndex,:] = self.dataIn.data.copy()
                         self.profIndex += 1
                         if self.profIndex == self.nProfiles:
                             self.__updateObjFromInput()
                             self.dataOut.data_pre = self.buffer.copy()
                             self.dataOut.paramInterval = nSeconds
                             self.dataOut.flagNoData = False
                             self.buffer = None
                             self.firstdatatime = None
                             self.profIndex = 0
                             return
                     #----------------------    Spectra Data    ---------------------------
                     if self.dataIn.type == "Spectra":
                         self.dataOut.data_pre = self.dataIn.data_spc.copy()
                         self.dataOut.abscissaList = self.dataIn.getVelRange(1)
                         self.dataOut.noise = self.dataIn.getNoise()
                         self.dataOut.normFactor = self.dataIn.normFactor
                         self.dataOut.groupList = self.dataIn.pairsList
                         self.dataOut.flagNoData = False
                     #----------------------    Correlation Data    ---------------------------
                     if self.dataIn.type == "Correlation":
                         lagRRange = self.dataIn.lagR
                         indR = numpy.where(lagRRange == 0)[0][0]
                         self.dataOut.data_pre = self.dataIn.data_corr.copy()[:,:,indR,:]
                         self.dataOut.abscissaList = self.dataIn.getLagTRange(1)
                         self.dataOut.noise = self.dataIn.noise
                         self.dataOut.normFactor = self.dataIn.normFactor
                         self.dataOut.data_SNR = self.dataIn.SNR
                         self.dataOut.groupList = self.dataIn.pairsList
                         self.dataOut.flagNoData = False
                             #----------------------    Correlation Data    ---------------------------
                     if self.dataIn.type == "Parameters":
                         self.dataOut.copy(self.dataIn)
                         self.dataOut.flagNoData = False
                         return True
                     self.__updateObjFromInput()
                     self.firstdatatime = None
                     self.dataOut.utctimeInit = self.dataIn.utctime
                     self.dataOut.outputInterval = self.dataIn.timeInterval
                 #-------------------    Get Moments    ----------------------------------
                 def GetMoments(self, channelList = None):
                     '''
                     Function GetMoments()
                     Input:
                         channelList    :    simple channel list to select e.g. [2,3,7]
                         self.dataOut.data_pre
                         self.dataOut.abscissaList
                         self.dataOut.noise
                     Affected:
                         self.dataOut.data_param
                         self.dataOut.data_SNR
                     '''
                     data = self.dataOut.data_pre
                     absc = self.dataOut.abscissaList[:-1]
                     noise = self.dataOut.noise
                     data_param = numpy.zeros((data.shape[0], 4, data.shape[2]))
                     if channelList== None:
                         channelList = self.dataIn.channelList
                     self.dataOut.channelList = channelList
                     for ind in channelList:
                         data_param[ind,:,:] = self.__calculateMoments(data[ind,:,:], absc, noise[ind])
                     self.dataOut.data_param = data_param[:,1:,:]
                     self.dataOut.data_SNR = data_param[:,0]
                     return
                 def __calculateMoments(self, oldspec, oldfreq, n0, nicoh = None, graph = None, smooth = None, type1 = None, fwindow = None, snrth = None, dc = None, aliasing = None, oldfd = None, wwauto = None):
                     if (nicoh == None): nicoh = 1
                     if (graph == None): graph = 0
                     if (smooth == None): smooth = 0
                     elif (self.smooth < 3): smooth = 0
                     if (type1 == None): type1 = 0
                     if (fwindow == None): fwindow = numpy.zeros(oldfreq.size) + 1
                     if (snrth == None): snrth = -3
                     if (dc == None): dc = 0
                     if (aliasing == None): aliasing = 0
                     if (oldfd == None): oldfd = 0
                     if (wwauto == None): wwauto = 0
                     if (n0 < 1.e-20):   n0 = 1.e-20
                     freq = oldfreq
                     vec_power = numpy.zeros(oldspec.shape[1])
                     vec_fd = numpy.zeros(oldspec.shape[1])
                     vec_w = numpy.zeros(oldspec.shape[1])
                     vec_snr = numpy.zeros(oldspec.shape[1])
                     for ind in range(oldspec.shape[1]):
                         spec = oldspec[:,ind]
                         aux = spec*fwindow
                         max_spec = aux.max()
                         m = list(aux).index(max_spec)
                         #Smooth
                         if (smooth == 0):   spec2 = spec
                         else:   spec2 = scipy.ndimage.filters.uniform_filter1d(spec,size=smooth)
                         #    Calculo de Momentos
                         bb = spec2[range(m,spec2.size)]
                         bb = (bb<n0).nonzero()
                         bb = bb[0]
                         ss = spec2[range(0,m + 1)]
                         ss = (ss<n0).nonzero()
                         ss = ss[0]
                         if (bb.size == 0):
                             bb0 = spec.size - 1 - m
                         else:
                             bb0 = bb[0] - 1
                             if (bb0 < 0):
                                 bb0 = 0
                         if (ss.size == 0):   ss1 = 1
                         else: ss1 = max(ss) + 1
                         if (ss1 > m):   ss1 = m
                         valid = numpy.asarray(range(int(m + bb0 - ss1 + 1))) + ss1
                         power = ((spec2[valid] - n0)*fwindow[valid]).sum()
                         fd = ((spec2[valid]- n0)*freq[valid]*fwindow[valid]).sum()/power
                         w = math.sqrt(((spec2[valid] - n0)*fwindow[valid]*(freq[valid]- fd)**2).sum()/power)
                         snr = (spec2.mean()-n0)/n0
                         if (snr < 1.e-20) :
                             snr = 1.e-20
                         vec_power[ind] = power
                         vec_fd[ind] = fd
                         vec_w[ind] = w
                         vec_snr[ind] = snr
                     moments = numpy.vstack((vec_snr, vec_power, vec_fd, vec_w))
                     return moments
                 #------------------    Get SA Parameters    --------------------------
                 def GetSAParameters(self):
-                    data = self.dataOut.data_pre
+                    pairslist = self.dataOut.groupList
-                    crossdata = self.dataIn.data_cspc
+                    num_pairs = len(pairslist)
-                    a = 1
+                    vel = self.dataOut.abscissaList
+                    spectra = self.dataOut.data_pre
+                    cspectra = self.dataIn.data_cspc
+                    delta_v = vel[1] - vel[0]
+                    #Calculating the power spectrum
+                    spc_pow = numpy.sum(spectra, 3)*delta_v
+                    #Normalizing Spectra
+                    norm_spectra = spectra/spc_pow
+                    #Calculating the norm_spectra at peak
+                    max_spectra = numpy.max(norm_spectra, 3)
+                    #Normalizing Cross Spectra
+                    norm_cspectra = numpy.zeros(cspectra.shape)
+                    for i in range(num_chan):
+                        norm_cspectra[i,:,:] = cspectra[i,:,:]/numpy.sqrt(spc_pow[pairslist[i][0],:]*spc_pow[pairslist[i][1],:])
+                    max_cspectra = numpy.max(norm_cspectra,2)
+                    max_cspectra_index = numpy.argmax(norm_cspectra, 2)
+                    for i in range(num_pairs):
+                        cspc_par[i,:,:] = __calculateMoments(norm_cspectra)
                 #-------------------    Get Lags    ----------------------------------
                 def GetLags(self):
                     '''
                     Function GetMoments()
                     Input:
                         self.dataOut.data_pre
                         self.dataOut.abscissaList
                         self.dataOut.noise
                         self.dataOut.normFactor
                         self.dataOut.data_SNR
                         self.dataOut.groupList
                         self.dataOut.nChannels
                     Affected:
                         self.dataOut.data_param
                     '''
                     data = self.dataOut.data_pre
                     normFactor = self.dataOut.normFactor
                     nHeights = self.dataOut.nHeights
                     absc = self.dataOut.abscissaList[:-1]
                     noise = self.dataOut.noise
                     SNR = self.dataOut.data_SNR
                     pairsList = self.dataOut.groupList
                     nChannels = self.dataOut.nChannels
                     pairsAutoCorr, pairsCrossCorr = self.__getPairsAutoCorr(pairsList, nChannels)
                     self.dataOut.data_param = numpy.zeros((len(pairsCrossCorr)*2 + 1, nHeights))
                     dataNorm = numpy.abs(data)
                     for l in range(len(pairsList)):
                         dataNorm[l,:,:] = dataNorm[l,:,:]/normFactor[l,:]
                     self.dataOut.data_param[:-1,:] = self.__calculateTaus(dataNorm, pairsCrossCorr, pairsAutoCorr, absc)
                     self.dataOut.data_param[-1,:] = self.__calculateLag1Phase(data, pairsAutoCorr, absc)
                     return
                 def __getPairsAutoCorr(self, pairsList, nChannels):
                     pairsAutoCorr = numpy.zeros(nChannels, dtype = 'int')*numpy.nan
                     for l in range(len(pairsList)):
                         firstChannel = pairsList[l][0]
                         secondChannel = pairsList[l][1]
                         #Obteniendo pares de Autocorrelacion
                         if firstChannel == secondChannel:
                             pairsAutoCorr[firstChannel] = int(l)
                     pairsAutoCorr = pairsAutoCorr.astype(int)
                     pairsCrossCorr = range(len(pairsList))
                     pairsCrossCorr = numpy.delete(pairsCrossCorr,pairsAutoCorr)
                     return pairsAutoCorr, pairsCrossCorr
                 def __calculateTaus(self, data, pairsCrossCorr, pairsAutoCorr, lagTRange):
                     Pt0 = data.shape[1]/2
                     #Funcion de Autocorrelacion
                     dataAutoCorr = stats.nanmean(data[pairsAutoCorr,:,:], axis = 0)
                     #Obtencion Indice de TauCross
                     indCross = data[pairsCrossCorr,:,:].argmax(axis = 1)
                     #Obtencion Indice de TauAuto
                     indAuto = numpy.zeros(indCross.shape,dtype = 'int')
                     CCValue = data[pairsCrossCorr,Pt0,:]
                     for i in range(pairsCrossCorr.size):
                         indAuto[i,:] = numpy.abs(dataAutoCorr - CCValue[i,:]).argmin(axis = 0)
                     #Obtencion de TauCross y TauAuto
                     tauCross = lagTRange[indCross]
                     tauAuto  = lagTRange[indAuto]
                     Nan1, Nan2 = numpy.where(tauCross == lagTRange[0])
                     tauCross[Nan1,Nan2] = numpy.nan
                     tauAuto[Nan1,Nan2] = numpy.nan
                     tau = numpy.vstack((tauCross,tauAuto))
                     return tau
                 def __calculateLag1Phase(self, data, pairs, lagTRange):
                     data1 = stats.nanmean(data[pairs,:,:], axis = 0)
                     lag1 = numpy.where(lagTRange == 0)[0][0] + 1
                     phase = numpy.angle(data1[lag1,:])
                     return phase
                     #-------------------    Detect Meteors  ------------------------------
-                def DetectMeteors(self, hei_ref = None, tauindex = 0,
+                def MeteorDetection(self, hei_ref = None, tauindex = 0,
-                                  predefinedPhaseShifts = None, centerReceiverIndex = 2,
+                                  predefinedPhaseShifts = None, centerReceiverIndex = 2, saveAll = False,
                                   cohDetection = False, cohDet_timeStep = 1, cohDet_thresh = 25,
                                   noise_timeStep = 4, noise_multiple = 4,
                                   multDet_timeLimit = 1, multDet_rangeLimit = 3,
                                   phaseThresh = 20, SNRThresh = 8,
                                   hmin = 70, hmax=110, azimuth = 0) :
                     '''
                     Function DetectMeteors()
                         Project developed with paper:
                         HOLDSWORTH ET AL. 2004
                     Input:
                         self.dataOut.data_pre
                         centerReceiverIndex:      From the channels, which is the center receiver
                         hei_ref:                  Height reference for the Beacon signal extraction
                         tauindex:
                         predefinedPhaseShifts:    Predefined phase offset for the voltge signals
                         cohDetection:             Whether to user Coherent detection or not
                         cohDet_timeStep:          Coherent Detection calculation time step
                         cohDet_thresh:            Coherent Detection phase threshold to correct phases
                         noise_timeStep:           Noise calculation time step
                         noise_multiple:           Noise multiple to define signal threshold
                         multDet_timeLimit:        Multiple Detection Removal time limit in seconds
                         multDet_rangeLimit:       Multiple Detection Removal range limit in km
                         phaseThresh:              Maximum phase difference between receiver to be consider a meteor
                         SNRThresh:                Minimum SNR threshold of the meteor signal to be consider a meteor
                         hmin:                     Minimum Height of the meteor to use it in the further wind estimations
                         hmax:                     Maximum Height of the meteor to use it in the further wind estimations
                         azimuth:                  Azimuth angle correction
                     Affected:
                         self.dataOut.data_param
                     Rejection Criteria (Errors):
 : No error; analysis OK
 : SNR < SNR threshold
 : angle of arrival (AOA) ambiguously determined
 : AOA estimate not feasible
 : Large difference in AOAs obtained from different antenna baselines
 : echo at start or end of time series
 : echo less than 5 examples long; too short for analysis
 : echo rise exceeds 0.3s
 : echo decay time less than twice rise time
 : large power level before echo
 : large power level after echo
 : poor fit to amplitude for estimation of decay time
 : poor fit to CCF phase variation for estimation of radial drift velocity
 : height unresolvable echo: not valid height within 70 to 110 km
 : height ambiguous echo: more then one possible height within 70 to 110 km
 : radial drift velocity or projected horizontal velocity exceeds 200 m/s
 : oscilatory echo, indicating event most likely not an underdense echo
 : phase difference in meteor Reestimation
                     Data Storage:
                         Meteors for Wind Estimation   (8):
                         Day    Hour    |    Range    Height
                         Azimuth    Zenith    errorCosDir
                         VelRad    errorVelRad
                         TypeError
                      '''
                     #Get Beacon signal
                     newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
                     if hei_ref != None:
                         newheis = numpy.where(self.dataOut.heightList>hei_ref)
                     heiRang = self.dataOut.getHeiRange()
                     #Pairs List
                     pairslist = []
                     nChannel = self.dataOut.nChannels
                     for i in range(nChannel):
                         if i != centerReceiverIndex:
                             pairslist.append((centerReceiverIndex,i))
                     #****************REMOVING HARDWARE PHASE DIFFERENCES***************
                     # see if the user put in pre defined phase shifts
                     voltsPShift = self.dataOut.data_pre.copy()
                     if predefinedPhaseShifts != None:
                         hardwarePhaseShifts = numpy.array(predefinedPhaseShifts)*numpy.pi/180
-                    else:
+                    elif beaconPhaseShifts:
                         #get hardware phase shifts using beacon signal
                         hardwarePhaseShifts = self.__getHardwarePhaseDiff(self.dataOut.data_pre, pairslist, newheis, 10)
                         hardwarePhaseShifts = numpy.insert(hardwarePhaseShifts,centerReceiverIndex,0)
+                    else:
+                        hardwarePhaseShifts = numpy.zeros(5)
                     voltsPShift = numpy.zeros((self.dataOut.data_pre.shape[0],self.dataOut.data_pre.shape[1],self.dataOut.data_pre.shape[2]), dtype = 'complex')
                     for i in range(self.dataOut.data_pre.shape[0]):
                         voltsPShift[i,:,:] = self.__shiftPhase(self.dataOut.data_pre[i,:,:], hardwarePhaseShifts[i])
                     #******************END OF REMOVING HARDWARE PHASE DIFFERENCES*********
                     #Remove DC
                     voltsDC = numpy.mean(voltsPShift,1)
                     voltsDC = numpy.mean(voltsDC,1)
                     for i in range(voltsDC.shape[0]):
                         voltsPShift[i] = voltsPShift[i] - voltsDC[i]
                     #Don't considerate last heights, theyre used to calculate Hardware Phase Shift
                     voltsPShift = voltsPShift[:,:,:newheis[0][0]]
                     #************ FIND POWER OF DATA W/COH OR NON COH DETECTION (3.4) **********
                     #Coherent Detection
                     if cohDetection:
                         #use coherent detection to get the net power
                         cohDet_thresh = cohDet_thresh*numpy.pi/180
                         voltsPShift = self.__coherentDetection(voltsPShift, cohDet_timeStep, self.dataOut.timeInterval, pairslist, cohDet_thresh)
                     #Non-coherent detection!
                     powerNet = numpy.nansum(numpy.abs(voltsPShift[:,:,:])**2,0)
                     #********** END OF COH/NON-COH POWER CALCULATION**********************
                     #********** FIND THE NOISE LEVEL AND POSSIBLE METEORS ****************
                     #Get noise
                     noise, noise1 = self.__getNoise(powerNet, noise_timeStep, self.dataOut.timeInterval)
             #         noise = self.getNoise1(powerNet, noise_timeStep, self.dataOut.timeInterval)
                     #Get signal threshold
                     signalThresh = noise_multiple*noise
                     #Meteor echoes detection
                     listMeteors = self.__findMeteors(powerNet, signalThresh)
                     #******* END OF NOISE LEVEL AND POSSIBLE METEORS CACULATION **********
                     #************** REMOVE MULTIPLE DETECTIONS (3.5) ***************************
                     #Parameters
                     heiRange = self.dataOut.getHeiRange()
                     rangeInterval = heiRange[1] - heiRange[0]
                     rangeLimit = multDet_rangeLimit/rangeInterval
                     timeLimit = multDet_timeLimit/self.dataOut.timeInterval
                     #Multiple detection removals
                     listMeteors1 = self.__removeMultipleDetections(listMeteors, rangeLimit, timeLimit)
                     #************ END OF REMOVE MULTIPLE DETECTIONS **********************
                     #*********************     METEOR REESTIMATION  (3.7, 3.8, 3.9, 3.10)   ********************
                     #Parameters
                     phaseThresh = phaseThresh*numpy.pi/180
                     thresh = [phaseThresh, noise_multiple, SNRThresh]
                     #Meteor reestimation  (Errors N 1, 6, 12, 17)
                     listMeteors2, listMeteorsPower, listMeteorsVolts = self.__meteorReestimation(listMeteors1, voltsPShift, pairslist, thresh, noise, self.dataOut.timeInterval, self.dataOut.frequency)
             #         listMeteors2, listMeteorsPower, listMeteorsVolts = self.meteorReestimation3(listMeteors2, listMeteorsPower, listMeteorsVolts, voltsPShift, pairslist, thresh, noise)
                     #Estimation of decay times (Errors N 7, 8, 11)
                     listMeteors3 = self.__estimateDecayTime(listMeteors2, listMeteorsPower, self.dataOut.timeInterval, self.dataOut.frequency)
                     #*******************     END OF METEOR REESTIMATION    *******************
                     #*********************    METEOR PARAMETERS CALCULATION (3.11, 3.12, 3.13)    **************************
                     #Calculating Radial Velocity (Error N 15)
                     radialStdThresh = 10
                     listMeteors4 = self.__getRadialVelocity(listMeteors3, listMeteorsVolts, radialStdThresh, pairslist, self.dataOut.timeInterval)
                     if len(listMeteors4) > 0:
-                        #Setting New Array
-                        date = repr(self.dataOut.datatime)
-                        arrayMeteors4, arrayParameters = self.__setNewArrays(listMeteors4, date, heiRang)
-                        #Calculate AOA (Error N 3, 4)
-                        #JONES ET AL. 1998
-                        AOAthresh = numpy.pi/8
-                        error = arrayParameters[:,-1]
-                        phases = -arrayMeteors4[:,9:13]
                         pairsList = []
-                        pairsList.append((0,3))
+                        pairx = (0,3)
-                        pairsList.append((1,2))
+                        pairy = (1,2)
-                        arrayParameters[:,4:7], arrayParameters[:,-1] = self.__getAOA(phases, pairsList, error, AOAthresh, azimuth)
+                        pairsList.append(pairx)
+                        pairsList.append(pairy)
-                        #Calculate Heights (Error N 13 and 14)
+                        #Setting New Array
-                        error = arrayParameters[:,-1]
+                        date = repr(self.dataOut.datatime)
-                        Ranges = arrayParameters[:,2]
+                        arrayParameters = self.__setNewArrays(listMeteors4, date, heiRang)
-                        zenith = arrayParameters[:,5]
-                        arrayParameters[:,3], arrayParameters[:,-1] = self.__getHeights(Ranges, zenith, error, hmin, hmax)
+                        meteorOps = MeteorOperations()
+                        jph = numpy.array([0,0,0,0])
+                        h = (hmin,hmax)
+                        arrayParameters = meteorOps.getMeteorParams(arrayParameters, azimuth, h, pairsList, jph)
+            #             #Calculate AOA (Error N 3, 4)
+            #             #JONES ET AL. 1998
+            #             error = arrayParameters[:,-1]
+            #             AOAthresh = numpy.pi/8
+            #             phases = -arrayParameters[:,9:13]
+            #             arrayParameters[:,4:7], arrayParameters[:,-1] = meteorOps.getAOA(phases, pairsList, error, AOAthresh, azimuth)
+            #
+            #             #Calculate Heights (Error N 13 and 14)
+            #             error = arrayParameters[:,-1]
+            #             Ranges = arrayParameters[:,2]
+            #             zenith = arrayParameters[:,5]
+            #             arrayParameters[:,3], arrayParameters[:,-1] = meteorOps.getHeights(Ranges, zenith, error, hmin, hmax)
+            #             error = arrayParameters[:,-1]
                     #*********************    END OF PARAMETERS CALCULATION    **************************
-                    #***************************+     SAVE DATA IN HDF5 FORMAT    **********************
+                    #***************************+     PASS DATA TO NEXT STEP    **********************
-                        self.dataOut.data_param = arrayParameters
+                        arrayFinal = arrayParameters.reshape((1,arrayParameters.shape[0],arrayParameters.shape[1]))
+                        self.dataOut.data_param = arrayFinal
                     return
                 def __getHardwarePhaseDiff(self, voltage0, pairslist, newheis, n):
                     minIndex = min(newheis[0])
                     maxIndex = max(newheis[0])
                     voltage = voltage0[:,:,minIndex:maxIndex+1]
                     nLength = voltage.shape[1]/n
                     nMin = 0
                     nMax = 0
                     phaseOffset = numpy.zeros((len(pairslist),n))
                     for i in range(n):
                         nMax += nLength
                         phaseCCF = -numpy.angle(self.__calculateCCF(voltage[:,nMin:nMax,:], pairslist, [0]))
                         phaseCCF = numpy.mean(phaseCCF, axis = 2)
                         phaseOffset[:,i] = phaseCCF.transpose()
                         nMin = nMax
             #         phaseDiff, phaseArrival = self.estimatePhaseDifference(voltage, pairslist)
                     #Remove Outliers
                     factor = 2
                     wt = phaseOffset - signal.medfilt(phaseOffset,(1,5))
                     dw = numpy.std(wt,axis = 1)
                     dw = dw.reshape((dw.size,1))
                     ind = numpy.where(numpy.logical_or(wt>dw*factor,wt<-dw*factor))
                     phaseOffset[ind] = numpy.nan
                     phaseOffset = stats.nanmean(phaseOffset, axis=1)
                     return phaseOffset
                 def __shiftPhase(self, data, phaseShift):
                     #this will shift the phase of a complex number
                     dataShifted = numpy.abs(data) * numpy.exp((numpy.angle(data)+phaseShift)*1j)
                     return dataShifted
                 def __estimatePhaseDifference(self, array, pairslist):
                     nChannel = array.shape[0]
                     nHeights = array.shape[2]
                     numPairs = len(pairslist)
             #         phaseCCF = numpy.zeros((nChannel, 5, nHeights))
                     phaseCCF = numpy.angle(self.__calculateCCF(array, pairslist, [-2,-1,0,1,2]))
                     #Correct phases
                     derPhaseCCF = phaseCCF[:,1:,:] - phaseCCF[:,0:-1,:]
                     indDer = numpy.where(numpy.abs(derPhaseCCF) > numpy.pi)
                     if indDer[0].shape[0] > 0:
                         for i in range(indDer[0].shape[0]):
                             signo = -numpy.sign(derPhaseCCF[indDer[0][i],indDer[1][i],indDer[2][i]])
                             phaseCCF[indDer[0][i],indDer[1][i]+1:,:] += signo*2*numpy.pi
             #         for j in range(numSides):
             #             phaseCCFAux = self.calculateCCF(arrayCenter, arraySides[j,:,:], [-2,1,0,1,2])
             #             phaseCCF[j,:,:] = numpy.angle(phaseCCFAux)
             #
                     #Linear
                     phaseInt = numpy.zeros((numPairs,1))
                     angAllCCF = phaseCCF[:,[0,1,3,4],0]
                     for j in range(numPairs):
                         fit = stats.linregress([-2,-1,1,2],angAllCCF[j,:])
                         phaseInt[j] = fit[1]
                     #Phase Differences
                     phaseDiff = phaseInt - phaseCCF[:,2,:]
                     phaseArrival = phaseInt.reshape(phaseInt.size)
                     #Dealias
                     indAlias = numpy.where(phaseArrival > numpy.pi)
                     phaseArrival[indAlias] -= 2*numpy.pi
                     indAlias = numpy.where(phaseArrival < -numpy.pi)
                     phaseArrival[indAlias] += 2*numpy.pi
                     return phaseDiff, phaseArrival
                 def __coherentDetection(self, volts, timeSegment, timeInterval, pairslist, thresh):
                     #this function will run the coherent detection used in Holdworth et al. 2004 and return the net power
                     #find the phase shifts of each channel over 1 second intervals
                     #only look at ranges below the beacon signal
                     numProfPerBlock = numpy.ceil(timeSegment/timeInterval)
                     numBlocks = int(volts.shape[1]/numProfPerBlock)
                     numHeights = volts.shape[2]
                     nChannel = volts.shape[0]
                     voltsCohDet = volts.copy()
                     pairsarray = numpy.array(pairslist)
                     indSides = pairsarray[:,1]
             #         indSides = numpy.array(range(nChannel))
             #         indSides = numpy.delete(indSides, indCenter)
             #
             #         listCenter = numpy.array_split(volts[indCenter,:,:], numBlocks, 0)
                     listBlocks = numpy.array_split(volts, numBlocks, 1)
                     startInd = 0
                     endInd = 0
                     for i in range(numBlocks):
                         startInd = endInd
                         endInd = endInd + listBlocks[i].shape[1]
                         arrayBlock = listBlocks[i]
             #             arrayBlockCenter = listCenter[i]
                         #Estimate the Phase Difference
                         phaseDiff, aux = self.__estimatePhaseDifference(arrayBlock, pairslist)
                         #Phase Difference RMS
                         arrayPhaseRMS = numpy.abs(phaseDiff)
                         phaseRMSaux = numpy.sum(arrayPhaseRMS < thresh,0)
                         indPhase = numpy.where(phaseRMSaux==4)
                         #Shifting
                         if indPhase[0].shape[0] > 0:
                             for j in range(indSides.size):
                                 arrayBlock[indSides[j],:,indPhase] = self.__shiftPhase(arrayBlock[indSides[j],:,indPhase], phaseDiff[j,indPhase].transpose())
                             voltsCohDet[:,startInd:endInd,:] = arrayBlock
                     return voltsCohDet
                 def __calculateCCF(self, volts, pairslist ,laglist):
                     nHeights = volts.shape[2]
                     nPoints = volts.shape[1]
                     voltsCCF = numpy.zeros((len(pairslist), len(laglist), nHeights),dtype = 'complex')
                     for i in range(len(pairslist)):
                         volts1 = volts[pairslist[i][0]]
                         volts2 = volts[pairslist[i][1]]
                         for t in range(len(laglist)):
                             idxT = laglist[t]
                             if idxT >= 0:
                                 vStacked = numpy.vstack((volts2[idxT:,:],
                                                        numpy.zeros((idxT, nHeights),dtype='complex')))
                             else:
                                 vStacked = numpy.vstack((numpy.zeros((-idxT, nHeights),dtype='complex'),
                                                        volts2[:(nPoints + idxT),:]))
                             voltsCCF[i,t,:] = numpy.sum((numpy.conjugate(volts1)*vStacked),axis=0)
                             vStacked = None
                     return voltsCCF
                 def __getNoise(self, power, timeSegment, timeInterval):
                     numProfPerBlock = numpy.ceil(timeSegment/timeInterval)
                     numBlocks = int(power.shape[0]/numProfPerBlock)
                     numHeights = power.shape[1]
                     listPower = numpy.array_split(power, numBlocks, 0)
                     noise = numpy.zeros((power.shape[0], power.shape[1]))
                     noise1 = numpy.zeros((power.shape[0], power.shape[1]))
                     startInd = 0
                     endInd = 0
                     for i in range(numBlocks):             #split por canal
                         startInd = endInd
                         endInd = endInd + listPower[i].shape[0]
                         arrayBlock = listPower[i]
                         noiseAux = numpy.mean(arrayBlock, 0)
             #             noiseAux = numpy.median(noiseAux)
             #             noiseAux = numpy.mean(arrayBlock)
                         noise[startInd:endInd,:] = noise[startInd:endInd,:] + noiseAux
                         noiseAux1 = numpy.mean(arrayBlock)
                         noise1[startInd:endInd,:] = noise1[startInd:endInd,:] + noiseAux1
                     return noise, noise1
                 def __findMeteors(self, power, thresh):
                     nProf = power.shape[0]
                     nHeights = power.shape[1]
                     listMeteors = []
                     for i in range(nHeights):
                         powerAux = power[:,i]
                         threshAux = thresh[:,i]
                         indUPthresh = numpy.where(powerAux > threshAux)[0]
                         indDNthresh = numpy.where(powerAux <= threshAux)[0]
                         j = 0
                         while (j < indUPthresh.size - 2):
                             if (indUPthresh[j + 2] == indUPthresh[j] + 2):
                                 indDNAux = numpy.where(indDNthresh > indUPthresh[j])
                                 indDNthresh = indDNthresh[indDNAux]
                                 if (indDNthresh.size > 0):
                                     indEnd = indDNthresh[0] - 1
                                     indInit = indUPthresh[j]
                                     meteor = powerAux[indInit:indEnd + 1]
                                     indPeak = meteor.argmax() + indInit
                                     FLA = sum(numpy.conj(meteor)*numpy.hstack((meteor[1:],0)))
                                     listMeteors.append(numpy.array([i,indInit,indPeak,indEnd,FLA])) #CHEQUEAR!!!!!
                                     j = numpy.where(indUPthresh == indEnd)[0] + 1
                                 else: j+=1
                             else: j+=1
                     return listMeteors
                 def __removeMultipleDetections(self,listMeteors, rangeLimit, timeLimit):
                     arrayMeteors = numpy.asarray(listMeteors)
                     listMeteors1 = []
                     while arrayMeteors.shape[0] > 0:
                         FLAs = arrayMeteors[:,4]
                         maxFLA = FLAs.argmax()
                         listMeteors1.append(arrayMeteors[maxFLA,:])
                         MeteorInitTime = arrayMeteors[maxFLA,1]
                         MeteorEndTime = arrayMeteors[maxFLA,3]
                         MeteorHeight = arrayMeteors[maxFLA,0]
                         #Check neighborhood
                         maxHeightIndex = MeteorHeight + rangeLimit
                         minHeightIndex = MeteorHeight - rangeLimit
                         minTimeIndex = MeteorInitTime - timeLimit
                         maxTimeIndex = MeteorEndTime + timeLimit
                         #Check Heights
                         indHeight = numpy.logical_and(arrayMeteors[:,0] >= minHeightIndex, arrayMeteors[:,0] <= maxHeightIndex)
                         indTime = numpy.logical_and(arrayMeteors[:,3] >= minTimeIndex, arrayMeteors[:,1] <= maxTimeIndex)
                         indBoth = numpy.where(numpy.logical_and(indTime,indHeight))
                         arrayMeteors = numpy.delete(arrayMeteors, indBoth, axis = 0)
                     return listMeteors1
                 def __meteorReestimation(self, listMeteors, volts, pairslist, thresh, noise, timeInterval,frequency):
                     numHeights = volts.shape[2]
                     nChannel = volts.shape[0]
                     thresholdPhase = thresh[0]
                     thresholdNoise = thresh[1]
                     thresholdDB = float(thresh[2])
                     thresholdDB1 = 10**(thresholdDB/10)
                     pairsarray = numpy.array(pairslist)
                     indSides = pairsarray[:,1]
                     pairslist1 = list(pairslist)
                     pairslist1.append((0,1))
                     pairslist1.append((3,4))
                     listMeteors1 = []
                     listPowerSeries = []
                     listVoltageSeries = []
                     #volts has the war data
                     if frequency == 30e6:
                         timeLag = 45*10**-3
                     else:
                         timeLag = 15*10**-3
                     lag = numpy.ceil(timeLag/timeInterval)
                     for i in range(len(listMeteors)):
                         ######################   3.6 - 3.7 PARAMETERS REESTIMATION    #########################
                         meteorAux = numpy.zeros(16)
                         #Loading meteor Data (mHeight, mStart, mPeak, mEnd)
                         mHeight = listMeteors[i][0]
                         mStart = listMeteors[i][1]
                         mPeak = listMeteors[i][2]
                         mEnd = listMeteors[i][3]
                         #get the volt data between the start and end times of the meteor
                         meteorVolts = volts[:,mStart:mEnd+1,mHeight]
                         meteorVolts = meteorVolts.reshape(meteorVolts.shape[0], meteorVolts.shape[1], 1)
                         #3.6. Phase Difference estimation
                         phaseDiff, aux = self.__estimatePhaseDifference(meteorVolts, pairslist)
                         #3.7. Phase difference removal & meteor start, peak and end times reestimated
                         #meteorVolts0.- all Channels, all Profiles
                         meteorVolts0 = volts[:,:,mHeight]
                         meteorThresh = noise[:,mHeight]*thresholdNoise
                         meteorNoise = noise[:,mHeight]
                         meteorVolts0[indSides,:] = self.__shiftPhase(meteorVolts0[indSides,:], phaseDiff) #Phase Shifting
                         powerNet0 = numpy.nansum(numpy.abs(meteorVolts0)**2, axis = 0)  #Power
                         #Times reestimation
                         mStart1 = numpy.where(powerNet0[:mPeak] < meteorThresh[:mPeak])[0]
                         if mStart1.size > 0:
                             mStart1 = mStart1[-1] + 1
                         else:
                             mStart1 = mPeak
                         mEnd1 = numpy.where(powerNet0[mPeak:] < meteorThresh[mPeak:])[0][0] + mPeak - 1
                         mEndDecayTime1 = numpy.where(powerNet0[mPeak:] < meteorNoise[mPeak:])[0]
                         if mEndDecayTime1.size == 0:
                             mEndDecayTime1 = powerNet0.size
                         else:
                             mEndDecayTime1 = mEndDecayTime1[0] + mPeak - 1
             #             mPeak1 = meteorVolts0[mStart1:mEnd1 + 1].argmax()
                         #meteorVolts1.- all Channels, from start to end
                         meteorVolts1 = meteorVolts0[:,mStart1:mEnd1 + 1]
                         meteorVolts2 = meteorVolts0[:,mPeak + lag:mEnd1 + 1]
                         if meteorVolts2.shape[1] == 0:
                             meteorVolts2 = meteorVolts0[:,mPeak:mEnd1 + 1]
                         meteorVolts1 = meteorVolts1.reshape(meteorVolts1.shape[0], meteorVolts1.shape[1], 1)
                         meteorVolts2 = meteorVolts2.reshape(meteorVolts2.shape[0], meteorVolts2.shape[1], 1)
                         #####################    END PARAMETERS REESTIMATION    #########################
                         #####################   3.8 PHASE DIFFERENCE REESTIMATION  ########################
             #             if mEnd1 - mStart1 > 4:       #Error Number 6: echo less than 5 samples long; too short for analysis
                         if meteorVolts2.shape[1] > 0:
                             #Phase Difference re-estimation
                             phaseDiff1, phaseDiffint = self.__estimatePhaseDifference(meteorVolts2, pairslist1)   #Phase Difference Estimation
             #                 phaseDiff1, phaseDiffint = self.estimatePhaseDifference(meteorVolts2, pairslist)
                             meteorVolts2 = meteorVolts2.reshape(meteorVolts2.shape[0], meteorVolts2.shape[1])
                             phaseDiff11 = numpy.reshape(phaseDiff1, (phaseDiff1.shape[0],1))
                             meteorVolts2[indSides,:] = self.__shiftPhase(meteorVolts2[indSides,:], phaseDiff11[0:4])     #Phase Shifting
                             #Phase Difference RMS
                             phaseRMS1 = numpy.sqrt(numpy.mean(numpy.square(phaseDiff1)))
                             powerNet1 = numpy.nansum(numpy.abs(meteorVolts1[:,:])**2,0)
                             #Data from Meteor
                             mPeak1 = powerNet1.argmax() + mStart1
                             mPeakPower1 = powerNet1.max()
                             noiseAux = sum(noise[mStart1:mEnd1 + 1,mHeight])
                             mSNR1 = (sum(powerNet1)-noiseAux)/noiseAux
                             Meteor1 = numpy.array([mHeight, mStart1, mPeak1, mEnd1, mPeakPower1, mSNR1, phaseRMS1])
                             Meteor1 = numpy.hstack((Meteor1,phaseDiffint))
                             PowerSeries  = powerNet0[mStart1:mEndDecayTime1 + 1]
                             #Vectorize
                             meteorAux[0:7] = [mHeight, mStart1, mPeak1, mEnd1, mPeakPower1, mSNR1, phaseRMS1]
                             meteorAux[7:11] = phaseDiffint[0:4]
                             #Rejection Criterions
                             if phaseRMS1 > thresholdPhase:  #Error Number 17: Phase variation
                                 meteorAux[-1] = 17
                             elif mSNR1 < thresholdDB1:      #Error Number 1: SNR < threshold dB
                                 meteorAux[-1] = 1
                         else:
                             meteorAux[0:4] = [mHeight, mStart, mPeak, mEnd]
                             meteorAux[-1] = 6 #Error Number 6: echo less than 5 samples long; too short for analysis
                             PowerSeries = 0
                         listMeteors1.append(meteorAux)
                         listPowerSeries.append(PowerSeries)
                         listVoltageSeries.append(meteorVolts1)
                     return listMeteors1, listPowerSeries, listVoltageSeries
                 def __estimateDecayTime(self, listMeteors, listPower, timeInterval, frequency):
                     threshError = 10
                     #Depending if it is 30 or 50 MHz
                     if frequency == 30e6:
                         timeLag = 45*10**-3
                     else:
                         timeLag = 15*10**-3
                     lag = numpy.ceil(timeLag/timeInterval)
                     listMeteors1 = []
                     for i in range(len(listMeteors)):
                         meteorPower = listPower[i]
                         meteorAux = listMeteors[i]
                         if meteorAux[-1] == 0:
                             try:
                                 indmax = meteorPower.argmax()
                                 indlag = indmax + lag
                                 y = meteorPower[indlag:]
                                 x = numpy.arange(0, y.size)*timeLag
                                 #first guess
                                 a = y[0]
                                 tau = timeLag
                                 #exponential fit
                                 popt, pcov = optimize.curve_fit(self.__exponential_function, x, y, p0 = [a, tau])
                                 y1 = self.__exponential_function(x, *popt)
                                 #error estimation
                                 error = sum((y - y1)**2)/(numpy.var(y)*(y.size - popt.size))
                                 decayTime = popt[1]
                                 riseTime = indmax*timeInterval
                                 meteorAux[11:13] = [decayTime, error]
                                 #Table items 7, 8 and 11
                                 if (riseTime > 0.3):            #Number 7: Echo rise exceeds 0.3s
                                     meteorAux[-1] = 7
                                 elif (decayTime < 2*riseTime) : #Number 8: Echo decay time less than than twice rise time
                                     meteorAux[-1] = 8
                                 if (error > threshError):       #Number 11: Poor fit to amplitude for estimation of decay time
                                     meteorAux[-1] = 11
                             except:
                                 meteorAux[-1] = 11
                         listMeteors1.append(meteorAux)
                     return listMeteors1
                 #Exponential Function
                 def __exponential_function(self, x, a, tau):
                     y = a*numpy.exp(-x/tau)
                     return y
                 def __getRadialVelocity(self, listMeteors, listVolts, radialStdThresh, pairslist,  timeInterval):
                     pairslist1 = list(pairslist)
                     pairslist1.append((0,1))
                     pairslist1.append((3,4))
                     numPairs = len(pairslist1)
                     #Time Lag
                     timeLag = 45*10**-3
                     c = 3e8
                     lag = numpy.ceil(timeLag/timeInterval)
                     freq = 30e6
                     listMeteors1 = []
                     for i in range(len(listMeteors)):
-                        meteor = listMeteors[i]
+                        meteorAux = listMeteors[i]
-                        meteorAux = numpy.hstack((meteor[:-1], 0, 0, meteor[-1]))
+                        if meteorAux[-1] == 0:
-                        if meteor[-1] == 0:
                             mStart = listMeteors[i][1]
                             mPeak = listMeteors[i][2]
                             mLag = mPeak - mStart + lag
                             #get the volt data between the start and end times of the meteor
                             meteorVolts = listVolts[i]
                             meteorVolts = meteorVolts.reshape(meteorVolts.shape[0], meteorVolts.shape[1], 1)
                             #Get CCF
                             allCCFs = self.__calculateCCF(meteorVolts, pairslist1, [-2,-1,0,1,2])
                             #Method 2
                             slopes = numpy.zeros(numPairs)
                             time = numpy.array([-2,-1,1,2])*timeInterval
                             angAllCCF = numpy.angle(allCCFs[:,[0,1,3,4],0])
                             #Correct phases
                             derPhaseCCF = angAllCCF[:,1:] - angAllCCF[:,0:-1]
                             indDer = numpy.where(numpy.abs(derPhaseCCF) > numpy.pi)
                             if indDer[0].shape[0] > 0:
                                 for i in range(indDer[0].shape[0]):
                                     signo = -numpy.sign(derPhaseCCF[indDer[0][i],indDer[1][i]])
                                     angAllCCF[indDer[0][i],indDer[1][i]+1:] += signo*2*numpy.pi
             #                     fit = scipy.stats.linregress(numpy.array([-2,-1,1,2])*timeInterval, numpy.array([phaseLagN2s[i],phaseLagN1s[i],phaseLag1s[i],phaseLag2s[i]]))
                             for j in range(numPairs):
                                 fit = stats.linregress(time, angAllCCF[j,:])
                                 slopes[j] = fit[0]
                             #Remove Outlier
             #                 indOut = numpy.argmax(numpy.abs(slopes - numpy.mean(slopes)))
             #                 slopes = numpy.delete(slopes,indOut)
             #                 indOut = numpy.argmax(numpy.abs(slopes - numpy.mean(slopes)))
             #                 slopes = numpy.delete(slopes,indOut)
                             radialVelocity = -numpy.mean(slopes)*(0.25/numpy.pi)*(c/freq)
                             radialError = numpy.std(slopes)*(0.25/numpy.pi)*(c/freq)
                             meteorAux[-2] = radialError
                             meteorAux[-3] = radialVelocity
                             #Setting Error
                             #Number 15: Radial Drift velocity or projected horizontal velocity exceeds 200 m/s
                             if numpy.abs(radialVelocity) > 200:
                                 meteorAux[-1] = 15
                             #Number 12: Poor fit to CCF variation for estimation of radial drift velocity
                             elif radialError > radialStdThresh:
                                 meteorAux[-1] = 12
                         listMeteors1.append(meteorAux)
                     return listMeteors1
                 def __setNewArrays(self, listMeteors, date, heiRang):
                     #New arrays
                     arrayMeteors = numpy.array(listMeteors)
-                    arrayParameters = numpy.zeros((len(listMeteors),10))
+                    arrayParameters = numpy.zeros((len(listMeteors), 14))
                     #Date inclusion
                     date = re.findall(r'\((.*?)\)', date)
                     date = date[0].split(',')
                     date = map(int, date)
                     date = [date[0]*10000 + date[1]*100 + date[2], date[3]*10000 + date[4]*100 + date[5]]
                     arrayDate = numpy.tile(date, (len(listMeteors), 1))
                     #Meteor array
-                    arrayMeteors[:,0] = heiRang[arrayMeteors[:,0].astype(int)]
+            #         arrayMeteors[:,0] = heiRang[arrayMeteors[:,0].astype(int)]
-                    arrayMeteors = numpy.hstack((arrayDate, arrayMeteors))
+            #         arrayMeteors = numpy.hstack((arrayDate, arrayMeteors))
                     #Parameters Array
-                    arrayParameters[:,0:3] = arrayMeteors[:,0:3]
+                    arrayParameters[:,:2] = arrayDate #Date
-                    arrayParameters[:,-3:] = arrayMeteors[:,-3:]
+                    arrayParameters[:,2] = heiRang[arrayMeteors[:,0].astype(int)] #Range
+                    arrayParameters[:,7:9] = arrayMeteors[:,-3:-1] #Radial velocity and its error
+                    arrayParameters[:,9:13] = arrayMeteors[:,7:11]  #Phases
+                    arrayParameters[:,-1] = arrayMeteors[:,-1]  #Error
-                    return arrayMeteors, arrayParameters
+                    return arrayParameters
                 def __getAOA(self, phases, pairsList, error, AOAthresh, azimuth):
                     arrayAOA = numpy.zeros((phases.shape[0],3))
                     cosdir0, cosdir = self.__getDirectionCosines(phases, pairsList)
                     arrayAOA[:,:2] = self.__calculateAOA(cosdir, azimuth)
                     cosDirError = numpy.sum(numpy.abs(cosdir0 - cosdir), axis = 1)
                     arrayAOA[:,2] = cosDirError
                     azimuthAngle = arrayAOA[:,0]
                     zenithAngle = arrayAOA[:,1]
                     #Setting Error
                     #Number 3: AOA not fesible
                     indInvalid = numpy.where(numpy.logical_and((numpy.logical_or(numpy.isnan(zenithAngle), numpy.isnan(azimuthAngle))),error == 0))[0]
                     error[indInvalid] = 3
                     #Number 4: Large difference in AOAs obtained from different antenna baselines
                     indInvalid = numpy.where(numpy.logical_and(cosDirError > AOAthresh,error == 0))[0]
                     error[indInvalid] = 4
                     return arrayAOA, error
                 def __getDirectionCosines(self, arrayPhase, pairsList):
                     #Initializing some variables
                     ang_aux = numpy.array([-8,-7,-6,-5,-4,-3,-2,-1,0,1,2,3,4,5,6,7,8])*2*numpy.pi
                     ang_aux = ang_aux.reshape(1,ang_aux.size)
                     cosdir = numpy.zeros((arrayPhase.shape[0],2))
                     cosdir0 = numpy.zeros((arrayPhase.shape[0],2))
                     for i in range(2):
                         #First Estimation
                         phi0_aux = arrayPhase[:,pairsList[i][0]] + arrayPhase[:,pairsList[i][1]]
                         #Dealias
                         indcsi = numpy.where(phi0_aux > numpy.pi)
                         phi0_aux[indcsi] -= 2*numpy.pi
                         indcsi = numpy.where(phi0_aux < -numpy.pi)
                         phi0_aux[indcsi] += 2*numpy.pi
                         #Direction Cosine 0
                         cosdir0[:,i] = -(phi0_aux)/(2*numpy.pi*0.5)
                         #Most-Accurate Second Estimation
                         phi1_aux =  arrayPhase[:,pairsList[i][0]] - arrayPhase[:,pairsList[i][1]]
                         phi1_aux = phi1_aux.reshape(phi1_aux.size,1)
                         #Direction Cosine 1
                         cosdir1 = -(phi1_aux + ang_aux)/(2*numpy.pi*4.5)
                         #Searching the correct Direction Cosine
                         cosdir0_aux = cosdir0[:,i]
                         cosdir0_aux = cosdir0_aux.reshape(cosdir0_aux.size,1)
                         #Minimum Distance
                         cosDiff = (cosdir1 - cosdir0_aux)**2
                         indcos = cosDiff.argmin(axis = 1)
                         #Saving Value obtained
                         cosdir[:,i] = cosdir1[numpy.arange(len(indcos)),indcos]
                     return cosdir0, cosdir
                 def __calculateAOA(self, cosdir, azimuth):
                     cosdirX = cosdir[:,0]
                     cosdirY = cosdir[:,1]
                     zenithAngle = numpy.arccos(numpy.sqrt(1 - cosdirX**2 - cosdirY**2))*180/numpy.pi
                     azimuthAngle = numpy.arctan2(cosdirX,cosdirY)*180/numpy.pi + azimuth #0 deg north, 90 deg east
                     angles = numpy.vstack((azimuthAngle, zenithAngle)).transpose()
                     return angles
                 def __getHeights(self, Ranges, zenith, error, minHeight, maxHeight):
                     Ramb = 375  #Ramb = c/(2*PRF)
                     Re = 6371   #Earth Radius
                     heights = numpy.zeros(Ranges.shape)
                     R_aux = numpy.array([0,1,2])*Ramb
                     R_aux = R_aux.reshape(1,R_aux.size)
                     Ranges = Ranges.reshape(Ranges.size,1)
                     Ri = Ranges + R_aux
                     hi = numpy.sqrt(Re**2 + Ri**2 + (2*Re*numpy.cos(zenith*numpy.pi/180)*Ri.transpose()).transpose()) - Re
                     #Check if there is a height between 70 and 110 km
                     h_bool = numpy.sum(numpy.logical_and(hi > minHeight, hi < maxHeight), axis = 1)
                     ind_h = numpy.where(h_bool == 1)[0]
                     hCorr = hi[ind_h, :]
                     ind_hCorr = numpy.where(numpy.logical_and(hi > minHeight, hi < maxHeight))
                     hCorr = hi[ind_hCorr]
                     heights[ind_h] = hCorr
                     #Setting Error
                     #Number 13: Height unresolvable echo: not valid height within 70 to 110 km
                     #Number 14: Height ambiguous echo: more than one possible height within 70 to 110 km
                     indInvalid2 = numpy.where(numpy.logical_and(h_bool > 1, error == 0))[0]
                     error[indInvalid2] = 14
                     indInvalid1 = numpy.where(numpy.logical_and(h_bool == 0, error == 0))[0]
                     error[indInvalid1] = 13
                     return heights, error
                 def SpectralFitting(self, getSNR = True, path=None, file=None, groupList=None):
                     '''
                     Function GetMoments()
                     Input:
                     Output:
                     Variables modified:
                     '''
                     if path != None:
                         sys.path.append(path)
                     self.dataOut.library = importlib.import_module(file)
                     #To be inserted as a parameter
                     groupArray = numpy.array(groupList)
             #         groupArray = numpy.array([[0,1],[2,3]])
                     self.dataOut.groupList = groupArray
                     nGroups = groupArray.shape[0]
                     nChannels = self.dataIn.nChannels
                     nHeights=self.dataIn.heightList.size
                     #Parameters Array
                     self.dataOut.data_param = None
                     #Set constants
                     constants = self.dataOut.library.setConstants(self.dataIn)
                     self.dataOut.constants = constants
                     M = self.dataIn.normFactor
                     N = self.dataIn.nFFTPoints
                     ippSeconds = self.dataIn.ippSeconds
                     K = self.dataIn.nIncohInt
                     pairsArray = numpy.array(self.dataIn.pairsList)
                     #List of possible combinations
                     listComb = itertools.combinations(numpy.arange(groupArray.shape[1]),2)
                     indCross = numpy.zeros(len(list(listComb)), dtype = 'int')
                     if getSNR:
                         listChannels = groupArray.reshape((groupArray.size))
                         listChannels.sort()
                         noise = self.dataIn.getNoise()
                         self.dataOut.data_SNR = self.__getSNR(self.dataIn.data_spc[listChannels,:,:], noise[listChannels])
                     for i in range(nGroups):
                         coord = groupArray[i,:]
                         #Input data array
                         data = self.dataIn.data_spc[coord,:,:]/(M*N)
                         data = data.reshape((data.shape[0]*data.shape[1],data.shape[2]))
                         #Cross Spectra data array for Covariance Matrixes
                         ind = 0
                         for pairs in listComb:
                             pairsSel = numpy.array([coord[x],coord[y]])
                             indCross[ind] = int(numpy.where(numpy.all(pairsArray == pairsSel, axis = 1))[0][0])
                             ind += 1
                         dataCross = self.dataIn.data_cspc[indCross,:,:]/(M*N)
                         dataCross = dataCross**2/K
                         for h in range(nHeights):
             #                 print self.dataOut.heightList[h]
                             #Input
                             d = data[:,h]
                             #Covariance Matrix
                             D = numpy.diag(d**2/K)
                             ind = 0
                             for pairs in listComb:
                                 #Coordinates in Covariance Matrix
                                 x = pairs[0]
                                 y = pairs[1]
                                 #Channel Index
                                 S12 = dataCross[ind,:,h]
                                 D12 = numpy.diag(S12)
                                 #Completing Covariance Matrix with Cross Spectras
                                 D[x*N:(x+1)*N,y*N:(y+1)*N] = D12
                                 D[y*N:(y+1)*N,x*N:(x+1)*N] = D12
                                 ind += 1
                             Dinv=numpy.linalg.inv(D)
                             L=numpy.linalg.cholesky(Dinv)
                             LT=L.T
                             dp = numpy.dot(LT,d)
                             #Initial values
                             data_spc = self.dataIn.data_spc[coord,:,h]
                             if (h>0)and(error1[3]<5):
                                 p0 = self.dataOut.data_param[i,:,h-1]
                             else:
                                 p0 = numpy.array(self.dataOut.library.initialValuesFunction(data_spc, constants, i))
                             try:
                                 #Least Squares
                                 minp,covp,infodict,mesg,ier = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants),full_output=True)
             #                   minp,covp = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants))
                                 #Chi square error
                                 error0 = numpy.sum(infodict['fvec']**2)/(2*N)
                                 #Error with Jacobian
                                 error1 = self.dataOut.library.errorFunction(minp,constants,LT)
                             except:
                                 minp = p0*numpy.nan
                                 error0 = numpy.nan
                                 error1 = p0*numpy.nan
                             #Save
                             if self.dataOut.data_param == None:
                                 self.dataOut.data_param = numpy.zeros((nGroups, p0.size, nHeights))*numpy.nan
                                 self.dataOut.data_error = numpy.zeros((nGroups, p0.size + 1, nHeights))*numpy.nan
                             self.dataOut.data_error[i,:,h] = numpy.hstack((error0,error1))
                             self.dataOut.data_param[i,:,h] = minp
                     return
                 def __residFunction(self, p, dp, LT, constants):
                     fm = self.dataOut.library.modelFunction(p, constants)
                     fmp=numpy.dot(LT,fm)
                     return  dp-fmp
                 def __getSNR(self, z, noise):
                     avg = numpy.average(z, axis=1)
                     SNR = (avg.T-noise)/noise
                     SNR = SNR.T
                     return SNR
                 def __chisq(p,chindex,hindex):
                     #similar to Resid but calculates CHI**2
                     [LT,d,fm]=setupLTdfm(p,chindex,hindex)
                     dp=numpy.dot(LT,d)
                     fmp=numpy.dot(LT,fm)
                     chisq=numpy.dot((dp-fmp).T,(dp-fmp))
                     return chisq
             class WindProfiler(Operation):
                 __isConfig = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer = None
                 __dataReady = False
                 __firstdata = None
                 n = None
                 def __init__(self):
                     Operation.__init__(self)
                 def __calculateCosDir(self, elev, azim):
                     zen = (90 - elev)*numpy.pi/180
                     azim = azim*numpy.pi/180
                     cosDirX = numpy.sqrt((1-numpy.cos(zen)**2)/((1+numpy.tan(azim)**2)))
                     cosDirY = numpy.sqrt(1-numpy.cos(zen)**2-cosDirX**2)
                     signX = numpy.sign(numpy.cos(azim))
                     signY = numpy.sign(numpy.sin(azim))
                     cosDirX = numpy.copysign(cosDirX, signX)
                     cosDirY = numpy.copysign(cosDirY, signY)
                     return cosDirX, cosDirY
                 def __calculateAngles(self, theta_x, theta_y, azimuth):
                     dir_cosw = numpy.sqrt(1-theta_x**2-theta_y**2)
                     zenith_arr = numpy.arccos(dir_cosw)
                     azimuth_arr = numpy.arctan2(theta_x,theta_y) + azimuth*math.pi/180
                     dir_cosu = numpy.sin(azimuth_arr)*numpy.sin(zenith_arr)
                     dir_cosv = numpy.cos(azimuth_arr)*numpy.sin(zenith_arr)
                     return azimuth_arr, zenith_arr, dir_cosu, dir_cosv, dir_cosw
                 def __calculateMatA(self, dir_cosu, dir_cosv, dir_cosw, horOnly):
             #
                     if horOnly:
                         A = numpy.c_[dir_cosu,dir_cosv]
                     else:
                         A = numpy.c_[dir_cosu,dir_cosv,dir_cosw]
                     A = numpy.asmatrix(A)
                     A1 = numpy.linalg.inv(A.transpose()*A)*A.transpose()
                     return A1
                 def __correctValues(self, heiRang, phi, velRadial, SNR):
                     listPhi = phi.tolist()
                     maxid = listPhi.index(max(listPhi))
                     minid = listPhi.index(min(listPhi))
                     rango = range(len(phi))
                #     rango = numpy.delete(rango,maxid)
                     heiRang1 = heiRang*math.cos(phi[maxid])
                     heiRangAux = heiRang*math.cos(phi[minid])
                     indOut = (heiRang1 < heiRangAux[0]).nonzero()
                     heiRang1 = numpy.delete(heiRang1,indOut)
                     velRadial1 = numpy.zeros([len(phi),len(heiRang1)])
                     SNR1 = numpy.zeros([len(phi),len(heiRang1)])
                     for i in rango:
                         x = heiRang*math.cos(phi[i])
                         y1 = velRadial[i,:]
                         f1 = interpolate.interp1d(x,y1,kind = 'cubic')
                         x1 = heiRang1
                         y11 = f1(x1)
                         y2 = SNR[i,:]
                         f2 = interpolate.interp1d(x,y2,kind = 'cubic')
                         y21 = f2(x1)
                         velRadial1[i,:] = y11
                         SNR1[i,:] = y21
                     return heiRang1, velRadial1, SNR1
                 def __calculateVelUVW(self, A, velRadial):
                     #Operacion Matricial
             #         velUVW = numpy.zeros((velRadial.shape[1],3))
             #         for ind in range(velRadial.shape[1]):
             #             velUVW[ind,:] = numpy.dot(A,velRadial[:,ind])
             #         velUVW = velUVW.transpose()
                     velUVW = numpy.zeros((A.shape[0],velRadial.shape[1]))
                     velUVW[:,:] = numpy.dot(A,velRadial)
                     return velUVW
                 def techniqueDBS(self, velRadial0, dirCosx, disrCosy, azimuth, correct, horizontalOnly, heiRang, SNR0):
                     """
                     Function that implements Doppler Beam Swinging (DBS) technique.
                     Input:    Radial velocities, Direction cosines (x and y) of the Beam, Antenna azimuth,
                                 Direction correction (if necessary), Ranges and SNR
                     Output:    Winds estimation (Zonal, Meridional and Vertical)
                     Parameters affected:    Winds, height range, SNR
                     """
                     azimuth_arr, zenith_arr, dir_cosu, dir_cosv, dir_cosw = self.__calculateAngles(dirCosx, disrCosy, azimuth)
                     heiRang1, velRadial1, SNR1 = self.__correctValues(heiRang, zenith_arr, correct*velRadial0, SNR0)
                     A = self.__calculateMatA(dir_cosu, dir_cosv, dir_cosw, horizontalOnly)
                     #Calculo de Componentes de la velocidad con DBS
                     winds = self.__calculateVelUVW(A,velRadial1)
                     return winds, heiRang1, SNR1
                 def __calculateDistance(self, posx, posy, pairsCrossCorr, pairsList, pairs, azimuth = None):
                     posx = numpy.asarray(posx)
                     posy = numpy.asarray(posy)
                     #Rotacion Inversa para alinear con el azimuth
                     if azimuth!= None:
                         azimuth = azimuth*math.pi/180
                         posx1 = posx*math.cos(azimuth) + posy*math.sin(azimuth)
                         posy1 = -posx*math.sin(azimuth) + posy*math.cos(azimuth)
                     else:
                         posx1 = posx
                         posy1 = posy
                     #Calculo de Distancias
                     distx = numpy.zeros(pairsCrossCorr.size)
                     disty = numpy.zeros(pairsCrossCorr.size)
                     dist = numpy.zeros(pairsCrossCorr.size)
                     ang = numpy.zeros(pairsCrossCorr.size)
                     for i in range(pairsCrossCorr.size):
                         distx[i] = posx1[pairsList[pairsCrossCorr[i]][1]] - posx1[pairsList[pairsCrossCorr[i]][0]]
                         disty[i] = posy1[pairsList[pairsCrossCorr[i]][1]] - posy1[pairsList[pairsCrossCorr[i]][0]]
                         dist[i] = numpy.sqrt(distx[i]**2 + disty[i]**2)
                         ang[i] = numpy.arctan2(disty[i],distx[i])
                     #Calculo de Matrices
                     nPairs = len(pairs)
                     ang1 = numpy.zeros((nPairs, 2, 1))
                     dist1 = numpy.zeros((nPairs, 2, 1))
                     for j in range(nPairs):
                         dist1[j,0,0] = dist[pairs[j][0]]
                         dist1[j,1,0] = dist[pairs[j][1]]
                         ang1[j,0,0] = ang[pairs[j][0]]
                         ang1[j,1,0] = ang[pairs[j][1]]
                     return distx,disty, dist1,ang1
                 def __calculateVelVer(self, phase, lagTRange, _lambda):
                     Ts = lagTRange[1] - lagTRange[0]
                     velW = -_lambda*phase/(4*math.pi*Ts)
                     return velW
                 def __calculateVelHorDir(self, dist, tau1, tau2, ang):
                     nPairs = tau1.shape[0]
                     vel = numpy.zeros((nPairs,3,tau1.shape[2]))
                     angCos = numpy.cos(ang)
                     angSin = numpy.sin(ang)
                     vel0 = dist*tau1/(2*tau2**2)
                     vel[:,0,:] = (vel0*angCos).sum(axis = 1)
                     vel[:,1,:] = (vel0*angSin).sum(axis = 1)
                     ind = numpy.where(numpy.isinf(vel))
                     vel[ind] = numpy.nan
                     return vel
                 def __getPairsAutoCorr(self, pairsList, nChannels):
                     pairsAutoCorr = numpy.zeros(nChannels, dtype = 'int')*numpy.nan
                     for l in range(len(pairsList)):
                         firstChannel = pairsList[l][0]
                         secondChannel = pairsList[l][1]
                         #Obteniendo pares de Autocorrelacion
                         if firstChannel == secondChannel:
                             pairsAutoCorr[firstChannel] = int(l)
                     pairsAutoCorr = pairsAutoCorr.astype(int)
                     pairsCrossCorr = range(len(pairsList))
                     pairsCrossCorr = numpy.delete(pairsCrossCorr,pairsAutoCorr)
                     return pairsAutoCorr, pairsCrossCorr
                 def techniqueSA(self, pairsSelected, pairsList, nChannels, tau, azimuth, _lambda, position_x, position_y, lagTRange, correctFactor):
                     """
                     Function that implements Spaced Antenna (SA) technique.
                     Input:    Radial velocities, Direction cosines (x and y) of the Beam, Antenna azimuth,
                                 Direction correction (if necessary), Ranges and SNR
                     Output:    Winds estimation (Zonal, Meridional and Vertical)
                     Parameters affected:    Winds
                     """
                     #Cross Correlation pairs obtained
                     pairsAutoCorr, pairsCrossCorr = self.__getPairsAutoCorr(pairsList, nChannels)
                     pairsArray = numpy.array(pairsList)[pairsCrossCorr]
                     pairsSelArray = numpy.array(pairsSelected)
                     pairs = []
                     #Wind estimation pairs obtained
                     for i in range(pairsSelArray.shape[0]/2):
                         ind1 = numpy.where(numpy.all(pairsArray == pairsSelArray[2*i], axis = 1))[0][0]
                         ind2 = numpy.where(numpy.all(pairsArray == pairsSelArray[2*i + 1], axis = 1))[0][0]
                         pairs.append((ind1,ind2))
                     indtau = tau.shape[0]/2
                     tau1 = tau[:indtau,:]
                     tau2 = tau[indtau:-1,:]
                     tau1 = tau1[pairs,:]
                     tau2 = tau2[pairs,:]
                     phase1 = tau[-1,:]
                     #---------------------------------------------------------------------
                     #Metodo Directo
                     distx, disty, dist, ang = self.__calculateDistance(position_x, position_y, pairsCrossCorr, pairsList, pairs,azimuth)
                     winds = self.__calculateVelHorDir(dist, tau1, tau2, ang)
                     winds = stats.nanmean(winds, axis=0)
                     #---------------------------------------------------------------------
                     #Metodo General
             #         distx, disty, dist = self.calculateDistance(position_x,position_y,pairsCrossCorr, pairsList, azimuth)
             #         #Calculo Coeficientes de Funcion de Correlacion
             #         F,G,A,B,H = self.calculateCoef(tau1,tau2,distx,disty,n)
             #         #Calculo de Velocidades
             #         winds = self.calculateVelUV(F,G,A,B,H)
                     #---------------------------------------------------------------------
                     winds[2,:] = self.__calculateVelVer(phase1, lagTRange, _lambda)
                     winds = correctFactor*winds
                     return winds
                 def __checkTime(self, currentTime, paramInterval, outputInterval):
                     dataTime = currentTime + paramInterval
                     deltaTime = dataTime - self.__initime
                     if deltaTime >= outputInterval or deltaTime < 0:
                         self.__dataReady = True
                     return
                 def techniqueMeteors(self, arrayMeteor, meteorThresh, heightMin, heightMax):
                     '''
                     Function that implements winds estimation technique with detected meteors.
                     Input:    Detected meteors, Minimum meteor quantity to wind estimation
                     Output:    Winds estimation (Zonal and Meridional)
                     Parameters affected:    Winds
                     '''
+                    print arrayMeteor.shape
                     #Settings
                     nInt = (heightMax - heightMin)/2
                     winds = numpy.zeros((2,nInt))*numpy.nan
                     #Filter errors
-                    error = numpy.where(arrayMeteor[:,-1] == 0)[0]
+                    error = numpy.where(arrayMeteor[0,:,-1] == 0)[0]
-                    finalMeteor = arrayMeteor[error,:]
+                    finalMeteor = arrayMeteor[0,error,:]
                     #Meteor Histogram
                     finalHeights = finalMeteor[:,3]
                     hist = numpy.histogram(finalHeights, bins = nInt, range = (heightMin,heightMax))
                     nMeteorsPerI = hist[0]
                     heightPerI = hist[1]
                     #Sort of meteors
                     indSort = finalHeights.argsort()
                     finalMeteor2 = finalMeteor[indSort,:]
                     #    Calculating winds
                     ind1 = 0
                     ind2 = 0
                     for i in range(nInt):
                         nMet = nMeteorsPerI[i]
                         ind1 = ind2
                         ind2 = ind1 + nMet
                         meteorAux = finalMeteor2[ind1:ind2,:]
                         if meteorAux.shape[0] >= meteorThresh:
                             vel = meteorAux[:, 7]
                             zen = meteorAux[:, 5]*numpy.pi/180
                             azim = meteorAux[:, 4]*numpy.pi/180
                             n = numpy.cos(zen)
                     #         m = (1 - n**2)/(1 - numpy.tan(azim)**2)
                     #         l = m*numpy.tan(azim)
                             l = numpy.sin(zen)*numpy.sin(azim)
                             m = numpy.sin(zen)*numpy.cos(azim)
                             A = numpy.vstack((l, m)).transpose()
                             A1 = numpy.dot(numpy.linalg.inv( numpy.dot(A.transpose(),A) ),A.transpose())
                             windsAux = numpy.dot(A1, vel)
                             winds[0,i] = windsAux[0]
                             winds[1,i] = windsAux[1]
                     return winds, heightPerI[:-1]
                 def run(self, dataOut, technique, **kwargs):
                     param = dataOut.data_param
-                    if dataOut.abscissaList != None:
+            #         if dataOut.abscissaList != None:
-                        absc = dataOut.abscissaList[:-1]
+            #             absc = dataOut.abscissaList[:-1]
                     noise = dataOut.noise
-                    heightList = dataOut.getHeiRange()
+                    heightList = dataOut.heightList
                     SNR = dataOut.data_SNR
                     if technique == 'DBS':
                         if kwargs.has_key('dirCosx') and kwargs.has_key('dirCosy'):
                             theta_x = numpy.array(kwargs['dirCosx'])
                             theta_y = numpy.array(kwargs['dirCosy'])
                         else:
                             elev = numpy.array(kwargs['elevation'])
                             azim = numpy.array(kwargs['azimuth'])
                             theta_x, theta_y = self.__calculateCosDir(elev, azim)
                         azimuth = kwargs['correctAzimuth']
                         if kwargs.has_key('horizontalOnly'):
                             horizontalOnly = kwargs['horizontalOnly']
                         else:   horizontalOnly = False
                         if kwargs.has_key('correctFactor'):
                             correctFactor = kwargs['correctFactor']
                         else:   correctFactor = 1
                         if kwargs.has_key('channelList'):
                             channelList = kwargs['channelList']
+                            if len(channelList) == 2:
+                                horizontalOnly = True
                             arrayChannel = numpy.array(channelList)
                             param = param[arrayChannel,:,:]
                             theta_x = theta_x[arrayChannel]
                             theta_y = theta_y[arrayChannel]
                         velRadial0 = param[:,1,:] #Radial velocity
-                        if velRadial0.shape[0] != theta_x.shape[0] or velRadial0.shape[0] != theta_y.shape[0]:
-                            raise ValueError, "The max number of channels is %d, and the length of cosine director is %d. Please check: dirCosX, dirCosY, elevation or azimuth arguments" %(velRadial0.shape[0], theta_x.shape[0])
-                        if theta_x.shape[0] == 2:
-                                horizontalOnly = True
                         dataOut.data_output, dataOut.heightList, dataOut.data_SNR = self.techniqueDBS(velRadial0, theta_x, theta_y, azimuth, correctFactor, horizontalOnly, heightList, SNR) #DBS Function
                         dataOut.utctimeInit = dataOut.utctime
                         dataOut.outputInterval = dataOut.timeInterval
                     elif technique == 'SA':
                         #Parameters
                         position_x = kwargs['positionX']
                         position_y = kwargs['positionY']
                         azimuth = kwargs['azimuth']
                         if kwargs.has_key('crosspairsList'):
                             pairs = kwargs['crosspairsList']
                         else:
                             pairs = None
                         if kwargs.has_key('correctFactor'):
                             correctFactor = kwargs['correctFactor']
                         else:
                             correctFactor = 1
                         tau = dataOut.data_param
                         _lambda = dataOut.C/dataOut.frequency
                         pairsList = dataOut.groupList
                         nChannels = dataOut.nChannels
                         dataOut.data_output = self.techniqueSA(pairs, pairsList, nChannels, tau, azimuth, _lambda, position_x, position_y, absc, correctFactor)
                         dataOut.utctimeInit = dataOut.utctime
                         dataOut.outputInterval = dataOut.timeInterval
                     elif technique == 'Meteors':
                         dataOut.flagNoData = True
                         self.__dataReady = False
                         if kwargs.has_key('nHours'):
                             nHours = kwargs['nHours']
                         else:
                             nHours = 1
                         if kwargs.has_key('meteorsPerBin'):
                             meteorThresh = kwargs['meteorsPerBin']
                         else:
                             meteorThresh = 6
                         if kwargs.has_key('hmin'):
                             hmin = kwargs['hmin']
                         else:   hmin = 70
                         if kwargs.has_key('hmax'):
                             hmax = kwargs['hmax']
                         else:   hmax = 110
                         dataOut.outputInterval = nHours*3600
                         if self.__isConfig == False:
             #                 self.__initime = dataOut.datatime.replace(minute = 0, second = 0, microsecond = 03)
                             #Get Initial LTC time
-                            self.__initime = datetime.datetime.utcfromtimestamp(self.dataOut.utctime)
+                            self.__initime = datetime.datetime.utcfromtimestamp(dataOut.utctime)
                             self.__initime = (self.__initime.replace(minute = 0, second = 0, microsecond = 0) - datetime.datetime(1970, 1, 1)).total_seconds()
                             self.__isConfig = True
                         if self.__buffer == None:
                             self.__buffer = dataOut.data_param
                             self.__firstdata = copy.copy(dataOut)
                         else:
-                            self.__buffer = numpy.vstack((self.__buffer, dataOut.data_param))
+                            self.__buffer = numpy.hstack((self.__buffer, dataOut.data_param))
                         self.__checkTime(dataOut.utctime, dataOut.paramInterval, dataOut.outputInterval) #Check if the buffer is ready
                         if self.__dataReady:
                             dataOut.utctimeInit = self.__initime
                             self.__initime += dataOut.outputInterval #to erase time offset
                             dataOut.data_output, dataOut.heightList = self.techniqueMeteors(self.__buffer, meteorThresh, hmin, hmax)
                             dataOut.flagNoData = False
                             self.__buffer = None
                     return
             class EWDriftsEstimation(Operation):
                 def __init__(self):
                     Operation.__init__(self)
                 def __correctValues(self, heiRang, phi, velRadial, SNR):
                     listPhi = phi.tolist()
                     maxid = listPhi.index(max(listPhi))
                     minid = listPhi.index(min(listPhi))
                     rango = range(len(phi))
                #     rango = numpy.delete(rango,maxid)
                     heiRang1 = heiRang*math.cos(phi[maxid])
                     heiRangAux = heiRang*math.cos(phi[minid])
                     indOut = (heiRang1 < heiRangAux[0]).nonzero()
                     heiRang1 = numpy.delete(heiRang1,indOut)
                     velRadial1 = numpy.zeros([len(phi),len(heiRang1)])
                     SNR1 = numpy.zeros([len(phi),len(heiRang1)])
                     for i in rango:
                         x = heiRang*math.cos(phi[i])
                         y1 = velRadial[i,:]
                         f1 = interpolate.interp1d(x,y1,kind = 'cubic')
                         x1 = heiRang1
                         y11 = f1(x1)
                         y2 = SNR[i,:]
                         f2 = interpolate.interp1d(x,y2,kind = 'cubic')
                         y21 = f2(x1)
                         velRadial1[i,:] = y11
                         SNR1[i,:] = y21
                     return heiRang1, velRadial1, SNR1
                 def run(self, dataOut, zenith, zenithCorrection):
                     heiRang = dataOut.heightList
                     velRadial = dataOut.data_param[:,3,:]
                     SNR = dataOut.data_SNR
                     zenith = numpy.array(zenith)
                     zenith -= zenithCorrection
                     zenith *= numpy.pi/180
                     heiRang1, velRadial1, SNR1 = self.__correctValues(heiRang, numpy.abs(zenith), velRadial, SNR)
                     alp = zenith[0]
                     bet = zenith[1]
                     w_w = velRadial1[0,:]
                     w_e = velRadial1[1,:]
                     w = (w_w*numpy.sin(bet) - w_e*numpy.sin(alp))/(numpy.cos(alp)*numpy.sin(bet) - numpy.cos(bet)*numpy.sin(alp))
                     u = (w_w*numpy.cos(bet) - w_e*numpy.cos(alp))/(numpy.sin(alp)*numpy.cos(bet) - numpy.sin(bet)*numpy.cos(alp))
                     winds = numpy.vstack((u,w))
                     dataOut.heightList = heiRang1
                     dataOut.data_output = winds
                     dataOut.data_SNR = SNR1
                     dataOut.utctimeInit = dataOut.utctime
                     dataOut.outputInterval = dataOut.timeInterval
                     return
+            class PhaseCalibration(Operation):
+                __buffer = None
+                __initime = None
+                __dataReady = False
+                __isConfig = False
+                def __checkTime(self, currentTime, initTime, paramInterval, outputInterval):
+                    dataTime = currentTime + paramInterval
+                    deltaTime = dataTime - initTime
+                    if deltaTime >= outputInterval or deltaTime < 0:
+                        return True
+                    return False
+                def __getGammas(self, pairs, k, d, phases):
+                    gammas = numpy.zeros(2)
+                    for i in range(len(pairs)):
+                        pairi = pairs[i]
+                        #Calculating gamma
+                        jdcos = phases[:,pairi[1]]/(k*d[pairi[1]])
+                        jgamma = numpy.angle(numpy.exp(1j*(k*d[pairi[0]]*jdcos - phases[:,pairi[0]])))
+                        #Revised distribution
+                        jgammaArray = numpy.hstack((jgamma,jgamma+0.5*numpy.pi,jgamma-0.5*numpy.pi))
+                        #Histogram
+                        nBins = 64.0
+                        rmin = -0.5*numpy.pi
+                        rmax = 0.5*numpy.pi
+                        phaseHisto = numpy.histogram(jgammaArray, bins=nBins, range=(rmin,rmax))
+                        meteorsY = phaseHisto[0]
+                        phasesX = phaseHisto[1][:-1]
+                        width = phasesX[1] - phasesX[0]
+                        phasesX += width/2
+                        #Gaussian aproximation
+                        bpeak = meteorsY.argmax()
+                        peak = meteorsY.max()
+                        jmin = bpeak - 5
+                        jmax = bpeak + 5 + 1
+                        if jmin<0:
+                            jmin = 0
+                            jmax = 6
+                        elif jmax > meteorsY.size:
+                            jmin = meteorsY.size - 6
+                            jmax = meteorsY.size
+                        x0 = numpy.array([peak,bpeak,50])
+                        coeff = optimize.leastsq(self.__residualFunction, x0, args=(meteorsY[jmin:jmax], phasesX[jmin:jmax]))
+                        #Gammas
+                        gammas[i] = coeff[0][1]
+            #             gammas[i] = bpeak
+                    return gammas
+                def __residualFunction(self, coeffs, y, t):
+                    return y - self.__gauss_function(t, coeffs)
+                def __gauss_function(self, t, coeffs):
+                    return coeffs[0]*numpy.exp(-0.5*((t - coeffs[1]) / coeffs[2])**2)
+                def __getPhases(self, azimuth, h, pairsList, d, gammas, meteorsArray):
+                    meteorOps = MeteorOperations()
+                    nchan = 4
+                    pairx = pairsList[0]
+                    pairy = pairsList[1]
+                    center_xangle = 0
+                    center_yangle = 0
+                    range_angle = numpy.array([8*numpy.pi,numpy.pi,numpy.pi/2,numpy.pi/4])
+                    ntimes = len(range_angle)
+                    nstepsx = 20.0
+                    nstepsy = 20.0
+                    for iz in range(ntimes):
+                        min_xangle = -range_angle[iz]/2 + center_xangle
+                        max_xangle = range_angle[iz]/2 + center_xangle
+                        min_yangle = -range_angle[iz]/2 + center_yangle
+                        max_yangle = range_angle[iz]/2 + center_yangle
+                        inc_x = (max_xangle-min_xangle)/nstepsx
+                        inc_y = (max_yangle-min_yangle)/nstepsy
+                        alpha_y = numpy.arange(nstepsy)*inc_y + min_yangle
+                        alpha_x = numpy.arange(nstepsx)*inc_x + min_xangle
+                        penalty = numpy.zeros((nstepsx,nstepsy))
+                        jph_array = numpy.zeros((nchan,nstepsx,nstepsy))
+                        jph = numpy.zeros(nchan)
+                        # Iterations looking for the offset
+                        for iy in range(int(nstepsy)):
+                            for ix in range(int(nstepsx)):
+                                jph[pairy[1]] = alpha_y[iy]
+                                jph[pairy[0]] = -gammas[1] + alpha_y[iy]*d[pairy[0]]/d[pairy[1]]
+                                jph[pairx[1]] = alpha_x[ix]
+                                jph[pairx[0]] = -gammas[0] + alpha_x[ix]*d[pairx[0]]/d[pairx[1]]
+                                jph_array[:,ix,iy] = jph
+                                meteorsArray1 = meteorOps.getMeteorParams(meteorsArray, azimuth, h, pairsList, jph)
+                                error = meteorsArray1[:,-1]
+                                ind1 = numpy.where(error==0)[0]
+                                penalty[ix,iy] = ind1.size
+                        i,j = numpy.unravel_index(penalty.argmax(), penalty.shape)
+                        phOffset = jph_array[:,i,j]
+                        center_xangle = phOffset[pairx[1]]
+                        center_yangle = phOffset[pairy[1]]
+                    phOffset = numpy.angle(numpy.exp(1j*jph_array[:,i,j]))
+                    phOffset = phOffset*180/numpy.pi
+                    return phOffset
+                def run(self, dataOut, pairs, distances, hmin, hmax, nHours = None):
+                    dataOut.flagNoData = True
+                    self.__dataReady = False
+                    if nHours == None:
+                        nHours = 1
+                    dataOut.outputInterval = nHours*3600
+                    if self.__isConfig == False:
+            #                 self.__initime = dataOut.datatime.replace(minute = 0, second = 0, microsecond = 03)
+                        #Get Initial LTC time
+                        self.__initime = datetime.datetime.utcfromtimestamp(dataOut.utctime)
+                        self.__initime = (self.__initime.replace(minute = 0, second = 0, microsecond = 0) - datetime.datetime(1970, 1, 1)).total_seconds()
+                        self.__isConfig = True
+                    if self.__buffer == None:
+                        self.__buffer = dataOut.data_param.copy()
+                    else:
+                        self.__buffer = numpy.hstack((self.__buffer, dataOut.data_param))
+                    self.__dataReady = self.__checkTime(dataOut.utctime, self.__initime, dataOut.paramInterval, dataOut.outputInterval) #Check if the buffer is ready
+                    if self.__dataReady:
+                        dataOut.utctimeInit = self.__initime
+                        self.__initime += dataOut.outputInterval #to erase time offset
+                        freq = dataOut.frequency
+                        c = dataOut.C #m/s
+                        lamb = c/freq
+                        k = 2*numpy.pi/lamb
+                        azimuth = 0
+                        h = (hmin, hmax)
+                        pairsList = ((0,3),(1,2))
+                        meteorsArray = self.__buffer[0,:,:]
+                        error = meteorsArray[:,-1]
+                        boolError = (error==0)|(error==3)|(error==4)|(error==13)|(error==14)
+                        ind1 = numpy.where(boolError)[0]
+                        meteorsArray = meteorsArray[ind1,:]
+                        meteorsArray[:,-1] = 0
+                        phases = meteorsArray[:,9:13]
+                        #Calculate Gammas
+                        gammas = self.__getGammas(pairs, k, distances, phases)
+            #             gammas = numpy.array([-21.70409463,45.76935864])*numpy.pi/180
+                        #Calculate Phases
+                        phasesOff = self.__getPhases(azimuth, h, pairsList, distances, gammas, meteorsArray)
+                        phasesOff = phasesOff.reshape((1,phasesOff.size))
+                        dataOut.data_output = -phasesOff
+                        dataOut.flagNoData = False
+                        self.__buffer = None
+                    return
+            class MeteorOperations():
+                def __init__(self):
+                    return
+                def getMeteorParams(self, arrayParameters0, azimuth, h, pairsList, jph):
+                    arrayParameters = arrayParameters0.copy()
+                    hmin = h[0]
+                    hmax = h[1]
+                    #Calculate AOA (Error N 3, 4)
+                    #JONES ET AL. 1998
+                    AOAthresh = numpy.pi/8
+                    error = arrayParameters[:,-1]
+                    phases = -arrayParameters[:,9:13] + jph
+                    arrayParameters[:,4:7], arrayParameters[:,-1] = self.__getAOA(phases, pairsList, error, AOAthresh, azimuth)
+                    #Calculate Heights (Error N 13 and 14)
+                    error = arrayParameters[:,-1]
+                    Ranges = arrayParameters[:,2]
+                    zenith = arrayParameters[:,5]
+                    arrayParameters[:,3], arrayParameters[:,-1] = self.__getHeights(Ranges, zenith, error, hmin, hmax)
+                    #----------------------- Get Final data    ------------------------------------
+            #         error = arrayParameters[:,-1]
+            #         ind1 = numpy.where(error==0)[0]
+            #         arrayParameters = arrayParameters[ind1,:]
+                    return arrayParameters
+                def __getAOA(self, phases, pairsList, error, AOAthresh, azimuth):
+                    arrayAOA = numpy.zeros((phases.shape[0],3))
+                    cosdir0, cosdir = self.__getDirectionCosines(phases, pairsList)
+                    arrayAOA[:,:2] = self.__calculateAOA(cosdir, azimuth)
+                    cosDirError = numpy.sum(numpy.abs(cosdir0 - cosdir), axis = 1)
+                    arrayAOA[:,2] = cosDirError
+                    azimuthAngle = arrayAOA[:,0]
+                    zenithAngle = arrayAOA[:,1]
+                    #Setting Error
+                    #Number 3: AOA not fesible
+                    indInvalid = numpy.where(numpy.logical_and((numpy.logical_or(numpy.isnan(zenithAngle), numpy.isnan(azimuthAngle))),error == 0))[0]
+                    error[indInvalid] = 3
+                    #Number 4: Large difference in AOAs obtained from different antenna baselines
+                    indInvalid = numpy.where(numpy.logical_and(cosDirError > AOAthresh,error == 0))[0]
+                    error[indInvalid] = 4
+                    return arrayAOA, error
+                def __getDirectionCosines(self, arrayPhase, pairsList):
+                    #Initializing some variables
+                    ang_aux = numpy.array([-8,-7,-6,-5,-4,-3,-2,-1,0,1,2,3,4,5,6,7,8])*2*numpy.pi
+                    ang_aux = ang_aux.reshape(1,ang_aux.size)
+                    cosdir = numpy.zeros((arrayPhase.shape[0],2))
+                    cosdir0 = numpy.zeros((arrayPhase.shape[0],2))
+                    for i in range(2):
+                        #First Estimation
+                        phi0_aux = arrayPhase[:,pairsList[i][0]] + arrayPhase[:,pairsList[i][1]]
+                        #Dealias
+                        indcsi = numpy.where(phi0_aux > numpy.pi)
+                        phi0_aux[indcsi] -= 2*numpy.pi
+                        indcsi = numpy.where(phi0_aux < -numpy.pi)
+                        phi0_aux[indcsi] += 2*numpy.pi
+                        #Direction Cosine 0
+                        cosdir0[:,i] = -(phi0_aux)/(2*numpy.pi*0.5)
+                        #Most-Accurate Second Estimation
+                        phi1_aux =  arrayPhase[:,pairsList[i][0]] - arrayPhase[:,pairsList[i][1]]
+                        phi1_aux = phi1_aux.reshape(phi1_aux.size,1)
+                        #Direction Cosine 1
+                        cosdir1 = -(phi1_aux + ang_aux)/(2*numpy.pi*4.5)
+                        #Searching the correct Direction Cosine
+                        cosdir0_aux = cosdir0[:,i]
+                        cosdir0_aux = cosdir0_aux.reshape(cosdir0_aux.size,1)
+                        #Minimum Distance
+                        cosDiff = (cosdir1 - cosdir0_aux)**2
+                        indcos = cosDiff.argmin(axis = 1)
+                        #Saving Value obtained
+                        cosdir[:,i] = cosdir1[numpy.arange(len(indcos)),indcos]
+                    return cosdir0, cosdir
+                def __calculateAOA(self, cosdir, azimuth):
+                    cosdirX = cosdir[:,0]
+                    cosdirY = cosdir[:,1]
+                    zenithAngle = numpy.arccos(numpy.sqrt(1 - cosdirX**2 - cosdirY**2))*180/numpy.pi
+                    azimuthAngle = numpy.arctan2(cosdirX,cosdirY)*180/numpy.pi + azimuth #0 deg north, 90 deg east
+                    angles = numpy.vstack((azimuthAngle, zenithAngle)).transpose()
+                    return angles
+                def __getHeights(self, Ranges, zenith, error, minHeight, maxHeight):
+                    Ramb = 375  #Ramb = c/(2*PRF)
+                    Re = 6371   #Earth Radius
+                    heights = numpy.zeros(Ranges.shape)
+                    R_aux = numpy.array([0,1,2])*Ramb
+                    R_aux = R_aux.reshape(1,R_aux.size)
+                    Ranges = Ranges.reshape(Ranges.size,1)
+                    Ri = Ranges + R_aux
+                    hi = numpy.sqrt(Re**2 + Ri**2 + (2*Re*numpy.cos(zenith*numpy.pi/180)*Ri.transpose()).transpose()) - Re
+                    #Check if there is a height between 70 and 110 km
+                    h_bool = numpy.sum(numpy.logical_and(hi > minHeight, hi < maxHeight), axis = 1)
+                    ind_h = numpy.where(h_bool == 1)[0]
+                    hCorr = hi[ind_h, :]
+                    ind_hCorr = numpy.where(numpy.logical_and(hi > minHeight, hi < maxHeight))
+                    hCorr = hi[ind_hCorr]
+                    heights[ind_h] = hCorr
+                    #Setting Error
+                    #Number 13: Height unresolvable echo: not valid height within 70 to 110 km
+                    #Number 14: Height ambiguous echo: more than one possible height within 70 to 110 km
+                    indInvalid2 = numpy.where(numpy.logical_and(h_bool > 1, error == 0))[0]
+                    error[indInvalid2] = 14
+                    indInvalid1 = numpy.where(numpy.logical_and(h_bool == 0, error == 0))[0]
+                    error[indInvalid1] = 13
+                    return heights, error
   No newline at end of file

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