SpectraProcessor.py
619 lines
| 20.2 KiB
| text/x-python
|
PythonLexer
|
r10 | ''' | |
Created on Feb 7, 2012 | |||
|
r16 | @author $Author$ | |
@version $Id$ | |||
|
r10 | ''' | |
|
r72 | import os, sys | |
import numpy | |||
path = os.path.split(os.getcwd())[0] | |||
sys.path.append(path) | |||
from Model.Spectra import Spectra | |||
from IO.SpectraIO import SpectraWriter | |||
from Graphics.SpectraPlot import Spectrum | |||
|
r10 | ||
class SpectraProcessor: | |||
''' | |||
classdocs | |||
''' | |||
|
r72 | ||
|
r99 | dataInObj = None | |
dataOutObj = None | |||
integratorObjIndex = None | |||
decoderObjIndex = None | |||
writerObjIndex = None | |||
plotterObjIndex = None | |||
integratorObjList = [] | |||
decoderObjList = [] | |||
writerObjList = [] | |||
plotterObjList = [] | |||
buffer = None | |||
ptsId = 0 | |||
nFFTPoints = None | |||
pairList = None | |||
|
r103 | def __init__(self, dataInObj=None, dataOutObj=None): | |
|
r10 | ''' | |
Constructor | |||
''' | |||
|
r85 | self.dataInObj = dataInObj | |
|
r72 | ||
|
r85 | if dataOutObj == None: | |
self.dataOutObj = Spectra() | |||
|
r72 | else: | |
|
r85 | self.dataOutObj = dataOutObj | |
|
r72 | ||
|
r99 | self.integratorObjIndex = None | |
self.decoderObjIndex = None | |||
self.writerObjIndex = None | |||
self.plotterObjIndex = None | |||
|
r72 | ||
|
r99 | self.integratorObjList = [] | |
self.decoderObjList = [] | |||
self.writerObjList = [] | |||
self.plotterObjList = [] | |||
|
r72 | ||
self.buffer = None | |||
self.ptsId = 0 | |||
|
r103 | def setIO(self,inputObject, outputObject): | |
|
r85 | ||
|
r103 | # if not( isinstance(inputObject, Voltage) ): | |
# print 'InputObject must be an instance from Voltage()' | |||
# sys.exit(0) | |||
|
r72 | ||
|
r103 | if not( isinstance(outputObject, Spectra) ): | |
print 'OutputObject must be an instance from Spectra()' | |||
sys.exit(0) | |||
self.dataInObj = inputObject | |||
self.dataOutObj = outputObject | |||
def setup(self,nFFTPoints=None, pairList=None): | |||
|
r85 | if nFFTPoints == None: | |
|
r89 | nFFTPoints = self.dataOutObj.nFFTPoints | |
|
r85 | ||
self.nFFTPoints = nFFTPoints | |||
self.pairList = pairList | |||
|
r103 | ||
# def init(self, nFFTPoints, pairList=None): | |||
def init(self): | |||
self.integratorObjIndex = 0 | |||
self.decoderObjIndex = 0 | |||
self.writerObjIndex = 0 | |||
self.plotterObjIndex = 0 | |||
# if nFFTPoints == None: | |||
# nFFTPoints = self.dataOutObj.nFFTPoints | |||
# | |||
# self.nFFTPoints = nFFTPoints | |||
# self.pairList = pairList | |||
# | |||
|
r85 | if not( isinstance(self.dataInObj, Spectra) ): | |
self.__getFft() | |||
|
r72 | else: | |
|
r85 | self.dataOutObj.copy(self.dataInObj) | |
|
r72 | ||
|
r85 | def __getFft(self): | |
""" | |||
Convierte valores de Voltaje a Spectra | |||
Affected: | |||
self.dataOutObj.data_spc | |||
self.dataOutObj.data_cspc | |||
self.dataOutObj.data_dc | |||
self.dataOutObj.heightList | |||
self.dataOutObj.m_BasicHeader | |||
self.dataOutObj.m_ProcessingHeader | |||
self.dataOutObj.m_RadarControllerHeader | |||
self.dataOutObj.m_SystemHeader | |||
self.ptsId | |||
self.buffer | |||
self.dataOutObj.flagNoData | |||
self.dataOutObj.dataType | |||
self.dataOutObj.nPairs | |||
self.dataOutObj.nChannels | |||
self.dataOutObj.nProfiles | |||
self.dataOutObj.m_SystemHeader.numChannels | |||
self.dataOutObj.m_ProcessingHeader.totalSpectra | |||
self.dataOutObj.m_ProcessingHeader.profilesPerBlock | |||
self.dataOutObj.m_ProcessingHeader.numHeights | |||
self.dataOutObj.m_ProcessingHeader.spectraComb | |||
self.dataOutObj.m_ProcessingHeader.shif_fft | |||
""" | |||
|
r104 | if self.dataInObj.flagNoData: | |
return 0 | |||
|
r85 | blocksize = 0 | |
|
r89 | nFFTPoints = self.nFFTPoints | |
nChannels, nheis = self.dataInObj.data.shape | |||
|
r72 | ||
if self.buffer == None: | |||
|
r89 | self.buffer = numpy.zeros((nChannels, nFFTPoints, nheis), dtype='complex') | |
|
r73 | ||
|
r85 | self.buffer[:,self.ptsId,:] = self.dataInObj.data | |
|
r72 | self.ptsId += 1 | |
|
r104 | if self.ptsId < self.nFFTPoints: | |
|
r85 | self.dataOutObj.flagNoData = True | |
return | |||
|
r72 | ||
|
r85 | fft_volt = numpy.fft.fft(self.buffer,axis=1) | |
dc = fft_volt[:,0,:] | |||
#calculo de self-spectra | |||
fft_volt = numpy.fft.fftshift(fft_volt,axes=(1,)) | |||
spc = numpy.abs(fft_volt * numpy.conjugate(fft_volt)) | |||
blocksize += dc.size | |||
blocksize += spc.size | |||
cspc = None | |||
|
r89 | nPair = 0 | |
|
r85 | if self.pairList != None: | |
|
r72 | #calculo de cross-spectra | |
|
r89 | nPairs = len(self.pairList) | |
cspc = numpy.zeros((nPairs, nFFTPoints, nheis), dtype='complex') | |||
|
r85 | for pair in self.pairList: | |
|
r89 | cspc[nPair,:,:] = numpy.abs(fft_volt[pair[0],:,:] * numpy.conjugate(fft_volt[pair[1],:,:])) | |
nPair += 1 | |||
|
r85 | blocksize += cspc.size | |
self.dataOutObj.data_spc = spc | |||
self.dataOutObj.data_cspc = cspc | |||
self.dataOutObj.data_dc = dc | |||
self.ptsId = 0 | |||
self.buffer = None | |||
self.dataOutObj.flagNoData = False | |||
self.dataOutObj.heightList = self.dataInObj.heightList | |||
self.dataOutObj.channelList = self.dataInObj.channelList | |||
self.dataOutObj.m_BasicHeader = self.dataInObj.m_BasicHeader.copy() | |||
self.dataOutObj.m_ProcessingHeader = self.dataInObj.m_ProcessingHeader.copy() | |||
self.dataOutObj.m_RadarControllerHeader = self.dataInObj.m_RadarControllerHeader.copy() | |||
self.dataOutObj.m_SystemHeader = self.dataInObj.m_SystemHeader.copy() | |||
self.dataOutObj.dataType = self.dataInObj.dataType | |||
|
r89 | self.dataOutObj.nPairs = nPair | |
self.dataOutObj.nChannels = nChannels | |||
self.dataOutObj.nProfiles = nFFTPoints | |||
|
r85 | self.dataOutObj.nHeights = nheis | |
|
r89 | self.dataOutObj.nFFTPoints = nFFTPoints | |
|
r85 | #self.dataOutObj.data = None | |
|
r89 | self.dataOutObj.m_SystemHeader.numChannels = nChannels | |
self.dataOutObj.m_SystemHeader.nProfiles = nFFTPoints | |||
|
r85 | ||
self.dataOutObj.m_ProcessingHeader.blockSize = blocksize | |||
|
r89 | self.dataOutObj.m_ProcessingHeader.totalSpectra = nChannels + nPair | |
self.dataOutObj.m_ProcessingHeader.profilesPerBlock = nFFTPoints | |||
|
r85 | self.dataOutObj.m_ProcessingHeader.numHeights = nheis | |
self.dataOutObj.m_ProcessingHeader.shif_fft = True | |||
|
r89 | spectraComb = numpy.zeros( (nChannels+nPair)*2,numpy.dtype('u1')) | |
|
r85 | k = 0 | |
|
r89 | for i in range( 0,nChannels*2,2 ): | |
|
r85 | spectraComb[i] = k | |
spectraComb[i+1] = k | |||
k += 1 | |||
k *= 2 | |||
if self.pairList != None: | |||
|
r72 | ||
|
r85 | for pair in self.pairList: | |
spectraComb[k] = pair[0] | |||
spectraComb[k+1] = pair[1] | |||
k += 2 | |||
|
r72 | ||
|
r85 | self.dataOutObj.m_ProcessingHeader.spectraComb = spectraComb | |
|
r73 | ||
|
r72 | def addWriter(self,wrpath): | |
|
r85 | objWriter = SpectraWriter(self.dataOutObj) | |
|
r72 | objWriter.setup(wrpath) | |
|
r99 | self.writerObjList.append(objWriter) | |
|
r72 | ||
|
r73 | def addPlotter(self, index=None): | |
if index==None: | |||
|
r99 | index = self.plotterObjIndex | |
|
r72 | ||
|
r85 | plotObj = Spectrum(self.dataOutObj, index) | |
|
r99 | self.plotterObjList.append(plotObj) | |
|
r72 | ||
|
r104 | def addIntegrator(self,N,timeInterval): | |
|
r72 | ||
|
r104 | objIncohInt = IncoherentIntegration(N,timeInterval) | |
|
r99 | self.integratorObjList.append(objIncohInt) | |
|
r72 | ||
|
r85 | def writeData(self, wrpath): | |
if self.dataOutObj.flagNoData: | |||
|
r72 | return 0 | |
|
r99 | if len(self.writerObjList) <= self.writerObjIndex: | |
|
r72 | self.addWriter(wrpath) | |
|
r99 | self.writerObjList[self.writerObjIndex].putData() | |
|
r72 | ||
|
r99 | self.writerObjIndex += 1 | |
|
r72 | ||
|
r73 | def plotData(self,xmin=None, xmax=None, ymin=None, ymax=None, winTitle='', index=None): | |
|
r85 | if self.dataOutObj.flagNoData: | |
|
r72 | return 0 | |
|
r99 | if len(self.plotterObjList) <= self.plotterObjIndex: | |
|
r73 | self.addPlotter(index) | |
|
r72 | ||
|
r99 | self.plotterObjList[self.plotterObjIndex].plotData(xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,winTitle=winTitle) | |
|
r72 | ||
|
r99 | self.plotterObjIndex += 1 | |
|
r72 | ||
|
r104 | def integrator(self, N=None, timeInterval=None): | |
|
r85 | if self.dataOutObj.flagNoData: | |
|
r72 | return 0 | |
|
r99 | if len(self.integratorObjList) <= self.integratorObjIndex: | |
|
r104 | self.addIntegrator(N,timeInterval) | |
|
r72 | ||
|
r104 | myIncohIntObj = self.integratorObjList[self.integratorObjIndex] | |
myIncohIntObj.exe(data=self.dataOutObj.data_spc,timeOfData=self.dataOutObj.m_BasicHeader.utc) | |||
|
r72 | ||
|
r104 | if myIncohIntObj.isReady: | |
self.dataOutObj.data_spc = myIncohIntObj.data | |||
self.dataOutObj.nAvg = myIncohIntObj.navg | |||
self.dataOutObj.m_ProcessingHeader.incoherentInt *= myIncohIntObj.navg | |||
#print "myIncohIntObj.navg: ",myIncohIntObj.navg | |||
|
r85 | self.dataOutObj.flagNoData = False | |
|
r104 | ||
|
r72 | else: | |
|
r85 | self.dataOutObj.flagNoData = True | |
|
r72 | ||
|
r99 | self.integratorObjIndex += 1 | |
|
r85 | ||
def removeDC(self, type): | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
pass | |||
def removeInterference(self): | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
pass | |||
def removeSatellites(self): | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
pass | |||
def selectChannels(self, channelList, pairList=None): | |||
""" | |||
Selecciona un bloque de datos en base a canales y pares segun el channelList y el pairList | |||
Input: | |||
channelList : lista sencilla de canales a seleccionar por ej. (2,3,7) | |||
pairList : tupla de pares que se desea selecionar por ej. ( (0,1), (0,2) ) | |||
Affected: | |||
self.dataOutObj.data_spc | |||
self.dataOutObj.data_cspc | |||
self.dataOutObj.data_dc | |||
self.dataOutObj.nChannels | |||
self.dataOutObj.nPairs | |||
self.dataOutObj.m_ProcessingHeader.spectraComb | |||
self.dataOutObj.m_SystemHeader.numChannels | |||
Return: | |||
None | |||
""" | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
|
r96 | channelIndexList = [] | |
for channel in channelList: | |||
if channel in self.dataOutObj.channelList: | |||
index = self.dataOutObj.channelList.index(channel) | |||
channelIndexList.append(index) | |||
continue | |||
raise ValueError, "The value %d in channelList is not valid" %channel | |||
|
r89 | nProfiles = self.dataOutObj.nProfiles | |
|
r96 | #dataType = self.dataOutObj.dataType | |
nHeights = self.dataOutObj.nHeights #m_ProcessingHeader.numHeights | |||
|
r85 | blocksize = 0 | |
#self spectra | |||
|
r96 | nChannels = len(channelIndexList) | |
spc = numpy.zeros( (nChannels,nProfiles,nHeights), dtype='float' ) #dataType[0] ) | |||
|
r72 | ||
|
r96 | for index, channel in enumerate(channelIndexList): | |
spc[index,:,:] = self.dataOutObj.data_spc[index,:,:] | |||
|
r85 | ||
#DC channel | |||
|
r89 | dc = numpy.zeros( (nChannels,nHeights), dtype='complex' ) | |
|
r96 | for index, channel in enumerate(channelIndexList): | |
|
r85 | dc[index,:] = self.dataOutObj.data_dc[channel,:] | |
blocksize += dc.size | |||
blocksize += spc.size | |||
|
r89 | nPairs = 0 | |
|
r85 | cspc = None | |
if pairList == None: | |||
pairList = self.pairList | |||
|
r96 | if (pairList != None) and (self.dataOutObj.data_cspc != None): | |
|
r85 | #cross spectra | |
|
r89 | nPairs = len(pairList) | |
cspc = numpy.zeros( (nPairs,nProfiles,nHeights), dtype='complex' ) | |||
|
r85 | ||
spectraComb = self.dataOutObj.m_ProcessingHeader.spectraComb | |||
totalSpectra = len(spectraComb) | |||
nchan = self.dataOutObj.nChannels | |||
|
r96 | pairIndexList = [] | |
|
r85 | ||
for pair in pairList: #busco el par en la lista de pares del Spectra Combinations | |||
for index in range(0,totalSpectra,2): | |||
if pair[0] == spectraComb[index] and pair[1] == spectraComb[index+1]: | |||
|
r96 | pairIndexList.append( index/2 - nchan ) | |
|
r85 | ||
|
r96 | for index, pair in enumerate(pairIndexList): | |
|
r85 | cspc[index,:,:] = self.dataOutObj.data_cspc[pair,:,:] | |
blocksize += cspc.size | |||
else: | |||
pairList = self.pairList | |||
cspc = self.dataOutObj.data_cspc | |||
if cspc != None: | |||
blocksize += cspc.size | |||
|
r89 | spectraComb = numpy.zeros( (nChannels+nPairs)*2,numpy.dtype('u1')) | |
|
r85 | i = 0 | |
for val in channelList: | |||
spectraComb[i] = val | |||
spectraComb[i+1] = val | |||
i += 2 | |||
if pairList != None: | |||
for pair in pairList: | |||
spectraComb[i] = pair[0] | |||
spectraComb[i+1] = pair[1] | |||
i += 2 | |||
self.dataOutObj.data_spc = spc | |||
self.dataOutObj.data_cspc = cspc | |||
self.dataOutObj.data_dc = dc | |||
|
r89 | self.dataOutObj.nChannels = nChannels | |
self.dataOutObj.nPairs = nPairs | |||
|
r85 | ||
self.dataOutObj.channelList = channelList | |||
self.dataOutObj.m_ProcessingHeader.spectraComb = spectraComb | |||
|
r89 | self.dataOutObj.m_ProcessingHeader.totalSpectra = nChannels + nPairs | |
self.dataOutObj.m_SystemHeader.numChannels = nChannels | |||
self.dataOutObj.nChannels = nChannels | |||
|
r85 | self.dataOutObj.m_ProcessingHeader.blockSize = blocksize | |
def selectHeightsByValue(self, minHei, maxHei): | |||
""" | |||
Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango | |||
minHei <= height <= maxHei | |||
Input: | |||
minHei : valor minimo de altura a considerar | |||
maxHei : valor maximo de altura a considerar | |||
Affected: | |||
Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex | |||
Return: | |||
None | |||
""" | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
|
r96 | if (minHei < self.dataOutObj.heightList[0]) or (minHei > maxHei): | |
raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei) | |||
if (maxHei > self.dataOutObj.heightList[-1]): | |||
raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei) | |||
|
r85 | minIndex = 0 | |
maxIndex = 0 | |||
data = self.dataOutObj.heightList | |||
for i,val in enumerate(data): | |||
if val < minHei: | |||
continue | |||
else: | |||
minIndex = i; | |||
break | |||
for i,val in enumerate(data): | |||
if val <= maxHei: | |||
maxIndex = i; | |||
else: | |||
break | |||
self.selectHeightsByIndex(minIndex, maxIndex) | |||
|
r104 | ||
|
r85 | def selectHeightsByIndex(self, minIndex, maxIndex): | |
""" | |||
Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango | |||
minIndex <= index <= maxIndex | |||
Input: | |||
minIndex : valor minimo de altura a considerar | |||
maxIndex : valor maximo de altura a considerar | |||
Affected: | |||
self.dataOutObj.data_spc | |||
self.dataOutObj.data_cspc | |||
self.dataOutObj.data_dc | |||
self.dataOutObj.heightList | |||
self.dataOutObj.nHeights | |||
self.dataOutObj.m_ProcessingHeader.numHeights | |||
self.dataOutObj.m_ProcessingHeader.blockSize | |||
self.dataOutObj.m_ProcessingHeader.firstHeight | |||
self.dataOutObj.m_RadarControllerHeader.numHeights | |||
Return: | |||
None | |||
""" | |||
if self.dataOutObj.flagNoData: | |||
return 0 | |||
|
r96 | if (minIndex < 0) or (minIndex > maxIndex): | |
raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex) | |||
if (maxIndex >= self.dataOutObj.nHeights): | |||
raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex) | |||
|
r89 | nChannels = self.dataOutObj.nChannels | |
nPairs = self.dataOutObj.nPairs | |||
nProfiles = self.dataOutObj.nProfiles | |||
|
r85 | dataType = self.dataOutObj.dataType | |
|
r96 | nHeights = maxIndex - minIndex + 1 | |
|
r85 | blockSize = 0 | |
#self spectra | |||
|
r96 | spc = self.dataOutObj.data_spc[:,:,minIndex:maxIndex+1] | |
blockSize += spc.size | |||
|
r85 | ||
#cross spectra | |||
|
r96 | cspc = None | |
if self.dataOutObj.data_cspc != None: | |||
cspc = self.dataOutObj.data_cspc[:,:,minIndex:maxIndex+1] | |||
blockSize += cspc.size | |||
|
r85 | ||
#DC channel | |||
|
r96 | dc = self.dataOutObj.data_dc[:,minIndex:maxIndex+1] | |
blockSize += dc.size | |||
|
r85 | ||
self.dataOutObj.data_spc = spc | |||
|
r96 | if cspc != None: | |
self.dataOutObj.data_cspc = cspc | |||
|
r85 | self.dataOutObj.data_dc = dc | |
firstHeight = self.dataOutObj.heightList[minIndex] | |||
|
r96 | self.dataOutObj.nHeights = nHeights | |
self.dataOutObj.m_ProcessingHeader.blockSize = blockSize | |||
self.dataOutObj.m_ProcessingHeader.numHeights = nHeights | |||
|
r85 | self.dataOutObj.m_ProcessingHeader.firstHeight = firstHeight | |
|
r96 | self.dataOutObj.m_RadarControllerHeader.numHeights = nHeights | |
|
r85 | ||
|
r96 | self.dataOutObj.heightList = self.dataOutObj.heightList[minIndex:maxIndex+1] | |
|
r85 | ||
|
r72 | class IncoherentIntegration: | |
|
r99 | ||
|
r104 | integ_counter = None | |
|
r99 | data = None | |
|
r104 | navg = None | |
|
r99 | buffer = None | |
nIncohInt = None | |||
|
r104 | def __init__(self, N = None, timeInterval = None): | |
""" | |||
N | |||
timeInterval - interval time [min], integer value | |||
""" | |||
|
r99 | ||
|
r72 | self.data = None | |
|
r104 | self.navg = None | |
|
r72 | self.buffer = None | |
|
r104 | self.timeOut = None | |
self.exitCondition = False | |||
self.isReady = False | |||
|
r72 | self.nIncohInt = N | |
|
r104 | self.integ_counter = 0 | |
if timeInterval!=None: | |||
self.timeIntervalInSeconds = timeInterval * 60. #if (type(timeInterval)!=integer) -> change this line | |||
if ((timeInterval==None) and (N==None)): | |||
print 'N = None ; timeInterval = None' | |||
sys.exit(0) | |||
elif timeInterval == None: | |||
self.timeFlag = False | |||
else: | |||
self.timeFlag = True | |||
|
r72 | ||
|
r104 | def exe(self,data,timeOfData): | |
""" | |||
data | |||
timeOfData [seconds] | |||
""" | |||
|
r72 | ||
|
r104 | if self.timeFlag: | |
if self.timeOut == None: | |||
self.timeOut = timeOfData + self.timeIntervalInSeconds | |||
if timeOfData < self.timeOut: | |||
if self.buffer == None: | |||
self.buffer = data | |||
else: | |||
self.buffer = self.buffer + data | |||
self.integ_counter += 1 | |||
else: | |||
self.exitCondition = True | |||
|
r72 | else: | |
|
r104 | if self.integ_counter < self.nIncohInt: | |
if self.buffer == None: | |||
self.buffer = data | |||
else: | |||
self.buffer = self.buffer + data | |||
self.integ_counter += 1 | |||
if self.integ_counter == self.nIncohInt: | |||
self.exitCondition = True | |||
if self.exitCondition: | |||
|
r72 | self.data = self.buffer | |
|
r104 | self.navg = self.integ_counter | |
self.isReady = True | |||
|
r72 | self.buffer = None | |
|
r104 | self.timeOut = None | |
self.integ_counter = 0 | |||
self.exitCondition = False | |||
|
r72 | ||
|
r104 | if self.timeFlag: | |
self.buffer = data | |||
self.timeOut = timeOfData + self.timeIntervalInSeconds | |||
else: | |||
self.isReady = False | |||