##// END OF EJS Templates
AMISR modules to read hdf5 files and link to SignalChain Objects...
Daniel Valdez -
r491:3acb3a6f21eb
parent child
Show More
@@ -0,0 +1,59
1 import numpy
2
3 class AMISR:
4 def __init__(self):
5 self.flagNoData = True
6 self.data = None
7 self.utctime = None
8 self.type = "AMISR"
9
10 #propiedades para compatibilidad con Voltages
11 self.timeZone = 300#timezone like jroheader, difference in minutes between UTC and localtime
12 self.dstFlag = 0#self.dataIn.dstFlag
13 self.errorCount = 0#self.dataIn.errorCount
14 self.useLocalTime = True#self.dataIn.useLocalTime
15
16 self.radarControllerHeaderObj = None#self.dataIn.radarControllerHeaderObj.copy()
17 self.systemHeaderObj = None#self.dataIn.systemHeaderObj.copy()
18 self.channelList = [0]#self.dataIn.channelList esto solo aplica para el caso de AMISR
19 self.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
20
21 self.flagTimeBlock = None#self.dataIn.flagTimeBlock
22 #self.utctime = #self.firstdatatime
23 self.flagDecodeData = None#self.dataIn.flagDecodeData #asumo q la data esta decodificada
24 self.flagDeflipData = None#self.dataIn.flagDeflipData #asumo q la data esta sin flip
25
26 self.nCohInt = 1#self.dataIn.nCohInt
27 self.nIncohInt = 1
28 self.ippSeconds = None#self.dataIn.ippSeconds, segun el filename/Setup/Tufile
29 self.windowOfFilter = None#self.dataIn.windowOfFilter
30
31 self.timeInterval = None#self.dataIn.timeInterval*self.dataOut.nFFTPoints*self.dataOut.nIncohInt
32 self.frequency = None#self.dataIn.frequency
33 self.realtime = 0#self.dataIn.realtime
34
35 #actualizar en la lectura de datos
36 self.heightList = None#self.dataIn.heightList
37 self.nProfiles = None#Number of samples or nFFTPoints
38 self.nRecords = None
39 self.nBeams = None
40 self.nBaud = None#self.dataIn.nBaud
41 self.nCode = None#self.dataIn.nCode
42 self.code = None#self.dataIn.code
43
44 #consideracion para los Beams
45 self.beamCodeDict = None
46 self.beamRangeDict = None
47
48 def copy(self, inputObj=None):
49
50 if inputObj == None:
51 return copy.deepcopy(self)
52
53 for key in inputObj.__dict__.keys():
54 self.__dict__[key] = inputObj.__dict__[key]
55
56
57 def isEmpty(self):
58
59 return self.flagNoData No newline at end of file
This diff has been collapsed as it changes many lines, (542 lines changed) Show them Hide them
@@ -0,0 +1,542
1 '''
2 @author: Daniel Suarez
3 '''
4
5 import os
6 import sys
7 import glob
8 import fnmatch
9 import datetime
10 import re
11 import h5py
12 import numpy
13
14 from model.proc.jroproc_base import ProcessingUnit, Operation
15 from model.data.jroamisr import AMISR
16
17 class RadacHeader():
18 def __init__(self, fp):
19 header = 'Raw11/Data/RadacHeader'
20 self.beamCodeByPulse = fp.get(header+'/BeamCode')
21 self.beamCode = fp.get('Raw11/Data/Beamcodes')
22 self.code = fp.get(header+'/Code')
23 self.frameCount = fp.get(header+'/FrameCount')
24 self.modeGroup = fp.get(header+'/ModeGroup')
25 self.nsamplesPulse = fp.get(header+'/NSamplesPulse')
26 self.pulseCount = fp.get(header+'/PulseCount')
27 self.radacTime = fp.get(header+'/RadacTime')
28 self.timeCount = fp.get(header+'/TimeCount')
29 self.timeStatus = fp.get(header+'/TimeStatus')
30
31 self.nrecords = self.pulseCount.shape[0] #nblocks
32 self.npulses = self.pulseCount.shape[1] #nprofile
33 self.nsamples = self.nsamplesPulse[0,0] #ngates
34 self.nbeams = self.beamCode.shape[1]
35
36
37 def getIndexRangeToPulse(self, idrecord=0):
38 indexToZero = numpy.where(self.pulseCount.value[idrecord,:]==0)
39 startPulseCountId = indexToZero[0][0]
40 endPulseCountId = startPulseCountId - 1
41 range1 = numpy.arange(startPulseCountId,self.npulses,1)
42 range2 = numpy.arange(0,startPulseCountId,1)
43 return range1, range2
44
45
46 class AMISRReader(ProcessingUnit):
47
48 path = None
49 startDate = None
50 endDate = None
51 startTime = None
52 endTime = None
53 walk = None
54 isConfig = False
55
56 def __init__(self):
57 self.set = None
58 self.subset = None
59 self.extension_file = '.h5'
60 self.dtc_str = 'dtc'
61 self.dtc_id = 0
62 self.status = True
63 self.isConfig = False
64 self.dirnameList = []
65 self.filenameList = []
66 self.fileIndex = None
67 self.flagNoMoreFiles = False
68 self.flagIsNewFile = 0
69 self.filename = ''
70 self.amisrFilePointer = None
71 self.radacHeaderObj = None
72 self.dataOut = self.__createObjByDefault()
73 self.datablock = None
74 self.rest_datablock = None
75 self.range = None
76 self.idrecord_count = 0
77 self.profileIndex = 0
78 self.idpulse_range1 = None
79 self.idpulse_range2 = None
80 self.beamCodeByFrame = None
81 self.radacTimeByFrame = None
82 #atributos originales tal y como esta en el archivo de datos
83 self.beamCodesFromFile = None
84 self.radacTimeFromFile = None
85 self.rangeFromFile = None
86 self.dataByFrame = None
87 self.dataset = None
88
89 self.beamCodeDict = {}
90 self.beamRangeDict = {}
91
92 #experiment cgf file
93 self.npulsesint_fromfile = None
94 self.recordsperfile_fromfile = None
95 self.nbeamcodes_fromfile = None
96 self.ngates_fromfile = None
97 self.ippSeconds_fromfile = None
98 self.frequency_h5file = None
99
100
101 self.__firstFile = True
102 self.buffer_radactime = None
103
104 def __createObjByDefault(self):
105
106 dataObj = AMISR()
107
108 return dataObj
109
110 def __setParameters(self,path,startDate,endDate,startTime,endTime,walk):
111 self.path = path
112 self.startDate = startDate
113 self.endDate = endDate
114 self.startTime = startTime
115 self.endTime = endTime
116 self.walk = walk
117
118 def __checkPath(self):
119 if os.path.exists(self.path):
120 self.status = 1
121 else:
122 self.status = 0
123 print 'Path:%s does not exists'%self.path
124
125 return
126
127 def __selDates(self, amisr_dirname_format):
128 try:
129 year = int(amisr_dirname_format[0:4])
130 month = int(amisr_dirname_format[4:6])
131 dom = int(amisr_dirname_format[6:8])
132 thisDate = datetime.date(year,month,dom)
133
134 if (thisDate>=self.startDate and thisDate <= self.endDate):
135 return amisr_dirname_format
136 except:
137 return None
138
139 def __findDataForDates(self):
140
141
142
143 if not(self.status):
144 return None
145
146 pat = '\d+.\d+'
147 dirnameList = [re.search(pat,x) for x in os.listdir(self.path)]
148 dirnameList = filter(lambda x:x!=None,dirnameList)
149 dirnameList = [x.string for x in dirnameList]
150 dirnameList = [self.__selDates(x) for x in dirnameList]
151 dirnameList = filter(lambda x:x!=None,dirnameList)
152 if len(dirnameList)>0:
153 self.status = 1
154 self.dirnameList = dirnameList
155 self.dirnameList.sort()
156 else:
157 self.status = 0
158 return None
159
160 def __getTimeFromData(self):
161 pass
162
163 def __filterByGlob1(self, dirName):
164 filter_files = glob.glob1(dirName, '*.*%s'%self.extension_file)
165 filterDict = {}
166 filterDict.setdefault(dirName)
167 filterDict[dirName] = filter_files
168 return filterDict
169
170 def __getFilenameList(self, fileListInKeys, dirList):
171 for value in fileListInKeys:
172 dirName = value.keys()[0]
173 for file in value[dirName]:
174 filename = os.path.join(dirName, file)
175 self.filenameList.append(filename)
176
177
178 def __selectDataForTimes(self):
179 #aun no esta implementado el filtro for tiempo
180 if not(self.status):
181 return None
182
183 dirList = [os.path.join(self.path,x) for x in self.dirnameList]
184
185 fileListInKeys = [self.__filterByGlob1(x) for x in dirList]
186
187 self.__getFilenameList(fileListInKeys, dirList)
188
189 if len(self.filenameList)>0:
190 self.status = 1
191 self.filenameList.sort()
192 else:
193 self.status = 0
194 return None
195
196
197 def __searchFilesOffline(self,
198 path,
199 startDate,
200 endDate,
201 startTime=datetime.time(0,0,0),
202 endTime=datetime.time(23,59,59),
203 walk=True):
204
205 self.__setParameters(path, startDate, endDate, startTime, endTime, walk)
206
207 self.__checkPath()
208
209 self.__findDataForDates()
210
211 self.__selectDataForTimes()
212
213 for i in range(len(self.filenameList)):
214 print "%s" %(self.filenameList[i])
215
216 return
217
218 def __setNextFileOffline(self):
219 idFile = self.fileIndex
220
221 while (True):
222 idFile += 1
223 if not(idFile < len(self.filenameList)):
224 self.flagNoMoreFiles = 1
225 print "No more Files"
226 return 0
227
228 filename = self.filenameList[idFile]
229
230 amisrFilePointer = h5py.File(filename,'r')
231
232 break
233
234 self.flagIsNewFile = 1
235 self.fileIndex = idFile
236 self.filename = filename
237
238 self.amisrFilePointer = amisrFilePointer
239
240 print "Setting the file: %s"%self.filename
241
242 return 1
243
244 def __readHeader(self):
245 self.radacHeaderObj = RadacHeader(self.amisrFilePointer)
246
247 #update values from experiment cfg file
248 if self.radacHeaderObj.nrecords == self.recordsperfile_fromfile:
249 self.radacHeaderObj.nrecords = self.recordsperfile_fromfile
250 self.radacHeaderObj.nbeams = self.nbeamcodes_fromfile
251 self.radacHeaderObj.npulses = self.npulsesint_fromfile
252 self.radacHeaderObj.nsamples = self.ngates_fromfile
253
254 #get tuning frequency
255 frequency_h5file_dataset = self.amisrFilePointer.get('Rx'+'/TuningFrequency')
256 self.frequency_h5file = frequency_h5file_dataset[0,0]
257
258 self.flagIsNewFile = 1
259
260 def __getBeamCode(self):
261 self.beamCodeDict = {}
262 self.beamRangeDict = {}
263
264 for i in range(len(self.radacHeaderObj.beamCode[0,:])):
265 self.beamCodeDict.setdefault(i)
266 self.beamRangeDict.setdefault(i)
267 self.beamCodeDict[i] = self.radacHeaderObj.beamCode[0,i]
268
269
270 just4record0 = self.radacHeaderObj.beamCodeByPulse[0,:]
271
272 for i in range(len(self.beamCodeDict.values())):
273 xx = numpy.where(just4record0==self.beamCodeDict.values()[i])
274 self.beamRangeDict[i] = xx[0]
275
276 def __getExpParameters(self):
277 if not(self.status):
278 return None
279
280 experimentCfgPath = os.path.join(self.path, self.dirnameList[0], 'Setup')
281
282 expFinder = glob.glob1(experimentCfgPath,'*.exp')
283 if len(expFinder)== 0:
284 self.status = 0
285 return None
286
287 experimentFilename = os.path.join(experimentCfgPath,expFinder[0])
288
289 f = open(experimentFilename)
290 lines = f.readlines()
291 f.close()
292
293 parmsList = ['npulsesint*','recordsperfile*','nbeamcodes*','ngates*']
294 filterList = [fnmatch.filter(lines, x) for x in parmsList]
295
296
297 values = [re.sub(r'\D',"",x[0]) for x in filterList]
298
299 self.npulsesint_fromfile = int(values[0])
300 self.recordsperfile_fromfile = int(values[1])
301 self.nbeamcodes_fromfile = int(values[2])
302 self.ngates_fromfile = int(values[3])
303
304 tufileFinder = fnmatch.filter(lines, 'tufile=*')
305 tufile = tufileFinder[0].split('=')[1].split('\n')[0]
306 tufilename = os.path.join(experimentCfgPath,tufile)
307
308 f = open(tufilename)
309 lines = f.readlines()
310 f.close()
311 self.ippSeconds_fromfile = float(lines[1].split()[2])/1E6
312
313
314 self.status = 1
315
316 def __setIdsAndArrays(self):
317 self.dataByFrame = self.__setDataByFrame()
318 self.beamCodeByFrame = self.amisrFilePointer.get('Raw11/Data/RadacHeader/BeamCode').value[0, :]
319 self.readRanges()
320 self.idpulse_range1, self.idpulse_range2 = self.radacHeaderObj.getIndexRangeToPulse(0)
321 self.radacTimeByFrame = numpy.zeros(self.radacHeaderObj.npulses)
322 self.buffer_radactime = numpy.zeros_like(self.radacTimeByFrame)
323
324
325 def __setNextFile(self):
326
327 newFile = self.__setNextFileOffline()
328
329 if not(newFile):
330 return 0
331
332 self.__readHeader()
333
334 if self.__firstFile:
335 self.__setIdsAndArrays()
336 self.__firstFile = False
337
338 self.__getBeamCode()
339 self.readDataBlock()
340
341
342 def setup(self,path=None,
343 startDate=None,
344 endDate=None,
345 startTime=datetime.time(0,0,0),
346 endTime=datetime.time(23,59,59),
347 walk=True):
348
349 #Busqueda de archivos offline
350 self.__searchFilesOffline(path, startDate, endDate, startTime, endTime, walk)
351
352 if not(self.filenameList):
353 print "There is no files into the folder: %s"%(path)
354
355 sys.exit(-1)
356
357 self.__getExpParameters()
358
359 self.fileIndex = -1
360
361 self.__setNextFile()
362
363 def readRanges(self):
364 dataset = self.amisrFilePointer.get('Raw11/Data/Samples/Range')
365 #self.rangeFromFile = dataset.value
366 self.rangeFromFile = numpy.reshape(dataset.value,(-1))
367 return range
368
369
370 def readRadacTime(self,idrecord, range1, range2):
371 self.radacTimeFromFile = self.radacHeaderObj.radacTime.value
372
373 radacTimeByFrame = numpy.zeros((self.radacHeaderObj.npulses))
374 #radacTimeByFrame = dataset[idrecord - 1,range1]
375 #radacTimeByFrame = dataset[idrecord,range2]
376
377 return radacTimeByFrame
378
379 def readBeamCode(self, idrecord, range1, range2):
380 dataset = self.amisrFilePointer.get('Raw11/Data/RadacHeader/BeamCode')
381 beamcodeByFrame = numpy.zeros((self.radacHeaderObj.npulses))
382 self.beamCodesFromFile = dataset.value
383
384 #beamcodeByFrame[range1] = dataset[idrecord - 1, range1]
385 #beamcodeByFrame[range2] = dataset[idrecord, range2]
386 beamcodeByFrame[range1] = dataset[idrecord, range1]
387 beamcodeByFrame[range2] = dataset[idrecord, range2]
388
389 return beamcodeByFrame
390
391
392 def __setDataByFrame(self):
393 ndata = 2 # porque es complejo
394 dataByFrame = numpy.zeros((self.radacHeaderObj.npulses, self.radacHeaderObj.nsamples, ndata))
395 return dataByFrame
396
397 def __readDataSet(self):
398 dataset = self.amisrFilePointer.get('Raw11/Data/Samples/Data')
399 return dataset
400
401 def __setDataBlock(self,):
402 real = self.dataByFrame[:,:,0] #asumo que 0 es real
403 imag = self.dataByFrame[:,:,1] #asumo que 1 es imaginario
404 datablock = real + imag*1j #armo el complejo
405 return datablock
406
407 def readSamples_version1(self,idrecord):
408 #estas tres primeras lineas solo se deben ejecutar una vez
409 if self.flagIsNewFile:
410 #reading dataset
411 self.dataset = self.__readDataSet()
412 self.flagIsNewFile = 0
413
414 if idrecord == 0:
415 #if self.buffer_last_record == None:
416 selectorById = self.radacHeaderObj.pulseCount[0,self.idpulse_range2]
417
418 self.dataByFrame[selectorById,:,:] = self.dataset[0, self.idpulse_range2,:,:]
419
420 self.radacTimeByFrame[selectorById] = self.radacHeaderObj.radacTime[0, self.idpulse_range2]
421
422 selectorById = self.radacHeaderObj.pulseCount[0,self.idpulse_range1]
423
424 self.radacTimeByFrame[selectorById] = self.buffer_radactime[selectorById]
425
426 datablock = self.__setDataBlock()
427
428 return datablock
429
430 selectorById = self.radacHeaderObj.pulseCount[idrecord-1,self.idpulse_range1]
431 self.dataByFrame[selectorById,:,:] = self.dataset[idrecord-1, self.idpulse_range1, :, :]
432 self.radacTimeByFrame[selectorById] = self.radacHeaderObj.radacTime[idrecord-1, self.idpulse_range1]
433
434 selectorById = self.radacHeaderObj.pulseCount[idrecord,self.idpulse_range2]#data incompleta ultimo archivo de carpeta, verifica el record real segun la dimension del arreglo de datos
435 self.dataByFrame[selectorById,:,:] = self.dataset[idrecord, self.idpulse_range2, :, :]
436 self.radacTimeByFrame[selectorById] = self.radacHeaderObj.radacTime[idrecord, self.idpulse_range2]
437
438 datablock = self.__setDataBlock()
439
440 selectorById = self.radacHeaderObj.pulseCount[idrecord,self.idpulse_range1]
441 self.dataByFrame[selectorById,:,:] = self.dataset[idrecord, self.idpulse_range1, :, :]
442 self.buffer_radactime[selectorById] = self.radacHeaderObj.radacTime[idrecord, self.idpulse_range1]
443
444 return datablock
445
446
447 def readSamples(self,idrecord):
448 if self.flagIsNewFile:
449 self.dataByFrame = self.__setDataByFrame()
450 self.beamCodeByFrame = self.amisrFilePointer.get('Raw11/Data/RadacHeader/BeamCode').value[idrecord, :]
451
452 #reading ranges
453 self.readRanges()
454 #reading dataset
455 self.dataset = self.__readDataSet()
456
457 self.flagIsNewFile = 0
458 self.radacTimeByFrame = self.radacHeaderObj.radacTime.value[idrecord, :]
459 self.dataByFrame = self.dataset[idrecord, :, :, :]
460 datablock = self.__setDataBlock()
461 return datablock
462
463
464 def readDataBlock(self):
465
466 self.datablock = self.readSamples_version1(self.idrecord_count)
467 #self.datablock = self.readSamples(self.idrecord_count)
468 #print 'record:', self.idrecord_count
469
470 self.idrecord_count += 1
471 self.profileIndex = 0
472
473 if self.idrecord_count >= self.radacHeaderObj.nrecords:
474 self.idrecord_count = 0
475 self.flagIsNewFile = 1
476
477 def readNextBlock(self):
478
479 self.readDataBlock()
480
481 if self.flagIsNewFile:
482 self.__setNextFile()
483 pass
484
485 def __hasNotDataInBuffer(self):
486 #self.radacHeaderObj.npulses debe ser otra variable para considerar el numero de pulsos a tomar en el primer y ultimo record
487 if self.profileIndex >= self.radacHeaderObj.npulses:
488 return 1
489 return 0
490
491 def printUTC(self):
492 print self.dataOut.utctime
493 print ''
494
495 def setObjProperties(self):
496 self.dataOut.heightList = self.rangeFromFile/1000.0 #km
497 self.dataOut.nProfiles = self.radacHeaderObj.npulses
498 self.dataOut.nRecords = self.radacHeaderObj.nrecords
499 self.dataOut.nBeams = self.radacHeaderObj.nbeams
500 self.dataOut.ippSeconds = self.ippSeconds_fromfile
501 self.dataOut.timeInterval = self.dataOut.ippSeconds * self.dataOut.nCohInt
502 self.dataOut.frequency = self.frequency_h5file
503 self.dataOut.nBaud = None
504 self.dataOut.nCode = None
505 self.dataOut.code = None
506
507 self.dataOut.beamCodeDict = self.beamCodeDict
508 self.dataOut.beamRangeDict = self.beamRangeDict
509
510 def getData(self):
511
512 if self.flagNoMoreFiles:
513 self.dataOut.flagNoData = True
514 print 'Process finished'
515 return 0
516
517 if self.__hasNotDataInBuffer():
518 self.readNextBlock()
519
520
521 if self.datablock == None: # setear esta condicion cuando no hayan datos por leers
522 self.dataOut.flagNoData = True
523 return 0
524
525 self.dataOut.data = numpy.reshape(self.datablock[self.profileIndex,:],(1,-1))
526
527 self.dataOut.utctime = self.radacTimeByFrame[self.profileIndex]
528
529 self.dataOut.flagNoData = False
530
531 self.profileIndex += 1
532
533 return self.dataOut.data
534
535
536 def run(self, **kwargs):
537 if not(self.isConfig):
538 self.setup(**kwargs)
539 self.setObjProperties()
540 self.isConfig = True
541
542 self.getData()
@@ -0,0 +1,82
1 '''
2 @author: Daniel Suarez
3 '''
4
5 from jroproc_base import ProcessingUnit, Operation
6 from model.data.jroamisr import AMISR
7
8 class AMISRProc(ProcessingUnit):
9 def __init__(self):
10 ProcessingUnit.__init__(self)
11 self.objectDict = {}
12 self.dataOut = AMISR()
13
14 def run(self):
15 if self.dataIn.type == 'AMISR':
16 self.dataOut.copy(self.dataIn)
17
18
19 class PrintInfo(Operation):
20 def __init__(self):
21 pass
22
23 def run(self, dataOut):
24
25 print 'Number of Records by File: %d'%dataOut.nRecords
26 print 'Number of Pulses: %d'%dataOut.nProfiles
27 print 'Number of Samples by Pulse: %d'%len(dataOut.heightList)
28 print 'Ipp Seconds: %f'%dataOut.ippSeconds
29 print 'Number of Beams: %d'%dataOut.nBeams
30 print 'BeamCodes:'
31 beamStrList = ['Beam %d -> Code %d'%(k,v) for k,v in dataOut.beamCodeDict.items()]
32 for b in beamStrList:
33 print b
34
35
36 class BeamSelector(Operation):
37 profileIndex = None
38 nProfiles = None
39
40 def __init__(self):
41
42 self.profileIndex = 0
43
44 def incIndex(self):
45 self.profileIndex += 1
46
47 if self.profileIndex >= self.nProfiles:
48 self.profileIndex = 0
49
50 def isProfileInRange(self, minIndex, maxIndex):
51
52 if self.profileIndex < minIndex:
53 return False
54
55 if self.profileIndex > maxIndex:
56 return False
57
58 return True
59
60 def isProfileInList(self, profileList):
61
62 if self.profileIndex not in profileList:
63 return False
64
65 return True
66
67 def run(self, dataOut, beam=None):
68
69 dataOut.flagNoData = True
70 self.nProfiles = dataOut.nProfiles
71
72 if beam != None:
73 if self.isProfileInList(dataOut.beamRangeDict[beam]):
74 dataOut.flagNoData = False
75
76 self.incIndex()
77 return 1
78
79 else:
80 raise ValueError, "BeamSelector needs beam value"
81
82 return 0 No newline at end of file
@@ -1,782 +1,725
1 1 '''
2 2
3 3 $Author: murco $
4 4 $Id: JROData.py 173 2012-11-20 15:06:21Z murco $
5 5 '''
6 6
7 7 import copy
8 8 import numpy
9 9 import datetime
10 10
11 11 from jroheaderIO import SystemHeader, RadarControllerHeader
12 12
13 13 def getNumpyDtype(dataTypeCode):
14 14
15 15 if dataTypeCode == 0:
16 16 numpyDtype = numpy.dtype([('real','<i1'),('imag','<i1')])
17 17 elif dataTypeCode == 1:
18 18 numpyDtype = numpy.dtype([('real','<i2'),('imag','<i2')])
19 19 elif dataTypeCode == 2:
20 20 numpyDtype = numpy.dtype([('real','<i4'),('imag','<i4')])
21 21 elif dataTypeCode == 3:
22 22 numpyDtype = numpy.dtype([('real','<i8'),('imag','<i8')])
23 23 elif dataTypeCode == 4:
24 24 numpyDtype = numpy.dtype([('real','<f4'),('imag','<f4')])
25 25 elif dataTypeCode == 5:
26 26 numpyDtype = numpy.dtype([('real','<f8'),('imag','<f8')])
27 27 else:
28 28 raise ValueError, 'dataTypeCode was not defined'
29 29
30 30 return numpyDtype
31 31
32 32 def getDataTypeCode(numpyDtype):
33 33
34 34 if numpyDtype == numpy.dtype([('real','<i1'),('imag','<i1')]):
35 35 datatype = 0
36 36 elif numpyDtype == numpy.dtype([('real','<i2'),('imag','<i2')]):
37 37 datatype = 1
38 38 elif numpyDtype == numpy.dtype([('real','<i4'),('imag','<i4')]):
39 39 datatype = 2
40 40 elif numpyDtype == numpy.dtype([('real','<i8'),('imag','<i8')]):
41 41 datatype = 3
42 42 elif numpyDtype == numpy.dtype([('real','<f4'),('imag','<f4')]):
43 43 datatype = 4
44 44 elif numpyDtype == numpy.dtype([('real','<f8'),('imag','<f8')]):
45 45 datatype = 5
46 46 else:
47 47 datatype = None
48 48
49 49 return datatype
50 50
51 51 def hildebrand_sekhon(data, navg):
52 52
53 53 data = data.copy()
54 54
55 55 sortdata = numpy.sort(data,axis=None)
56 56 lenOfData = len(sortdata)
57 57 nums_min = lenOfData/10
58 58
59 59 if (lenOfData/10) > 2:
60 60 nums_min = lenOfData/10
61 61 else:
62 62 nums_min = 2
63 63
64 64 sump = 0.
65 65
66 66 sumq = 0.
67 67
68 68 j = 0
69 69
70 70 cont = 1
71 71
72 72 while((cont==1)and(j<lenOfData)):
73 73
74 74 sump += sortdata[j]
75 75
76 76 sumq += sortdata[j]**2
77 77
78 78 j += 1
79 79
80 80 if j > nums_min:
81 81 rtest = float(j)/(j-1) + 1.0/navg
82 82 if ((sumq*j) > (rtest*sump**2)):
83 83 j = j - 1
84 84 sump = sump - sortdata[j]
85 85 sumq = sumq - sortdata[j]**2
86 86 cont = 0
87 87
88 88 lnoise = sump /j
89 89 stdv = numpy.sqrt((sumq - lnoise**2)/(j - 1))
90 90 return lnoise
91 91
92 92 class GenericData(object):
93 93
94 94 flagNoData = True
95 95
96 96 def __init__(self):
97 97
98 98 raise ValueError, "This class has not been implemented"
99 99
100 100 def copy(self, inputObj=None):
101 101
102 102 if inputObj == None:
103 103 return copy.deepcopy(self)
104 104
105 105 for key in inputObj.__dict__.keys():
106 106 self.__dict__[key] = inputObj.__dict__[key]
107 107
108 108 def deepcopy(self):
109 109
110 110 return copy.deepcopy(self)
111 111
112 112 def isEmpty(self):
113 113
114 114 return self.flagNoData
115 115
116 116 class JROData(GenericData):
117 117
118 118 # m_BasicHeader = BasicHeader()
119 119 # m_ProcessingHeader = ProcessingHeader()
120 120
121 121 systemHeaderObj = SystemHeader()
122 122
123 123 radarControllerHeaderObj = RadarControllerHeader()
124 124
125 125 # data = None
126 126
127 127 type = None
128 128
129 129 datatype = None #dtype but in string
130 130
131 131 # dtype = None
132 132
133 133 # nChannels = None
134 134
135 135 # nHeights = None
136 136
137 137 nProfiles = None
138 138
139 139 heightList = None
140 140
141 141 channelList = None
142 142
143 143 flagTimeBlock = False
144 144
145 145 useLocalTime = False
146 146
147 147 utctime = None
148 148
149 149 timeZone = None
150 150
151 151 dstFlag = None
152 152
153 153 errorCount = None
154 154
155 155 blocksize = None
156 156
157 157 nCode = None
158 158
159 159 nBaud = None
160 160
161 161 code = None
162 162
163 163 flagDecodeData = False #asumo q la data no esta decodificada
164 164
165 165 flagDeflipData = False #asumo q la data no esta sin flip
166 166
167 167 flagShiftFFT = False
168 168
169 169 # ippSeconds = None
170 170
171 171 timeInterval = None
172 172
173 173 nCohInt = None
174 174
175 175 noise = None
176 176
177 177 windowOfFilter = 1
178 178
179 179 #Speed of ligth
180 180 C = 3e8
181 181
182 182 frequency = 49.92e6
183 183
184 184 realtime = False
185 185
186 186 beacon_heiIndexList = None
187 187
188 188 last_block = None
189 189
190 190 blocknow = None
191 191
192 192 def __init__(self):
193 193
194 194 raise ValueError, "This class has not been implemented"
195 195
196 196 def getNoise(self):
197 197
198 198 raise ValueError, "Not implemented"
199 199
200 200 def getNChannels(self):
201 201
202 202 return len(self.channelList)
203 203
204 204 def getChannelIndexList(self):
205 205
206 206 return range(self.nChannels)
207 207
208 208 def getNHeights(self):
209 209
210 210 return len(self.heightList)
211 211
212 212 def getHeiRange(self, extrapoints=0):
213 213
214 214 heis = self.heightList
215 215 # deltah = self.heightList[1] - self.heightList[0]
216 216 #
217 217 # heis.append(self.heightList[-1])
218 218
219 219 return heis
220 220
221 221 def getltctime(self):
222 222
223 223 if self.useLocalTime:
224 224 return self.utctime - self.timeZone*60
225 225
226 226 return self.utctime
227 227
228 228 def getDatatime(self):
229 229
230 230 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
231 231 return datatimeValue
232 232
233 233 def getTimeRange(self):
234 234
235 235 datatime = []
236 236
237 237 datatime.append(self.ltctime)
238 238 datatime.append(self.ltctime + self.timeInterval)
239 239
240 240 datatime = numpy.array(datatime)
241 241
242 242 return datatime
243 243
244 244 def getFmax(self):
245 245
246 246 PRF = 1./(self.ippSeconds * self.nCohInt)
247 247
248 248 fmax = PRF/2.
249 249
250 250 return fmax
251 251
252 252 def getVmax(self):
253 253
254 254 _lambda = self.C/self.frequency
255 255
256 256 vmax = self.getFmax() * _lambda
257 257
258 258 return vmax
259 259
260 260 def get_ippSeconds(self):
261 261 '''
262 262 '''
263 263 return self.radarControllerHeaderObj.ippSeconds
264 264
265 265 def set_ippSeconds(self, ippSeconds):
266 266 '''
267 267 '''
268 268
269 269 self.radarControllerHeaderObj.ippSeconds = ippSeconds
270 270
271 271 return
272 272
273 273 def get_dtype(self):
274 274 '''
275 275 '''
276 276 return getNumpyDtype(self.datatype)
277 277
278 278 def set_dtype(self, numpyDtype):
279 279 '''
280 280 '''
281 281
282 282 self.datatype = getDataTypeCode(numpyDtype)
283 283
284 284 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
285 285 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
286 286 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
287 287 #noise = property(getNoise, "I'm the 'nHeights' property.")
288 288 datatime = property(getDatatime, "I'm the 'datatime' property")
289 289 ltctime = property(getltctime, "I'm the 'ltctime' property")
290 290 ippSeconds = property(get_ippSeconds, set_ippSeconds)
291 291 dtype = property(get_dtype, set_dtype)
292 292
293 293 class Voltage(JROData):
294 294
295 295 #data es un numpy array de 2 dmensiones (canales, alturas)
296 296 data = None
297 297
298 298 def __init__(self):
299 299 '''
300 300 Constructor
301 301 '''
302 302
303 303 self.radarControllerHeaderObj = RadarControllerHeader()
304 304
305 305 self.systemHeaderObj = SystemHeader()
306 306
307 307 self.type = "Voltage"
308 308
309 309 self.data = None
310 310
311 311 # self.dtype = None
312 312
313 313 # self.nChannels = 0
314 314
315 315 # self.nHeights = 0
316 316
317 317 self.nProfiles = None
318 318
319 319 self.heightList = None
320 320
321 321 self.channelList = None
322 322
323 323 # self.channelIndexList = None
324 324
325 325 self.flagNoData = True
326 326
327 327 self.flagTimeBlock = False
328 328
329 329 self.utctime = None
330 330
331 331 self.timeZone = None
332 332
333 333 self.dstFlag = None
334 334
335 335 self.errorCount = None
336 336
337 337 self.nCohInt = None
338 338
339 339 self.blocksize = None
340 340
341 341 self.flagDecodeData = False #asumo q la data no esta decodificada
342 342
343 343 self.flagDeflipData = False #asumo q la data no esta sin flip
344 344
345 345 self.flagShiftFFT = False
346 346
347 347
348 348 def getNoisebyHildebrand(self):
349 349 """
350 350 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
351 351
352 352 Return:
353 353 noiselevel
354 354 """
355 355
356 356 for channel in range(self.nChannels):
357 357 daux = self.data_spc[channel,:,:]
358 358 self.noise[channel] = hildebrand_sekhon(daux, self.nCohInt)
359 359
360 360 return self.noise
361 361
362 362 def getNoise(self, type = 1):
363 363
364 364 self.noise = numpy.zeros(self.nChannels)
365 365
366 366 if type == 1:
367 367 noise = self.getNoisebyHildebrand()
368 368
369 369 return 10*numpy.log10(noise)
370 370
371 371 noise = property(getNoise, "I'm the 'nHeights' property.")
372 372
373 373 class Spectra(JROData):
374 374
375 375 #data es un numpy array de 2 dmensiones (canales, perfiles, alturas)
376 376 data_spc = None
377 377
378 378 #data es un numpy array de 2 dmensiones (canales, pares, alturas)
379 379 data_cspc = None
380 380
381 381 #data es un numpy array de 2 dmensiones (canales, alturas)
382 382 data_dc = None
383 383
384 384 nFFTPoints = None
385 385
386 386 # nPairs = None
387 387
388 388 pairsList = None
389 389
390 390 nIncohInt = None
391 391
392 392 wavelength = None #Necesario para cacular el rango de velocidad desde la frecuencia
393 393
394 394 nCohInt = None #se requiere para determinar el valor de timeInterval
395 395
396 396 ippFactor = None
397 397
398 398 def __init__(self):
399 399 '''
400 400 Constructor
401 401 '''
402 402
403 403 self.radarControllerHeaderObj = RadarControllerHeader()
404 404
405 405 self.systemHeaderObj = SystemHeader()
406 406
407 407 self.type = "Spectra"
408 408
409 409 # self.data = None
410 410
411 411 # self.dtype = None
412 412
413 413 # self.nChannels = 0
414 414
415 415 # self.nHeights = 0
416 416
417 417 self.nProfiles = None
418 418
419 419 self.heightList = None
420 420
421 421 self.channelList = None
422 422
423 423 # self.channelIndexList = None
424 424
425 425 self.pairsList = None
426 426
427 427 self.flagNoData = True
428 428
429 429 self.flagTimeBlock = False
430 430
431 431 self.utctime = None
432 432
433 433 self.nCohInt = None
434 434
435 435 self.nIncohInt = None
436 436
437 437 self.blocksize = None
438 438
439 439 self.nFFTPoints = None
440 440
441 441 self.wavelength = None
442 442
443 443 self.flagDecodeData = False #asumo q la data no esta decodificada
444 444
445 445 self.flagDeflipData = False #asumo q la data no esta sin flip
446 446
447 447 self.flagShiftFFT = False
448 448
449 449 self.ippFactor = 1
450 450
451 451 #self.noise = None
452 452
453 453 self.beacon_heiIndexList = []
454 454
455 455 self.noise_estimation = None
456 456
457 457
458 458 def getNoisebyHildebrand(self):
459 459 """
460 460 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
461 461
462 462 Return:
463 463 noiselevel
464 464 """
465 465
466 466 noise = numpy.zeros(self.nChannels)
467 467 for channel in range(self.nChannels):
468 468 daux = self.data_spc[channel,:,:]
469 469 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
470 470
471 471 return noise
472 472
473 473 def getNoise(self):
474 474 if self.noise_estimation != None:
475 475 return self.noise_estimation #this was estimated by getNoise Operation defined in jroproc_spectra.py
476 476 else:
477 477 noise = self.getNoisebyHildebrand()
478 478 return noise
479 479
480 480
481 481 def getFreqRange(self, extrapoints=0):
482 482
483 483 deltafreq = self.getFmax() / (self.nFFTPoints*self.ippFactor)
484 484 freqrange = deltafreq*(numpy.arange(self.nFFTPoints+extrapoints)-self.nFFTPoints/2.) - deltafreq/2
485 485
486 486 return freqrange
487 487
488 488 def getVelRange(self, extrapoints=0):
489 489
490 490 deltav = self.getVmax() / (self.nFFTPoints*self.ippFactor)
491 491 velrange = deltav*(numpy.arange(self.nFFTPoints+extrapoints)-self.nFFTPoints/2.) - deltav/2
492 492
493 493 return velrange
494 494
495 495 def getNPairs(self):
496 496
497 497 return len(self.pairsList)
498 498
499 499 def getPairsIndexList(self):
500 500
501 501 return range(self.nPairs)
502 502
503 503 def getNormFactor(self):
504 504 pwcode = 1
505 505 if self.flagDecodeData:
506 506 pwcode = numpy.sum(self.code[0]**2)
507 507 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
508 508 normFactor = self.nProfiles*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
509 509
510 510 return normFactor
511 511
512 512 def getFlagCspc(self):
513 513
514 514 if self.data_cspc == None:
515 515 return True
516 516
517 517 return False
518 518
519 519 def getFlagDc(self):
520 520
521 521 if self.data_dc == None:
522 522 return True
523 523
524 524 return False
525 525
526 526 nPairs = property(getNPairs, "I'm the 'nPairs' property.")
527 527 pairsIndexList = property(getPairsIndexList, "I'm the 'pairsIndexList' property.")
528 528 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
529 529 flag_cspc = property(getFlagCspc)
530 530 flag_dc = property(getFlagDc)
531 531 noise = property(getNoise, "I'm the 'nHeights' property.")
532 532
533 533 class SpectraHeis(Spectra):
534 534
535 535 data_spc = None
536 536
537 537 data_cspc = None
538 538
539 539 data_dc = None
540 540
541 541 nFFTPoints = None
542 542
543 543 # nPairs = None
544 544
545 545 pairsList = None
546 546
547 547 nIncohInt = None
548 548
549 549 def __init__(self):
550 550
551 551 self.radarControllerHeaderObj = RadarControllerHeader()
552 552
553 553 self.systemHeaderObj = SystemHeader()
554 554
555 555 self.type = "SpectraHeis"
556 556
557 557 # self.dtype = None
558 558
559 559 # self.nChannels = 0
560 560
561 561 # self.nHeights = 0
562 562
563 563 self.nProfiles = None
564 564
565 565 self.heightList = None
566 566
567 567 self.channelList = None
568 568
569 569 # self.channelIndexList = None
570 570
571 571 self.flagNoData = True
572 572
573 573 self.flagTimeBlock = False
574 574
575 575 # self.nPairs = 0
576 576
577 577 self.utctime = None
578 578
579 579 self.blocksize = None
580 580
581 581 class Fits:
582 582
583 583 heightList = None
584 584
585 585 channelList = None
586 586
587 587 flagNoData = True
588 588
589 589 flagTimeBlock = False
590 590
591 591 useLocalTime = False
592 592
593 593 utctime = None
594 594
595 595 timeZone = None
596 596
597 597 # ippSeconds = None
598 598
599 599 timeInterval = None
600 600
601 601 nCohInt = None
602 602
603 603 nIncohInt = None
604 604
605 605 noise = None
606 606
607 607 windowOfFilter = 1
608 608
609 609 #Speed of ligth
610 610 C = 3e8
611 611
612 612 frequency = 49.92e6
613 613
614 614 realtime = False
615 615
616 616
617 617 def __init__(self):
618 618
619 619 self.type = "Fits"
620 620
621 621 self.nProfiles = None
622 622
623 623 self.heightList = None
624 624
625 625 self.channelList = None
626 626
627 627 # self.channelIndexList = None
628 628
629 629 self.flagNoData = True
630 630
631 631 self.utctime = None
632 632
633 633 self.nCohInt = None
634 634
635 635 self.nIncohInt = None
636 636
637 637 self.useLocalTime = True
638 638
639 639 # self.utctime = None
640 640 # self.timeZone = None
641 641 # self.ltctime = None
642 642 # self.timeInterval = None
643 643 # self.header = None
644 644 # self.data_header = None
645 645 # self.data = None
646 646 # self.datatime = None
647 647 # self.flagNoData = False
648 648 # self.expName = ''
649 649 # self.nChannels = None
650 650 # self.nSamples = None
651 651 # self.dataBlocksPerFile = None
652 652 # self.comments = ''
653 653 #
654 654
655 655
656 656 def getltctime(self):
657 657
658 658 if self.useLocalTime:
659 659 return self.utctime - self.timeZone*60
660 660
661 661 return self.utctime
662 662
663 663 def getDatatime(self):
664 664
665 665 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
666 666 return datatime
667 667
668 668 def getTimeRange(self):
669 669
670 670 datatime = []
671 671
672 672 datatime.append(self.ltctime)
673 673 datatime.append(self.ltctime + self.timeInterval)
674 674
675 675 datatime = numpy.array(datatime)
676 676
677 677 return datatime
678 678
679 679 def getHeiRange(self):
680 680
681 681 heis = self.heightList
682 682
683 683 return heis
684 684
685 685 def isEmpty(self):
686 686
687 687 return self.flagNoData
688 688
689 689 def getNHeights(self):
690 690
691 691 return len(self.heightList)
692 692
693 693 def getNChannels(self):
694 694
695 695 return len(self.channelList)
696 696
697 697 def getChannelIndexList(self):
698 698
699 699 return range(self.nChannels)
700 700
701 701 def getNoise(self, type = 1):
702 702
703 703 self.noise = numpy.zeros(self.nChannels)
704 704
705 705 if type == 1:
706 706 noise = self.getNoisebyHildebrand()
707 707
708 708 if type == 2:
709 709 noise = self.getNoisebySort()
710 710
711 711 if type == 3:
712 712 noise = self.getNoisebyWindow()
713 713
714 714 return noise
715 715
716 716 datatime = property(getDatatime, "I'm the 'datatime' property")
717 717 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
718 718 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
719 719 channelIndexList = property(getChannelIndexList, "I'm the 'channelIndexList' property.")
720 720 noise = property(getNoise, "I'm the 'nHeights' property.")
721 721 datatime = property(getDatatime, "I'm the 'datatime' property")
722 722 ltctime = property(getltctime, "I'm the 'ltctime' property")
723 723
724 724 ltctime = property(getltctime, "I'm the 'ltctime' property")
725 725
726 class AMISR:
727 def __init__(self):
728 self.flagNoData = True
729 self.data = None
730 self.utctime = None
731 self.type = "AMISR"
732
733 #propiedades para compatibilidad con Voltages
734 self.timeZone = 300#timezone like jroheader, difference in minutes between UTC and localtime
735 self.dstFlag = 0#self.dataIn.dstFlag
736 self.errorCount = 0#self.dataIn.errorCount
737 self.useLocalTime = True#self.dataIn.useLocalTime
738
739 self.radarControllerHeaderObj = None#self.dataIn.radarControllerHeaderObj.copy()
740 self.systemHeaderObj = None#self.dataIn.systemHeaderObj.copy()
741 self.channelList = [0]#self.dataIn.channelList esto solo aplica para el caso de AMISR
742 self.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
743
744 self.flagTimeBlock = None#self.dataIn.flagTimeBlock
745 #self.utctime = #self.firstdatatime
746 self.flagDecodeData = None#self.dataIn.flagDecodeData #asumo q la data esta decodificada
747 self.flagDeflipData = None#self.dataIn.flagDeflipData #asumo q la data esta sin flip
748
749 self.nCohInt = 1#self.dataIn.nCohInt
750 self.nIncohInt = 1
751 self.ippSeconds = None#self.dataIn.ippSeconds, segun el filename/Setup/Tufile
752 self.windowOfFilter = None#self.dataIn.windowOfFilter
753
754 self.timeInterval = None#self.dataIn.timeInterval*self.dataOut.nFFTPoints*self.dataOut.nIncohInt
755 self.frequency = None#self.dataIn.frequency
756 self.realtime = 0#self.dataIn.realtime
757
758 #actualizar en la lectura de datos
759 self.heightList = None#self.dataIn.heightList
760 self.nProfiles = None#Number of samples or nFFTPoints
761 self.nRecords = None
762 self.nBeams = None
763 self.nBaud = None#self.dataIn.nBaud
764 self.nCode = None#self.dataIn.nCode
765 self.code = None#self.dataIn.code
766
767 #consideracion para los Beams
768 self.beamCodeDict = None
769 self.beamRangeDict = None
770
771 def copy(self, inputObj=None):
772
773 if inputObj == None:
774 return copy.deepcopy(self)
775
776 for key in inputObj.__dict__.keys():
777 self.__dict__[key] = inputObj.__dict__[key]
778
779
780 def isEmpty(self):
781
782 return self.flagNoData No newline at end of file
@@ -1,1334 +1,1334
1 1 '''
2 2
3 3 '''
4 4 import os
5 5 import sys
6 6 import glob
7 7 import time
8 8 import numpy
9 9 import fnmatch
10 10 import time, datetime
11 import h5py
11 #import h5py
12 12 import traceback
13 13
14 14 #try:
15 15 # import pyfits
16 16 #except:
17 17 # print "pyfits module has not been imported, it should be installed to save files in fits format"
18 18
19 19 #from jrodata import *
20 20 #from jroheaderIO import *
21 21 #from jroprocessing import *
22 22
23 23 #import re
24 24 #from xml.etree.ElementTree import Element, SubElement, ElementTree
25 25
26 26
27 27 LOCALTIME = True #-18000
28 28
29 29 from model.data.jroheaderIO import PROCFLAG, BasicHeader, SystemHeader, RadarControllerHeader, ProcessingHeader
30 30
31 31 def isNumber(str):
32 32 """
33 33 Chequea si el conjunto de caracteres que componen un string puede ser convertidos a un numero.
34 34
35 35 Excepciones:
36 36 Si un determinado string no puede ser convertido a numero
37 37 Input:
38 38 str, string al cual se le analiza para determinar si convertible a un numero o no
39 39
40 40 Return:
41 41 True : si el string es uno numerico
42 42 False : no es un string numerico
43 43 """
44 44 try:
45 45 float( str )
46 46 return True
47 47 except:
48 48 return False
49 49
50 50 def isThisFileinRange(filename, startUTSeconds, endUTSeconds):
51 51 """
52 52 Esta funcion determina si un archivo de datos se encuentra o no dentro del rango de fecha especificado.
53 53
54 54 Inputs:
55 55 filename : nombre completo del archivo de datos en formato Jicamarca (.r)
56 56
57 57 startUTSeconds : fecha inicial del rango seleccionado. La fecha esta dada en
58 58 segundos contados desde 01/01/1970.
59 59 endUTSeconds : fecha final del rango seleccionado. La fecha esta dada en
60 60 segundos contados desde 01/01/1970.
61 61
62 62 Return:
63 63 Boolean : Retorna True si el archivo de datos contiene datos en el rango de
64 64 fecha especificado, de lo contrario retorna False.
65 65
66 66 Excepciones:
67 67 Si el archivo no existe o no puede ser abierto
68 68 Si la cabecera no puede ser leida.
69 69
70 70 """
71 71 basicHeaderObj = BasicHeader(LOCALTIME)
72 72
73 73 try:
74 74 fp = open(filename,'rb')
75 75 except IOError:
76 76 traceback.print_exc()
77 77 raise IOError, "The file %s can't be opened" %(filename)
78 78
79 79 sts = basicHeaderObj.read(fp)
80 80 fp.close()
81 81
82 82 if not(sts):
83 83 print "Skipping the file %s because it has not a valid header" %(filename)
84 84 return 0
85 85
86 86 if not ((startUTSeconds <= basicHeaderObj.utc) and (endUTSeconds > basicHeaderObj.utc)):
87 87 return 0
88 88
89 89 return 1
90 90
91 91 def isFileinThisTime(filename, startTime, endTime):
92 92 """
93 93 Retorna 1 si el archivo de datos se encuentra dentro del rango de horas especificado.
94 94
95 95 Inputs:
96 96 filename : nombre completo del archivo de datos en formato Jicamarca (.r)
97 97
98 98 startTime : tiempo inicial del rango seleccionado en formato datetime.time
99 99
100 100 endTime : tiempo final del rango seleccionado en formato datetime.time
101 101
102 102 Return:
103 103 Boolean : Retorna True si el archivo de datos contiene datos en el rango de
104 104 fecha especificado, de lo contrario retorna False.
105 105
106 106 Excepciones:
107 107 Si el archivo no existe o no puede ser abierto
108 108 Si la cabecera no puede ser leida.
109 109
110 110 """
111 111
112 112
113 113 try:
114 114 fp = open(filename,'rb')
115 115 except IOError:
116 116 traceback.print_exc()
117 117 raise IOError, "The file %s can't be opened" %(filename)
118 118
119 119 basicHeaderObj = BasicHeader(LOCALTIME)
120 120 sts = basicHeaderObj.read(fp)
121 121 fp.close()
122 122
123 123 thisDatetime = basicHeaderObj.datatime
124 124 thisTime = thisDatetime.time()
125 125
126 126 if not(sts):
127 127 print "Skipping the file %s because it has not a valid header" %(filename)
128 128 return None
129 129
130 130 if not ((startTime <= thisTime) and (endTime > thisTime)):
131 131 return None
132 132
133 133 return thisDatetime
134 134
135 135 def getFileFromSet(path, ext, set):
136 136 validFilelist = []
137 137 fileList = os.listdir(path)
138 138
139 139 # 0 1234 567 89A BCDE
140 140 # H YYYY DDD SSS .ext
141 141
142 142 for thisFile in fileList:
143 143 try:
144 144 year = int(thisFile[1:5])
145 145 doy = int(thisFile[5:8])
146 146 except:
147 147 continue
148 148
149 149 if (os.path.splitext(thisFile)[-1].lower() != ext.lower()):
150 150 continue
151 151
152 152 validFilelist.append(thisFile)
153 153
154 154 myfile = fnmatch.filter(validFilelist,'*%4.4d%3.3d%3.3d*'%(year,doy,set))
155 155
156 156 if len(myfile)!= 0:
157 157 return myfile[0]
158 158 else:
159 159 filename = '*%4.4d%3.3d%3.3d%s'%(year,doy,set,ext.lower())
160 160 print 'the filename %s does not exist'%filename
161 161 print '...going to the last file: '
162 162
163 163 if validFilelist:
164 164 validFilelist = sorted( validFilelist, key=str.lower )
165 165 return validFilelist[-1]
166 166
167 167 return None
168 168
169 169 def getlastFileFromPath(path, ext):
170 170 """
171 171 Depura el fileList dejando solo los que cumplan el formato de "PYYYYDDDSSS.ext"
172 172 al final de la depuracion devuelve el ultimo file de la lista que quedo.
173 173
174 174 Input:
175 175 fileList : lista conteniendo todos los files (sin path) que componen una determinada carpeta
176 176 ext : extension de los files contenidos en una carpeta
177 177
178 178 Return:
179 179 El ultimo file de una determinada carpeta, no se considera el path.
180 180 """
181 181 validFilelist = []
182 182 fileList = os.listdir(path)
183 183
184 184 # 0 1234 567 89A BCDE
185 185 # H YYYY DDD SSS .ext
186 186
187 187 for thisFile in fileList:
188 188
189 189 year = thisFile[1:5]
190 190 if not isNumber(year):
191 191 continue
192 192
193 193 doy = thisFile[5:8]
194 194 if not isNumber(doy):
195 195 continue
196 196
197 197 year = int(year)
198 198 doy = int(doy)
199 199
200 200 if (os.path.splitext(thisFile)[-1].lower() != ext.lower()):
201 201 continue
202 202
203 203 validFilelist.append(thisFile)
204 204
205 205 if validFilelist:
206 206 validFilelist = sorted( validFilelist, key=str.lower )
207 207 return validFilelist[-1]
208 208
209 209 return None
210 210
211 211 def checkForRealPath(path, foldercounter, year, doy, set, ext):
212 212 """
213 213 Por ser Linux Case Sensitive entonces checkForRealPath encuentra el nombre correcto de un path,
214 214 Prueba por varias combinaciones de nombres entre mayusculas y minusculas para determinar
215 215 el path exacto de un determinado file.
216 216
217 217 Example :
218 218 nombre correcto del file es .../.../D2009307/P2009307367.ext
219 219
220 220 Entonces la funcion prueba con las siguientes combinaciones
221 221 .../.../y2009307367.ext
222 222 .../.../Y2009307367.ext
223 223 .../.../x2009307/y2009307367.ext
224 224 .../.../x2009307/Y2009307367.ext
225 225 .../.../X2009307/y2009307367.ext
226 226 .../.../X2009307/Y2009307367.ext
227 227 siendo para este caso, la ultima combinacion de letras, identica al file buscado
228 228
229 229 Return:
230 230 Si encuentra la cobinacion adecuada devuelve el path completo y el nombre del file
231 231 caso contrario devuelve None como path y el la ultima combinacion de nombre en mayusculas
232 232 para el filename
233 233 """
234 234 fullfilename = None
235 235 find_flag = False
236 236 filename = None
237 237
238 238 prefixDirList = [None,'d','D']
239 239 if ext.lower() == ".r": #voltage
240 240 prefixFileList = ['d','D']
241 241 elif ext.lower() == ".pdata": #spectra
242 242 prefixFileList = ['p','P']
243 243 else:
244 244 return None, filename
245 245
246 246 #barrido por las combinaciones posibles
247 247 for prefixDir in prefixDirList:
248 248 thispath = path
249 249 if prefixDir != None:
250 250 #formo el nombre del directorio xYYYYDDD (x=d o x=D)
251 251 if foldercounter == 0:
252 252 thispath = os.path.join(path, "%s%04d%03d" % ( prefixDir, year, doy ))
253 253 else:
254 254 thispath = os.path.join(path, "%s%04d%03d_%02d" % ( prefixDir, year, doy , foldercounter))
255 255 for prefixFile in prefixFileList: #barrido por las dos combinaciones posibles de "D"
256 256 filename = "%s%04d%03d%03d%s" % ( prefixFile, year, doy, set, ext ) #formo el nombre del file xYYYYDDDSSS.ext
257 257 fullfilename = os.path.join( thispath, filename ) #formo el path completo
258 258
259 259 if os.path.exists( fullfilename ): #verifico que exista
260 260 find_flag = True
261 261 break
262 262 if find_flag:
263 263 break
264 264
265 265 if not(find_flag):
266 266 return None, filename
267 267
268 268 return fullfilename, filename
269 269
270 270 def isDoyFolder(folder):
271 271 try:
272 272 year = int(folder[1:5])
273 273 except:
274 274 return 0
275 275
276 276 try:
277 277 doy = int(folder[5:8])
278 278 except:
279 279 return 0
280 280
281 281 return 1
282 282
283 283 class JRODataIO:
284 284
285 285 c = 3E8
286 286
287 287 isConfig = False
288 288
289 289 basicHeaderObj = None
290 290
291 291 systemHeaderObj = None
292 292
293 293 radarControllerHeaderObj = None
294 294
295 295 processingHeaderObj = None
296 296
297 297 online = 0
298 298
299 299 dtype = None
300 300
301 301 pathList = []
302 302
303 303 filenameList = []
304 304
305 305 filename = None
306 306
307 307 ext = None
308 308
309 309 flagIsNewFile = 1
310 310
311 311 flagTimeBlock = 0
312 312
313 313 flagIsNewBlock = 0
314 314
315 315 fp = None
316 316
317 317 firstHeaderSize = 0
318 318
319 319 basicHeaderSize = 24
320 320
321 321 versionFile = 1103
322 322
323 323 fileSize = None
324 324
325 325 # ippSeconds = None
326 326
327 327 fileSizeByHeader = None
328 328
329 329 fileIndex = None
330 330
331 331 profileIndex = None
332 332
333 333 blockIndex = None
334 334
335 335 nTotalBlocks = None
336 336
337 337 maxTimeStep = 30
338 338
339 339 lastUTTime = None
340 340
341 341 datablock = None
342 342
343 343 dataOut = None
344 344
345 345 blocksize = None
346 346
347 347 def __init__(self):
348 348
349 349 raise ValueError, "Not implemented"
350 350
351 351 def run(self):
352 352
353 353 raise ValueError, "Not implemented"
354 354
355 355 class JRODataReader(JRODataIO):
356 356
357 357 nReadBlocks = 0
358 358
359 359 delay = 10 #number of seconds waiting a new file
360 360
361 361 nTries = 3 #quantity tries
362 362
363 363 nFiles = 3 #number of files for searching
364 364
365 365 path = None
366 366
367 367 foldercounter = 0
368 368
369 369 flagNoMoreFiles = 0
370 370
371 371 datetimeList = []
372 372
373 373 __isFirstTimeOnline = 1
374 374
375 375 __printInfo = True
376 376
377 377 profileIndex = None
378 378
379 379 def __init__(self):
380 380
381 381 """
382 382
383 383 """
384 384
385 385 raise ValueError, "This method has not been implemented"
386 386
387 387
388 388 def createObjByDefault(self):
389 389 """
390 390
391 391 """
392 392 raise ValueError, "This method has not been implemented"
393 393
394 394 def getBlockDimension(self):
395 395
396 396 raise ValueError, "No implemented"
397 397
398 398 def __searchFilesOffLine(self,
399 399 path,
400 400 startDate,
401 401 endDate,
402 402 startTime=datetime.time(0,0,0),
403 403 endTime=datetime.time(23,59,59),
404 404 set=None,
405 405 expLabel='',
406 406 ext='.r',
407 407 walk=True):
408 408
409 409 pathList = []
410 410
411 411 if not walk:
412 412 #pathList.append(path)
413 413 multi_path = path.split(',')
414 414 for single_path in multi_path:
415 415 pathList.append(single_path)
416 416
417 417 else:
418 418 #dirList = []
419 419 multi_path = path.split(',')
420 420 for single_path in multi_path:
421 421 dirList = []
422 422 for thisPath in os.listdir(single_path):
423 423 if not os.path.isdir(os.path.join(single_path,thisPath)):
424 424 continue
425 425 if not isDoyFolder(thisPath):
426 426 continue
427 427
428 428 dirList.append(thisPath)
429 429
430 430 if not(dirList):
431 431 return None, None
432 432
433 433 thisDate = startDate
434 434
435 435 while(thisDate <= endDate):
436 436 year = thisDate.timetuple().tm_year
437 437 doy = thisDate.timetuple().tm_yday
438 438
439 439 matchlist = fnmatch.filter(dirList, '?' + '%4.4d%3.3d' % (year,doy) + '*')
440 440 if len(matchlist) == 0:
441 441 thisDate += datetime.timedelta(1)
442 442 continue
443 443 for match in matchlist:
444 444 pathList.append(os.path.join(single_path,match,expLabel))
445 445
446 446 thisDate += datetime.timedelta(1)
447 447
448 448 if pathList == []:
449 449 print "Any folder was found for the date range: %s-%s" %(startDate, endDate)
450 450 return None, None
451 451
452 452 print "%d folder(s) was(were) found for the date range: %s - %s" %(len(pathList), startDate, endDate)
453 453
454 454 filenameList = []
455 455 datetimeList = []
456 456 pathDict = {}
457 457 filenameList_to_sort = []
458 458
459 459 for i in range(len(pathList)):
460 460
461 461 thisPath = pathList[i]
462 462
463 463 fileList = glob.glob1(thisPath, "*%s" %ext)
464 464 fileList.sort()
465 465 pathDict.setdefault(fileList[0])
466 466 pathDict[fileList[0]] = i
467 467 filenameList_to_sort.append(fileList[0])
468 468
469 469 filenameList_to_sort.sort()
470 470
471 471 for file in filenameList_to_sort:
472 472 thisPath = pathList[pathDict[file]]
473 473
474 474 fileList = glob.glob1(thisPath, "*%s" %ext)
475 475 fileList.sort()
476 476
477 477 for file in fileList:
478 478
479 479 filename = os.path.join(thisPath,file)
480 480 thisDatetime = isFileinThisTime(filename, startTime, endTime)
481 481
482 482 if not(thisDatetime):
483 483 continue
484 484
485 485 filenameList.append(filename)
486 486 datetimeList.append(thisDatetime)
487 487
488 488 if not(filenameList):
489 489 print "Any file was found for the time range %s - %s" %(startTime, endTime)
490 490 return None, None
491 491
492 492 print "%d file(s) was(were) found for the time range: %s - %s" %(len(filenameList), startTime, endTime)
493 493 print
494 494
495 495 for i in range(len(filenameList)):
496 496 print "%s -> [%s]" %(filenameList[i], datetimeList[i].ctime())
497 497
498 498 self.filenameList = filenameList
499 499 self.datetimeList = datetimeList
500 500
501 501 return pathList, filenameList
502 502
503 503 def __searchFilesOnLine(self, path, expLabel = "", ext = None, walk=True, set=None):
504 504
505 505 """
506 506 Busca el ultimo archivo de la ultima carpeta (determinada o no por startDateTime) y
507 507 devuelve el archivo encontrado ademas de otros datos.
508 508
509 509 Input:
510 510 path : carpeta donde estan contenidos los files que contiene data
511 511
512 512 expLabel : Nombre del subexperimento (subfolder)
513 513
514 514 ext : extension de los files
515 515
516 516 walk : Si es habilitado no realiza busquedas dentro de los ubdirectorios (doypath)
517 517
518 518 Return:
519 519 directory : eL directorio donde esta el file encontrado
520 520 filename : el ultimo file de una determinada carpeta
521 521 year : el anho
522 522 doy : el numero de dia del anho
523 523 set : el set del archivo
524 524
525 525
526 526 """
527 527 dirList = []
528 528
529 529 if not walk:
530 530 fullpath = path
531 531 foldercounter = 0
532 532 else:
533 533 #Filtra solo los directorios
534 534 for thisPath in os.listdir(path):
535 535 if not os.path.isdir(os.path.join(path,thisPath)):
536 536 continue
537 537 if not isDoyFolder(thisPath):
538 538 continue
539 539
540 540 dirList.append(thisPath)
541 541
542 542 if not(dirList):
543 543 return None, None, None, None, None, None
544 544
545 545 dirList = sorted( dirList, key=str.lower )
546 546
547 547 doypath = dirList[-1]
548 548 foldercounter = int(doypath.split('_')[1]) if len(doypath.split('_'))>1 else 0
549 549 fullpath = os.path.join(path, doypath, expLabel)
550 550
551 551
552 552 print "%s folder was found: " %(fullpath )
553 553
554 554 if set == None:
555 555 filename = getlastFileFromPath(fullpath, ext)
556 556 else:
557 557 filename = getFileFromSet(fullpath, ext, set)
558 558
559 559 if not(filename):
560 560 return None, None, None, None, None, None
561 561
562 562 print "%s file was found" %(filename)
563 563
564 564 if not(self.__verifyFile(os.path.join(fullpath, filename))):
565 565 return None, None, None, None, None, None
566 566
567 567 year = int( filename[1:5] )
568 568 doy = int( filename[5:8] )
569 569 set = int( filename[8:11] )
570 570
571 571 return fullpath, foldercounter, filename, year, doy, set
572 572
573 573 def __setNextFileOffline(self):
574 574
575 575 idFile = self.fileIndex
576 576
577 577 while (True):
578 578 idFile += 1
579 579 if not(idFile < len(self.filenameList)):
580 580 self.flagNoMoreFiles = 1
581 581 print "No more Files"
582 582 return 0
583 583
584 584 filename = self.filenameList[idFile]
585 585
586 586 if not(self.__verifyFile(filename)):
587 587 continue
588 588
589 589 fileSize = os.path.getsize(filename)
590 590 fp = open(filename,'rb')
591 591 break
592 592
593 593 self.flagIsNewFile = 1
594 594 self.fileIndex = idFile
595 595 self.filename = filename
596 596 self.fileSize = fileSize
597 597 self.fp = fp
598 598
599 599 print "Setting the file: %s"%self.filename
600 600
601 601 return 1
602 602
603 603 def __setNextFileOnline(self):
604 604 """
605 605 Busca el siguiente file que tenga suficiente data para ser leida, dentro de un folder especifico, si
606 606 no encuentra un file valido espera un tiempo determinado y luego busca en los posibles n files
607 607 siguientes.
608 608
609 609 Affected:
610 610 self.flagIsNewFile
611 611 self.filename
612 612 self.fileSize
613 613 self.fp
614 614 self.set
615 615 self.flagNoMoreFiles
616 616
617 617 Return:
618 618 0 : si luego de una busqueda del siguiente file valido este no pudo ser encontrado
619 619 1 : si el file fue abierto con exito y esta listo a ser leido
620 620
621 621 Excepciones:
622 622 Si un determinado file no puede ser abierto
623 623 """
624 624 nFiles = 0
625 625 fileOk_flag = False
626 626 firstTime_flag = True
627 627
628 628 self.set += 1
629 629
630 630 if self.set > 999:
631 631 self.set = 0
632 632 self.foldercounter += 1
633 633
634 634 #busca el 1er file disponible
635 635 fullfilename, filename = checkForRealPath( self.path, self.foldercounter, self.year, self.doy, self.set, self.ext )
636 636 if fullfilename:
637 637 if self.__verifyFile(fullfilename, False):
638 638 fileOk_flag = True
639 639
640 640 #si no encuentra un file entonces espera y vuelve a buscar
641 641 if not(fileOk_flag):
642 642 for nFiles in range(self.nFiles+1): #busco en los siguientes self.nFiles+1 files posibles
643 643
644 644 if firstTime_flag: #si es la 1era vez entonces hace el for self.nTries veces
645 645 tries = self.nTries
646 646 else:
647 647 tries = 1 #si no es la 1era vez entonces solo lo hace una vez
648 648
649 649 for nTries in range( tries ):
650 650 if firstTime_flag:
651 651 print "\tWaiting %0.2f sec for the file \"%s\" , try %03d ..." % ( self.delay, filename, nTries+1 )
652 652 time.sleep( self.delay )
653 653 else:
654 654 print "\tSearching next \"%s%04d%03d%03d%s\" file ..." % (self.optchar, self.year, self.doy, self.set, self.ext)
655 655
656 656 fullfilename, filename = checkForRealPath( self.path, self.foldercounter, self.year, self.doy, self.set, self.ext )
657 657 if fullfilename:
658 658 if self.__verifyFile(fullfilename):
659 659 fileOk_flag = True
660 660 break
661 661
662 662 if fileOk_flag:
663 663 break
664 664
665 665 firstTime_flag = False
666 666
667 667 print "\tSkipping the file \"%s\" due to this file doesn't exist" % filename
668 668 self.set += 1
669 669
670 670 if nFiles == (self.nFiles-1): #si no encuentro el file buscado cambio de carpeta y busco en la siguiente carpeta
671 671 self.set = 0
672 672 self.doy += 1
673 673 self.foldercounter = 0
674 674
675 675 if fileOk_flag:
676 676 self.fileSize = os.path.getsize( fullfilename )
677 677 self.filename = fullfilename
678 678 self.flagIsNewFile = 1
679 679 if self.fp != None: self.fp.close()
680 680 self.fp = open(fullfilename, 'rb')
681 681 self.flagNoMoreFiles = 0
682 682 print 'Setting the file: %s' % fullfilename
683 683 else:
684 684 self.fileSize = 0
685 685 self.filename = None
686 686 self.flagIsNewFile = 0
687 687 self.fp = None
688 688 self.flagNoMoreFiles = 1
689 689 print 'No more Files'
690 690
691 691 return fileOk_flag
692 692
693 693 def setNextFile(self):
694 694 if self.fp != None:
695 695 self.fp.close()
696 696
697 697 if self.online:
698 698 newFile = self.__setNextFileOnline()
699 699 else:
700 700 newFile = self.__setNextFileOffline()
701 701
702 702 if not(newFile):
703 703 return 0
704 704
705 705 self.__readFirstHeader()
706 706 self.nReadBlocks = 0
707 707 return 1
708 708
709 709 def __waitNewBlock(self):
710 710 """
711 711 Return 1 si se encontro un nuevo bloque de datos, 0 de otra forma.
712 712
713 713 Si el modo de lectura es OffLine siempre retorn 0
714 714 """
715 715 if not self.online:
716 716 return 0
717 717
718 718 if (self.nReadBlocks >= self.processingHeaderObj.dataBlocksPerFile):
719 719 return 0
720 720
721 721 currentPointer = self.fp.tell()
722 722
723 723 neededSize = self.processingHeaderObj.blockSize + self.basicHeaderSize
724 724
725 725 for nTries in range( self.nTries ):
726 726
727 727 self.fp.close()
728 728 self.fp = open( self.filename, 'rb' )
729 729 self.fp.seek( currentPointer )
730 730
731 731 self.fileSize = os.path.getsize( self.filename )
732 732 currentSize = self.fileSize - currentPointer
733 733
734 734 if ( currentSize >= neededSize ):
735 735 self.basicHeaderObj.read(self.fp)
736 736 return 1
737 737
738 738 if self.fileSize == self.fileSizeByHeader:
739 739 # self.flagEoF = True
740 740 return 0
741 741
742 742 print "\tWaiting %0.2f seconds for the next block, try %03d ..." % (self.delay, nTries+1)
743 743 time.sleep( self.delay )
744 744
745 745
746 746 return 0
747 747
748 748 def waitDataBlock(self,pointer_location):
749 749
750 750 currentPointer = pointer_location
751 751
752 752 neededSize = self.processingHeaderObj.blockSize #+ self.basicHeaderSize
753 753
754 754 for nTries in range( self.nTries ):
755 755 self.fp.close()
756 756 self.fp = open( self.filename, 'rb' )
757 757 self.fp.seek( currentPointer )
758 758
759 759 self.fileSize = os.path.getsize( self.filename )
760 760 currentSize = self.fileSize - currentPointer
761 761
762 762 if ( currentSize >= neededSize ):
763 763 return 1
764 764
765 765 print "\tWaiting %0.2f seconds for the next block, try %03d ..." % (self.delay, nTries+1)
766 766 time.sleep( self.delay )
767 767
768 768 return 0
769 769
770 770 def __jumpToLastBlock(self):
771 771
772 772 if not(self.__isFirstTimeOnline):
773 773 return
774 774
775 775 csize = self.fileSize - self.fp.tell()
776 776 blocksize = self.processingHeaderObj.blockSize
777 777
778 778 #salta el primer bloque de datos
779 779 if csize > self.processingHeaderObj.blockSize:
780 780 self.fp.seek(self.fp.tell() + blocksize)
781 781 else:
782 782 return
783 783
784 784 csize = self.fileSize - self.fp.tell()
785 785 neededsize = self.processingHeaderObj.blockSize + self.basicHeaderSize
786 786 while True:
787 787
788 788 if self.fp.tell()<self.fileSize:
789 789 self.fp.seek(self.fp.tell() + neededsize)
790 790 else:
791 791 self.fp.seek(self.fp.tell() - neededsize)
792 792 break
793 793
794 794 # csize = self.fileSize - self.fp.tell()
795 795 # neededsize = self.processingHeaderObj.blockSize + self.basicHeaderSize
796 796 # factor = int(csize/neededsize)
797 797 # if factor > 0:
798 798 # self.fp.seek(self.fp.tell() + factor*neededsize)
799 799
800 800 self.flagIsNewFile = 0
801 801 self.__isFirstTimeOnline = 0
802 802
803 803 def __setNewBlock(self):
804 804
805 805 if self.fp == None:
806 806 return 0
807 807
808 808 if self.online:
809 809 self.__jumpToLastBlock()
810 810
811 811 if self.flagIsNewFile:
812 812 return 1
813 813
814 814 self.lastUTTime = self.basicHeaderObj.utc
815 815 currentSize = self.fileSize - self.fp.tell()
816 816 neededSize = self.processingHeaderObj.blockSize + self.basicHeaderSize
817 817
818 818 if (currentSize >= neededSize):
819 819 self.basicHeaderObj.read(self.fp)
820 820 return 1
821 821
822 822 if self.__waitNewBlock():
823 823 return 1
824 824
825 825 if not(self.setNextFile()):
826 826 return 0
827 827
828 828 deltaTime = self.basicHeaderObj.utc - self.lastUTTime #
829 829
830 830 self.flagTimeBlock = 0
831 831
832 832 if deltaTime > self.maxTimeStep:
833 833 self.flagTimeBlock = 1
834 834
835 835 return 1
836 836
837 837 def readNextBlock(self):
838 838 if not(self.__setNewBlock()):
839 839 return 0
840 840
841 841 if not(self.readBlock()):
842 842 return 0
843 843
844 844 return 1
845 845
846 846 def __readFirstHeader(self):
847 847
848 848 self.basicHeaderObj.read(self.fp)
849 849 self.systemHeaderObj.read(self.fp)
850 850 self.radarControllerHeaderObj.read(self.fp)
851 851 self.processingHeaderObj.read(self.fp)
852 852
853 853 self.firstHeaderSize = self.basicHeaderObj.size
854 854
855 855 datatype = int(numpy.log2((self.processingHeaderObj.processFlags & PROCFLAG.DATATYPE_MASK))-numpy.log2(PROCFLAG.DATATYPE_CHAR))
856 856 if datatype == 0:
857 857 datatype_str = numpy.dtype([('real','<i1'),('imag','<i1')])
858 858 elif datatype == 1:
859 859 datatype_str = numpy.dtype([('real','<i2'),('imag','<i2')])
860 860 elif datatype == 2:
861 861 datatype_str = numpy.dtype([('real','<i4'),('imag','<i4')])
862 862 elif datatype == 3:
863 863 datatype_str = numpy.dtype([('real','<i8'),('imag','<i8')])
864 864 elif datatype == 4:
865 865 datatype_str = numpy.dtype([('real','<f4'),('imag','<f4')])
866 866 elif datatype == 5:
867 867 datatype_str = numpy.dtype([('real','<f8'),('imag','<f8')])
868 868 else:
869 869 raise ValueError, 'Data type was not defined'
870 870
871 871 self.dtype = datatype_str
872 872 #self.ippSeconds = 2 * 1000 * self.radarControllerHeaderObj.ipp / self.c
873 873 self.fileSizeByHeader = self.processingHeaderObj.dataBlocksPerFile * self.processingHeaderObj.blockSize + self.firstHeaderSize + self.basicHeaderSize*(self.processingHeaderObj.dataBlocksPerFile - 1)
874 874 # self.dataOut.channelList = numpy.arange(self.systemHeaderObj.numChannels)
875 875 # self.dataOut.channelIndexList = numpy.arange(self.systemHeaderObj.numChannels)
876 876 self.getBlockDimension()
877 877
878 878 def __verifyFile(self, filename, msgFlag=True):
879 879 msg = None
880 880 try:
881 881 fp = open(filename, 'rb')
882 882 currentPosition = fp.tell()
883 883 except IOError:
884 884 traceback.print_exc()
885 885 if msgFlag:
886 886 print "The file %s can't be opened" % (filename)
887 887 return False
888 888
889 889 neededSize = self.processingHeaderObj.blockSize + self.firstHeaderSize
890 890
891 891 if neededSize == 0:
892 892 basicHeaderObj = BasicHeader(LOCALTIME)
893 893 systemHeaderObj = SystemHeader()
894 894 radarControllerHeaderObj = RadarControllerHeader()
895 895 processingHeaderObj = ProcessingHeader()
896 896
897 897 try:
898 898 if not( basicHeaderObj.read(fp) ): raise IOError
899 899 if not( systemHeaderObj.read(fp) ): raise IOError
900 900 if not( radarControllerHeaderObj.read(fp) ): raise IOError
901 901 if not( processingHeaderObj.read(fp) ): raise IOError
902 902 # data_type = int(numpy.log2((processingHeaderObj.processFlags & PROCFLAG.DATATYPE_MASK))-numpy.log2(PROCFLAG.DATATYPE_CHAR))
903 903
904 904 neededSize = processingHeaderObj.blockSize + basicHeaderObj.size
905 905
906 906 except IOError:
907 907 traceback.print_exc()
908 908 if msgFlag:
909 909 print "\tThe file %s is empty or it hasn't enough data" % filename
910 910
911 911 fp.close()
912 912 return False
913 913 else:
914 914 msg = "\tSkipping the file %s due to it hasn't enough data" %filename
915 915
916 916 fp.close()
917 917 fileSize = os.path.getsize(filename)
918 918 currentSize = fileSize - currentPosition
919 919 if currentSize < neededSize:
920 920 if msgFlag and (msg != None):
921 921 print msg #print"\tSkipping the file %s due to it hasn't enough data" %filename
922 922 return False
923 923
924 924 return True
925 925
926 926 def setup(self,
927 927 path=None,
928 928 startDate=None,
929 929 endDate=None,
930 930 startTime=datetime.time(0,0,0),
931 931 endTime=datetime.time(23,59,59),
932 932 set=None,
933 933 expLabel = "",
934 934 ext = None,
935 935 online = False,
936 936 delay = 60,
937 937 walk = True):
938 938
939 939 if path == None:
940 940 raise ValueError, "The path is not valid"
941 941
942 942 if ext == None:
943 943 ext = self.ext
944 944
945 945 if online:
946 946 print "Searching files in online mode..."
947 947
948 948 for nTries in range( self.nTries ):
949 949 fullpath, foldercounter, file, year, doy, set = self.__searchFilesOnLine(path=path, expLabel=expLabel, ext=ext, walk=walk, set=set)
950 950
951 951 if fullpath:
952 952 break
953 953
954 954 print '\tWaiting %0.2f sec for an valid file in %s: try %02d ...' % (self.delay, path, nTries+1)
955 955 time.sleep( self.delay )
956 956
957 957 if not(fullpath):
958 958 print "There 'isn't valied files in %s" % path
959 959 return None
960 960
961 961 self.year = year
962 962 self.doy = doy
963 963 self.set = set - 1
964 964 self.path = path
965 965 self.foldercounter = foldercounter
966 966 last_set = None
967 967
968 968 else:
969 969 print "Searching files in offline mode ..."
970 970 pathList, filenameList = self.__searchFilesOffLine(path, startDate=startDate, endDate=endDate,
971 971 startTime=startTime, endTime=endTime,
972 972 set=set, expLabel=expLabel, ext=ext,
973 973 walk=walk)
974 974
975 975 if not(pathList):
976 976 print "No *%s files into the folder %s \nfor the range: %s - %s"%(ext, path,
977 977 datetime.datetime.combine(startDate,startTime).ctime(),
978 978 datetime.datetime.combine(endDate,endTime).ctime())
979 979
980 980 sys.exit(-1)
981 981
982 982
983 983 self.fileIndex = -1
984 984 self.pathList = pathList
985 985 self.filenameList = filenameList
986 986 file_name = os.path.basename(filenameList[-1])
987 987 basename, ext = os.path.splitext(file_name)
988 988 last_set = int(basename[-3:])
989 989
990 990 self.online = online
991 991 self.delay = delay
992 992 ext = ext.lower()
993 993 self.ext = ext
994 994
995 995 if not(self.setNextFile()):
996 996 if (startDate!=None) and (endDate!=None):
997 997 print "No files in range: %s - %s" %(datetime.datetime.combine(startDate,startTime).ctime(), datetime.datetime.combine(endDate,endTime).ctime())
998 998 elif startDate != None:
999 999 print "No files in range: %s" %(datetime.datetime.combine(startDate,startTime).ctime())
1000 1000 else:
1001 1001 print "No files"
1002 1002
1003 1003 sys.exit(-1)
1004 1004
1005 1005 # self.updateDataHeader()
1006 1006 if last_set != None:
1007 1007 self.dataOut.last_block = last_set * self.processingHeaderObj.dataBlocksPerFile + self.basicHeaderObj.dataBlock
1008 1008 return
1009 1009
1010 1010 def getBasicHeader(self):
1011 1011
1012 1012 self.dataOut.utctime = self.basicHeaderObj.utc + self.basicHeaderObj.miliSecond/1000. + self.profileIndex * self.radarControllerHeaderObj.ippSeconds
1013 1013
1014 1014 self.dataOut.flagTimeBlock = self.flagTimeBlock
1015 1015
1016 1016 self.dataOut.timeZone = self.basicHeaderObj.timeZone
1017 1017
1018 1018 self.dataOut.dstFlag = self.basicHeaderObj.dstFlag
1019 1019
1020 1020 self.dataOut.errorCount = self.basicHeaderObj.errorCount
1021 1021
1022 1022 self.dataOut.useLocalTime = self.basicHeaderObj.useLocalTime
1023 1023
1024 1024 def getFirstHeader(self):
1025 1025
1026 1026 raise ValueError, "This method has not been implemented"
1027 1027
1028 1028 def getData(self):
1029 1029
1030 1030 raise ValueError, "This method has not been implemented"
1031 1031
1032 1032 def hasNotDataInBuffer(self):
1033 1033
1034 1034 raise ValueError, "This method has not been implemented"
1035 1035
1036 1036 def readBlock(self):
1037 1037
1038 1038 raise ValueError, "This method has not been implemented"
1039 1039
1040 1040 def isEndProcess(self):
1041 1041
1042 1042 return self.flagNoMoreFiles
1043 1043
1044 1044 def printReadBlocks(self):
1045 1045
1046 1046 print "Number of read blocks per file %04d" %self.nReadBlocks
1047 1047
1048 1048 def printTotalBlocks(self):
1049 1049
1050 1050 print "Number of read blocks %04d" %self.nTotalBlocks
1051 1051
1052 1052 def printNumberOfBlock(self):
1053 1053
1054 1054 if self.flagIsNewBlock:
1055 1055 print "Block No. %04d, Total blocks %04d -> %s" %(self.basicHeaderObj.dataBlock, self.nTotalBlocks, self.dataOut.datatime.ctime())
1056 1056 self.dataOut.blocknow = self.basicHeaderObj.dataBlock
1057 1057
1058 1058 def printInfo(self):
1059 1059
1060 1060 if self.__printInfo == False:
1061 1061 return
1062 1062
1063 1063 self.basicHeaderObj.printInfo()
1064 1064 self.systemHeaderObj.printInfo()
1065 1065 self.radarControllerHeaderObj.printInfo()
1066 1066 self.processingHeaderObj.printInfo()
1067 1067
1068 1068 self.__printInfo = False
1069 1069
1070 1070
1071 1071 def run(self, **kwargs):
1072 1072
1073 1073 if not(self.isConfig):
1074 1074
1075 1075 # self.dataOut = dataOut
1076 1076 self.setup(**kwargs)
1077 1077 self.isConfig = True
1078 1078
1079 1079 self.getData()
1080 1080
1081 1081 class JRODataWriter(JRODataIO):
1082 1082
1083 1083 """
1084 1084 Esta clase permite escribir datos a archivos procesados (.r o ,pdata). La escritura
1085 1085 de los datos siempre se realiza por bloques.
1086 1086 """
1087 1087
1088 1088 blockIndex = 0
1089 1089
1090 1090 path = None
1091 1091
1092 1092 setFile = None
1093 1093
1094 1094 profilesPerBlock = None
1095 1095
1096 1096 blocksPerFile = None
1097 1097
1098 1098 nWriteBlocks = 0
1099 1099
1100 1100 def __init__(self, dataOut=None):
1101 1101 raise ValueError, "Not implemented"
1102 1102
1103 1103
1104 1104 def hasAllDataInBuffer(self):
1105 1105 raise ValueError, "Not implemented"
1106 1106
1107 1107
1108 1108 def setBlockDimension(self):
1109 1109 raise ValueError, "Not implemented"
1110 1110
1111 1111
1112 1112 def writeBlock(self):
1113 1113 raise ValueError, "No implemented"
1114 1114
1115 1115
1116 1116 def putData(self):
1117 1117 raise ValueError, "No implemented"
1118 1118
1119 1119
1120 1120 def setBasicHeader(self):
1121 1121
1122 1122 self.basicHeaderObj.size = self.basicHeaderSize #bytes
1123 1123 self.basicHeaderObj.version = self.versionFile
1124 1124 self.basicHeaderObj.dataBlock = self.nTotalBlocks
1125 1125
1126 1126 utc = numpy.floor(self.dataOut.utctime)
1127 1127 milisecond = (self.dataOut.utctime - utc)* 1000.0
1128 1128
1129 1129 self.basicHeaderObj.utc = utc
1130 1130 self.basicHeaderObj.miliSecond = milisecond
1131 1131 self.basicHeaderObj.timeZone = self.dataOut.timeZone
1132 1132 self.basicHeaderObj.dstFlag = self.dataOut.dstFlag
1133 1133 self.basicHeaderObj.errorCount = self.dataOut.errorCount
1134 1134
1135 1135 def setFirstHeader(self):
1136 1136 """
1137 1137 Obtiene una copia del First Header
1138 1138
1139 1139 Affected:
1140 1140
1141 1141 self.basicHeaderObj
1142 1142 self.systemHeaderObj
1143 1143 self.radarControllerHeaderObj
1144 1144 self.processingHeaderObj self.
1145 1145
1146 1146 Return:
1147 1147 None
1148 1148 """
1149 1149
1150 1150 raise ValueError, "No implemented"
1151 1151
1152 1152 def __writeFirstHeader(self):
1153 1153 """
1154 1154 Escribe el primer header del file es decir el Basic header y el Long header (SystemHeader, RadarControllerHeader, ProcessingHeader)
1155 1155
1156 1156 Affected:
1157 1157 __dataType
1158 1158
1159 1159 Return:
1160 1160 None
1161 1161 """
1162 1162
1163 1163 # CALCULAR PARAMETROS
1164 1164
1165 1165 sizeLongHeader = self.systemHeaderObj.size + self.radarControllerHeaderObj.size + self.processingHeaderObj.size
1166 1166 self.basicHeaderObj.size = self.basicHeaderSize + sizeLongHeader
1167 1167
1168 1168 self.basicHeaderObj.write(self.fp)
1169 1169 self.systemHeaderObj.write(self.fp)
1170 1170 self.radarControllerHeaderObj.write(self.fp)
1171 1171 self.processingHeaderObj.write(self.fp)
1172 1172
1173 1173 self.dtype = self.dataOut.dtype
1174 1174
1175 1175 def __setNewBlock(self):
1176 1176 """
1177 1177 Si es un nuevo file escribe el First Header caso contrario escribe solo el Basic Header
1178 1178
1179 1179 Return:
1180 1180 0 : si no pudo escribir nada
1181 1181 1 : Si escribio el Basic el First Header
1182 1182 """
1183 1183 if self.fp == None:
1184 1184 self.setNextFile()
1185 1185
1186 1186 if self.flagIsNewFile:
1187 1187 return 1
1188 1188
1189 1189 if self.blockIndex < self.processingHeaderObj.dataBlocksPerFile:
1190 1190 self.basicHeaderObj.write(self.fp)
1191 1191 return 1
1192 1192
1193 1193 if not( self.setNextFile() ):
1194 1194 return 0
1195 1195
1196 1196 return 1
1197 1197
1198 1198
1199 1199 def writeNextBlock(self):
1200 1200 """
1201 1201 Selecciona el bloque siguiente de datos y los escribe en un file
1202 1202
1203 1203 Return:
1204 1204 0 : Si no hizo pudo escribir el bloque de datos
1205 1205 1 : Si no pudo escribir el bloque de datos
1206 1206 """
1207 1207 if not( self.__setNewBlock() ):
1208 1208 return 0
1209 1209
1210 1210 self.writeBlock()
1211 1211
1212 1212 return 1
1213 1213
1214 1214 def setNextFile(self):
1215 1215 """
1216 1216 Determina el siguiente file que sera escrito
1217 1217
1218 1218 Affected:
1219 1219 self.filename
1220 1220 self.subfolder
1221 1221 self.fp
1222 1222 self.setFile
1223 1223 self.flagIsNewFile
1224 1224
1225 1225 Return:
1226 1226 0 : Si el archivo no puede ser escrito
1227 1227 1 : Si el archivo esta listo para ser escrito
1228 1228 """
1229 1229 ext = self.ext
1230 1230 path = self.path
1231 1231
1232 1232 if self.fp != None:
1233 1233 self.fp.close()
1234 1234
1235 1235 timeTuple = time.localtime( self.dataOut.utctime)
1236 1236 subfolder = 'd%4.4d%3.3d' % (timeTuple.tm_year,timeTuple.tm_yday)
1237 1237
1238 1238 fullpath = os.path.join( path, subfolder )
1239 1239 if not( os.path.exists(fullpath) ):
1240 1240 os.mkdir(fullpath)
1241 1241 self.setFile = -1 #inicializo mi contador de seteo
1242 1242 else:
1243 1243 filesList = os.listdir( fullpath )
1244 1244 if len( filesList ) > 0:
1245 1245 filesList = sorted( filesList, key=str.lower )
1246 1246 filen = filesList[-1]
1247 1247 # el filename debera tener el siguiente formato
1248 1248 # 0 1234 567 89A BCDE (hex)
1249 1249 # x YYYY DDD SSS .ext
1250 1250 if isNumber( filen[8:11] ):
1251 1251 self.setFile = int( filen[8:11] ) #inicializo mi contador de seteo al seteo del ultimo file
1252 1252 else:
1253 1253 self.setFile = -1
1254 1254 else:
1255 1255 self.setFile = -1 #inicializo mi contador de seteo
1256 1256
1257 1257 setFile = self.setFile
1258 1258 setFile += 1
1259 1259
1260 1260 file = '%s%4.4d%3.3d%3.3d%s' % (self.optchar,
1261 1261 timeTuple.tm_year,
1262 1262 timeTuple.tm_yday,
1263 1263 setFile,
1264 1264 ext )
1265 1265
1266 1266 filename = os.path.join( path, subfolder, file )
1267 1267
1268 1268 fp = open( filename,'wb' )
1269 1269
1270 1270 self.blockIndex = 0
1271 1271
1272 1272 #guardando atributos
1273 1273 self.filename = filename
1274 1274 self.subfolder = subfolder
1275 1275 self.fp = fp
1276 1276 self.setFile = setFile
1277 1277 self.flagIsNewFile = 1
1278 1278
1279 1279 self.setFirstHeader()
1280 1280
1281 1281 print 'Writing the file: %s'%self.filename
1282 1282
1283 1283 self.__writeFirstHeader()
1284 1284
1285 1285 return 1
1286 1286
1287 1287 def setup(self, dataOut, path, blocksPerFile, profilesPerBlock=64, set=0, ext=None):
1288 1288 """
1289 1289 Setea el tipo de formato en la cual sera guardada la data y escribe el First Header
1290 1290
1291 1291 Inputs:
1292 1292 path : el path destino en el cual se escribiran los files a crear
1293 1293 format : formato en el cual sera salvado un file
1294 1294 set : el setebo del file
1295 1295
1296 1296 Return:
1297 1297 0 : Si no realizo un buen seteo
1298 1298 1 : Si realizo un buen seteo
1299 1299 """
1300 1300
1301 1301 if ext == None:
1302 1302 ext = self.ext
1303 1303
1304 1304 ext = ext.lower()
1305 1305
1306 1306 self.ext = ext
1307 1307
1308 1308 self.path = path
1309 1309
1310 1310 self.setFile = set - 1
1311 1311
1312 1312 self.blocksPerFile = blocksPerFile
1313 1313
1314 1314 self.profilesPerBlock = profilesPerBlock
1315 1315
1316 1316 self.dataOut = dataOut
1317 1317
1318 1318 if not(self.setNextFile()):
1319 1319 print "There isn't a next file"
1320 1320 return 0
1321 1321
1322 1322 self.setBlockDimension()
1323 1323
1324 1324 return 1
1325 1325
1326 1326 def run(self, dataOut, **kwargs):
1327 1327
1328 1328 if not(self.isConfig):
1329 1329
1330 1330 self.setup(dataOut, **kwargs)
1331 1331 self.isConfig = True
1332 1332
1333 1333 self.putData()
1334 1334
@@ -1,3 +1,4
1 1 from jroIO_voltage import *
2 2 from jroIO_spectra import *
3 3 from jroIO_heispectra import *
4 from jroIO_amisr import * No newline at end of file
@@ -1,1024 +1,1025
1 1 import numpy
2 2
3 3 from jroproc_base import ProcessingUnit, Operation
4 4 from model.data.jrodata import Spectra
5 from model.data.jrodata import hildebrand_sekhon
5 6
6 7 class SpectraProc(ProcessingUnit):
7 8
8 9 def __init__(self):
9 10
10 11 ProcessingUnit.__init__(self)
11 12
12 13 self.buffer = None
13 14 self.firstdatatime = None
14 15 self.profIndex = 0
15 16 self.dataOut = Spectra()
16 17
17 18 def __updateObjFromInput(self):
18 19
19 20 self.dataOut.timeZone = self.dataIn.timeZone
20 21 self.dataOut.dstFlag = self.dataIn.dstFlag
21 22 self.dataOut.errorCount = self.dataIn.errorCount
22 23 self.dataOut.useLocalTime = self.dataIn.useLocalTime
23 24
24 25 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
25 26 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
26 27 self.dataOut.channelList = self.dataIn.channelList
27 28 self.dataOut.heightList = self.dataIn.heightList
28 29 self.dataOut.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
29 30 # self.dataOut.nHeights = self.dataIn.nHeights
30 31 # self.dataOut.nChannels = self.dataIn.nChannels
31 32 self.dataOut.nBaud = self.dataIn.nBaud
32 33 self.dataOut.nCode = self.dataIn.nCode
33 34 self.dataOut.code = self.dataIn.code
34 35 self.dataOut.nProfiles = self.dataOut.nFFTPoints
35 36 # self.dataOut.channelIndexList = self.dataIn.channelIndexList
36 37 self.dataOut.flagTimeBlock = self.dataIn.flagTimeBlock
37 38 self.dataOut.utctime = self.firstdatatime
38 39 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
39 40 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
40 41 # self.dataOut.flagShiftFFT = self.dataIn.flagShiftFFT
41 42 self.dataOut.nCohInt = self.dataIn.nCohInt
42 43 self.dataOut.nIncohInt = 1
43 44 # self.dataOut.ippSeconds = self.dataIn.ippSeconds
44 45 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
45 46
46 47 self.dataOut.timeInterval = self.dataIn.timeInterval*self.dataOut.nFFTPoints*self.dataOut.nIncohInt
47 48 self.dataOut.frequency = self.dataIn.frequency
48 49 self.dataOut.realtime = self.dataIn.realtime
49 50
50 51 def __getFft(self):
51 52 """
52 53 Convierte valores de Voltaje a Spectra
53 54
54 55 Affected:
55 56 self.dataOut.data_spc
56 57 self.dataOut.data_cspc
57 58 self.dataOut.data_dc
58 59 self.dataOut.heightList
59 60 self.profIndex
60 61 self.buffer
61 62 self.dataOut.flagNoData
62 63 """
63 64 fft_volt = numpy.fft.fft(self.buffer,n=self.dataOut.nFFTPoints,axis=1)
64 65 fft_volt = fft_volt.astype(numpy.dtype('complex'))
65 66 dc = fft_volt[:,0,:]
66 67
67 68 #calculo de self-spectra
68 69 fft_volt = numpy.fft.fftshift(fft_volt,axes=(1,))
69 70 spc = fft_volt * numpy.conjugate(fft_volt)
70 71 spc = spc.real
71 72
72 73 blocksize = 0
73 74 blocksize += dc.size
74 75 blocksize += spc.size
75 76
76 77 cspc = None
77 78 pairIndex = 0
78 79 if self.dataOut.pairsList != None:
79 80 #calculo de cross-spectra
80 81 cspc = numpy.zeros((self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
81 82 for pair in self.dataOut.pairsList:
82 83 cspc[pairIndex,:,:] = fft_volt[pair[0],:,:] * numpy.conjugate(fft_volt[pair[1],:,:])
83 84 pairIndex += 1
84 85 blocksize += cspc.size
85 86
86 87 self.dataOut.data_spc = spc
87 88 self.dataOut.data_cspc = cspc
88 89 self.dataOut.data_dc = dc
89 90 self.dataOut.blockSize = blocksize
90 91 self.dataOut.flagShiftFFT = False
91 92
92 93 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], ippFactor=None):
93 94
94 95 self.dataOut.flagNoData = True
95 96
96 97 if self.dataIn.type == "Spectra":
97 98 self.dataOut.copy(self.dataIn)
98 99 return True
99 100
100 101 if self.dataIn.type == "Voltage":
101 102
102 103 if nFFTPoints == None:
103 104 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
104 105
105 106 nProfiles = nFFTPoints
106 107
107 108 # if pairsList == None:
108 109 # nPairs = 0
109 110 # else:
110 111 # nPairs = len(pairsList)
111 112
112 113 if ippFactor == None:
113 114 ippFactor = 1
114 115 self.dataOut.ippFactor = ippFactor
115 116
116 117 self.dataOut.nFFTPoints = nFFTPoints
117 118 self.dataOut.pairsList = pairsList
118 119 # self.dataOut.nPairs = nPairs
119 120
120 121 if self.buffer == None:
121 122 self.buffer = numpy.zeros((self.dataIn.nChannels,
122 123 nProfiles,
123 124 self.dataIn.nHeights),
124 125 dtype='complex')
125 126
126 127
127 128 self.buffer[:,self.profIndex,:] = self.dataIn.data.copy()
128 129 self.profIndex += 1
129 130
130 131 if self.firstdatatime == None:
131 132 self.firstdatatime = self.dataIn.utctime
132 133
133 134 if self.profIndex == nProfiles:
134 135 self.__updateObjFromInput()
135 136 self.__getFft()
136 137
137 138 self.dataOut.flagNoData = False
138 139
139 140 self.buffer = None
140 141 self.firstdatatime = None
141 142 self.profIndex = 0
142 143
143 144 return True
144 145
145 146 raise ValueError, "The type object %s is not valid"%(self.dataIn.type)
146 147
147 148 def selectChannels(self, channelList):
148 149
149 150 channelIndexList = []
150 151
151 152 for channel in channelList:
152 153 index = self.dataOut.channelList.index(channel)
153 154 channelIndexList.append(index)
154 155
155 156 self.selectChannelsByIndex(channelIndexList)
156 157
157 158 def selectChannelsByIndex(self, channelIndexList):
158 159 """
159 160 Selecciona un bloque de datos en base a canales segun el channelIndexList
160 161
161 162 Input:
162 163 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
163 164
164 165 Affected:
165 166 self.dataOut.data_spc
166 167 self.dataOut.channelIndexList
167 168 self.dataOut.nChannels
168 169
169 170 Return:
170 171 None
171 172 """
172 173
173 174 for channelIndex in channelIndexList:
174 175 if channelIndex not in self.dataOut.channelIndexList:
175 176 print channelIndexList
176 177 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
177 178
178 179 # nChannels = len(channelIndexList)
179 180
180 181 data_spc = self.dataOut.data_spc[channelIndexList,:]
181 182
182 183 self.dataOut.data_spc = data_spc
183 184 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
184 185 # self.dataOut.nChannels = nChannels
185 186
186 187 return 1
187 188
188 189 def selectHeights(self, minHei, maxHei):
189 190 """
190 191 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
191 192 minHei <= height <= maxHei
192 193
193 194 Input:
194 195 minHei : valor minimo de altura a considerar
195 196 maxHei : valor maximo de altura a considerar
196 197
197 198 Affected:
198 199 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
199 200
200 201 Return:
201 202 1 si el metodo se ejecuto con exito caso contrario devuelve 0
202 203 """
203 204 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
204 205 raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei)
205 206
206 207 if (maxHei > self.dataOut.heightList[-1]):
207 208 maxHei = self.dataOut.heightList[-1]
208 209 # raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei)
209 210
210 211 minIndex = 0
211 212 maxIndex = 0
212 213 heights = self.dataOut.heightList
213 214
214 215 inda = numpy.where(heights >= minHei)
215 216 indb = numpy.where(heights <= maxHei)
216 217
217 218 try:
218 219 minIndex = inda[0][0]
219 220 except:
220 221 minIndex = 0
221 222
222 223 try:
223 224 maxIndex = indb[0][-1]
224 225 except:
225 226 maxIndex = len(heights)
226 227
227 228 self.selectHeightsByIndex(minIndex, maxIndex)
228 229
229 230 return 1
230 231
231 232 def getBeaconSignal(self, tauindex = 0, channelindex = 0, hei_ref=None):
232 233 newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
233 234
234 235 if hei_ref != None:
235 236 newheis = numpy.where(self.dataOut.heightList>hei_ref)
236 237
237 238 minIndex = min(newheis[0])
238 239 maxIndex = max(newheis[0])
239 240 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
240 241 heightList = self.dataOut.heightList[minIndex:maxIndex+1]
241 242
242 243 # determina indices
243 244 nheis = int(self.dataOut.radarControllerHeaderObj.txB/(self.dataOut.heightList[1]-self.dataOut.heightList[0]))
244 245 avg_dB = 10*numpy.log10(numpy.sum(data_spc[channelindex,:,:],axis=0))
245 246 beacon_dB = numpy.sort(avg_dB)[-nheis:]
246 247 beacon_heiIndexList = []
247 248 for val in avg_dB.tolist():
248 249 if val >= beacon_dB[0]:
249 250 beacon_heiIndexList.append(avg_dB.tolist().index(val))
250 251
251 252 #data_spc = data_spc[:,:,beacon_heiIndexList]
252 253 data_cspc = None
253 254 if self.dataOut.data_cspc != None:
254 255 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
255 256 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
256 257
257 258 data_dc = None
258 259 if self.dataOut.data_dc != None:
259 260 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
260 261 #data_dc = data_dc[:,beacon_heiIndexList]
261 262
262 263 self.dataOut.data_spc = data_spc
263 264 self.dataOut.data_cspc = data_cspc
264 265 self.dataOut.data_dc = data_dc
265 266 self.dataOut.heightList = heightList
266 267 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
267 268
268 269 return 1
269 270
270 271
271 272 def selectHeightsByIndex(self, minIndex, maxIndex):
272 273 """
273 274 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
274 275 minIndex <= index <= maxIndex
275 276
276 277 Input:
277 278 minIndex : valor de indice minimo de altura a considerar
278 279 maxIndex : valor de indice maximo de altura a considerar
279 280
280 281 Affected:
281 282 self.dataOut.data_spc
282 283 self.dataOut.data_cspc
283 284 self.dataOut.data_dc
284 285 self.dataOut.heightList
285 286
286 287 Return:
287 288 1 si el metodo se ejecuto con exito caso contrario devuelve 0
288 289 """
289 290
290 291 if (minIndex < 0) or (minIndex > maxIndex):
291 292 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
292 293
293 294 if (maxIndex >= self.dataOut.nHeights):
294 295 maxIndex = self.dataOut.nHeights-1
295 296 # raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
296 297
297 298 # nHeights = maxIndex - minIndex + 1
298 299
299 300 #Spectra
300 301 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
301 302
302 303 data_cspc = None
303 304 if self.dataOut.data_cspc != None:
304 305 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
305 306
306 307 data_dc = None
307 308 if self.dataOut.data_dc != None:
308 309 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
309 310
310 311 self.dataOut.data_spc = data_spc
311 312 self.dataOut.data_cspc = data_cspc
312 313 self.dataOut.data_dc = data_dc
313 314
314 315 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex+1]
315 316
316 317 return 1
317 318
318 319 def removeDC(self, mode = 2):
319 320 jspectra = self.dataOut.data_spc
320 321 jcspectra = self.dataOut.data_cspc
321 322
322 323
323 324 num_chan = jspectra.shape[0]
324 325 num_hei = jspectra.shape[2]
325 326
326 327 if jcspectra != None:
327 328 jcspectraExist = True
328 329 num_pairs = jcspectra.shape[0]
329 330 else: jcspectraExist = False
330 331
331 332 freq_dc = jspectra.shape[1]/2
332 333 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
333 334
334 335 if ind_vel[0]<0:
335 336 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
336 337
337 338 if mode == 1:
338 339 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
339 340
340 341 if jcspectraExist:
341 342 jcspectra[:,freq_dc,:] = (jcspectra[:,ind_vel[1],:] + jcspectra[:,ind_vel[2],:])/2
342 343
343 344 if mode == 2:
344 345
345 346 vel = numpy.array([-2,-1,1,2])
346 347 xx = numpy.zeros([4,4])
347 348
348 349 for fil in range(4):
349 350 xx[fil,:] = vel[fil]**numpy.asarray(range(4))
350 351
351 352 xx_inv = numpy.linalg.inv(xx)
352 353 xx_aux = xx_inv[0,:]
353 354
354 355 for ich in range(num_chan):
355 356 yy = jspectra[ich,ind_vel,:]
356 357 jspectra[ich,freq_dc,:] = numpy.dot(xx_aux,yy)
357 358
358 359 junkid = jspectra[ich,freq_dc,:]<=0
359 360 cjunkid = sum(junkid)
360 361
361 362 if cjunkid.any():
362 363 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
363 364
364 365 if jcspectraExist:
365 366 for ip in range(num_pairs):
366 367 yy = jcspectra[ip,ind_vel,:]
367 368 jcspectra[ip,freq_dc,:] = numpy.dot(xx_aux,yy)
368 369
369 370
370 371 self.dataOut.data_spc = jspectra
371 372 self.dataOut.data_cspc = jcspectra
372 373
373 374 return 1
374 375
375 376 def removeInterference(self, interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None):
376 377
377 378 jspectra = self.dataOut.data_spc
378 379 jcspectra = self.dataOut.data_cspc
379 380 jnoise = self.dataOut.getNoise()
380 381 num_incoh = self.dataOut.nIncohInt
381 382
382 383 num_channel = jspectra.shape[0]
383 384 num_prof = jspectra.shape[1]
384 385 num_hei = jspectra.shape[2]
385 386
386 387 #hei_interf
387 388 if hei_interf == None:
388 389 count_hei = num_hei/2 #Como es entero no importa
389 390 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
390 391 hei_interf = numpy.asarray(hei_interf)[0]
391 392 #nhei_interf
392 393 if (nhei_interf == None):
393 394 nhei_interf = 5
394 395 if (nhei_interf < 1):
395 396 nhei_interf = 1
396 397 if (nhei_interf > count_hei):
397 398 nhei_interf = count_hei
398 399 if (offhei_interf == None):
399 400 offhei_interf = 0
400 401
401 402 ind_hei = range(num_hei)
402 403 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
403 404 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
404 405 mask_prof = numpy.asarray(range(num_prof))
405 406 num_mask_prof = mask_prof.size
406 407 comp_mask_prof = [0, num_prof/2]
407 408
408 409
409 410 #noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
410 411 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
411 412 jnoise = numpy.nan
412 413 noise_exist = jnoise[0] < numpy.Inf
413 414
414 415 #Subrutina de Remocion de la Interferencia
415 416 for ich in range(num_channel):
416 417 #Se ordena los espectros segun su potencia (menor a mayor)
417 418 power = jspectra[ich,mask_prof,:]
418 419 power = power[:,hei_interf]
419 420 power = power.sum(axis = 0)
420 421 psort = power.ravel().argsort()
421 422
422 423 #Se estima la interferencia promedio en los Espectros de Potencia empleando
423 424 junkspc_interf = jspectra[ich,:,hei_interf[psort[range(offhei_interf, nhei_interf + offhei_interf)]]]
424 425
425 426 if noise_exist:
426 427 # tmp_noise = jnoise[ich] / num_prof
427 428 tmp_noise = jnoise[ich]
428 429 junkspc_interf = junkspc_interf - tmp_noise
429 430 #junkspc_interf[:,comp_mask_prof] = 0
430 431
431 432 jspc_interf = junkspc_interf.sum(axis = 0) / nhei_interf
432 433 jspc_interf = jspc_interf.transpose()
433 434 #Calculando el espectro de interferencia promedio
434 435 noiseid = numpy.where(jspc_interf <= tmp_noise/ math.sqrt(num_incoh))
435 436 noiseid = noiseid[0]
436 437 cnoiseid = noiseid.size
437 438 interfid = numpy.where(jspc_interf > tmp_noise/ math.sqrt(num_incoh))
438 439 interfid = interfid[0]
439 440 cinterfid = interfid.size
440 441
441 442 if (cnoiseid > 0): jspc_interf[noiseid] = 0
442 443
443 444 #Expandiendo los perfiles a limpiar
444 445 if (cinterfid > 0):
445 446 new_interfid = (numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof)%num_prof
446 447 new_interfid = numpy.asarray(new_interfid)
447 448 new_interfid = {x for x in new_interfid}
448 449 new_interfid = numpy.array(list(new_interfid))
449 450 new_cinterfid = new_interfid.size
450 451 else: new_cinterfid = 0
451 452
452 453 for ip in range(new_cinterfid):
453 454 ind = junkspc_interf[:,new_interfid[ip]].ravel().argsort()
454 455 jspc_interf[new_interfid[ip]] = junkspc_interf[ind[nhei_interf/2],new_interfid[ip]]
455 456
456 457
457 458 jspectra[ich,:,ind_hei] = jspectra[ich,:,ind_hei] - jspc_interf #Corregir indices
458 459
459 460 #Removiendo la interferencia del punto de mayor interferencia
460 461 ListAux = jspc_interf[mask_prof].tolist()
461 462 maxid = ListAux.index(max(ListAux))
462 463
463 464
464 465 if cinterfid > 0:
465 466 for ip in range(cinterfid*(interf == 2) - 1):
466 467 ind = (jspectra[ich,interfid[ip],:] < tmp_noise*(1 + 1/math.sqrt(num_incoh))).nonzero()
467 468 cind = len(ind)
468 469
469 470 if (cind > 0):
470 471 jspectra[ich,interfid[ip],ind] = tmp_noise*(1 + (numpy.random.uniform(cind) - 0.5)/math.sqrt(num_incoh))
471 472
472 473 ind = numpy.array([-2,-1,1,2])
473 474 xx = numpy.zeros([4,4])
474 475
475 476 for id1 in range(4):
476 477 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
477 478
478 479 xx_inv = numpy.linalg.inv(xx)
479 480 xx = xx_inv[:,0]
480 481 ind = (ind + maxid + num_mask_prof)%num_mask_prof
481 482 yy = jspectra[ich,mask_prof[ind],:]
482 483 jspectra[ich,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
483 484
484 485
485 486 indAux = (jspectra[ich,:,:] < tmp_noise*(1-1/math.sqrt(num_incoh))).nonzero()
486 487 jspectra[ich,indAux[0],indAux[1]] = tmp_noise * (1 - 1/math.sqrt(num_incoh))
487 488
488 489 #Remocion de Interferencia en el Cross Spectra
489 490 if jcspectra == None: return jspectra, jcspectra
490 491 num_pairs = jcspectra.size/(num_prof*num_hei)
491 492 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
492 493
493 494 for ip in range(num_pairs):
494 495
495 496 #-------------------------------------------
496 497
497 498 cspower = numpy.abs(jcspectra[ip,mask_prof,:])
498 499 cspower = cspower[:,hei_interf]
499 500 cspower = cspower.sum(axis = 0)
500 501
501 502 cspsort = cspower.ravel().argsort()
502 503 junkcspc_interf = jcspectra[ip,:,hei_interf[cspsort[range(offhei_interf, nhei_interf + offhei_interf)]]]
503 504 junkcspc_interf = junkcspc_interf.transpose()
504 505 jcspc_interf = junkcspc_interf.sum(axis = 1)/nhei_interf
505 506
506 507 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
507 508
508 509 median_real = numpy.median(numpy.real(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
509 510 median_imag = numpy.median(numpy.imag(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
510 511 junkcspc_interf[comp_mask_prof,:] = numpy.complex(median_real, median_imag)
511 512
512 513 for iprof in range(num_prof):
513 514 ind = numpy.abs(junkcspc_interf[iprof,:]).ravel().argsort()
514 515 jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf/2]]
515 516
516 517 #Removiendo la Interferencia
517 518 jcspectra[ip,:,ind_hei] = jcspectra[ip,:,ind_hei] - jcspc_interf
518 519
519 520 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
520 521 maxid = ListAux.index(max(ListAux))
521 522
522 523 ind = numpy.array([-2,-1,1,2])
523 524 xx = numpy.zeros([4,4])
524 525
525 526 for id1 in range(4):
526 527 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
527 528
528 529 xx_inv = numpy.linalg.inv(xx)
529 530 xx = xx_inv[:,0]
530 531
531 532 ind = (ind + maxid + num_mask_prof)%num_mask_prof
532 533 yy = jcspectra[ip,mask_prof[ind],:]
533 534 jcspectra[ip,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
534 535
535 536 #Guardar Resultados
536 537 self.dataOut.data_spc = jspectra
537 538 self.dataOut.data_cspc = jcspectra
538 539
539 540 return 1
540 541
541 542 def setRadarFrequency(self, frequency=None):
542 543 if frequency != None:
543 544 self.dataOut.frequency = frequency
544 545
545 546 return 1
546 547
547 548 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
548 549 #validacion de rango
549 550 if minHei == None:
550 551 minHei = self.dataOut.heightList[0]
551 552
552 553 if maxHei == None:
553 554 maxHei = self.dataOut.heightList[-1]
554 555
555 556 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
556 557 print 'minHei: %.2f is out of the heights range'%(minHei)
557 558 print 'minHei is setting to %.2f'%(self.dataOut.heightList[0])
558 559 minHei = self.dataOut.heightList[0]
559 560
560 561 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
561 562 print 'maxHei: %.2f is out of the heights range'%(maxHei)
562 563 print 'maxHei is setting to %.2f'%(self.dataOut.heightList[-1])
563 564 maxHei = self.dataOut.heightList[-1]
564 565
565 566 # validacion de velocidades
566 567 velrange = self.dataOut.getVelRange(1)
567 568
568 569 if minVel == None:
569 570 minVel = velrange[0]
570 571
571 572 if maxVel == None:
572 573 maxVel = velrange[-1]
573 574
574 575 if (minVel < velrange[0]) or (minVel > maxVel):
575 576 print 'minVel: %.2f is out of the velocity range'%(minVel)
576 577 print 'minVel is setting to %.2f'%(velrange[0])
577 578 minVel = velrange[0]
578 579
579 580 if (maxVel > velrange[-1]) or (maxVel < minVel):
580 581 print 'maxVel: %.2f is out of the velocity range'%(maxVel)
581 582 print 'maxVel is setting to %.2f'%(velrange[-1])
582 583 maxVel = velrange[-1]
583 584
584 585 # seleccion de indices para rango
585 586 minIndex = 0
586 587 maxIndex = 0
587 588 heights = self.dataOut.heightList
588 589
589 590 inda = numpy.where(heights >= minHei)
590 591 indb = numpy.where(heights <= maxHei)
591 592
592 593 try:
593 594 minIndex = inda[0][0]
594 595 except:
595 596 minIndex = 0
596 597
597 598 try:
598 599 maxIndex = indb[0][-1]
599 600 except:
600 601 maxIndex = len(heights)
601 602
602 603 if (minIndex < 0) or (minIndex > maxIndex):
603 604 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
604 605
605 606 if (maxIndex >= self.dataOut.nHeights):
606 607 maxIndex = self.dataOut.nHeights-1
607 608
608 609 # seleccion de indices para velocidades
609 610 indminvel = numpy.where(velrange >= minVel)
610 611 indmaxvel = numpy.where(velrange <= maxVel)
611 612 try:
612 613 minIndexVel = indminvel[0][0]
613 614 except:
614 615 minIndexVel = 0
615 616
616 617 try:
617 618 maxIndexVel = indmaxvel[0][-1]
618 619 except:
619 620 maxIndexVel = len(velrange)
620 621
621 622 #seleccion del espectro
622 623 data_spc = self.dataOut.data_spc[:,minIndexVel:maxIndexVel+1,minIndex:maxIndex+1]
623 624 #estimacion de ruido
624 625 noise = numpy.zeros(self.dataOut.nChannels)
625 626
626 627 for channel in range(self.dataOut.nChannels):
627 628 daux = data_spc[channel,:,:]
628 629 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
629 630
630 631 self.dataOut.noise_estimation = noise.copy()
631 632
632 633 return 1
633 634
634 635 class IncohInt(Operation):
635 636
636 637
637 638 __profIndex = 0
638 639 __withOverapping = False
639 640
640 641 __byTime = False
641 642 __initime = None
642 643 __lastdatatime = None
643 644 __integrationtime = None
644 645
645 646 __buffer_spc = None
646 647 __buffer_cspc = None
647 648 __buffer_dc = None
648 649
649 650 __dataReady = False
650 651
651 652 __timeInterval = None
652 653
653 654 n = None
654 655
655 656
656 657
657 658 def __init__(self):
658 659
659 660 Operation.__init__(self)
660 661 # self.isConfig = False
661 662
662 663 def setup(self, n=None, timeInterval=None, overlapping=False):
663 664 """
664 665 Set the parameters of the integration class.
665 666
666 667 Inputs:
667 668
668 669 n : Number of coherent integrations
669 670 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
670 671 overlapping :
671 672
672 673 """
673 674
674 675 self.__initime = None
675 676 self.__lastdatatime = 0
676 677 self.__buffer_spc = None
677 678 self.__buffer_cspc = None
678 679 self.__buffer_dc = None
679 680 self.__dataReady = False
680 681
681 682
682 683 if n == None and timeInterval == None:
683 684 raise ValueError, "n or timeInterval should be specified ..."
684 685
685 686 if n != None:
686 687 self.n = n
687 688 self.__byTime = False
688 689 else:
689 690 self.__integrationtime = timeInterval #if (type(timeInterval)!=integer) -> change this line
690 691 self.n = 9999
691 692 self.__byTime = True
692 693
693 694 if overlapping:
694 695 self.__withOverapping = True
695 696 else:
696 697 self.__withOverapping = False
697 698 self.__buffer_spc = 0
698 699 self.__buffer_cspc = 0
699 700 self.__buffer_dc = 0
700 701
701 702 self.__profIndex = 0
702 703
703 704 def putData(self, data_spc, data_cspc, data_dc):
704 705
705 706 """
706 707 Add a profile to the __buffer_spc and increase in one the __profileIndex
707 708
708 709 """
709 710
710 711 if not self.__withOverapping:
711 712 self.__buffer_spc += data_spc
712 713
713 714 if data_cspc == None:
714 715 self.__buffer_cspc = None
715 716 else:
716 717 self.__buffer_cspc += data_cspc
717 718
718 719 if data_dc == None:
719 720 self.__buffer_dc = None
720 721 else:
721 722 self.__buffer_dc += data_dc
722 723
723 724 self.__profIndex += 1
724 725 return
725 726
726 727 #Overlapping data
727 728 nChannels, nFFTPoints, nHeis = data_spc.shape
728 729 data_spc = numpy.reshape(data_spc, (1, nChannels, nFFTPoints, nHeis))
729 730 if data_cspc != None:
730 731 data_cspc = numpy.reshape(data_cspc, (1, -1, nFFTPoints, nHeis))
731 732 if data_dc != None:
732 733 data_dc = numpy.reshape(data_dc, (1, -1, nHeis))
733 734
734 735 #If the buffer is empty then it takes the data value
735 736 if self.__buffer_spc == None:
736 737 self.__buffer_spc = data_spc
737 738
738 739 if data_cspc == None:
739 740 self.__buffer_cspc = None
740 741 else:
741 742 self.__buffer_cspc += data_cspc
742 743
743 744 if data_dc == None:
744 745 self.__buffer_dc = None
745 746 else:
746 747 self.__buffer_dc += data_dc
747 748
748 749 self.__profIndex += 1
749 750 return
750 751
751 752 #If the buffer length is lower than n then stakcing the data value
752 753 if self.__profIndex < self.n:
753 754 self.__buffer_spc = numpy.vstack((self.__buffer_spc, data_spc))
754 755
755 756 if data_cspc != None:
756 757 self.__buffer_cspc = numpy.vstack((self.__buffer_cspc, data_cspc))
757 758
758 759 if data_dc != None:
759 760 self.__buffer_dc = numpy.vstack((self.__buffer_dc, data_dc))
760 761
761 762 self.__profIndex += 1
762 763 return
763 764
764 765 #If the buffer length is equal to n then replacing the last buffer value with the data value
765 766 self.__buffer_spc = numpy.roll(self.__buffer_spc, -1, axis=0)
766 767 self.__buffer_spc[self.n-1] = data_spc
767 768
768 769 if data_cspc != None:
769 770 self.__buffer_cspc = numpy.roll(self.__buffer_cspc, -1, axis=0)
770 771 self.__buffer_cspc[self.n-1] = data_cspc
771 772
772 773 if data_dc != None:
773 774 self.__buffer_dc = numpy.roll(self.__buffer_dc, -1, axis=0)
774 775 self.__buffer_dc[self.n-1] = data_dc
775 776
776 777 self.__profIndex = self.n
777 778 return
778 779
779 780
780 781 def pushData(self):
781 782 """
782 783 Return the sum of the last profiles and the profiles used in the sum.
783 784
784 785 Affected:
785 786
786 787 self.__profileIndex
787 788
788 789 """
789 790 data_spc = None
790 791 data_cspc = None
791 792 data_dc = None
792 793
793 794 if not self.__withOverapping:
794 795 data_spc = self.__buffer_spc
795 796 data_cspc = self.__buffer_cspc
796 797 data_dc = self.__buffer_dc
797 798
798 799 n = self.__profIndex
799 800
800 801 self.__buffer_spc = 0
801 802 self.__buffer_cspc = 0
802 803 self.__buffer_dc = 0
803 804 self.__profIndex = 0
804 805
805 806 return data_spc, data_cspc, data_dc, n
806 807
807 808 #Integration with Overlapping
808 809 data_spc = numpy.sum(self.__buffer_spc, axis=0)
809 810
810 811 if self.__buffer_cspc != None:
811 812 data_cspc = numpy.sum(self.__buffer_cspc, axis=0)
812 813
813 814 if self.__buffer_dc != None:
814 815 data_dc = numpy.sum(self.__buffer_dc, axis=0)
815 816
816 817 n = self.__profIndex
817 818
818 819 return data_spc, data_cspc, data_dc, n
819 820
820 821 def byProfiles(self, *args):
821 822
822 823 self.__dataReady = False
823 824 avgdata_spc = None
824 825 avgdata_cspc = None
825 826 avgdata_dc = None
826 827 # n = None
827 828
828 829 self.putData(*args)
829 830
830 831 if self.__profIndex == self.n:
831 832
832 833 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
833 834 self.__dataReady = True
834 835
835 836 return avgdata_spc, avgdata_cspc, avgdata_dc
836 837
837 838 def byTime(self, datatime, *args):
838 839
839 840 self.__dataReady = False
840 841 avgdata_spc = None
841 842 avgdata_cspc = None
842 843 avgdata_dc = None
843 844 n = None
844 845
845 846 self.putData(*args)
846 847
847 848 if (datatime - self.__initime) >= self.__integrationtime:
848 849 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
849 850 self.n = n
850 851 self.__dataReady = True
851 852
852 853 return avgdata_spc, avgdata_cspc, avgdata_dc
853 854
854 855 def integrate(self, datatime, *args):
855 856
856 857 if self.__initime == None:
857 858 self.__initime = datatime
858 859
859 860 if self.__byTime:
860 861 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(datatime, *args)
861 862 else:
862 863 avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
863 864
864 865 self.__lastdatatime = datatime
865 866
866 867 if avgdata_spc == None:
867 868 return None, None, None, None
868 869
869 870 avgdatatime = self.__initime
870 871 try:
871 872 self.__timeInterval = (self.__lastdatatime - self.__initime)/(self.n - 1)
872 873 except:
873 874 self.__timeInterval = self.__lastdatatime - self.__initime
874 875
875 876 deltatime = datatime -self.__lastdatatime
876 877
877 878 if not self.__withOverapping:
878 879 self.__initime = datatime
879 880 else:
880 881 self.__initime += deltatime
881 882
882 883 return avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc
883 884
884 885 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
885 886
886 887 if n==1:
887 888 dataOut.flagNoData = False
888 889 return
889 890
890 891 if not self.isConfig:
891 892 self.setup(n, timeInterval, overlapping)
892 893 self.isConfig = True
893 894
894 895 avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
895 896 dataOut.data_spc,
896 897 dataOut.data_cspc,
897 898 dataOut.data_dc)
898 899
899 900 # dataOut.timeInterval *= n
900 901 dataOut.flagNoData = True
901 902
902 903 if self.__dataReady:
903 904
904 905 dataOut.data_spc = avgdata_spc
905 906 dataOut.data_cspc = avgdata_cspc
906 907 dataOut.data_dc = avgdata_dc
907 908
908 909 dataOut.nIncohInt *= self.n
909 910 dataOut.utctime = avgdatatime
910 911 #dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt * dataOut.nIncohInt * dataOut.nFFTPoints
911 912 dataOut.timeInterval = self.__timeInterval*self.n
912 913 dataOut.flagNoData = False
913 914
914 915 class ProfileConcat(Operation):
915 916
916 917 isConfig = False
917 918 buffer = None
918 919
919 920 def __init__(self):
920 921
921 922 Operation.__init__(self)
922 923 self.profileIndex = 0
923 924
924 925 def reset(self):
925 926 self.buffer = numpy.zeros_like(self.buffer)
926 927 self.start_index = 0
927 928 self.times = 1
928 929
929 930 def setup(self, data, m, n=1):
930 931 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
931 932 self.profiles = data.shape[1]
932 933 self.start_index = 0
933 934 self.times = 1
934 935
935 936 def concat(self, data):
936 937
937 938 self.buffer[:,self.start_index:self.profiles*self.times] = data.copy()
938 939 self.start_index = self.start_index + self.profiles
939 940
940 941 def run(self, dataOut, m):
941 942
942 943 dataOut.flagNoData = True
943 944
944 945 if not self.isConfig:
945 946 self.setup(dataOut.data, m, 1)
946 947 self.isConfig = True
947 948
948 949 self.concat(dataOut.data)
949 950 self.times += 1
950 951 if self.times > m:
951 952 dataOut.data = self.buffer
952 953 self.reset()
953 954 dataOut.flagNoData = False
954 955 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
955 956 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
956 957 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * 5
957 958 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
958 959
959 960 class ProfileSelector(Operation):
960 961
961 962 profileIndex = None
962 963 # Tamanho total de los perfiles
963 964 nProfiles = None
964 965
965 966 def __init__(self):
966 967
967 968 Operation.__init__(self)
968 969 self.profileIndex = 0
969 970
970 971 def incIndex(self):
971 972 self.profileIndex += 1
972 973
973 974 if self.profileIndex >= self.nProfiles:
974 975 self.profileIndex = 0
975 976
976 977 def isProfileInRange(self, minIndex, maxIndex):
977 978
978 979 if self.profileIndex < minIndex:
979 980 return False
980 981
981 982 if self.profileIndex > maxIndex:
982 983 return False
983 984
984 985 return True
985 986
986 987 def isProfileInList(self, profileList):
987 988
988 989 if self.profileIndex not in profileList:
989 990 return False
990 991
991 992 return True
992 993
993 994 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None):
994 995
995 996 dataOut.flagNoData = True
996 997 self.nProfiles = dataOut.nProfiles
997 998
998 999 if profileList != None:
999 1000 if self.isProfileInList(profileList):
1000 1001 dataOut.flagNoData = False
1001 1002
1002 1003 self.incIndex()
1003 1004 return 1
1004 1005
1005 1006
1006 1007 elif profileRangeList != None:
1007 1008 minIndex = profileRangeList[0]
1008 1009 maxIndex = profileRangeList[1]
1009 1010 if self.isProfileInRange(minIndex, maxIndex):
1010 1011 dataOut.flagNoData = False
1011 1012
1012 1013 self.incIndex()
1013 1014 return 1
1014 1015 elif beam != None:
1015 1016 if self.isProfileInList(dataOut.beamRangeDict[beam]):
1016 1017 dataOut.flagNoData = False
1017 1018
1018 1019 self.incIndex()
1019 1020 return 1
1020 1021
1021 1022 else:
1022 1023 raise ValueError, "ProfileSelector needs profileList or profileRangeList"
1023 1024
1024 1025 return 0
@@ -1,518 +1,524
1 1 import numpy
2 2
3 3 from jroproc_base import ProcessingUnit, Operation
4 4 from model.data.jrodata import Voltage
5 5
6 6 class VoltageProc(ProcessingUnit):
7 7
8 8
9 9 def __init__(self):
10 10
11 11 ProcessingUnit.__init__(self)
12 12
13 13 # self.objectDict = {}
14 14 self.dataOut = Voltage()
15 15 self.flip = 1
16 16
17 17 def run(self):
18 if self.dataIn.type == 'AMISR':
19 self.__updateObjFromAmisrInput()
20
21 if self.dataIn.type == 'Voltage':
18 22 self.dataOut.copy(self.dataIn)
19 23
20 # def __updateObjFromAmisrInput(self):
21 #
22 # self.dataOut.timeZone = self.dataIn.timeZone
23 # self.dataOut.dstFlag = self.dataIn.dstFlag
24 # self.dataOut.errorCount = self.dataIn.errorCount
25 # self.dataOut.useLocalTime = self.dataIn.useLocalTime
26 #
27 # self.dataOut.flagNoData = self.dataIn.flagNoData
28 # self.dataOut.data = self.dataIn.data
29 # self.dataOut.utctime = self.dataIn.utctime
30 # self.dataOut.channelList = self.dataIn.channelList
31 # self.dataOut.timeInterval = self.dataIn.timeInterval
32 # self.dataOut.heightList = self.dataIn.heightList
33 # self.dataOut.nProfiles = self.dataIn.nProfiles
34 #
35 # self.dataOut.nCohInt = self.dataIn.nCohInt
36 # self.dataOut.ippSeconds = self.dataIn.ippSeconds
37 # self.dataOut.frequency = self.dataIn.frequency
24 # self.dataOut.copy(self.dataIn)
25
26 def __updateObjFromAmisrInput(self):
27
28 self.dataOut.timeZone = self.dataIn.timeZone
29 self.dataOut.dstFlag = self.dataIn.dstFlag
30 self.dataOut.errorCount = self.dataIn.errorCount
31 self.dataOut.useLocalTime = self.dataIn.useLocalTime
32
33 self.dataOut.flagNoData = self.dataIn.flagNoData
34 self.dataOut.data = self.dataIn.data
35 self.dataOut.utctime = self.dataIn.utctime
36 self.dataOut.channelList = self.dataIn.channelList
37 self.dataOut.timeInterval = self.dataIn.timeInterval
38 self.dataOut.heightList = self.dataIn.heightList
39 self.dataOut.nProfiles = self.dataIn.nProfiles
40
41 self.dataOut.nCohInt = self.dataIn.nCohInt
42 self.dataOut.ippSeconds = self.dataIn.ippSeconds
43 self.dataOut.frequency = self.dataIn.frequency
38 44 #
39 45 # pass#
40 46 #
41 47 # def init(self):
42 48 #
43 49 #
44 50 # if self.dataIn.type == 'AMISR':
45 51 # self.__updateObjFromAmisrInput()
46 52 #
47 53 # if self.dataIn.type == 'Voltage':
48 54 # self.dataOut.copy(self.dataIn)
49 55 # # No necesita copiar en cada init() los atributos de dataIn
50 56 # # la copia deberia hacerse por cada nuevo bloque de datos
51 57
52 58 def selectChannels(self, channelList):
53 59
54 60 channelIndexList = []
55 61
56 62 for channel in channelList:
57 63 index = self.dataOut.channelList.index(channel)
58 64 channelIndexList.append(index)
59 65
60 66 self.selectChannelsByIndex(channelIndexList)
61 67
62 68 def selectChannelsByIndex(self, channelIndexList):
63 69 """
64 70 Selecciona un bloque de datos en base a canales segun el channelIndexList
65 71
66 72 Input:
67 73 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
68 74
69 75 Affected:
70 76 self.dataOut.data
71 77 self.dataOut.channelIndexList
72 78 self.dataOut.nChannels
73 79 self.dataOut.m_ProcessingHeader.totalSpectra
74 80 self.dataOut.systemHeaderObj.numChannels
75 81 self.dataOut.m_ProcessingHeader.blockSize
76 82
77 83 Return:
78 84 None
79 85 """
80 86
81 87 for channelIndex in channelIndexList:
82 88 if channelIndex not in self.dataOut.channelIndexList:
83 89 print channelIndexList
84 90 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
85 91
86 92 # nChannels = len(channelIndexList)
87 93
88 94 data = self.dataOut.data[channelIndexList,:]
89 95
90 96 self.dataOut.data = data
91 97 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
92 98 # self.dataOut.nChannels = nChannels
93 99
94 100 return 1
95 101
96 102 def selectHeights(self, minHei=None, maxHei=None):
97 103 """
98 104 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
99 105 minHei <= height <= maxHei
100 106
101 107 Input:
102 108 minHei : valor minimo de altura a considerar
103 109 maxHei : valor maximo de altura a considerar
104 110
105 111 Affected:
106 112 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
107 113
108 114 Return:
109 115 1 si el metodo se ejecuto con exito caso contrario devuelve 0
110 116 """
111 117
112 118 if minHei == None:
113 119 minHei = self.dataOut.heightList[0]
114 120
115 121 if maxHei == None:
116 122 maxHei = self.dataOut.heightList[-1]
117 123
118 124 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
119 125 raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei)
120 126
121 127
122 128 if (maxHei > self.dataOut.heightList[-1]):
123 129 maxHei = self.dataOut.heightList[-1]
124 130 # raise ValueError, "some value in (%d,%d) is not valid" % (minHei, maxHei)
125 131
126 132 minIndex = 0
127 133 maxIndex = 0
128 134 heights = self.dataOut.heightList
129 135
130 136 inda = numpy.where(heights >= minHei)
131 137 indb = numpy.where(heights <= maxHei)
132 138
133 139 try:
134 140 minIndex = inda[0][0]
135 141 except:
136 142 minIndex = 0
137 143
138 144 try:
139 145 maxIndex = indb[0][-1]
140 146 except:
141 147 maxIndex = len(heights)
142 148
143 149 self.selectHeightsByIndex(minIndex, maxIndex)
144 150
145 151 return 1
146 152
147 153
148 154 def selectHeightsByIndex(self, minIndex, maxIndex):
149 155 """
150 156 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
151 157 minIndex <= index <= maxIndex
152 158
153 159 Input:
154 160 minIndex : valor de indice minimo de altura a considerar
155 161 maxIndex : valor de indice maximo de altura a considerar
156 162
157 163 Affected:
158 164 self.dataOut.data
159 165 self.dataOut.heightList
160 166
161 167 Return:
162 168 1 si el metodo se ejecuto con exito caso contrario devuelve 0
163 169 """
164 170
165 171 if (minIndex < 0) or (minIndex > maxIndex):
166 172 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
167 173
168 174 if (maxIndex >= self.dataOut.nHeights):
169 175 maxIndex = self.dataOut.nHeights-1
170 176 # raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
171 177
172 178 # nHeights = maxIndex - minIndex + 1
173 179
174 180 #voltage
175 181 data = self.dataOut.data[:,minIndex:maxIndex+1]
176 182
177 183 # firstHeight = self.dataOut.heightList[minIndex]
178 184
179 185 self.dataOut.data = data
180 186 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex+1]
181 187
182 188 return 1
183 189
184 190
185 191 def filterByHeights(self, window):
186 192 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
187 193
188 194 if window == None:
189 195 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
190 196
191 197 newdelta = deltaHeight * window
192 198 r = self.dataOut.data.shape[1] % window
193 199 buffer = self.dataOut.data[:,0:self.dataOut.data.shape[1]-r]
194 200 buffer = buffer.reshape(self.dataOut.data.shape[0],self.dataOut.data.shape[1]/window,window)
195 201 buffer = numpy.sum(buffer,2)
196 202 self.dataOut.data = buffer
197 203 self.dataOut.heightList = numpy.arange(self.dataOut.heightList[0],newdelta*(self.dataOut.nHeights-r)/window,newdelta)
198 204 self.dataOut.windowOfFilter = window
199 205
200 206 def deFlip(self):
201 207 self.dataOut.data *= self.flip
202 208 self.flip *= -1.
203 209
204 210 def setRadarFrequency(self, frequency=None):
205 211 if frequency != None:
206 212 self.dataOut.frequency = frequency
207 213
208 214 return 1
209 215
210 216 class CohInt(Operation):
211 217
212 218 isConfig = False
213 219
214 220 __profIndex = 0
215 221 __withOverapping = False
216 222
217 223 __byTime = False
218 224 __initime = None
219 225 __lastdatatime = None
220 226 __integrationtime = None
221 227
222 228 __buffer = None
223 229
224 230 __dataReady = False
225 231
226 232 n = None
227 233
228 234
229 235 def __init__(self):
230 236
231 237 Operation.__init__(self)
232 238
233 239 # self.isConfig = False
234 240
235 241 def setup(self, n=None, timeInterval=None, overlapping=False):
236 242 """
237 243 Set the parameters of the integration class.
238 244
239 245 Inputs:
240 246
241 247 n : Number of coherent integrations
242 248 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
243 249 overlapping :
244 250
245 251 """
246 252
247 253 self.__initime = None
248 254 self.__lastdatatime = 0
249 255 self.__buffer = None
250 256 self.__dataReady = False
251 257
252 258
253 259 if n == None and timeInterval == None:
254 260 raise ValueError, "n or timeInterval should be specified ..."
255 261
256 262 if n != None:
257 263 self.n = n
258 264 self.__byTime = False
259 265 else:
260 266 self.__integrationtime = timeInterval * 60. #if (type(timeInterval)!=integer) -> change this line
261 267 self.n = 9999
262 268 self.__byTime = True
263 269
264 270 if overlapping:
265 271 self.__withOverapping = True
266 272 self.__buffer = None
267 273 else:
268 274 self.__withOverapping = False
269 275 self.__buffer = 0
270 276
271 277 self.__profIndex = 0
272 278
273 279 def putData(self, data):
274 280
275 281 """
276 282 Add a profile to the __buffer and increase in one the __profileIndex
277 283
278 284 """
279 285
280 286 if not self.__withOverapping:
281 287 self.__buffer += data.copy()
282 288 self.__profIndex += 1
283 289 return
284 290
285 291 #Overlapping data
286 292 nChannels, nHeis = data.shape
287 293 data = numpy.reshape(data, (1, nChannels, nHeis))
288 294
289 295 #If the buffer is empty then it takes the data value
290 296 if self.__buffer == None:
291 297 self.__buffer = data
292 298 self.__profIndex += 1
293 299 return
294 300
295 301 #If the buffer length is lower than n then stakcing the data value
296 302 if self.__profIndex < self.n:
297 303 self.__buffer = numpy.vstack((self.__buffer, data))
298 304 self.__profIndex += 1
299 305 return
300 306
301 307 #If the buffer length is equal to n then replacing the last buffer value with the data value
302 308 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
303 309 self.__buffer[self.n-1] = data
304 310 self.__profIndex = self.n
305 311 return
306 312
307 313
308 314 def pushData(self):
309 315 """
310 316 Return the sum of the last profiles and the profiles used in the sum.
311 317
312 318 Affected:
313 319
314 320 self.__profileIndex
315 321
316 322 """
317 323
318 324 if not self.__withOverapping:
319 325 data = self.__buffer
320 326 n = self.__profIndex
321 327
322 328 self.__buffer = 0
323 329 self.__profIndex = 0
324 330
325 331 return data, n
326 332
327 333 #Integration with Overlapping
328 334 data = numpy.sum(self.__buffer, axis=0)
329 335 n = self.__profIndex
330 336
331 337 return data, n
332 338
333 339 def byProfiles(self, data):
334 340
335 341 self.__dataReady = False
336 342 avgdata = None
337 343 # n = None
338 344
339 345 self.putData(data)
340 346
341 347 if self.__profIndex == self.n:
342 348
343 349 avgdata, n = self.pushData()
344 350 self.__dataReady = True
345 351
346 352 return avgdata
347 353
348 354 def byTime(self, data, datatime):
349 355
350 356 self.__dataReady = False
351 357 avgdata = None
352 358 n = None
353 359
354 360 self.putData(data)
355 361
356 362 if (datatime - self.__initime) >= self.__integrationtime:
357 363 avgdata, n = self.pushData()
358 364 self.n = n
359 365 self.__dataReady = True
360 366
361 367 return avgdata
362 368
363 369 def integrate(self, data, datatime=None):
364 370
365 371 if self.__initime == None:
366 372 self.__initime = datatime
367 373
368 374 if self.__byTime:
369 375 avgdata = self.byTime(data, datatime)
370 376 else:
371 377 avgdata = self.byProfiles(data)
372 378
373 379
374 380 self.__lastdatatime = datatime
375 381
376 382 if avgdata == None:
377 383 return None, None
378 384
379 385 avgdatatime = self.__initime
380 386
381 387 deltatime = datatime -self.__lastdatatime
382 388
383 389 if not self.__withOverapping:
384 390 self.__initime = datatime
385 391 else:
386 392 self.__initime += deltatime
387 393
388 394 return avgdata, avgdatatime
389 395
390 396 def run(self, dataOut, **kwargs):
391 397
392 398 if not self.isConfig:
393 399 self.setup(**kwargs)
394 400 self.isConfig = True
395 401
396 402 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
397 403
398 404 # dataOut.timeInterval *= n
399 405 dataOut.flagNoData = True
400 406
401 407 if self.__dataReady:
402 408 dataOut.data = avgdata
403 409 dataOut.nCohInt *= self.n
404 410 dataOut.utctime = avgdatatime
405 411 dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
406 412 dataOut.flagNoData = False
407 413
408 414 class Decoder(Operation):
409 415
410 416 isConfig = False
411 417 __profIndex = 0
412 418
413 419 code = None
414 420
415 421 nCode = None
416 422 nBaud = None
417 423
418 424 def __init__(self):
419 425
420 426 Operation.__init__(self)
421 427 # self.isConfig = False
422 428
423 429 def setup(self, code, shape):
424 430
425 431 self.__profIndex = 0
426 432
427 433 self.code = code
428 434
429 435 self.nCode = len(code)
430 436 self.nBaud = len(code[0])
431 437
432 438 self.__nChannels, self.__nHeis = shape
433 439
434 440 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
435 441
436 442 __codeBuffer[:,0:self.nBaud] = self.code
437 443
438 444 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
439 445
440 446 self.ndatadec = self.__nHeis - self.nBaud + 1
441 447
442 448 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
443 449
444 450 def convolutionInFreq(self, data):
445 451
446 452 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
447 453
448 454 fft_data = numpy.fft.fft(data, axis=1)
449 455
450 456 conv = fft_data*fft_code
451 457
452 458 data = numpy.fft.ifft(conv,axis=1)
453 459
454 460 datadec = data[:,:-self.nBaud+1]
455 461
456 462 return datadec
457 463
458 464 def convolutionInFreqOpt(self, data):
459 465
460 466 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
461 467
462 468 data = cfunctions.decoder(fft_code, data)
463 469
464 470 datadec = data[:,:-self.nBaud+1]
465 471
466 472 return datadec
467 473
468 474 def convolutionInTime(self, data):
469 475
470 476 code = self.code[self.__profIndex]
471 477
472 478 for i in range(self.__nChannels):
473 479 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='valid')
474 480
475 481 return self.datadecTime
476 482
477 483 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0):
478 484
479 485 if code == None:
480 486 code = dataOut.code
481 487 else:
482 488 code = numpy.array(code).reshape(nCode,nBaud)
483 489 dataOut.code = code
484 490 dataOut.nCode = nCode
485 491 dataOut.nBaud = nBaud
486 492 dataOut.radarControllerHeaderObj.code = code
487 493 dataOut.radarControllerHeaderObj.nCode = nCode
488 494 dataOut.radarControllerHeaderObj.nBaud = nBaud
489 495
490 496
491 497 if not self.isConfig:
492 498
493 499 self.setup(code, dataOut.data.shape)
494 500 self.isConfig = True
495 501
496 502 if mode == 0:
497 503 datadec = self.convolutionInTime(dataOut.data)
498 504
499 505 if mode == 1:
500 506 datadec = self.convolutionInFreq(dataOut.data)
501 507
502 508 if mode == 2:
503 509 datadec = self.convolutionInFreqOpt(dataOut.data)
504 510
505 511 dataOut.data = datadec
506 512
507 513 dataOut.heightList = dataOut.heightList[0:self.ndatadec]
508 514
509 515 dataOut.flagDecodeData = True #asumo q la data no esta decodificada
510 516
511 517 if self.__profIndex == self.nCode-1:
512 518 self.__profIndex = 0
513 519 return 1
514 520
515 521 self.__profIndex += 1
516 522
517 523 return 1
518 524 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
@@ -1,3 +1,4
1 1 from jroproc_voltage import *
2 2 from jroproc_spectra import *
3 3 from jroproc_heispectra import *
4 from jroproc_amisr import * No newline at end of file
@@ -1,79 +1,79
1 1 import os, sys
2 2
3 3 path = os.path.split(os.getcwd())[0]
4 4 sys.path.append(path)
5 5
6 6 from controller import *
7 7
8 8 desc = "AMISR Experiment"
9 9
10 10 filename = "amisr_reader.xml"
11 11
12 12 controllerObj = Project()
13 13
14 14 controllerObj.setup(id = '191', name='test01', description=desc)
15 15
16 16 path = os.path.join(os.environ['HOME'],'Documents/amisr')
17 17
18 18 figpath = os.path.join(os.environ['HOME'],'Pictures/amisr')
19 19
20 readUnitConfObj = controllerObj.addReadUnit(datatype='AMISR',
20 readUnitConfObj = controllerObj.addReadUnit(datatype='AMISRReader',
21 21 path=path,
22 22 startDate='2014/08/18',
23 23 endDate='2014/08/18',
24 24 startTime='00:00:00',
25 25 endTime='23:59:59',
26 26 walk=1)
27 27
28 28 #AMISR Processing Unit
29 procUnitAMISRBeam0 = controllerObj.addProcUnit(datatype='AMISR', inputId=readUnitConfObj.getId())
29 procUnitAMISRBeam0 = controllerObj.addProcUnit(datatype='AMISRProc', inputId=readUnitConfObj.getId())
30 30
31 31 #Beam Selector
32 32 opObj11 = procUnitAMISRBeam0.addOperation(name='BeamSelector', optype='other')
33 33 opObj11.addParameter(name='beam', value='0', format='int')
34 34
35 35 #Voltage Processing Unit
36 procUnitConfObjBeam0 = controllerObj.addProcUnit(datatype='Voltage', inputId=procUnitAMISRBeam0.getId())
36 procUnitConfObjBeam0 = controllerObj.addProcUnit(datatype='VoltageProc', inputId=procUnitAMISRBeam0.getId())
37 37 #Coherent Integration
38 38 opObj11 = procUnitConfObjBeam0.addOperation(name='CohInt', optype='other')
39 39 opObj11.addParameter(name='n', value='128', format='int')
40 40 #Spectra Unit Processing, getting spectras with nProfiles and nFFTPoints
41 procUnitConfObjSpectraBeam0 = controllerObj.addProcUnit(datatype='Spectra', inputId=procUnitConfObjBeam0.getId())
41 procUnitConfObjSpectraBeam0 = controllerObj.addProcUnit(datatype='SpectraProc', inputId=procUnitConfObjBeam0.getId())
42 42 procUnitConfObjSpectraBeam0.addParameter(name='nFFTPoints', value=32, format='int')
43 43 procUnitConfObjSpectraBeam0.addParameter(name='nProfiles', value=32, format='int')
44 44 #Noise Estimation
45 45 opObj11 = procUnitConfObjSpectraBeam0.addOperation(name='getNoise')
46 46 opObj11.addParameter(name='minHei', value='100', format='float')
47 47 opObj11.addParameter(name='maxHei', value='450', format='float')
48 48
49 49 #SpectraPlot
50 50 opObj11 = procUnitConfObjSpectraBeam0.addOperation(name='SpectraPlot', optype='other')
51 51 opObj11.addParameter(name='id', value='100', format='int')
52 52 opObj11.addParameter(name='wintitle', value='AMISR Beam 0', format='str')
53 53
54 54 #RTIPlot
55 55 title0 = 'RTI AMISR Beam 0'
56 56 opObj11 = procUnitConfObjSpectraBeam0.addOperation(name='RTIPlot', optype='other')
57 57 opObj11.addParameter(name='id', value='200', format='int')
58 58 opObj11.addParameter(name='wintitle', value=title0, format='str')
59 59 opObj11.addParameter(name='showprofile', value='0', format='int')
60 60 #Setting RTI time using xmin,xmax
61 61 opObj11.addParameter(name='xmin', value='0', format='int')
62 62 opObj11.addParameter(name='xmax', value='18', format='int')
63 63 #Setting dB range with zmin, zmax
64 64 opObj11.addParameter(name='zmin', value='45', format='int')
65 65 opObj11.addParameter(name='zmax', value='70', format='int')
66 66 #Save RTI
67 67 figfile0 = 'amisr_rti_beam0.png'
68 68 opObj11.addParameter(name='figpath', value=figpath, format='str')
69 opObj11.addParameter(name='figfile', value=figfile0, format='str')
69 #opObj11.addParameter(name='figfile', value=figfile0, format='str')
70 70
71 71
72 72 print "Escribiendo el archivo XML"
73 73 controllerObj.writeXml(filename)
74 74 print "Leyendo el archivo XML"
75 75 controllerObj.readXml(filename)
76 76
77 77 controllerObj.createObjects()
78 78 controllerObj.connectObjects()
79 79 controllerObj.run()
@@ -1,39 +1,39
1 1 import os, sys
2 2
3 3 path = os.path.split(os.getcwd())[0]
4 4 sys.path.append(path)
5 5
6 6 from controller import *
7 7
8 8 desc = "AMISR Experiment"
9 9
10 10 filename = "amisr_reader.xml"
11 11
12 12 controllerObj = Project()
13 13
14 14 controllerObj.setup(id = '191', name='test01', description=desc)
15 15
16 16 path = os.path.join(os.environ['HOME'],'Documents/amisr') #'/home/signalchain/Documents/amisr'
17 17
18 18 figpath = os.path.join(os.environ['HOME'],'Pictures/amisr')
19 19
20 readUnitConfObj = controllerObj.addReadUnit(datatype='AMISR',
20 readUnitConfObj = controllerObj.addReadUnit(datatype='AMISRReader',
21 21 path=path,
22 22 startDate='2014/08/18',
23 23 endDate='2014/08/18',
24 24 startTime='00:00:00',
25 25 endTime='23:59:59',
26 26 walk=1)
27 27
28 procUnitAMISR = controllerObj.addProcUnit(datatype='AMISR', inputId=readUnitConfObj.getId())
28 procUnitAMISR = controllerObj.addProcUnit(datatype='AMISRProc', inputId=readUnitConfObj.getId())
29 29
30 30 opObj11 = procUnitAMISR.addOperation(name='PrintInfo', optype='other')
31 31
32 32 print "Escribiendo el archivo XML"
33 33 controllerObj.writeXml(filename)
34 34 print "Leyendo el archivo XML"
35 35 controllerObj.readXml(filename)
36 36
37 37 controllerObj.createObjects()
38 38 controllerObj.connectObjects()
39 39 controllerObj.run()
General Comments 0
You need to be logged in to leave comments. Login now