##// END OF EJS Templates
Fix saving shifted pdata, add shift_fft parameter in SpectraProc to correct shifted pdata
jespinoza -
r1132:b323193ace39
parent child
Show More
@@ -1,679 +1,679
1 '''
1 '''
2 Created on Jul 2, 2014
2 Created on Jul 2, 2014
3
3
4 @author: roj-idl71
4 @author: roj-idl71
5 '''
5 '''
6 import numpy
6 import numpy
7
7
8 from jroIO_base import LOCALTIME, JRODataReader, JRODataWriter
8 from jroIO_base import LOCALTIME, JRODataReader, JRODataWriter
9 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
9 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
10 from schainpy.model.data.jroheaderIO import PROCFLAG, BasicHeader, SystemHeader, RadarControllerHeader, ProcessingHeader
10 from schainpy.model.data.jroheaderIO import PROCFLAG, BasicHeader, SystemHeader, RadarControllerHeader, ProcessingHeader
11 from schainpy.model.data.jrodata import Spectra
11 from schainpy.model.data.jrodata import Spectra
12
12
13 class SpectraReader(JRODataReader, ProcessingUnit):
13 class SpectraReader(JRODataReader, ProcessingUnit):
14 """
14 """
15 Esta clase permite leer datos de espectros desde archivos procesados (.pdata). La lectura
15 Esta clase permite leer datos de espectros desde archivos procesados (.pdata). La lectura
16 de los datos siempre se realiza por bloques. Los datos leidos (array de 3 dimensiones)
16 de los datos siempre se realiza por bloques. Los datos leidos (array de 3 dimensiones)
17 son almacenados en tres buffer's para el Self Spectra, el Cross Spectra y el DC Channel.
17 son almacenados en tres buffer's para el Self Spectra, el Cross Spectra y el DC Channel.
18
18
19 paresCanalesIguales * alturas * perfiles (Self Spectra)
19 paresCanalesIguales * alturas * perfiles (Self Spectra)
20 paresCanalesDiferentes * alturas * perfiles (Cross Spectra)
20 paresCanalesDiferentes * alturas * perfiles (Cross Spectra)
21 canales * alturas (DC Channels)
21 canales * alturas (DC Channels)
22
22
23 Esta clase contiene instancias (objetos) de las clases BasicHeader, SystemHeader,
23 Esta clase contiene instancias (objetos) de las clases BasicHeader, SystemHeader,
24 RadarControllerHeader y Spectra. Los tres primeros se usan para almacenar informacion de la
24 RadarControllerHeader y Spectra. Los tres primeros se usan para almacenar informacion de la
25 cabecera de datos (metadata), y el cuarto (Spectra) para obtener y almacenar un bloque de
25 cabecera de datos (metadata), y el cuarto (Spectra) para obtener y almacenar un bloque de
26 datos desde el "buffer" cada vez que se ejecute el metodo "getData".
26 datos desde el "buffer" cada vez que se ejecute el metodo "getData".
27
27
28 Example:
28 Example:
29 dpath = "/home/myuser/data"
29 dpath = "/home/myuser/data"
30
30
31 startTime = datetime.datetime(2010,1,20,0,0,0,0,0,0)
31 startTime = datetime.datetime(2010,1,20,0,0,0,0,0,0)
32
32
33 endTime = datetime.datetime(2010,1,21,23,59,59,0,0,0)
33 endTime = datetime.datetime(2010,1,21,23,59,59,0,0,0)
34
34
35 readerObj = SpectraReader()
35 readerObj = SpectraReader()
36
36
37 readerObj.setup(dpath, startTime, endTime)
37 readerObj.setup(dpath, startTime, endTime)
38
38
39 while(True):
39 while(True):
40
40
41 readerObj.getData()
41 readerObj.getData()
42
42
43 print readerObj.data_spc
43 print readerObj.data_spc
44
44
45 print readerObj.data_cspc
45 print readerObj.data_cspc
46
46
47 print readerObj.data_dc
47 print readerObj.data_dc
48
48
49 if readerObj.flagNoMoreFiles:
49 if readerObj.flagNoMoreFiles:
50 break
50 break
51
51
52 """
52 """
53
53
54 pts2read_SelfSpectra = 0
54 pts2read_SelfSpectra = 0
55
55
56 pts2read_CrossSpectra = 0
56 pts2read_CrossSpectra = 0
57
57
58 pts2read_DCchannels = 0
58 pts2read_DCchannels = 0
59
59
60 ext = ".pdata"
60 ext = ".pdata"
61
61
62 optchar = "P"
62 optchar = "P"
63
63
64 dataOut = None
64 dataOut = None
65
65
66 nRdChannels = None
66 nRdChannels = None
67
67
68 nRdPairs = None
68 nRdPairs = None
69
69
70 rdPairList = []
70 rdPairList = []
71
71
72 def __init__(self, **kwargs):
72 def __init__(self, **kwargs):
73 """
73 """
74 Inicializador de la clase SpectraReader para la lectura de datos de espectros.
74 Inicializador de la clase SpectraReader para la lectura de datos de espectros.
75
75
76 Inputs:
76 Inputs:
77 dataOut : Objeto de la clase Spectra. Este objeto sera utilizado para
77 dataOut : Objeto de la clase Spectra. Este objeto sera utilizado para
78 almacenar un perfil de datos cada vez que se haga un requerimiento
78 almacenar un perfil de datos cada vez que se haga un requerimiento
79 (getData). El perfil sera obtenido a partir del buffer de datos,
79 (getData). El perfil sera obtenido a partir del buffer de datos,
80 si el buffer esta vacio se hara un nuevo proceso de lectura de un
80 si el buffer esta vacio se hara un nuevo proceso de lectura de un
81 bloque de datos.
81 bloque de datos.
82 Si este parametro no es pasado se creara uno internamente.
82 Si este parametro no es pasado se creara uno internamente.
83
83
84 Affected:
84 Affected:
85 self.dataOut
85 self.dataOut
86
86
87 Return : None
87 Return : None
88 """
88 """
89
89
90 #Eliminar de la base la herencia
90 #Eliminar de la base la herencia
91 ProcessingUnit.__init__(self, **kwargs)
91 ProcessingUnit.__init__(self, **kwargs)
92
92
93 # self.isConfig = False
93 # self.isConfig = False
94
94
95 self.pts2read_SelfSpectra = 0
95 self.pts2read_SelfSpectra = 0
96
96
97 self.pts2read_CrossSpectra = 0
97 self.pts2read_CrossSpectra = 0
98
98
99 self.pts2read_DCchannels = 0
99 self.pts2read_DCchannels = 0
100
100
101 self.datablock = None
101 self.datablock = None
102
102
103 self.utc = None
103 self.utc = None
104
104
105 self.ext = ".pdata"
105 self.ext = ".pdata"
106
106
107 self.optchar = "P"
107 self.optchar = "P"
108
108
109 self.basicHeaderObj = BasicHeader(LOCALTIME)
109 self.basicHeaderObj = BasicHeader(LOCALTIME)
110
110
111 self.systemHeaderObj = SystemHeader()
111 self.systemHeaderObj = SystemHeader()
112
112
113 self.radarControllerHeaderObj = RadarControllerHeader()
113 self.radarControllerHeaderObj = RadarControllerHeader()
114
114
115 self.processingHeaderObj = ProcessingHeader()
115 self.processingHeaderObj = ProcessingHeader()
116
116
117 self.online = 0
117 self.online = 0
118
118
119 self.fp = None
119 self.fp = None
120
120
121 self.idFile = None
121 self.idFile = None
122
122
123 self.dtype = None
123 self.dtype = None
124
124
125 self.fileSizeByHeader = None
125 self.fileSizeByHeader = None
126
126
127 self.filenameList = []
127 self.filenameList = []
128
128
129 self.filename = None
129 self.filename = None
130
130
131 self.fileSize = None
131 self.fileSize = None
132
132
133 self.firstHeaderSize = 0
133 self.firstHeaderSize = 0
134
134
135 self.basicHeaderSize = 24
135 self.basicHeaderSize = 24
136
136
137 self.pathList = []
137 self.pathList = []
138
138
139 self.lastUTTime = 0
139 self.lastUTTime = 0
140
140
141 self.maxTimeStep = 30
141 self.maxTimeStep = 30
142
142
143 self.flagNoMoreFiles = 0
143 self.flagNoMoreFiles = 0
144
144
145 self.set = 0
145 self.set = 0
146
146
147 self.path = None
147 self.path = None
148
148
149 self.delay = 60 #seconds
149 self.delay = 60 #seconds
150
150
151 self.nTries = 3 #quantity tries
151 self.nTries = 3 #quantity tries
152
152
153 self.nFiles = 3 #number of files for searching
153 self.nFiles = 3 #number of files for searching
154
154
155 self.nReadBlocks = 0
155 self.nReadBlocks = 0
156
156
157 self.flagIsNewFile = 1
157 self.flagIsNewFile = 1
158
158
159 self.__isFirstTimeOnline = 1
159 self.__isFirstTimeOnline = 1
160
160
161 # self.ippSeconds = 0
161 # self.ippSeconds = 0
162
162
163 self.flagDiscontinuousBlock = 0
163 self.flagDiscontinuousBlock = 0
164
164
165 self.flagIsNewBlock = 0
165 self.flagIsNewBlock = 0
166
166
167 self.nTotalBlocks = 0
167 self.nTotalBlocks = 0
168
168
169 self.blocksize = 0
169 self.blocksize = 0
170
170
171 self.dataOut = self.createObjByDefault()
171 self.dataOut = self.createObjByDefault()
172
172
173 self.profileIndex = 1 #Always
173 self.profileIndex = 1 #Always
174
174
175
175
176 def createObjByDefault(self):
176 def createObjByDefault(self):
177
177
178 dataObj = Spectra()
178 dataObj = Spectra()
179
179
180 return dataObj
180 return dataObj
181
181
182 def __hasNotDataInBuffer(self):
182 def __hasNotDataInBuffer(self):
183 return 1
183 return 1
184
184
185
185
186 def getBlockDimension(self):
186 def getBlockDimension(self):
187 """
187 """
188 Obtiene la cantidad de puntos a leer por cada bloque de datos
188 Obtiene la cantidad de puntos a leer por cada bloque de datos
189
189
190 Affected:
190 Affected:
191 self.nRdChannels
191 self.nRdChannels
192 self.nRdPairs
192 self.nRdPairs
193 self.pts2read_SelfSpectra
193 self.pts2read_SelfSpectra
194 self.pts2read_CrossSpectra
194 self.pts2read_CrossSpectra
195 self.pts2read_DCchannels
195 self.pts2read_DCchannels
196 self.blocksize
196 self.blocksize
197 self.dataOut.nChannels
197 self.dataOut.nChannels
198 self.dataOut.nPairs
198 self.dataOut.nPairs
199
199
200 Return:
200 Return:
201 None
201 None
202 """
202 """
203 self.nRdChannels = 0
203 self.nRdChannels = 0
204 self.nRdPairs = 0
204 self.nRdPairs = 0
205 self.rdPairList = []
205 self.rdPairList = []
206
206
207 for i in range(0, self.processingHeaderObj.totalSpectra*2, 2):
207 for i in range(0, self.processingHeaderObj.totalSpectra*2, 2):
208 if self.processingHeaderObj.spectraComb[i] == self.processingHeaderObj.spectraComb[i+1]:
208 if self.processingHeaderObj.spectraComb[i] == self.processingHeaderObj.spectraComb[i+1]:
209 self.nRdChannels = self.nRdChannels + 1 #par de canales iguales
209 self.nRdChannels = self.nRdChannels + 1 #par de canales iguales
210 else:
210 else:
211 self.nRdPairs = self.nRdPairs + 1 #par de canales diferentes
211 self.nRdPairs = self.nRdPairs + 1 #par de canales diferentes
212 self.rdPairList.append((self.processingHeaderObj.spectraComb[i], self.processingHeaderObj.spectraComb[i+1]))
212 self.rdPairList.append((self.processingHeaderObj.spectraComb[i], self.processingHeaderObj.spectraComb[i+1]))
213
213
214 pts2read = self.processingHeaderObj.nHeights * self.processingHeaderObj.profilesPerBlock
214 pts2read = self.processingHeaderObj.nHeights * self.processingHeaderObj.profilesPerBlock
215
215
216 self.pts2read_SelfSpectra = int(self.nRdChannels * pts2read)
216 self.pts2read_SelfSpectra = int(self.nRdChannels * pts2read)
217 self.blocksize = self.pts2read_SelfSpectra
217 self.blocksize = self.pts2read_SelfSpectra
218
218
219 if self.processingHeaderObj.flag_cspc:
219 if self.processingHeaderObj.flag_cspc:
220 self.pts2read_CrossSpectra = int(self.nRdPairs * pts2read)
220 self.pts2read_CrossSpectra = int(self.nRdPairs * pts2read)
221 self.blocksize += self.pts2read_CrossSpectra
221 self.blocksize += self.pts2read_CrossSpectra
222
222
223 if self.processingHeaderObj.flag_dc:
223 if self.processingHeaderObj.flag_dc:
224 self.pts2read_DCchannels = int(self.systemHeaderObj.nChannels * self.processingHeaderObj.nHeights)
224 self.pts2read_DCchannels = int(self.systemHeaderObj.nChannels * self.processingHeaderObj.nHeights)
225 self.blocksize += self.pts2read_DCchannels
225 self.blocksize += self.pts2read_DCchannels
226
226
227 # self.blocksize = self.pts2read_SelfSpectra + self.pts2read_CrossSpectra + self.pts2read_DCchannels
227 # self.blocksize = self.pts2read_SelfSpectra + self.pts2read_CrossSpectra + self.pts2read_DCchannels
228
228
229
229
230 def readBlock(self):
230 def readBlock(self):
231 """
231 """
232 Lee el bloque de datos desde la posicion actual del puntero del archivo
232 Lee el bloque de datos desde la posicion actual del puntero del archivo
233 (self.fp) y actualiza todos los parametros relacionados al bloque de datos
233 (self.fp) y actualiza todos los parametros relacionados al bloque de datos
234 (metadata + data). La data leida es almacenada en el buffer y el contador del buffer
234 (metadata + data). La data leida es almacenada en el buffer y el contador del buffer
235 es seteado a 0
235 es seteado a 0
236
236
237 Return: None
237 Return: None
238
238
239 Variables afectadas:
239 Variables afectadas:
240
240
241 self.flagIsNewFile
241 self.flagIsNewFile
242 self.flagIsNewBlock
242 self.flagIsNewBlock
243 self.nTotalBlocks
243 self.nTotalBlocks
244 self.data_spc
244 self.data_spc
245 self.data_cspc
245 self.data_cspc
246 self.data_dc
246 self.data_dc
247
247
248 Exceptions:
248 Exceptions:
249 Si un bloque leido no es un bloque valido
249 Si un bloque leido no es un bloque valido
250 """
250 """
251 blockOk_flag = False
251 blockOk_flag = False
252 fpointer = self.fp.tell()
252 fpointer = self.fp.tell()
253
253
254 spc = numpy.fromfile( self.fp, self.dtype[0], self.pts2read_SelfSpectra )
254 spc = numpy.fromfile( self.fp, self.dtype[0], self.pts2read_SelfSpectra )
255 spc = spc.reshape( (self.nRdChannels, self.processingHeaderObj.nHeights, self.processingHeaderObj.profilesPerBlock) ) #transforma a un arreglo 3D
255 spc = spc.reshape( (self.nRdChannels, self.processingHeaderObj.nHeights, self.processingHeaderObj.profilesPerBlock) ) #transforma a un arreglo 3D
256
256
257 if self.processingHeaderObj.flag_cspc:
257 if self.processingHeaderObj.flag_cspc:
258 cspc = numpy.fromfile( self.fp, self.dtype, self.pts2read_CrossSpectra )
258 cspc = numpy.fromfile( self.fp, self.dtype, self.pts2read_CrossSpectra )
259 cspc = cspc.reshape( (self.nRdPairs, self.processingHeaderObj.nHeights, self.processingHeaderObj.profilesPerBlock) ) #transforma a un arreglo 3D
259 cspc = cspc.reshape( (self.nRdPairs, self.processingHeaderObj.nHeights, self.processingHeaderObj.profilesPerBlock) ) #transforma a un arreglo 3D
260
260
261 if self.processingHeaderObj.flag_dc:
261 if self.processingHeaderObj.flag_dc:
262 dc = numpy.fromfile( self.fp, self.dtype, self.pts2read_DCchannels ) #int(self.processingHeaderObj.nHeights*self.systemHeaderObj.nChannels) )
262 dc = numpy.fromfile( self.fp, self.dtype, self.pts2read_DCchannels ) #int(self.processingHeaderObj.nHeights*self.systemHeaderObj.nChannels) )
263 dc = dc.reshape( (self.systemHeaderObj.nChannels, self.processingHeaderObj.nHeights) ) #transforma a un arreglo 2D
263 dc = dc.reshape( (self.systemHeaderObj.nChannels, self.processingHeaderObj.nHeights) ) #transforma a un arreglo 2D
264
264
265
265
266 if not(self.processingHeaderObj.shif_fft):
266 if self.processingHeaderObj.shif_fft:
267 #desplaza a la derecha en el eje 2 determinadas posiciones
267 #desplaza a la derecha en el eje 2 determinadas posiciones
268 shift = int(self.processingHeaderObj.profilesPerBlock/2)
268 shift = int(self.processingHeaderObj.profilesPerBlock/2)
269 spc = numpy.roll( spc, shift , axis=2 )
269 spc = numpy.roll( spc, shift , axis=2 )
270
270
271 if self.processingHeaderObj.flag_cspc:
271 if self.processingHeaderObj.flag_cspc:
272 #desplaza a la derecha en el eje 2 determinadas posiciones
272 #desplaza a la derecha en el eje 2 determinadas posiciones
273 cspc = numpy.roll( cspc, shift, axis=2 )
273 cspc = numpy.roll( cspc, shift, axis=2 )
274
274
275 #Dimensions : nChannels, nProfiles, nSamples
275 #Dimensions : nChannels, nProfiles, nSamples
276 spc = numpy.transpose( spc, (0,2,1) )
276 spc = numpy.transpose( spc, (0,2,1) )
277 self.data_spc = spc
277 self.data_spc = spc
278
278
279 if self.processingHeaderObj.flag_cspc:
279 if self.processingHeaderObj.flag_cspc:
280 cspc = numpy.transpose( cspc, (0,2,1) )
280 cspc = numpy.transpose( cspc, (0,2,1) )
281 self.data_cspc = cspc['real'] + cspc['imag']*1j
281 self.data_cspc = cspc['real'] + cspc['imag']*1j
282 else:
282 else:
283 self.data_cspc = None
283 self.data_cspc = None
284
284
285 if self.processingHeaderObj.flag_dc:
285 if self.processingHeaderObj.flag_dc:
286 self.data_dc = dc['real'] + dc['imag']*1j
286 self.data_dc = dc['real'] + dc['imag']*1j
287 else:
287 else:
288 self.data_dc = None
288 self.data_dc = None
289
289
290 self.flagIsNewFile = 0
290 self.flagIsNewFile = 0
291 self.flagIsNewBlock = 1
291 self.flagIsNewBlock = 1
292
292
293 self.nTotalBlocks += 1
293 self.nTotalBlocks += 1
294 self.nReadBlocks += 1
294 self.nReadBlocks += 1
295
295
296 return 1
296 return 1
297
297
298 def getFirstHeader(self):
298 def getFirstHeader(self):
299
299
300 self.getBasicHeader()
300 self.getBasicHeader()
301
301
302 self.dataOut.systemHeaderObj = self.systemHeaderObj.copy()
302 self.dataOut.systemHeaderObj = self.systemHeaderObj.copy()
303
303
304 self.dataOut.radarControllerHeaderObj = self.radarControllerHeaderObj.copy()
304 self.dataOut.radarControllerHeaderObj = self.radarControllerHeaderObj.copy()
305
305
306 # self.dataOut.ippSeconds = self.ippSeconds
306 # self.dataOut.ippSeconds = self.ippSeconds
307
307
308 # self.dataOut.timeInterval = self.radarControllerHeaderObj.ippSeconds * self.processingHeaderObj.nCohInt * self.processingHeaderObj.nIncohInt * self.processingHeaderObj.profilesPerBlock
308 # self.dataOut.timeInterval = self.radarControllerHeaderObj.ippSeconds * self.processingHeaderObj.nCohInt * self.processingHeaderObj.nIncohInt * self.processingHeaderObj.profilesPerBlock
309
309
310 self.dataOut.dtype = self.dtype
310 self.dataOut.dtype = self.dtype
311
311
312 # self.dataOut.nPairs = self.nPairs
312 # self.dataOut.nPairs = self.nPairs
313
313
314 self.dataOut.pairsList = self.rdPairList
314 self.dataOut.pairsList = self.rdPairList
315
315
316 self.dataOut.nProfiles = self.processingHeaderObj.profilesPerBlock
316 self.dataOut.nProfiles = self.processingHeaderObj.profilesPerBlock
317
317
318 self.dataOut.nFFTPoints = self.processingHeaderObj.profilesPerBlock
318 self.dataOut.nFFTPoints = self.processingHeaderObj.profilesPerBlock
319
319
320 self.dataOut.nCohInt = self.processingHeaderObj.nCohInt
320 self.dataOut.nCohInt = self.processingHeaderObj.nCohInt
321
321
322 self.dataOut.nIncohInt = self.processingHeaderObj.nIncohInt
322 self.dataOut.nIncohInt = self.processingHeaderObj.nIncohInt
323
323
324 xf = self.processingHeaderObj.firstHeight + self.processingHeaderObj.nHeights*self.processingHeaderObj.deltaHeight
324 xf = self.processingHeaderObj.firstHeight + self.processingHeaderObj.nHeights*self.processingHeaderObj.deltaHeight
325
325
326 self.dataOut.heightList = numpy.arange(self.processingHeaderObj.firstHeight, xf, self.processingHeaderObj.deltaHeight)
326 self.dataOut.heightList = numpy.arange(self.processingHeaderObj.firstHeight, xf, self.processingHeaderObj.deltaHeight)
327
327
328 self.dataOut.channelList = range(self.systemHeaderObj.nChannels)
328 self.dataOut.channelList = range(self.systemHeaderObj.nChannels)
329
329
330 self.dataOut.flagShiftFFT = True #Data is always shifted
330 self.dataOut.flagShiftFFT = True #Data is always shifted
331
331
332 self.dataOut.flagDecodeData = self.processingHeaderObj.flag_decode #asumo q la data no esta decodificada
332 self.dataOut.flagDecodeData = self.processingHeaderObj.flag_decode #asumo q la data no esta decodificada
333
333
334 self.dataOut.flagDeflipData = self.processingHeaderObj.flag_deflip #asumo q la data esta sin flip
334 self.dataOut.flagDeflipData = self.processingHeaderObj.flag_deflip #asumo q la data esta sin flip
335
335
336 def getData(self):
336 def getData(self):
337 """
337 """
338 First method to execute before "RUN" is called.
338 First method to execute before "RUN" is called.
339
339
340 Copia el buffer de lectura a la clase "Spectra",
340 Copia el buffer de lectura a la clase "Spectra",
341 con todos los parametros asociados a este (metadata). cuando no hay datos en el buffer de
341 con todos los parametros asociados a este (metadata). cuando no hay datos en el buffer de
342 lectura es necesario hacer una nueva lectura de los bloques de datos usando "readNextBlock"
342 lectura es necesario hacer una nueva lectura de los bloques de datos usando "readNextBlock"
343
343
344 Return:
344 Return:
345 0 : Si no hay mas archivos disponibles
345 0 : Si no hay mas archivos disponibles
346 1 : Si hizo una buena copia del buffer
346 1 : Si hizo una buena copia del buffer
347
347
348 Affected:
348 Affected:
349 self.dataOut
349 self.dataOut
350
350
351 self.flagDiscontinuousBlock
351 self.flagDiscontinuousBlock
352 self.flagIsNewBlock
352 self.flagIsNewBlock
353 """
353 """
354
354
355 if self.flagNoMoreFiles:
355 if self.flagNoMoreFiles:
356 self.dataOut.flagNoData = True
356 self.dataOut.flagNoData = True
357 print 'Process finished'
357 print 'Process finished'
358 return 0
358 return 0
359
359
360 self.flagDiscontinuousBlock = 0
360 self.flagDiscontinuousBlock = 0
361 self.flagIsNewBlock = 0
361 self.flagIsNewBlock = 0
362
362
363 if self.__hasNotDataInBuffer():
363 if self.__hasNotDataInBuffer():
364
364
365 if not( self.readNextBlock() ):
365 if not( self.readNextBlock() ):
366 self.dataOut.flagNoData = True
366 self.dataOut.flagNoData = True
367 return 0
367 return 0
368
368
369 #data es un numpy array de 3 dmensiones (perfiles, alturas y canales)
369 #data es un numpy array de 3 dmensiones (perfiles, alturas y canales)
370
370
371 if self.data_spc is None:
371 if self.data_spc is None:
372 self.dataOut.flagNoData = True
372 self.dataOut.flagNoData = True
373 return 0
373 return 0
374
374
375 self.getBasicHeader()
375 self.getBasicHeader()
376
376
377 self.getFirstHeader()
377 self.getFirstHeader()
378
378
379 self.dataOut.data_spc = self.data_spc
379 self.dataOut.data_spc = self.data_spc
380
380
381 self.dataOut.data_cspc = self.data_cspc
381 self.dataOut.data_cspc = self.data_cspc
382
382
383 self.dataOut.data_dc = self.data_dc
383 self.dataOut.data_dc = self.data_dc
384
384
385 self.dataOut.flagNoData = False
385 self.dataOut.flagNoData = False
386
386
387 self.dataOut.realtime = self.online
387 self.dataOut.realtime = self.online
388
388
389 return self.dataOut.data_spc
389 return self.dataOut.data_spc
390
390
391 class SpectraWriter(JRODataWriter, Operation):
391 class SpectraWriter(JRODataWriter, Operation):
392
392
393 """
393 """
394 Esta clase permite escribir datos de espectros a archivos procesados (.pdata). La escritura
394 Esta clase permite escribir datos de espectros a archivos procesados (.pdata). La escritura
395 de los datos siempre se realiza por bloques.
395 de los datos siempre se realiza por bloques.
396 """
396 """
397
397
398 ext = ".pdata"
398 ext = ".pdata"
399
399
400 optchar = "P"
400 optchar = "P"
401
401
402 shape_spc_Buffer = None
402 shape_spc_Buffer = None
403
403
404 shape_cspc_Buffer = None
404 shape_cspc_Buffer = None
405
405
406 shape_dc_Buffer = None
406 shape_dc_Buffer = None
407
407
408 data_spc = None
408 data_spc = None
409
409
410 data_cspc = None
410 data_cspc = None
411
411
412 data_dc = None
412 data_dc = None
413
413
414 # dataOut = None
414 # dataOut = None
415
415
416 def __init__(self, **kwargs):
416 def __init__(self, **kwargs):
417 """
417 """
418 Inicializador de la clase SpectraWriter para la escritura de datos de espectros.
418 Inicializador de la clase SpectraWriter para la escritura de datos de espectros.
419
419
420 Affected:
420 Affected:
421 self.dataOut
421 self.dataOut
422 self.basicHeaderObj
422 self.basicHeaderObj
423 self.systemHeaderObj
423 self.systemHeaderObj
424 self.radarControllerHeaderObj
424 self.radarControllerHeaderObj
425 self.processingHeaderObj
425 self.processingHeaderObj
426
426
427 Return: None
427 Return: None
428 """
428 """
429
429
430 Operation.__init__(self, **kwargs)
430 Operation.__init__(self, **kwargs)
431
431
432 self.isConfig = False
432 self.isConfig = False
433
433
434 self.nTotalBlocks = 0
434 self.nTotalBlocks = 0
435
435
436 self.data_spc = None
436 self.data_spc = None
437
437
438 self.data_cspc = None
438 self.data_cspc = None
439
439
440 self.data_dc = None
440 self.data_dc = None
441
441
442 self.fp = None
442 self.fp = None
443
443
444 self.flagIsNewFile = 1
444 self.flagIsNewFile = 1
445
445
446 self.nTotalBlocks = 0
446 self.nTotalBlocks = 0
447
447
448 self.flagIsNewBlock = 0
448 self.flagIsNewBlock = 0
449
449
450 self.setFile = None
450 self.setFile = None
451
451
452 self.dtype = None
452 self.dtype = None
453
453
454 self.path = None
454 self.path = None
455
455
456 self.noMoreFiles = 0
456 self.noMoreFiles = 0
457
457
458 self.filename = None
458 self.filename = None
459
459
460 self.basicHeaderObj = BasicHeader(LOCALTIME)
460 self.basicHeaderObj = BasicHeader(LOCALTIME)
461
461
462 self.systemHeaderObj = SystemHeader()
462 self.systemHeaderObj = SystemHeader()
463
463
464 self.radarControllerHeaderObj = RadarControllerHeader()
464 self.radarControllerHeaderObj = RadarControllerHeader()
465
465
466 self.processingHeaderObj = ProcessingHeader()
466 self.processingHeaderObj = ProcessingHeader()
467
467
468
468
469 def hasAllDataInBuffer(self):
469 def hasAllDataInBuffer(self):
470 return 1
470 return 1
471
471
472
472
473 def setBlockDimension(self):
473 def setBlockDimension(self):
474 """
474 """
475 Obtiene las formas dimensionales del los subbloques de datos que componen un bloque
475 Obtiene las formas dimensionales del los subbloques de datos que componen un bloque
476
476
477 Affected:
477 Affected:
478 self.shape_spc_Buffer
478 self.shape_spc_Buffer
479 self.shape_cspc_Buffer
479 self.shape_cspc_Buffer
480 self.shape_dc_Buffer
480 self.shape_dc_Buffer
481
481
482 Return: None
482 Return: None
483 """
483 """
484 self.shape_spc_Buffer = (self.dataOut.nChannels,
484 self.shape_spc_Buffer = (self.dataOut.nChannels,
485 self.processingHeaderObj.nHeights,
485 self.processingHeaderObj.nHeights,
486 self.processingHeaderObj.profilesPerBlock)
486 self.processingHeaderObj.profilesPerBlock)
487
487
488 self.shape_cspc_Buffer = (self.dataOut.nPairs,
488 self.shape_cspc_Buffer = (self.dataOut.nPairs,
489 self.processingHeaderObj.nHeights,
489 self.processingHeaderObj.nHeights,
490 self.processingHeaderObj.profilesPerBlock)
490 self.processingHeaderObj.profilesPerBlock)
491
491
492 self.shape_dc_Buffer = (self.dataOut.nChannels,
492 self.shape_dc_Buffer = (self.dataOut.nChannels,
493 self.processingHeaderObj.nHeights)
493 self.processingHeaderObj.nHeights)
494
494
495
495
496 def writeBlock(self):
496 def writeBlock(self):
497 """
497 """
498 Escribe el buffer en el file designado
498 Escribe el buffer en el file designado
499
499
500 Affected:
500 Affected:
501 self.data_spc
501 self.data_spc
502 self.data_cspc
502 self.data_cspc
503 self.data_dc
503 self.data_dc
504 self.flagIsNewFile
504 self.flagIsNewFile
505 self.flagIsNewBlock
505 self.flagIsNewBlock
506 self.nTotalBlocks
506 self.nTotalBlocks
507 self.nWriteBlocks
507 self.nWriteBlocks
508
508
509 Return: None
509 Return: None
510 """
510 """
511
511
512 spc = numpy.transpose( self.data_spc, (0,2,1) )
512 spc = numpy.transpose( self.data_spc, (0,2,1) )
513 if not( self.processingHeaderObj.shif_fft ):
513 if self.processingHeaderObj.shif_fft:
514 spc = numpy.roll( spc, self.processingHeaderObj.profilesPerBlock/2, axis=2 ) #desplaza a la derecha en el eje 2 determinadas posiciones
514 spc = numpy.roll( spc, self.processingHeaderObj.profilesPerBlock/2, axis=2 ) #desplaza a la derecha en el eje 2 determinadas posiciones
515 data = spc.reshape((-1))
515 data = spc.reshape((-1))
516 data = data.astype(self.dtype[0])
516 data = data.astype(self.dtype[0])
517 data.tofile(self.fp)
517 data.tofile(self.fp)
518
518
519 if self.data_cspc is not None:
519 if self.data_cspc is not None:
520 data = numpy.zeros( self.shape_cspc_Buffer, self.dtype )
520 data = numpy.zeros( self.shape_cspc_Buffer, self.dtype )
521 cspc = numpy.transpose( self.data_cspc, (0,2,1) )
521 cspc = numpy.transpose( self.data_cspc, (0,2,1) )
522 if not( self.processingHeaderObj.shif_fft ):
522 if self.processingHeaderObj.shif_fft:
523 cspc = numpy.roll( cspc, self.processingHeaderObj.profilesPerBlock/2, axis=2 ) #desplaza a la derecha en el eje 2 determinadas posiciones
523 cspc = numpy.roll( cspc, self.processingHeaderObj.profilesPerBlock/2, axis=2 ) #desplaza a la derecha en el eje 2 determinadas posiciones
524 data['real'] = cspc.real
524 data['real'] = cspc.real
525 data['imag'] = cspc.imag
525 data['imag'] = cspc.imag
526 data = data.reshape((-1))
526 data = data.reshape((-1))
527 data.tofile(self.fp)
527 data.tofile(self.fp)
528
528
529 if self.data_dc is not None:
529 if self.data_dc is not None:
530 data = numpy.zeros( self.shape_dc_Buffer, self.dtype )
530 data = numpy.zeros( self.shape_dc_Buffer, self.dtype )
531 dc = self.data_dc
531 dc = self.data_dc
532 data['real'] = dc.real
532 data['real'] = dc.real
533 data['imag'] = dc.imag
533 data['imag'] = dc.imag
534 data = data.reshape((-1))
534 data = data.reshape((-1))
535 data.tofile(self.fp)
535 data.tofile(self.fp)
536
536
537 # self.data_spc.fill(0)
537 # self.data_spc.fill(0)
538 #
538 #
539 # if self.data_dc is not None:
539 # if self.data_dc is not None:
540 # self.data_dc.fill(0)
540 # self.data_dc.fill(0)
541 #
541 #
542 # if self.data_cspc is not None:
542 # if self.data_cspc is not None:
543 # self.data_cspc.fill(0)
543 # self.data_cspc.fill(0)
544
544
545 self.flagIsNewFile = 0
545 self.flagIsNewFile = 0
546 self.flagIsNewBlock = 1
546 self.flagIsNewBlock = 1
547 self.nTotalBlocks += 1
547 self.nTotalBlocks += 1
548 self.nWriteBlocks += 1
548 self.nWriteBlocks += 1
549 self.blockIndex += 1
549 self.blockIndex += 1
550
550
551 # print "[Writing] Block = %d04" %self.blockIndex
551 # print "[Writing] Block = %d04" %self.blockIndex
552
552
553 def putData(self):
553 def putData(self):
554 """
554 """
555 Setea un bloque de datos y luego los escribe en un file
555 Setea un bloque de datos y luego los escribe en un file
556
556
557 Affected:
557 Affected:
558 self.data_spc
558 self.data_spc
559 self.data_cspc
559 self.data_cspc
560 self.data_dc
560 self.data_dc
561
561
562 Return:
562 Return:
563 0 : Si no hay data o no hay mas files que puedan escribirse
563 0 : Si no hay data o no hay mas files que puedan escribirse
564 1 : Si se escribio la data de un bloque en un file
564 1 : Si se escribio la data de un bloque en un file
565 """
565 """
566
566
567 if self.dataOut.flagNoData:
567 if self.dataOut.flagNoData:
568 return 0
568 return 0
569
569
570 self.flagIsNewBlock = 0
570 self.flagIsNewBlock = 0
571
571
572 if self.dataOut.flagDiscontinuousBlock:
572 if self.dataOut.flagDiscontinuousBlock:
573 self.data_spc.fill(0)
573 self.data_spc.fill(0)
574 if self.dataOut.data_cspc is not None:
574 if self.dataOut.data_cspc is not None:
575 self.data_cspc.fill(0)
575 self.data_cspc.fill(0)
576 if self.dataOut.data_dc is not None:
576 if self.dataOut.data_dc is not None:
577 self.data_dc.fill(0)
577 self.data_dc.fill(0)
578 self.setNextFile()
578 self.setNextFile()
579
579
580 if self.flagIsNewFile == 0:
580 if self.flagIsNewFile == 0:
581 self.setBasicHeader()
581 self.setBasicHeader()
582
582
583 self.data_spc = self.dataOut.data_spc.copy()
583 self.data_spc = self.dataOut.data_spc.copy()
584
584
585 if self.dataOut.data_cspc is not None:
585 if self.dataOut.data_cspc is not None:
586 self.data_cspc = self.dataOut.data_cspc.copy()
586 self.data_cspc = self.dataOut.data_cspc.copy()
587
587
588 if self.dataOut.data_dc is not None:
588 if self.dataOut.data_dc is not None:
589 self.data_dc = self.dataOut.data_dc.copy()
589 self.data_dc = self.dataOut.data_dc.copy()
590
590
591 # #self.processingHeaderObj.dataBlocksPerFile)
591 # #self.processingHeaderObj.dataBlocksPerFile)
592 if self.hasAllDataInBuffer():
592 if self.hasAllDataInBuffer():
593 # self.setFirstHeader()
593 # self.setFirstHeader()
594 self.writeNextBlock()
594 self.writeNextBlock()
595
595
596 return 1
596 return 1
597
597
598 def __getBlockSize(self):
598 def __getBlockSize(self):
599 '''
599 '''
600 Este metodos determina el cantidad de bytes para un bloque de datos de tipo Spectra
600 Este metodos determina el cantidad de bytes para un bloque de datos de tipo Spectra
601 '''
601 '''
602
602
603 dtype_width = self.getDtypeWidth()
603 dtype_width = self.getDtypeWidth()
604
604
605 pts2write = self.dataOut.nHeights * self.dataOut.nFFTPoints
605 pts2write = self.dataOut.nHeights * self.dataOut.nFFTPoints
606
606
607 pts2write_SelfSpectra = int(self.dataOut.nChannels * pts2write)
607 pts2write_SelfSpectra = int(self.dataOut.nChannels * pts2write)
608 blocksize = (pts2write_SelfSpectra*dtype_width)
608 blocksize = (pts2write_SelfSpectra*dtype_width)
609
609
610 if self.dataOut.data_cspc is not None:
610 if self.dataOut.data_cspc is not None:
611 pts2write_CrossSpectra = int(self.dataOut.nPairs * pts2write)
611 pts2write_CrossSpectra = int(self.dataOut.nPairs * pts2write)
612 blocksize += (pts2write_CrossSpectra*dtype_width*2)
612 blocksize += (pts2write_CrossSpectra*dtype_width*2)
613
613
614 if self.dataOut.data_dc is not None:
614 if self.dataOut.data_dc is not None:
615 pts2write_DCchannels = int(self.dataOut.nChannels * self.dataOut.nHeights)
615 pts2write_DCchannels = int(self.dataOut.nChannels * self.dataOut.nHeights)
616 blocksize += (pts2write_DCchannels*dtype_width*2)
616 blocksize += (pts2write_DCchannels*dtype_width*2)
617
617
618 # blocksize = blocksize #* datatypeValue * 2 #CORREGIR ESTO
618 # blocksize = blocksize #* datatypeValue * 2 #CORREGIR ESTO
619
619
620 return blocksize
620 return blocksize
621
621
622 def setFirstHeader(self):
622 def setFirstHeader(self):
623
623
624 """
624 """
625 Obtiene una copia del First Header
625 Obtiene una copia del First Header
626
626
627 Affected:
627 Affected:
628 self.systemHeaderObj
628 self.systemHeaderObj
629 self.radarControllerHeaderObj
629 self.radarControllerHeaderObj
630 self.dtype
630 self.dtype
631
631
632 Return:
632 Return:
633 None
633 None
634 """
634 """
635
635
636 self.systemHeaderObj = self.dataOut.systemHeaderObj.copy()
636 self.systemHeaderObj = self.dataOut.systemHeaderObj.copy()
637 self.systemHeaderObj.nChannels = self.dataOut.nChannels
637 self.systemHeaderObj.nChannels = self.dataOut.nChannels
638 self.radarControllerHeaderObj = self.dataOut.radarControllerHeaderObj.copy()
638 self.radarControllerHeaderObj = self.dataOut.radarControllerHeaderObj.copy()
639
639
640 self.processingHeaderObj.dtype = 1 # Spectra
640 self.processingHeaderObj.dtype = 1 # Spectra
641 self.processingHeaderObj.blockSize = self.__getBlockSize()
641 self.processingHeaderObj.blockSize = self.__getBlockSize()
642 self.processingHeaderObj.profilesPerBlock = self.dataOut.nFFTPoints
642 self.processingHeaderObj.profilesPerBlock = self.dataOut.nFFTPoints
643 self.processingHeaderObj.dataBlocksPerFile = self.blocksPerFile
643 self.processingHeaderObj.dataBlocksPerFile = self.blocksPerFile
644 self.processingHeaderObj.nWindows = 1 #podria ser 1 o self.dataOut.processingHeaderObj.nWindows
644 self.processingHeaderObj.nWindows = 1 #podria ser 1 o self.dataOut.processingHeaderObj.nWindows
645 self.processingHeaderObj.nCohInt = self.dataOut.nCohInt# Se requiere para determinar el valor de timeInterval
645 self.processingHeaderObj.nCohInt = self.dataOut.nCohInt# Se requiere para determinar el valor de timeInterval
646 self.processingHeaderObj.nIncohInt = self.dataOut.nIncohInt
646 self.processingHeaderObj.nIncohInt = self.dataOut.nIncohInt
647 self.processingHeaderObj.totalSpectra = self.dataOut.nPairs + self.dataOut.nChannels
647 self.processingHeaderObj.totalSpectra = self.dataOut.nPairs + self.dataOut.nChannels
648 self.processingHeaderObj.shif_fft = self.dataOut.flagShiftFFT
648 self.processingHeaderObj.shif_fft = self.dataOut.flagShiftFFT
649
649
650 if self.processingHeaderObj.totalSpectra > 0:
650 if self.processingHeaderObj.totalSpectra > 0:
651 channelList = []
651 channelList = []
652 for channel in range(self.dataOut.nChannels):
652 for channel in range(self.dataOut.nChannels):
653 channelList.append(channel)
653 channelList.append(channel)
654 channelList.append(channel)
654 channelList.append(channel)
655
655
656 pairsList = []
656 pairsList = []
657 if self.dataOut.nPairs > 0:
657 if self.dataOut.nPairs > 0:
658 for pair in self.dataOut.pairsList:
658 for pair in self.dataOut.pairsList:
659 pairsList.append(pair[0])
659 pairsList.append(pair[0])
660 pairsList.append(pair[1])
660 pairsList.append(pair[1])
661
661
662 spectraComb = channelList + pairsList
662 spectraComb = channelList + pairsList
663 spectraComb = numpy.array(spectraComb, dtype="u1")
663 spectraComb = numpy.array(spectraComb, dtype="u1")
664 self.processingHeaderObj.spectraComb = spectraComb
664 self.processingHeaderObj.spectraComb = spectraComb
665
665
666 if self.dataOut.code is not None:
666 if self.dataOut.code is not None:
667 self.processingHeaderObj.code = self.dataOut.code
667 self.processingHeaderObj.code = self.dataOut.code
668 self.processingHeaderObj.nCode = self.dataOut.nCode
668 self.processingHeaderObj.nCode = self.dataOut.nCode
669 self.processingHeaderObj.nBaud = self.dataOut.nBaud
669 self.processingHeaderObj.nBaud = self.dataOut.nBaud
670
670
671 if self.processingHeaderObj.nWindows != 0:
671 if self.processingHeaderObj.nWindows != 0:
672 self.processingHeaderObj.firstHeight = self.dataOut.heightList[0]
672 self.processingHeaderObj.firstHeight = self.dataOut.heightList[0]
673 self.processingHeaderObj.deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
673 self.processingHeaderObj.deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
674 self.processingHeaderObj.nHeights = self.dataOut.nHeights
674 self.processingHeaderObj.nHeights = self.dataOut.nHeights
675 self.processingHeaderObj.samplesWin = self.dataOut.nHeights
675 self.processingHeaderObj.samplesWin = self.dataOut.nHeights
676
676
677 self.processingHeaderObj.processFlags = self.getProcessFlags()
677 self.processingHeaderObj.processFlags = self.getProcessFlags()
678
678
679 self.setBasicHeader()
679 self.setBasicHeader()
@@ -1,943 +1,952
1 import itertools
1 import itertools
2
2
3 import numpy
3 import numpy
4
4
5 from jroproc_base import ProcessingUnit, Operation
5 from jroproc_base import ProcessingUnit, Operation
6 from schainpy.model.data.jrodata import Spectra
6 from schainpy.model.data.jrodata import Spectra
7 from schainpy.model.data.jrodata import hildebrand_sekhon
7 from schainpy.model.data.jrodata import hildebrand_sekhon
8
8
9
9
10 class SpectraProc(ProcessingUnit):
10 class SpectraProc(ProcessingUnit):
11
11
12 def __init__(self, **kwargs):
12 def __init__(self, **kwargs):
13
13
14 ProcessingUnit.__init__(self, **kwargs)
14 ProcessingUnit.__init__(self, **kwargs)
15
15
16 self.buffer = None
16 self.buffer = None
17 self.firstdatatime = None
17 self.firstdatatime = None
18 self.profIndex = 0
18 self.profIndex = 0
19 self.dataOut = Spectra()
19 self.dataOut = Spectra()
20 self.id_min = None
20 self.id_min = None
21 self.id_max = None
21 self.id_max = None
22
22
23 def __updateSpecFromVoltage(self):
23 def __updateSpecFromVoltage(self):
24
24
25 self.dataOut.timeZone = self.dataIn.timeZone
25 self.dataOut.timeZone = self.dataIn.timeZone
26 self.dataOut.dstFlag = self.dataIn.dstFlag
26 self.dataOut.dstFlag = self.dataIn.dstFlag
27 self.dataOut.errorCount = self.dataIn.errorCount
27 self.dataOut.errorCount = self.dataIn.errorCount
28 self.dataOut.useLocalTime = self.dataIn.useLocalTime
28 self.dataOut.useLocalTime = self.dataIn.useLocalTime
29 try:
29 try:
30 self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
30 self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
31 except:
31 except:
32 pass
32 pass
33 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
33 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
34 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
34 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
35 self.dataOut.channelList = self.dataIn.channelList
35 self.dataOut.channelList = self.dataIn.channelList
36 self.dataOut.heightList = self.dataIn.heightList
36 self.dataOut.heightList = self.dataIn.heightList
37 self.dataOut.dtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
37 self.dataOut.dtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
38
38
39 self.dataOut.nBaud = self.dataIn.nBaud
39 self.dataOut.nBaud = self.dataIn.nBaud
40 self.dataOut.nCode = self.dataIn.nCode
40 self.dataOut.nCode = self.dataIn.nCode
41 self.dataOut.code = self.dataIn.code
41 self.dataOut.code = self.dataIn.code
42 self.dataOut.nProfiles = self.dataOut.nFFTPoints
42 self.dataOut.nProfiles = self.dataOut.nFFTPoints
43
43
44 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
44 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
45 self.dataOut.utctime = self.firstdatatime
45 self.dataOut.utctime = self.firstdatatime
46 # asumo q la data esta decodificada
46 # asumo q la data esta decodificada
47 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData
47 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData
48 # asumo q la data esta sin flip
48 # asumo q la data esta sin flip
49 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData
49 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData
50 self.dataOut.flagShiftFFT = False
50 self.dataOut.flagShiftFFT = False
51
51
52 self.dataOut.nCohInt = self.dataIn.nCohInt
52 self.dataOut.nCohInt = self.dataIn.nCohInt
53 self.dataOut.nIncohInt = 1
53 self.dataOut.nIncohInt = 1
54
54
55 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
55 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
56
56
57 self.dataOut.frequency = self.dataIn.frequency
57 self.dataOut.frequency = self.dataIn.frequency
58 self.dataOut.realtime = self.dataIn.realtime
58 self.dataOut.realtime = self.dataIn.realtime
59
59
60 self.dataOut.azimuth = self.dataIn.azimuth
60 self.dataOut.azimuth = self.dataIn.azimuth
61 self.dataOut.zenith = self.dataIn.zenith
61 self.dataOut.zenith = self.dataIn.zenith
62
62
63 self.dataOut.beam.codeList = self.dataIn.beam.codeList
63 self.dataOut.beam.codeList = self.dataIn.beam.codeList
64 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
64 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
65 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
65 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
66
66
67 def __getFft(self):
67 def __getFft(self):
68 """
68 """
69 Convierte valores de Voltaje a Spectra
69 Convierte valores de Voltaje a Spectra
70
70
71 Affected:
71 Affected:
72 self.dataOut.data_spc
72 self.dataOut.data_spc
73 self.dataOut.data_cspc
73 self.dataOut.data_cspc
74 self.dataOut.data_dc
74 self.dataOut.data_dc
75 self.dataOut.heightList
75 self.dataOut.heightList
76 self.profIndex
76 self.profIndex
77 self.buffer
77 self.buffer
78 self.dataOut.flagNoData
78 self.dataOut.flagNoData
79 """
79 """
80 fft_volt = numpy.fft.fft(
80 fft_volt = numpy.fft.fft(
81 self.buffer, n=self.dataOut.nFFTPoints, axis=1)
81 self.buffer, n=self.dataOut.nFFTPoints, axis=1)
82 fft_volt = fft_volt.astype(numpy.dtype('complex'))
82 fft_volt = fft_volt.astype(numpy.dtype('complex'))
83 dc = fft_volt[:, 0, :]
83 dc = fft_volt[:, 0, :]
84
84
85 # calculo de self-spectra
85 # calculo de self-spectra
86 fft_volt = numpy.fft.fftshift(fft_volt, axes=(1,))
86 fft_volt = numpy.fft.fftshift(fft_volt, axes=(1,))
87 spc = fft_volt * numpy.conjugate(fft_volt)
87 spc = fft_volt * numpy.conjugate(fft_volt)
88 spc = spc.real
88 spc = spc.real
89
89
90 blocksize = 0
90 blocksize = 0
91 blocksize += dc.size
91 blocksize += dc.size
92 blocksize += spc.size
92 blocksize += spc.size
93
93
94 cspc = None
94 cspc = None
95 pairIndex = 0
95 pairIndex = 0
96 if self.dataOut.pairsList != None:
96 if self.dataOut.pairsList != None:
97 # calculo de cross-spectra
97 # calculo de cross-spectra
98 cspc = numpy.zeros(
98 cspc = numpy.zeros(
99 (self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
99 (self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
100 for pair in self.dataOut.pairsList:
100 for pair in self.dataOut.pairsList:
101 if pair[0] not in self.dataOut.channelList:
101 if pair[0] not in self.dataOut.channelList:
102 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" % (
102 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" % (
103 str(pair), str(self.dataOut.channelList))
103 str(pair), str(self.dataOut.channelList))
104 if pair[1] not in self.dataOut.channelList:
104 if pair[1] not in self.dataOut.channelList:
105 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" % (
105 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" % (
106 str(pair), str(self.dataOut.channelList))
106 str(pair), str(self.dataOut.channelList))
107
107
108 cspc[pairIndex, :, :] = fft_volt[pair[0], :, :] * \
108 cspc[pairIndex, :, :] = fft_volt[pair[0], :, :] * \
109 numpy.conjugate(fft_volt[pair[1], :, :])
109 numpy.conjugate(fft_volt[pair[1], :, :])
110 pairIndex += 1
110 pairIndex += 1
111 blocksize += cspc.size
111 blocksize += cspc.size
112
112
113 self.dataOut.data_spc = spc
113 self.dataOut.data_spc = spc
114 self.dataOut.data_cspc = cspc
114 self.dataOut.data_cspc = cspc
115 self.dataOut.data_dc = dc
115 self.dataOut.data_dc = dc
116 self.dataOut.blockSize = blocksize
116 self.dataOut.blockSize = blocksize
117 self.dataOut.flagShiftFFT = True
117 self.dataOut.flagShiftFFT = True
118
118
119 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], ippFactor=None):
119 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], ippFactor=None, shift_fft=False):
120
120
121 self.dataOut.flagNoData = True
121 self.dataOut.flagNoData = True
122
122
123 if self.dataIn.type == "Spectra":
123 if self.dataIn.type == "Spectra":
124 self.dataOut.copy(self.dataIn)
124 self.dataOut.copy(self.dataIn)
125 if not pairsList:
125 if not pairsList:
126 pairsList = itertools.combinations(self.dataOut.channelList, 2)
126 pairsList = itertools.combinations(self.dataOut.channelList, 2)
127 if self.dataOut.data_cspc is not None:
127 if self.dataOut.data_cspc is not None:
128 self.__selectPairs(pairsList)
128 self.__selectPairs(pairsList)
129 if shift_fft:
130 #desplaza a la derecha en el eje 2 determinadas posiciones
131 shift = int(self.dataOut.nFFTPoints/2)
132 self.dataOut.data_spc = numpy.roll(self.dataOut.data_spc, shift , axis=1)
133
134 if self.dataOut.data_cspc is not None:
135 #desplaza a la derecha en el eje 2 determinadas posiciones
136 self.dataOut.data_cspc = numpy.roll(self.dataOut.cspc, shift, axis=1)
137
129 return True
138 return True
130
139
131 if self.dataIn.type == "Voltage":
140 if self.dataIn.type == "Voltage":
132
141
133 if nFFTPoints == None:
142 if nFFTPoints == None:
134 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
143 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
135
144
136 if nProfiles == None:
145 if nProfiles == None:
137 nProfiles = nFFTPoints
146 nProfiles = nFFTPoints
138
147
139 if ippFactor == None:
148 if ippFactor == None:
140 ippFactor = 1
149 ippFactor = 1
141
150
142 self.dataOut.ippFactor = ippFactor
151 self.dataOut.ippFactor = ippFactor
143
152
144 self.dataOut.nFFTPoints = nFFTPoints
153 self.dataOut.nFFTPoints = nFFTPoints
145 self.dataOut.pairsList = pairsList
154 self.dataOut.pairsList = pairsList
146
155
147 if self.buffer is None:
156 if self.buffer is None:
148 self.buffer = numpy.zeros((self.dataIn.nChannels,
157 self.buffer = numpy.zeros((self.dataIn.nChannels,
149 nProfiles,
158 nProfiles,
150 self.dataIn.nHeights),
159 self.dataIn.nHeights),
151 dtype='complex')
160 dtype='complex')
152
161
153 if self.dataIn.flagDataAsBlock:
162 if self.dataIn.flagDataAsBlock:
154 # data dimension: [nChannels, nProfiles, nSamples]
163 # data dimension: [nChannels, nProfiles, nSamples]
155 nVoltProfiles = self.dataIn.data.shape[1]
164 nVoltProfiles = self.dataIn.data.shape[1]
156 # nVoltProfiles = self.dataIn.nProfiles
165 # nVoltProfiles = self.dataIn.nProfiles
157
166
158 if nVoltProfiles == nProfiles:
167 if nVoltProfiles == nProfiles:
159 self.buffer = self.dataIn.data.copy()
168 self.buffer = self.dataIn.data.copy()
160 self.profIndex = nVoltProfiles
169 self.profIndex = nVoltProfiles
161
170
162 elif nVoltProfiles < nProfiles:
171 elif nVoltProfiles < nProfiles:
163
172
164 if self.profIndex == 0:
173 if self.profIndex == 0:
165 self.id_min = 0
174 self.id_min = 0
166 self.id_max = nVoltProfiles
175 self.id_max = nVoltProfiles
167
176
168 self.buffer[:, self.id_min:self.id_max,
177 self.buffer[:, self.id_min:self.id_max,
169 :] = self.dataIn.data
178 :] = self.dataIn.data
170 self.profIndex += nVoltProfiles
179 self.profIndex += nVoltProfiles
171 self.id_min += nVoltProfiles
180 self.id_min += nVoltProfiles
172 self.id_max += nVoltProfiles
181 self.id_max += nVoltProfiles
173 else:
182 else:
174 raise ValueError, "The type object %s has %d profiles, it should just has %d profiles" % (
183 raise ValueError, "The type object %s has %d profiles, it should just has %d profiles" % (
175 self.dataIn.type, self.dataIn.data.shape[1], nProfiles)
184 self.dataIn.type, self.dataIn.data.shape[1], nProfiles)
176 self.dataOut.flagNoData = True
185 self.dataOut.flagNoData = True
177 return 0
186 return 0
178 else:
187 else:
179 self.buffer[:, self.profIndex, :] = self.dataIn.data.copy()
188 self.buffer[:, self.profIndex, :] = self.dataIn.data.copy()
180 self.profIndex += 1
189 self.profIndex += 1
181
190
182 if self.firstdatatime == None:
191 if self.firstdatatime == None:
183 self.firstdatatime = self.dataIn.utctime
192 self.firstdatatime = self.dataIn.utctime
184
193
185 if self.profIndex == nProfiles:
194 if self.profIndex == nProfiles:
186 self.__updateSpecFromVoltage()
195 self.__updateSpecFromVoltage()
187 self.__getFft()
196 self.__getFft()
188
197
189 self.dataOut.flagNoData = False
198 self.dataOut.flagNoData = False
190 self.firstdatatime = None
199 self.firstdatatime = None
191 self.profIndex = 0
200 self.profIndex = 0
192
201
193 return True
202 return True
194
203
195 raise ValueError, "The type of input object '%s' is not valid" % (
204 raise ValueError, "The type of input object '%s' is not valid" % (
196 self.dataIn.type)
205 self.dataIn.type)
197
206
198 def __selectPairs(self, pairsList):
207 def __selectPairs(self, pairsList):
199
208
200 if not pairsList:
209 if not pairsList:
201 return
210 return
202
211
203 pairs = []
212 pairs = []
204 pairsIndex = []
213 pairsIndex = []
205
214
206 for pair in pairsList:
215 for pair in pairsList:
207 if pair[0] not in self.dataOut.channelList or pair[1] not in self.dataOut.channelList:
216 if pair[0] not in self.dataOut.channelList or pair[1] not in self.dataOut.channelList:
208 continue
217 continue
209 pairs.append(pair)
218 pairs.append(pair)
210 pairsIndex.append(pairs.index(pair))
219 pairsIndex.append(pairs.index(pair))
211
220
212 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndex]
221 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndex]
213 self.dataOut.pairsList = pairs
222 self.dataOut.pairsList = pairs
214
223
215 return
224 return
216
225
217 def __selectPairsByChannel(self, channelList=None):
226 def __selectPairsByChannel(self, channelList=None):
218
227
219 if channelList == None:
228 if channelList == None:
220 return
229 return
221
230
222 pairsIndexListSelected = []
231 pairsIndexListSelected = []
223 for pairIndex in self.dataOut.pairsIndexList:
232 for pairIndex in self.dataOut.pairsIndexList:
224 # First pair
233 # First pair
225 if self.dataOut.pairsList[pairIndex][0] not in channelList:
234 if self.dataOut.pairsList[pairIndex][0] not in channelList:
226 continue
235 continue
227 # Second pair
236 # Second pair
228 if self.dataOut.pairsList[pairIndex][1] not in channelList:
237 if self.dataOut.pairsList[pairIndex][1] not in channelList:
229 continue
238 continue
230
239
231 pairsIndexListSelected.append(pairIndex)
240 pairsIndexListSelected.append(pairIndex)
232
241
233 if not pairsIndexListSelected:
242 if not pairsIndexListSelected:
234 self.dataOut.data_cspc = None
243 self.dataOut.data_cspc = None
235 self.dataOut.pairsList = []
244 self.dataOut.pairsList = []
236 return
245 return
237
246
238 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
247 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
239 self.dataOut.pairsList = [self.dataOut.pairsList[i]
248 self.dataOut.pairsList = [self.dataOut.pairsList[i]
240 for i in pairsIndexListSelected]
249 for i in pairsIndexListSelected]
241
250
242 return
251 return
243
252
244 def selectChannels(self, channelList):
253 def selectChannels(self, channelList):
245
254
246 channelIndexList = []
255 channelIndexList = []
247
256
248 for channel in channelList:
257 for channel in channelList:
249 if channel not in self.dataOut.channelList:
258 if channel not in self.dataOut.channelList:
250 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" % (
259 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" % (
251 channel, str(self.dataOut.channelList))
260 channel, str(self.dataOut.channelList))
252
261
253 index = self.dataOut.channelList.index(channel)
262 index = self.dataOut.channelList.index(channel)
254 channelIndexList.append(index)
263 channelIndexList.append(index)
255
264
256 self.selectChannelsByIndex(channelIndexList)
265 self.selectChannelsByIndex(channelIndexList)
257
266
258 def selectChannelsByIndex(self, channelIndexList):
267 def selectChannelsByIndex(self, channelIndexList):
259 """
268 """
260 Selecciona un bloque de datos en base a canales segun el channelIndexList
269 Selecciona un bloque de datos en base a canales segun el channelIndexList
261
270
262 Input:
271 Input:
263 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
272 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
264
273
265 Affected:
274 Affected:
266 self.dataOut.data_spc
275 self.dataOut.data_spc
267 self.dataOut.channelIndexList
276 self.dataOut.channelIndexList
268 self.dataOut.nChannels
277 self.dataOut.nChannels
269
278
270 Return:
279 Return:
271 None
280 None
272 """
281 """
273
282
274 for channelIndex in channelIndexList:
283 for channelIndex in channelIndexList:
275 if channelIndex not in self.dataOut.channelIndexList:
284 if channelIndex not in self.dataOut.channelIndexList:
276 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " % (
285 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " % (
277 channelIndex, self.dataOut.channelIndexList)
286 channelIndex, self.dataOut.channelIndexList)
278
287
279 # nChannels = len(channelIndexList)
288 # nChannels = len(channelIndexList)
280
289
281 data_spc = self.dataOut.data_spc[channelIndexList, :]
290 data_spc = self.dataOut.data_spc[channelIndexList, :]
282 data_dc = self.dataOut.data_dc[channelIndexList, :]
291 data_dc = self.dataOut.data_dc[channelIndexList, :]
283
292
284 self.dataOut.data_spc = data_spc
293 self.dataOut.data_spc = data_spc
285 self.dataOut.data_dc = data_dc
294 self.dataOut.data_dc = data_dc
286
295
287 self.dataOut.channelList = [
296 self.dataOut.channelList = [
288 self.dataOut.channelList[i] for i in channelIndexList]
297 self.dataOut.channelList[i] for i in channelIndexList]
289 # self.dataOut.nChannels = nChannels
298 # self.dataOut.nChannels = nChannels
290
299
291 self.__selectPairsByChannel(self.dataOut.channelList)
300 self.__selectPairsByChannel(self.dataOut.channelList)
292
301
293 return 1
302 return 1
294
303
295 def selectHeights(self, minHei, maxHei):
304 def selectHeights(self, minHei, maxHei):
296 """
305 """
297 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
306 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
298 minHei <= height <= maxHei
307 minHei <= height <= maxHei
299
308
300 Input:
309 Input:
301 minHei : valor minimo de altura a considerar
310 minHei : valor minimo de altura a considerar
302 maxHei : valor maximo de altura a considerar
311 maxHei : valor maximo de altura a considerar
303
312
304 Affected:
313 Affected:
305 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
314 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
306
315
307 Return:
316 Return:
308 1 si el metodo se ejecuto con exito caso contrario devuelve 0
317 1 si el metodo se ejecuto con exito caso contrario devuelve 0
309 """
318 """
310
319
311 if (minHei > maxHei):
320 if (minHei > maxHei):
312 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (
321 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (
313 minHei, maxHei)
322 minHei, maxHei)
314
323
315 if (minHei < self.dataOut.heightList[0]):
324 if (minHei < self.dataOut.heightList[0]):
316 minHei = self.dataOut.heightList[0]
325 minHei = self.dataOut.heightList[0]
317
326
318 if (maxHei > self.dataOut.heightList[-1]):
327 if (maxHei > self.dataOut.heightList[-1]):
319 maxHei = self.dataOut.heightList[-1]
328 maxHei = self.dataOut.heightList[-1]
320
329
321 minIndex = 0
330 minIndex = 0
322 maxIndex = 0
331 maxIndex = 0
323 heights = self.dataOut.heightList
332 heights = self.dataOut.heightList
324
333
325 inda = numpy.where(heights >= minHei)
334 inda = numpy.where(heights >= minHei)
326 indb = numpy.where(heights <= maxHei)
335 indb = numpy.where(heights <= maxHei)
327
336
328 try:
337 try:
329 minIndex = inda[0][0]
338 minIndex = inda[0][0]
330 except:
339 except:
331 minIndex = 0
340 minIndex = 0
332
341
333 try:
342 try:
334 maxIndex = indb[0][-1]
343 maxIndex = indb[0][-1]
335 except:
344 except:
336 maxIndex = len(heights)
345 maxIndex = len(heights)
337
346
338 self.selectHeightsByIndex(minIndex, maxIndex)
347 self.selectHeightsByIndex(minIndex, maxIndex)
339
348
340 return 1
349 return 1
341
350
342 def getBeaconSignal(self, tauindex=0, channelindex=0, hei_ref=None):
351 def getBeaconSignal(self, tauindex=0, channelindex=0, hei_ref=None):
343 newheis = numpy.where(
352 newheis = numpy.where(
344 self.dataOut.heightList > self.dataOut.radarControllerHeaderObj.Taus[tauindex])
353 self.dataOut.heightList > self.dataOut.radarControllerHeaderObj.Taus[tauindex])
345
354
346 if hei_ref != None:
355 if hei_ref != None:
347 newheis = numpy.where(self.dataOut.heightList > hei_ref)
356 newheis = numpy.where(self.dataOut.heightList > hei_ref)
348
357
349 minIndex = min(newheis[0])
358 minIndex = min(newheis[0])
350 maxIndex = max(newheis[0])
359 maxIndex = max(newheis[0])
351 data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
360 data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
352 heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
361 heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
353
362
354 # determina indices
363 # determina indices
355 nheis = int(self.dataOut.radarControllerHeaderObj.txB /
364 nheis = int(self.dataOut.radarControllerHeaderObj.txB /
356 (self.dataOut.heightList[1] - self.dataOut.heightList[0]))
365 (self.dataOut.heightList[1] - self.dataOut.heightList[0]))
357 avg_dB = 10 * \
366 avg_dB = 10 * \
358 numpy.log10(numpy.sum(data_spc[channelindex, :, :], axis=0))
367 numpy.log10(numpy.sum(data_spc[channelindex, :, :], axis=0))
359 beacon_dB = numpy.sort(avg_dB)[-nheis:]
368 beacon_dB = numpy.sort(avg_dB)[-nheis:]
360 beacon_heiIndexList = []
369 beacon_heiIndexList = []
361 for val in avg_dB.tolist():
370 for val in avg_dB.tolist():
362 if val >= beacon_dB[0]:
371 if val >= beacon_dB[0]:
363 beacon_heiIndexList.append(avg_dB.tolist().index(val))
372 beacon_heiIndexList.append(avg_dB.tolist().index(val))
364
373
365 #data_spc = data_spc[:,:,beacon_heiIndexList]
374 #data_spc = data_spc[:,:,beacon_heiIndexList]
366 data_cspc = None
375 data_cspc = None
367 if self.dataOut.data_cspc is not None:
376 if self.dataOut.data_cspc is not None:
368 data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
377 data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
369 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
378 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
370
379
371 data_dc = None
380 data_dc = None
372 if self.dataOut.data_dc is not None:
381 if self.dataOut.data_dc is not None:
373 data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
382 data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
374 #data_dc = data_dc[:,beacon_heiIndexList]
383 #data_dc = data_dc[:,beacon_heiIndexList]
375
384
376 self.dataOut.data_spc = data_spc
385 self.dataOut.data_spc = data_spc
377 self.dataOut.data_cspc = data_cspc
386 self.dataOut.data_cspc = data_cspc
378 self.dataOut.data_dc = data_dc
387 self.dataOut.data_dc = data_dc
379 self.dataOut.heightList = heightList
388 self.dataOut.heightList = heightList
380 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
389 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
381
390
382 return 1
391 return 1
383
392
384 def selectHeightsByIndex(self, minIndex, maxIndex):
393 def selectHeightsByIndex(self, minIndex, maxIndex):
385 """
394 """
386 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
395 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
387 minIndex <= index <= maxIndex
396 minIndex <= index <= maxIndex
388
397
389 Input:
398 Input:
390 minIndex : valor de indice minimo de altura a considerar
399 minIndex : valor de indice minimo de altura a considerar
391 maxIndex : valor de indice maximo de altura a considerar
400 maxIndex : valor de indice maximo de altura a considerar
392
401
393 Affected:
402 Affected:
394 self.dataOut.data_spc
403 self.dataOut.data_spc
395 self.dataOut.data_cspc
404 self.dataOut.data_cspc
396 self.dataOut.data_dc
405 self.dataOut.data_dc
397 self.dataOut.heightList
406 self.dataOut.heightList
398
407
399 Return:
408 Return:
400 1 si el metodo se ejecuto con exito caso contrario devuelve 0
409 1 si el metodo se ejecuto con exito caso contrario devuelve 0
401 """
410 """
402
411
403 if (minIndex < 0) or (minIndex > maxIndex):
412 if (minIndex < 0) or (minIndex > maxIndex):
404 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (
413 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (
405 minIndex, maxIndex)
414 minIndex, maxIndex)
406
415
407 if (maxIndex >= self.dataOut.nHeights):
416 if (maxIndex >= self.dataOut.nHeights):
408 maxIndex = self.dataOut.nHeights - 1
417 maxIndex = self.dataOut.nHeights - 1
409
418
410 # Spectra
419 # Spectra
411 data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
420 data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
412
421
413 data_cspc = None
422 data_cspc = None
414 if self.dataOut.data_cspc is not None:
423 if self.dataOut.data_cspc is not None:
415 data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
424 data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
416
425
417 data_dc = None
426 data_dc = None
418 if self.dataOut.data_dc is not None:
427 if self.dataOut.data_dc is not None:
419 data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
428 data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
420
429
421 self.dataOut.data_spc = data_spc
430 self.dataOut.data_spc = data_spc
422 self.dataOut.data_cspc = data_cspc
431 self.dataOut.data_cspc = data_cspc
423 self.dataOut.data_dc = data_dc
432 self.dataOut.data_dc = data_dc
424
433
425 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
434 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
426
435
427 return 1
436 return 1
428
437
429 def removeDC(self, mode=2):
438 def removeDC(self, mode=2):
430 jspectra = self.dataOut.data_spc
439 jspectra = self.dataOut.data_spc
431 jcspectra = self.dataOut.data_cspc
440 jcspectra = self.dataOut.data_cspc
432
441
433 num_chan = jspectra.shape[0]
442 num_chan = jspectra.shape[0]
434 num_hei = jspectra.shape[2]
443 num_hei = jspectra.shape[2]
435
444
436 if jcspectra is not None:
445 if jcspectra is not None:
437 jcspectraExist = True
446 jcspectraExist = True
438 num_pairs = jcspectra.shape[0]
447 num_pairs = jcspectra.shape[0]
439 else:
448 else:
440 jcspectraExist = False
449 jcspectraExist = False
441
450
442 freq_dc = jspectra.shape[1] / 2
451 freq_dc = jspectra.shape[1] / 2
443 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
452 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
444
453
445 if ind_vel[0] < 0:
454 if ind_vel[0] < 0:
446 ind_vel[range(0, 1)] = ind_vel[range(0, 1)] + self.num_prof
455 ind_vel[range(0, 1)] = ind_vel[range(0, 1)] + self.num_prof
447
456
448 if mode == 1:
457 if mode == 1:
449 jspectra[:, freq_dc, :] = (
458 jspectra[:, freq_dc, :] = (
450 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
459 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
451
460
452 if jcspectraExist:
461 if jcspectraExist:
453 jcspectra[:, freq_dc, :] = (
462 jcspectra[:, freq_dc, :] = (
454 jcspectra[:, ind_vel[1], :] + jcspectra[:, ind_vel[2], :]) / 2
463 jcspectra[:, ind_vel[1], :] + jcspectra[:, ind_vel[2], :]) / 2
455
464
456 if mode == 2:
465 if mode == 2:
457
466
458 vel = numpy.array([-2, -1, 1, 2])
467 vel = numpy.array([-2, -1, 1, 2])
459 xx = numpy.zeros([4, 4])
468 xx = numpy.zeros([4, 4])
460
469
461 for fil in range(4):
470 for fil in range(4):
462 xx[fil, :] = vel[fil]**numpy.asarray(range(4))
471 xx[fil, :] = vel[fil]**numpy.asarray(range(4))
463
472
464 xx_inv = numpy.linalg.inv(xx)
473 xx_inv = numpy.linalg.inv(xx)
465 xx_aux = xx_inv[0, :]
474 xx_aux = xx_inv[0, :]
466
475
467 for ich in range(num_chan):
476 for ich in range(num_chan):
468 yy = jspectra[ich, ind_vel, :]
477 yy = jspectra[ich, ind_vel, :]
469 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
478 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
470
479
471 junkid = jspectra[ich, freq_dc, :] <= 0
480 junkid = jspectra[ich, freq_dc, :] <= 0
472 cjunkid = sum(junkid)
481 cjunkid = sum(junkid)
473
482
474 if cjunkid.any():
483 if cjunkid.any():
475 jspectra[ich, freq_dc, junkid.nonzero()] = (
484 jspectra[ich, freq_dc, junkid.nonzero()] = (
476 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
485 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
477
486
478 if jcspectraExist:
487 if jcspectraExist:
479 for ip in range(num_pairs):
488 for ip in range(num_pairs):
480 yy = jcspectra[ip, ind_vel, :]
489 yy = jcspectra[ip, ind_vel, :]
481 jcspectra[ip, freq_dc, :] = numpy.dot(xx_aux, yy)
490 jcspectra[ip, freq_dc, :] = numpy.dot(xx_aux, yy)
482
491
483 self.dataOut.data_spc = jspectra
492 self.dataOut.data_spc = jspectra
484 self.dataOut.data_cspc = jcspectra
493 self.dataOut.data_cspc = jcspectra
485
494
486 return 1
495 return 1
487
496
488 def removeInterference(self, interf=2, hei_interf=None, nhei_interf=None, offhei_interf=None):
497 def removeInterference(self, interf=2, hei_interf=None, nhei_interf=None, offhei_interf=None):
489
498
490 jspectra = self.dataOut.data_spc
499 jspectra = self.dataOut.data_spc
491 jcspectra = self.dataOut.data_cspc
500 jcspectra = self.dataOut.data_cspc
492 jnoise = self.dataOut.getNoise()
501 jnoise = self.dataOut.getNoise()
493 num_incoh = self.dataOut.nIncohInt
502 num_incoh = self.dataOut.nIncohInt
494
503
495 num_channel = jspectra.shape[0]
504 num_channel = jspectra.shape[0]
496 num_prof = jspectra.shape[1]
505 num_prof = jspectra.shape[1]
497 num_hei = jspectra.shape[2]
506 num_hei = jspectra.shape[2]
498
507
499 # hei_interf
508 # hei_interf
500 if hei_interf is None:
509 if hei_interf is None:
501 count_hei = num_hei / 2 # Como es entero no importa
510 count_hei = num_hei / 2 # Como es entero no importa
502 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
511 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
503 hei_interf = numpy.asarray(hei_interf)[0]
512 hei_interf = numpy.asarray(hei_interf)[0]
504 # nhei_interf
513 # nhei_interf
505 if (nhei_interf == None):
514 if (nhei_interf == None):
506 nhei_interf = 5
515 nhei_interf = 5
507 if (nhei_interf < 1):
516 if (nhei_interf < 1):
508 nhei_interf = 1
517 nhei_interf = 1
509 if (nhei_interf > count_hei):
518 if (nhei_interf > count_hei):
510 nhei_interf = count_hei
519 nhei_interf = count_hei
511 if (offhei_interf == None):
520 if (offhei_interf == None):
512 offhei_interf = 0
521 offhei_interf = 0
513
522
514 ind_hei = range(num_hei)
523 ind_hei = range(num_hei)
515 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
524 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
516 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
525 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
517 mask_prof = numpy.asarray(range(num_prof))
526 mask_prof = numpy.asarray(range(num_prof))
518 num_mask_prof = mask_prof.size
527 num_mask_prof = mask_prof.size
519 comp_mask_prof = [0, num_prof / 2]
528 comp_mask_prof = [0, num_prof / 2]
520
529
521 # noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
530 # noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
522 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
531 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
523 jnoise = numpy.nan
532 jnoise = numpy.nan
524 noise_exist = jnoise[0] < numpy.Inf
533 noise_exist = jnoise[0] < numpy.Inf
525
534
526 # Subrutina de Remocion de la Interferencia
535 # Subrutina de Remocion de la Interferencia
527 for ich in range(num_channel):
536 for ich in range(num_channel):
528 # Se ordena los espectros segun su potencia (menor a mayor)
537 # Se ordena los espectros segun su potencia (menor a mayor)
529 power = jspectra[ich, mask_prof, :]
538 power = jspectra[ich, mask_prof, :]
530 power = power[:, hei_interf]
539 power = power[:, hei_interf]
531 power = power.sum(axis=0)
540 power = power.sum(axis=0)
532 psort = power.ravel().argsort()
541 psort = power.ravel().argsort()
533
542
534 # Se estima la interferencia promedio en los Espectros de Potencia empleando
543 # Se estima la interferencia promedio en los Espectros de Potencia empleando
535 junkspc_interf = jspectra[ich, :, hei_interf[psort[range(
544 junkspc_interf = jspectra[ich, :, hei_interf[psort[range(
536 offhei_interf, nhei_interf + offhei_interf)]]]
545 offhei_interf, nhei_interf + offhei_interf)]]]
537
546
538 if noise_exist:
547 if noise_exist:
539 # tmp_noise = jnoise[ich] / num_prof
548 # tmp_noise = jnoise[ich] / num_prof
540 tmp_noise = jnoise[ich]
549 tmp_noise = jnoise[ich]
541 junkspc_interf = junkspc_interf - tmp_noise
550 junkspc_interf = junkspc_interf - tmp_noise
542 #junkspc_interf[:,comp_mask_prof] = 0
551 #junkspc_interf[:,comp_mask_prof] = 0
543
552
544 jspc_interf = junkspc_interf.sum(axis=0) / nhei_interf
553 jspc_interf = junkspc_interf.sum(axis=0) / nhei_interf
545 jspc_interf = jspc_interf.transpose()
554 jspc_interf = jspc_interf.transpose()
546 # Calculando el espectro de interferencia promedio
555 # Calculando el espectro de interferencia promedio
547 noiseid = numpy.where(
556 noiseid = numpy.where(
548 jspc_interf <= tmp_noise / numpy.sqrt(num_incoh))
557 jspc_interf <= tmp_noise / numpy.sqrt(num_incoh))
549 noiseid = noiseid[0]
558 noiseid = noiseid[0]
550 cnoiseid = noiseid.size
559 cnoiseid = noiseid.size
551 interfid = numpy.where(
560 interfid = numpy.where(
552 jspc_interf > tmp_noise / numpy.sqrt(num_incoh))
561 jspc_interf > tmp_noise / numpy.sqrt(num_incoh))
553 interfid = interfid[0]
562 interfid = interfid[0]
554 cinterfid = interfid.size
563 cinterfid = interfid.size
555
564
556 if (cnoiseid > 0):
565 if (cnoiseid > 0):
557 jspc_interf[noiseid] = 0
566 jspc_interf[noiseid] = 0
558
567
559 # Expandiendo los perfiles a limpiar
568 # Expandiendo los perfiles a limpiar
560 if (cinterfid > 0):
569 if (cinterfid > 0):
561 new_interfid = (
570 new_interfid = (
562 numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof) % num_prof
571 numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof) % num_prof
563 new_interfid = numpy.asarray(new_interfid)
572 new_interfid = numpy.asarray(new_interfid)
564 new_interfid = {x for x in new_interfid}
573 new_interfid = {x for x in new_interfid}
565 new_interfid = numpy.array(list(new_interfid))
574 new_interfid = numpy.array(list(new_interfid))
566 new_cinterfid = new_interfid.size
575 new_cinterfid = new_interfid.size
567 else:
576 else:
568 new_cinterfid = 0
577 new_cinterfid = 0
569
578
570 for ip in range(new_cinterfid):
579 for ip in range(new_cinterfid):
571 ind = junkspc_interf[:, new_interfid[ip]].ravel().argsort()
580 ind = junkspc_interf[:, new_interfid[ip]].ravel().argsort()
572 jspc_interf[new_interfid[ip]
581 jspc_interf[new_interfid[ip]
573 ] = junkspc_interf[ind[nhei_interf / 2], new_interfid[ip]]
582 ] = junkspc_interf[ind[nhei_interf / 2], new_interfid[ip]]
574
583
575 jspectra[ich, :, ind_hei] = jspectra[ich, :,
584 jspectra[ich, :, ind_hei] = jspectra[ich, :,
576 ind_hei] - jspc_interf # Corregir indices
585 ind_hei] - jspc_interf # Corregir indices
577
586
578 # Removiendo la interferencia del punto de mayor interferencia
587 # Removiendo la interferencia del punto de mayor interferencia
579 ListAux = jspc_interf[mask_prof].tolist()
588 ListAux = jspc_interf[mask_prof].tolist()
580 maxid = ListAux.index(max(ListAux))
589 maxid = ListAux.index(max(ListAux))
581
590
582 if cinterfid > 0:
591 if cinterfid > 0:
583 for ip in range(cinterfid * (interf == 2) - 1):
592 for ip in range(cinterfid * (interf == 2) - 1):
584 ind = (jspectra[ich, interfid[ip], :] < tmp_noise *
593 ind = (jspectra[ich, interfid[ip], :] < tmp_noise *
585 (1 + 1 / numpy.sqrt(num_incoh))).nonzero()
594 (1 + 1 / numpy.sqrt(num_incoh))).nonzero()
586 cind = len(ind)
595 cind = len(ind)
587
596
588 if (cind > 0):
597 if (cind > 0):
589 jspectra[ich, interfid[ip], ind] = tmp_noise * \
598 jspectra[ich, interfid[ip], ind] = tmp_noise * \
590 (1 + (numpy.random.uniform(cind) - 0.5) /
599 (1 + (numpy.random.uniform(cind) - 0.5) /
591 numpy.sqrt(num_incoh))
600 numpy.sqrt(num_incoh))
592
601
593 ind = numpy.array([-2, -1, 1, 2])
602 ind = numpy.array([-2, -1, 1, 2])
594 xx = numpy.zeros([4, 4])
603 xx = numpy.zeros([4, 4])
595
604
596 for id1 in range(4):
605 for id1 in range(4):
597 xx[:, id1] = ind[id1]**numpy.asarray(range(4))
606 xx[:, id1] = ind[id1]**numpy.asarray(range(4))
598
607
599 xx_inv = numpy.linalg.inv(xx)
608 xx_inv = numpy.linalg.inv(xx)
600 xx = xx_inv[:, 0]
609 xx = xx_inv[:, 0]
601 ind = (ind + maxid + num_mask_prof) % num_mask_prof
610 ind = (ind + maxid + num_mask_prof) % num_mask_prof
602 yy = jspectra[ich, mask_prof[ind], :]
611 yy = jspectra[ich, mask_prof[ind], :]
603 jspectra[ich, mask_prof[maxid], :] = numpy.dot(
612 jspectra[ich, mask_prof[maxid], :] = numpy.dot(
604 yy.transpose(), xx)
613 yy.transpose(), xx)
605
614
606 indAux = (jspectra[ich, :, :] < tmp_noise *
615 indAux = (jspectra[ich, :, :] < tmp_noise *
607 (1 - 1 / numpy.sqrt(num_incoh))).nonzero()
616 (1 - 1 / numpy.sqrt(num_incoh))).nonzero()
608 jspectra[ich, indAux[0], indAux[1]] = tmp_noise * \
617 jspectra[ich, indAux[0], indAux[1]] = tmp_noise * \
609 (1 - 1 / numpy.sqrt(num_incoh))
618 (1 - 1 / numpy.sqrt(num_incoh))
610
619
611 # Remocion de Interferencia en el Cross Spectra
620 # Remocion de Interferencia en el Cross Spectra
612 if jcspectra is None:
621 if jcspectra is None:
613 return jspectra, jcspectra
622 return jspectra, jcspectra
614 num_pairs = jcspectra.size / (num_prof * num_hei)
623 num_pairs = jcspectra.size / (num_prof * num_hei)
615 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
624 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
616
625
617 for ip in range(num_pairs):
626 for ip in range(num_pairs):
618
627
619 #-------------------------------------------
628 #-------------------------------------------
620
629
621 cspower = numpy.abs(jcspectra[ip, mask_prof, :])
630 cspower = numpy.abs(jcspectra[ip, mask_prof, :])
622 cspower = cspower[:, hei_interf]
631 cspower = cspower[:, hei_interf]
623 cspower = cspower.sum(axis=0)
632 cspower = cspower.sum(axis=0)
624
633
625 cspsort = cspower.ravel().argsort()
634 cspsort = cspower.ravel().argsort()
626 junkcspc_interf = jcspectra[ip, :, hei_interf[cspsort[range(
635 junkcspc_interf = jcspectra[ip, :, hei_interf[cspsort[range(
627 offhei_interf, nhei_interf + offhei_interf)]]]
636 offhei_interf, nhei_interf + offhei_interf)]]]
628 junkcspc_interf = junkcspc_interf.transpose()
637 junkcspc_interf = junkcspc_interf.transpose()
629 jcspc_interf = junkcspc_interf.sum(axis=1) / nhei_interf
638 jcspc_interf = junkcspc_interf.sum(axis=1) / nhei_interf
630
639
631 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
640 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
632
641
633 median_real = numpy.median(numpy.real(
642 median_real = numpy.median(numpy.real(
634 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
643 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
635 median_imag = numpy.median(numpy.imag(
644 median_imag = numpy.median(numpy.imag(
636 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
645 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
637 junkcspc_interf[comp_mask_prof, :] = numpy.complex(
646 junkcspc_interf[comp_mask_prof, :] = numpy.complex(
638 median_real, median_imag)
647 median_real, median_imag)
639
648
640 for iprof in range(num_prof):
649 for iprof in range(num_prof):
641 ind = numpy.abs(junkcspc_interf[iprof, :]).ravel().argsort()
650 ind = numpy.abs(junkcspc_interf[iprof, :]).ravel().argsort()
642 jcspc_interf[iprof] = junkcspc_interf[iprof,
651 jcspc_interf[iprof] = junkcspc_interf[iprof,
643 ind[nhei_interf / 2]]
652 ind[nhei_interf / 2]]
644
653
645 # Removiendo la Interferencia
654 # Removiendo la Interferencia
646 jcspectra[ip, :, ind_hei] = jcspectra[ip,
655 jcspectra[ip, :, ind_hei] = jcspectra[ip,
647 :, ind_hei] - jcspc_interf
656 :, ind_hei] - jcspc_interf
648
657
649 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
658 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
650 maxid = ListAux.index(max(ListAux))
659 maxid = ListAux.index(max(ListAux))
651
660
652 ind = numpy.array([-2, -1, 1, 2])
661 ind = numpy.array([-2, -1, 1, 2])
653 xx = numpy.zeros([4, 4])
662 xx = numpy.zeros([4, 4])
654
663
655 for id1 in range(4):
664 for id1 in range(4):
656 xx[:, id1] = ind[id1]**numpy.asarray(range(4))
665 xx[:, id1] = ind[id1]**numpy.asarray(range(4))
657
666
658 xx_inv = numpy.linalg.inv(xx)
667 xx_inv = numpy.linalg.inv(xx)
659 xx = xx_inv[:, 0]
668 xx = xx_inv[:, 0]
660
669
661 ind = (ind + maxid + num_mask_prof) % num_mask_prof
670 ind = (ind + maxid + num_mask_prof) % num_mask_prof
662 yy = jcspectra[ip, mask_prof[ind], :]
671 yy = jcspectra[ip, mask_prof[ind], :]
663 jcspectra[ip, mask_prof[maxid], :] = numpy.dot(yy.transpose(), xx)
672 jcspectra[ip, mask_prof[maxid], :] = numpy.dot(yy.transpose(), xx)
664
673
665 # Guardar Resultados
674 # Guardar Resultados
666 self.dataOut.data_spc = jspectra
675 self.dataOut.data_spc = jspectra
667 self.dataOut.data_cspc = jcspectra
676 self.dataOut.data_cspc = jcspectra
668
677
669 return 1
678 return 1
670
679
671 def setRadarFrequency(self, frequency=None):
680 def setRadarFrequency(self, frequency=None):
672
681
673 if frequency != None:
682 if frequency != None:
674 self.dataOut.frequency = frequency
683 self.dataOut.frequency = frequency
675
684
676 return 1
685 return 1
677
686
678 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
687 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
679 # validacion de rango
688 # validacion de rango
680 if minHei == None:
689 if minHei == None:
681 minHei = self.dataOut.heightList[0]
690 minHei = self.dataOut.heightList[0]
682
691
683 if maxHei == None:
692 if maxHei == None:
684 maxHei = self.dataOut.heightList[-1]
693 maxHei = self.dataOut.heightList[-1]
685
694
686 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
695 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
687 print 'minHei: %.2f is out of the heights range' % (minHei)
696 print 'minHei: %.2f is out of the heights range' % (minHei)
688 print 'minHei is setting to %.2f' % (self.dataOut.heightList[0])
697 print 'minHei is setting to %.2f' % (self.dataOut.heightList[0])
689 minHei = self.dataOut.heightList[0]
698 minHei = self.dataOut.heightList[0]
690
699
691 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
700 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
692 print 'maxHei: %.2f is out of the heights range' % (maxHei)
701 print 'maxHei: %.2f is out of the heights range' % (maxHei)
693 print 'maxHei is setting to %.2f' % (self.dataOut.heightList[-1])
702 print 'maxHei is setting to %.2f' % (self.dataOut.heightList[-1])
694 maxHei = self.dataOut.heightList[-1]
703 maxHei = self.dataOut.heightList[-1]
695
704
696 # validacion de velocidades
705 # validacion de velocidades
697 velrange = self.dataOut.getVelRange(1)
706 velrange = self.dataOut.getVelRange(1)
698
707
699 if minVel == None:
708 if minVel == None:
700 minVel = velrange[0]
709 minVel = velrange[0]
701
710
702 if maxVel == None:
711 if maxVel == None:
703 maxVel = velrange[-1]
712 maxVel = velrange[-1]
704
713
705 if (minVel < velrange[0]) or (minVel > maxVel):
714 if (minVel < velrange[0]) or (minVel > maxVel):
706 print 'minVel: %.2f is out of the velocity range' % (minVel)
715 print 'minVel: %.2f is out of the velocity range' % (minVel)
707 print 'minVel is setting to %.2f' % (velrange[0])
716 print 'minVel is setting to %.2f' % (velrange[0])
708 minVel = velrange[0]
717 minVel = velrange[0]
709
718
710 if (maxVel > velrange[-1]) or (maxVel < minVel):
719 if (maxVel > velrange[-1]) or (maxVel < minVel):
711 print 'maxVel: %.2f is out of the velocity range' % (maxVel)
720 print 'maxVel: %.2f is out of the velocity range' % (maxVel)
712 print 'maxVel is setting to %.2f' % (velrange[-1])
721 print 'maxVel is setting to %.2f' % (velrange[-1])
713 maxVel = velrange[-1]
722 maxVel = velrange[-1]
714
723
715 # seleccion de indices para rango
724 # seleccion de indices para rango
716 minIndex = 0
725 minIndex = 0
717 maxIndex = 0
726 maxIndex = 0
718 heights = self.dataOut.heightList
727 heights = self.dataOut.heightList
719
728
720 inda = numpy.where(heights >= minHei)
729 inda = numpy.where(heights >= minHei)
721 indb = numpy.where(heights <= maxHei)
730 indb = numpy.where(heights <= maxHei)
722
731
723 try:
732 try:
724 minIndex = inda[0][0]
733 minIndex = inda[0][0]
725 except:
734 except:
726 minIndex = 0
735 minIndex = 0
727
736
728 try:
737 try:
729 maxIndex = indb[0][-1]
738 maxIndex = indb[0][-1]
730 except:
739 except:
731 maxIndex = len(heights)
740 maxIndex = len(heights)
732
741
733 if (minIndex < 0) or (minIndex > maxIndex):
742 if (minIndex < 0) or (minIndex > maxIndex):
734 raise ValueError, "some value in (%d,%d) is not valid" % (
743 raise ValueError, "some value in (%d,%d) is not valid" % (
735 minIndex, maxIndex)
744 minIndex, maxIndex)
736
745
737 if (maxIndex >= self.dataOut.nHeights):
746 if (maxIndex >= self.dataOut.nHeights):
738 maxIndex = self.dataOut.nHeights - 1
747 maxIndex = self.dataOut.nHeights - 1
739
748
740 # seleccion de indices para velocidades
749 # seleccion de indices para velocidades
741 indminvel = numpy.where(velrange >= minVel)
750 indminvel = numpy.where(velrange >= minVel)
742 indmaxvel = numpy.where(velrange <= maxVel)
751 indmaxvel = numpy.where(velrange <= maxVel)
743 try:
752 try:
744 minIndexVel = indminvel[0][0]
753 minIndexVel = indminvel[0][0]
745 except:
754 except:
746 minIndexVel = 0
755 minIndexVel = 0
747
756
748 try:
757 try:
749 maxIndexVel = indmaxvel[0][-1]
758 maxIndexVel = indmaxvel[0][-1]
750 except:
759 except:
751 maxIndexVel = len(velrange)
760 maxIndexVel = len(velrange)
752
761
753 # seleccion del espectro
762 # seleccion del espectro
754 data_spc = self.dataOut.data_spc[:,
763 data_spc = self.dataOut.data_spc[:,
755 minIndexVel:maxIndexVel + 1, minIndex:maxIndex + 1]
764 minIndexVel:maxIndexVel + 1, minIndex:maxIndex + 1]
756 # estimacion de ruido
765 # estimacion de ruido
757 noise = numpy.zeros(self.dataOut.nChannels)
766 noise = numpy.zeros(self.dataOut.nChannels)
758
767
759 for channel in range(self.dataOut.nChannels):
768 for channel in range(self.dataOut.nChannels):
760 daux = data_spc[channel, :, :]
769 daux = data_spc[channel, :, :]
761 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
770 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
762
771
763 self.dataOut.noise_estimation = noise.copy()
772 self.dataOut.noise_estimation = noise.copy()
764
773
765 return 1
774 return 1
766
775
767
776
768 class IncohInt(Operation):
777 class IncohInt(Operation):
769
778
770 __profIndex = 0
779 __profIndex = 0
771 __withOverapping = False
780 __withOverapping = False
772
781
773 __byTime = False
782 __byTime = False
774 __initime = None
783 __initime = None
775 __lastdatatime = None
784 __lastdatatime = None
776 __integrationtime = None
785 __integrationtime = None
777
786
778 __buffer_spc = None
787 __buffer_spc = None
779 __buffer_cspc = None
788 __buffer_cspc = None
780 __buffer_dc = None
789 __buffer_dc = None
781
790
782 __dataReady = False
791 __dataReady = False
783
792
784 __timeInterval = None
793 __timeInterval = None
785
794
786 n = None
795 n = None
787
796
788 def __init__(self, **kwargs):
797 def __init__(self, **kwargs):
789
798
790 Operation.__init__(self, **kwargs)
799 Operation.__init__(self, **kwargs)
791 # self.isConfig = False
800 # self.isConfig = False
792
801
793 def setup(self, n=None, timeInterval=None, overlapping=False):
802 def setup(self, n=None, timeInterval=None, overlapping=False):
794 """
803 """
795 Set the parameters of the integration class.
804 Set the parameters of the integration class.
796
805
797 Inputs:
806 Inputs:
798
807
799 n : Number of coherent integrations
808 n : Number of coherent integrations
800 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
809 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
801 overlapping :
810 overlapping :
802
811
803 """
812 """
804
813
805 self.__initime = None
814 self.__initime = None
806 self.__lastdatatime = 0
815 self.__lastdatatime = 0
807
816
808 self.__buffer_spc = 0
817 self.__buffer_spc = 0
809 self.__buffer_cspc = 0
818 self.__buffer_cspc = 0
810 self.__buffer_dc = 0
819 self.__buffer_dc = 0
811
820
812 self.__profIndex = 0
821 self.__profIndex = 0
813 self.__dataReady = False
822 self.__dataReady = False
814 self.__byTime = False
823 self.__byTime = False
815
824
816 if n is None and timeInterval is None:
825 if n is None and timeInterval is None:
817 raise ValueError, "n or timeInterval should be specified ..."
826 raise ValueError, "n or timeInterval should be specified ..."
818
827
819 if n is not None:
828 if n is not None:
820 self.n = int(n)
829 self.n = int(n)
821 else:
830 else:
822 # if (type(timeInterval)!=integer) -> change this line
831 # if (type(timeInterval)!=integer) -> change this line
823 self.__integrationtime = int(timeInterval)
832 self.__integrationtime = int(timeInterval)
824 self.n = None
833 self.n = None
825 self.__byTime = True
834 self.__byTime = True
826
835
827 def putData(self, data_spc, data_cspc, data_dc):
836 def putData(self, data_spc, data_cspc, data_dc):
828 """
837 """
829 Add a profile to the __buffer_spc and increase in one the __profileIndex
838 Add a profile to the __buffer_spc and increase in one the __profileIndex
830
839
831 """
840 """
832
841
833 self.__buffer_spc += data_spc
842 self.__buffer_spc += data_spc
834
843
835 if data_cspc is None:
844 if data_cspc is None:
836 self.__buffer_cspc = None
845 self.__buffer_cspc = None
837 else:
846 else:
838 self.__buffer_cspc += data_cspc
847 self.__buffer_cspc += data_cspc
839
848
840 if data_dc is None:
849 if data_dc is None:
841 self.__buffer_dc = None
850 self.__buffer_dc = None
842 else:
851 else:
843 self.__buffer_dc += data_dc
852 self.__buffer_dc += data_dc
844
853
845 self.__profIndex += 1
854 self.__profIndex += 1
846
855
847 return
856 return
848
857
849 def pushData(self):
858 def pushData(self):
850 """
859 """
851 Return the sum of the last profiles and the profiles used in the sum.
860 Return the sum of the last profiles and the profiles used in the sum.
852
861
853 Affected:
862 Affected:
854
863
855 self.__profileIndex
864 self.__profileIndex
856
865
857 """
866 """
858
867
859 data_spc = self.__buffer_spc
868 data_spc = self.__buffer_spc
860 data_cspc = self.__buffer_cspc
869 data_cspc = self.__buffer_cspc
861 data_dc = self.__buffer_dc
870 data_dc = self.__buffer_dc
862 n = self.__profIndex
871 n = self.__profIndex
863
872
864 self.__buffer_spc = 0
873 self.__buffer_spc = 0
865 self.__buffer_cspc = 0
874 self.__buffer_cspc = 0
866 self.__buffer_dc = 0
875 self.__buffer_dc = 0
867 self.__profIndex = 0
876 self.__profIndex = 0
868
877
869 return data_spc, data_cspc, data_dc, n
878 return data_spc, data_cspc, data_dc, n
870
879
871 def byProfiles(self, *args):
880 def byProfiles(self, *args):
872
881
873 self.__dataReady = False
882 self.__dataReady = False
874 avgdata_spc = None
883 avgdata_spc = None
875 avgdata_cspc = None
884 avgdata_cspc = None
876 avgdata_dc = None
885 avgdata_dc = None
877
886
878 self.putData(*args)
887 self.putData(*args)
879
888
880 if self.__profIndex == self.n:
889 if self.__profIndex == self.n:
881
890
882 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
891 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
883 self.n = n
892 self.n = n
884 self.__dataReady = True
893 self.__dataReady = True
885
894
886 return avgdata_spc, avgdata_cspc, avgdata_dc
895 return avgdata_spc, avgdata_cspc, avgdata_dc
887
896
888 def byTime(self, datatime, *args):
897 def byTime(self, datatime, *args):
889
898
890 self.__dataReady = False
899 self.__dataReady = False
891 avgdata_spc = None
900 avgdata_spc = None
892 avgdata_cspc = None
901 avgdata_cspc = None
893 avgdata_dc = None
902 avgdata_dc = None
894
903
895 self.putData(*args)
904 self.putData(*args)
896
905
897 if (datatime - self.__initime) >= self.__integrationtime:
906 if (datatime - self.__initime) >= self.__integrationtime:
898 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
907 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
899 self.n = n
908 self.n = n
900 self.__dataReady = True
909 self.__dataReady = True
901
910
902 return avgdata_spc, avgdata_cspc, avgdata_dc
911 return avgdata_spc, avgdata_cspc, avgdata_dc
903
912
904 def integrate(self, datatime, *args):
913 def integrate(self, datatime, *args):
905
914
906 if self.__profIndex == 0:
915 if self.__profIndex == 0:
907 self.__initime = datatime
916 self.__initime = datatime
908
917
909 if self.__byTime:
918 if self.__byTime:
910 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
919 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
911 datatime, *args)
920 datatime, *args)
912 else:
921 else:
913 avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
922 avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
914
923
915 if not self.__dataReady:
924 if not self.__dataReady:
916 return None, None, None, None
925 return None, None, None, None
917
926
918 return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
927 return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
919
928
920 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
929 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
921 if n == 1:
930 if n == 1:
922 return
931 return
923
932
924 dataOut.flagNoData = True
933 dataOut.flagNoData = True
925
934
926 if not self.isConfig:
935 if not self.isConfig:
927 self.setup(n, timeInterval, overlapping)
936 self.setup(n, timeInterval, overlapping)
928 self.isConfig = True
937 self.isConfig = True
929
938
930 avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
939 avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
931 dataOut.data_spc,
940 dataOut.data_spc,
932 dataOut.data_cspc,
941 dataOut.data_cspc,
933 dataOut.data_dc)
942 dataOut.data_dc)
934
943
935 if self.__dataReady:
944 if self.__dataReady:
936
945
937 dataOut.data_spc = avgdata_spc
946 dataOut.data_spc = avgdata_spc
938 dataOut.data_cspc = avgdata_cspc
947 dataOut.data_cspc = avgdata_cspc
939 dataOut.data_dc = avgdata_dc
948 dataOut.data_dc = avgdata_dc
940
949
941 dataOut.nIncohInt *= self.n
950 dataOut.nIncohInt *= self.n
942 dataOut.utctime = avgdatatime
951 dataOut.utctime = avgdatatime
943 dataOut.flagNoData = False
952 dataOut.flagNoData = False
General Comments 0
You need to be logged in to leave comments. Login now