##// END OF EJS Templates
First Spectral Fitting and EW Drifts operative module inside Signal Chain TRUNK
First Spectral Fitting and EW Drifts operative module inside Signal Chain TRUNK

File last commit:

r502:8b1242a708e8
r513:19f921674eb5
Show More
jroproc_correlation.py
245 lines | 9.8 KiB | text/x-python | PythonLexer
/ schainpy / model / proc / jroproc_correlation.py
import numpy
from jroproc_base import ProcessingUnit, Operation
from model.data.jrodata import Correlation
class CorrelationProc(ProcessingUnit):
def __init__(self):
ProcessingUnit.__init__(self)
self.objectDict = {}
self.buffer = None
self.firstdatatime = None
self.profIndex = 0
self.dataOut = Correlation()
def __updateObjFromInput(self):
self.dataOut.timeZone = self.dataIn.timeZone
self.dataOut.dstFlag = self.dataIn.dstFlag
self.dataOut.errorCount = self.dataIn.errorCount
self.dataOut.useLocalTime = self.dataIn.useLocalTime
self.dataOut.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.flagTimeBlock = self.dataIn.flagTimeBlock
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.nPoints
def removeDC(self, jspectra):
nChannel = jspectra.shape[0]
for i in range(nChannel):
jspectra_tmp = jspectra[i,:,:]
jspectra_DC = numpy.mean(jspectra_tmp,axis = 0)
jspectra_tmp = jspectra_tmp - jspectra_DC
jspectra[i,:,:] = jspectra_tmp
return jspectra
def removeNoise(self, mode = 2):
indR = numpy.where(self.dataOut.lagR == 0)[0][0]
indT = numpy.where(self.dataOut.lagT == 0)[0][0]
jspectra = self.dataOut.data_corr[:,:,indR,:]
num_chan = jspectra.shape[0]
num_hei = jspectra.shape[2]
freq_dc = indT
ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
NPot = self.dataOut.getNoise(mode)
jspectra[:,freq_dc,:] = jspectra[:,freq_dc,:] - NPot
SPot = jspectra[:,freq_dc,:]
pairsAutoCorr = self.dataOut.getPairsAutoCorr()
# self.dataOut.signalPotency = SPot
self.dataOut.noise = NPot
self.dataOut.SNR = (SPot/NPot)[pairsAutoCorr]
self.dataOut.data_corr[:,:,indR,:] = jspectra
return 1
def calculateNormFactor(self):
pairsList = self.dataOut.pairsList
pairsAutoCorr = self.dataOut.pairsAutoCorr
nHeights = self.dataOut.nHeights
nPairs = len(pairsList)
normFactor = numpy.zeros((nPairs,nHeights))
indR = numpy.where(self.dataOut.lagR == 0)[0][0]
indT = numpy.where(self.dataOut.lagT == 0)[0][0]
for l in range(len(pairsList)):
firstChannel = pairsList[l][0]
secondChannel = pairsList[l][1]
AC1 = pairsAutoCorr[firstChannel]
AC2 = pairsAutoCorr[secondChannel]
if (AC1 >= 0 and AC2 >= 0):
data1 = numpy.abs(self.dataOut.data_corr[AC1,:,indR,:])
data2 = numpy.abs(self.dataOut.data_corr[AC2,:,indR,:])
maxim1 = data1.max(axis = 0)
maxim2 = data1.max(axis = 0)
maxim = numpy.sqrt(maxim1*maxim2)
else:
#In case there is no autocorrelation for the pair
data = numpy.abs(self.dataOut.data_corr[l,:,indR,:])
maxim = numpy.max(data, axis = 0)
normFactor[l,:] = maxim
self.dataOut.normFactor = normFactor
return 1
def run(self, lagT=None, lagR=None, pairsList=None,
nPoints=None, nAvg=None, bufferSize=None,
fullT = False, fullR = False, removeDC = False):
self.dataOut.flagNoData = True
if self.dataIn.type == "Correlation":
self.dataOut.copy(self.dataIn)
return
if self.dataIn.type == "Voltage":
if pairsList == None:
pairsList = [numpy.array([0,0])]
if nPoints == None:
nPoints = 128
#------------------------------------------------------------
#Condicionales para calcular Correlaciones en Tiempo y Rango
if fullT:
lagT = numpy.arange(nPoints*2 - 1) - nPoints + 1
elif lagT == None:
lagT = numpy.array([0])
else:
lagT = numpy.array(lagT)
if fullR:
lagR = numpy.arange(self.dataOut.nHeights)
elif lagR == None:
lagR = numpy.array([0])
#-------------------------------------------------------------
if nAvg == None:
nAvg = 1
if bufferSize == None:
bufferSize = 0
deltaH = self.dataIn.heightList[1] - self.dataIn.heightList[0]
self.dataOut.lagR = numpy.round(numpy.array(lagR)/deltaH)
self.dataOut.pairsList = pairsList
self.dataOut.nPoints = nPoints
# channels = numpy.sort(list(set(list(itertools.chain.from_iterable(pairsList)))))
if self.buffer == None:
self.buffer = numpy.zeros((self.dataIn.nChannels,self.dataIn.nProfiles,self.dataIn.nHeights),dtype='complex')
self.buffer[:,self.profIndex,:] = self.dataIn.data.copy()[:,:]
self.profIndex += 1
if self.firstdatatime == None:
self.firstdatatime = self.dataIn.utctime
if self.profIndex == nPoints:
tmp = self.buffer[:,0:nPoints,:]
self.buffer = None
self.buffer = tmp
#--------------- Remover DC ------------
if removeDC:
self.buffer = self.removeDC(self.buffer)
#---------------------------------------------
self.dataOut.data_volts = self.buffer
self.__updateObjFromInput()
self.dataOut.data_corr = numpy.zeros((len(pairsList),
len(lagT),len(lagR),
self.dataIn.nHeights),
dtype='complex')
for l in range(len(pairsList)):
firstChannel = pairsList[l][0]
secondChannel = pairsList[l][1]
tmp = None
tmp = numpy.zeros((len(lagT),len(lagR),self.dataIn.nHeights),dtype='complex')
for t in range(len(lagT)):
for r in range(len(lagR)):
idxT = lagT[t]
idxR = lagR[r]
if idxT >= 0:
vStacked = numpy.vstack((self.buffer[secondChannel,idxT:,:],
numpy.zeros((idxT,self.dataIn.nHeights),dtype='complex')))
else:
vStacked = numpy.vstack((numpy.zeros((-idxT,self.dataIn.nHeights),dtype='complex'),
self.buffer[secondChannel,:(nPoints + idxT),:]))
if idxR >= 0:
hStacked = numpy.hstack((vStacked[:,idxR:],numpy.zeros((nPoints,idxR),dtype='complex')))
else:
hStacked = numpy.hstack((numpy.zeros((nPoints,-idxR),dtype='complex'),vStacked[:,(self.dataOut.nHeights + idxR)]))
tmp[t,r,:] = numpy.sum((numpy.conjugate(self.buffer[firstChannel,:,:])*hStacked),axis=0)
hStacked = None
vStacked = None
self.dataOut.data_corr[l,:,:,:] = tmp[:,:,:]
#Se Calcula los factores de Normalizacion
self.dataOut.pairsAutoCorr = self.dataOut.getPairsAutoCorr()
self.dataOut.lagT = lagT*self.dataIn.ippSeconds*self.dataIn.nCohInt
self.dataOut.lagR = lagR
self.calculateNormFactor()
self.dataOut.flagNoData = False
self.buffer = None
self.firstdatatime = None
self.profIndex = 0
return