jroproc_correlation.py
178 lines
| 6.7 KiB
| text/x-python
|
PythonLexer
|
r502 | import numpy | ||
from jroproc_base import ProcessingUnit, Operation | ||||
|
r854 | from schainpy.model.data.jrodata import Correlation, hildebrand_sekhon | ||
|
r502 | |||
class CorrelationProc(ProcessingUnit): | ||||
|
r897 | |||
|
r854 | pairsList = None | ||
|
r897 | |||
|
r854 | data_cf = None | ||
|
r897 | |||
def __init__(self, **kwargs): | ||||
ProcessingUnit.__init__(self, **kwargs) | ||||
|
r502 | self.objectDict = {} | ||
self.buffer = None | ||||
self.firstdatatime = None | ||||
self.profIndex = 0 | ||||
self.dataOut = Correlation() | ||||
|
r897 | |||
|
r587 | def __updateObjFromVoltage(self): | ||
|
r897 | |||
|
r502 | self.dataOut.timeZone = self.dataIn.timeZone | ||
self.dataOut.dstFlag = self.dataIn.dstFlag | ||||
self.dataOut.errorCount = self.dataIn.errorCount | ||||
self.dataOut.useLocalTime = self.dataIn.useLocalTime | ||||
|
r897 | |||
|
r502 | 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 | ||||
|
r568 | self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock | ||
|
r502 | 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 | ||||
|
r762 | self.dataOut.nCohInt = self.dataIn.nCohInt | ||
|
r502 | # self.dataOut.nIncohInt = 1 | ||
self.dataOut.ippSeconds = self.dataIn.ippSeconds | ||||
|
r854 | self.dataOut.nProfiles = self.dataIn.nProfiles | ||
self.dataOut.utctime = self.dataIn.utctime | ||||
|
r502 | # self.dataOut.windowOfFilter = self.dataIn.windowOfFilter | ||
|
r897 | |||
|
r587 | # self.dataOut.timeInterval = self.dataIn.timeInterval*self.dataOut.nPoints | ||
|
r897 | |||
|
r502 | def removeDC(self, jspectra): | ||
|
r897 | |||
|
r502 | nChannel = jspectra.shape[0] | ||
|
r897 | |||
|
r502 | for i in range(nChannel): | ||
jspectra_tmp = jspectra[i,:,:] | ||||
jspectra_DC = numpy.mean(jspectra_tmp,axis = 0) | ||||
|
r897 | |||
|
r502 | jspectra_tmp = jspectra_tmp - jspectra_DC | ||
jspectra[i,:,:] = jspectra_tmp | ||||
|
r897 | |||
|
r502 | return jspectra | ||
|
r897 | |||
|
r502 | def removeNoise(self, mode = 2): | ||
indR = numpy.where(self.dataOut.lagR == 0)[0][0] | ||||
indT = numpy.where(self.dataOut.lagT == 0)[0][0] | ||||
|
r897 | |||
|
r502 | jspectra = self.dataOut.data_corr[:,:,indR,:] | ||
|
r897 | |||
|
r502 | num_chan = jspectra.shape[0] | ||
num_hei = jspectra.shape[2] | ||||
|
r897 | |||
|
r502 | freq_dc = indT | ||
ind_vel = numpy.array([-2,-1,1,2]) + freq_dc | ||||
|
r897 | |||
|
r502 | NPot = self.dataOut.getNoise(mode) | ||
|
r897 | jspectra[:,freq_dc,:] = jspectra[:,freq_dc,:] - NPot | ||
|
r502 | SPot = jspectra[:,freq_dc,:] | ||
pairsAutoCorr = self.dataOut.getPairsAutoCorr() | ||||
# self.dataOut.signalPotency = SPot | ||||
self.dataOut.noise = NPot | ||||
self.dataOut.SNR = (SPot/NPot)[pairsAutoCorr] | ||||
|
r897 | self.dataOut.data_corr[:,:,indR,:] = jspectra | ||
return 1 | ||||
|
r854 | def run(self, lags=None, mode = 'time', pairsList=None, fullBuffer=False, nAvg = 1, removeDC = False, splitCF=False): | ||
|
r897 | |||
|
r502 | self.dataOut.flagNoData = True | ||
|
r897 | |||
|
r502 | if self.dataIn.type == "Correlation": | ||
|
r897 | |||
|
r502 | self.dataOut.copy(self.dataIn) | ||
|
r897 | |||
|
r502 | return | ||
|
r897 | |||
|
r502 | if self.dataIn.type == "Voltage": | ||
|
r897 | |||
|
r854 | nChannels = self.dataIn.nChannels | ||
|
r897 | nProfiles = self.dataIn.nProfiles | ||
|
r854 | nHeights = self.dataIn.nHeights | ||
data_pre = self.dataIn.data | ||||
|
r897 | |||
|
r854 | #--------------- Remover DC ------------ | ||
if removeDC: | ||||
data_pre = self.removeDC(data_pre) | ||||
|
r897 | |||
|
r854 | #--------------------------------------------- | ||
# pairsList = list(ccfList) | ||||
# for i in acfList: | ||||
# pairsList.append((i,i)) | ||||
|
r897 | # | ||
|
r854 | # ccf_pairs = numpy.arange(len(ccfList)) | ||
# acf_pairs = numpy.arange(len(ccfList),len(pairsList)) | ||||
self.__updateObjFromVoltage() | ||||
#---------------------------------------------------------------------- | ||||
#Creating temporal buffers | ||||
if fullBuffer: | ||||
|
r897 | tmp = numpy.zeros((len(pairsList), len(lags), nProfiles, nHeights), dtype = 'complex')*numpy.nan | ||
|
r854 | elif mode == 'time': | ||
if lags == None: | ||||
lags = numpy.arange(-nProfiles+1, nProfiles) | ||||
tmp = numpy.zeros((len(pairsList), len(lags), nHeights),dtype='complex') | ||||
elif mode == 'height': | ||||
if lags == None: | ||||
lags = numpy.arange(-nHeights+1, nHeights) | ||||
tmp = numpy.zeros(len(pairsList), (len(lags), nProfiles),dtype='complex') | ||||
|
r897 | |||
|
r854 | #For loop | ||
for l in range(len(pairsList)): | ||||
|
r897 | |||
|
r854 | ch0 = pairsList[l][0] | ||
|
r897 | ch1 = pairsList[l][1] | ||
|
r854 | for i in range(len(lags)): | ||
idx = lags[i] | ||||
|
r897 | |||
|
r854 | if idx >= 0: | ||
if mode == 'time': | ||||
ccf0 = data_pre[ch0,:nProfiles-idx,:]*numpy.conj(data_pre[ch1,idx:,:]) #time | ||||
else: | ||||
ccf0 = data_pre[ch0,:,nHeights-idx]*numpy.conj(data_pre[ch1,:,idx:]) #heights | ||||
else: | ||||
if mode == 'time': | ||||
ccf0 = data_pre[ch0,-idx:,:]*numpy.conj(data_pre[ch1,:nProfiles+idx,:]) #time | ||||
else: | ||||
ccf0 = data_pre[ch0,:,-idx:]*numpy.conj(data_pre[ch1,:,:nHeights+idx]) #heights | ||||
|
r897 | |||
|
r854 | if fullBuffer: | ||
tmp[l,i,:ccf0.shape[0],:] = ccf0 | ||||
else: | ||||
tmp[l,i,:] = numpy.sum(ccf0, axis=0) | ||||
|
r897 | |||
|
r854 | #----------------------------------------------------------------- | ||
if fullBuffer: | ||||
tmp = numpy.sum(numpy.reshape(tmp,(tmp.shape[0],tmp.shape[1],tmp.shape[2]/nAvg,nAvg,tmp.shape[3])),axis=3) | ||||
self.dataOut.nAvg = nAvg | ||||
|
r897 | |||
|
r854 | self.dataOut.data_cf = tmp | ||
self.dataOut.mode = mode | ||||
self.dataOut.nLags = len(lags) | ||||
self.dataOut.pairsList = pairsList | ||||
self.dataOut.nPairs = len(pairsList) | ||||
|
r897 | |||
|
r854 | #Se Calcula los factores de Normalizacion | ||
if mode == 'time': | ||||
delta = self.dataIn.ippSeconds*self.dataIn.nCohInt | ||||
else: | ||||
delta = self.dataIn.heightList[1] - self.dataIn.heightList[0] | ||||
self.dataOut.lagRange = numpy.array(lags)*delta | ||||
|
r897 | # self.dataOut.nCohInt = self.dataIn.nCohInt*nAvg | ||
|
r854 | self.dataOut.flagNoData = False | ||
|
r989 | # a = self.dataOut.normFactor | ||
|
r897 | return | ||