##// END OF EJS Templates
-Function to interpolate corrupted profiles
Julio Valdez -
r836:0804a6804e7a
parent child
Show More
@@ -1,632 +1,632
1 1 '''
2 2 Created on Jul 2, 2014
3 3
4 4 @author: roj-idl71
5 5 '''
6 6
7 7 import numpy
8 8
9 9 from jroIO_base import LOCALTIME, JRODataReader, JRODataWriter
10 10 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
11 11 from schainpy.model.data.jroheaderIO import PROCFLAG, BasicHeader, SystemHeader, RadarControllerHeader, ProcessingHeader
12 12 from schainpy.model.data.jrodata import Voltage
13 13 # from _sha import blocksize
14 14
15 15 class VoltageReader(JRODataReader, ProcessingUnit):
16 16 """
17 17 Esta clase permite leer datos de voltage desde archivos en formato rawdata (.r). La lectura
18 18 de los datos siempre se realiza por bloques. Los datos leidos (array de 3 dimensiones:
19 19 perfiles*alturas*canales) son almacenados en la variable "buffer".
20 20
21 21 perfiles * alturas * canales
22 22
23 23 Esta clase contiene instancias (objetos) de las clases BasicHeader, SystemHeader,
24 24 RadarControllerHeader y Voltage. Los tres primeros se usan para almacenar informacion de la
25 25 cabecera de datos (metadata), y el cuarto (Voltage) para obtener y almacenar un perfil de
26 26 datos desde el "buffer" cada vez que se ejecute el metodo "getData".
27 27
28 28 Example:
29 29
30 30 dpath = "/home/myuser/data"
31 31
32 32 startTime = datetime.datetime(2010,1,20,0,0,0,0,0,0)
33 33
34 34 endTime = datetime.datetime(2010,1,21,23,59,59,0,0,0)
35 35
36 36 readerObj = VoltageReader()
37 37
38 38 readerObj.setup(dpath, startTime, endTime)
39 39
40 40 while(True):
41 41
42 42 #to get one profile
43 43 profile = readerObj.getData()
44 44
45 45 #print the profile
46 46 print profile
47 47
48 48 #If you want to see all datablock
49 49 print readerObj.datablock
50 50
51 51 if readerObj.flagNoMoreFiles:
52 52 break
53 53
54 54 """
55 55
56 56 ext = ".r"
57 57
58 58 optchar = "D"
59 59 dataOut = None
60 60
61 61 def __init__(self):
62 62 """
63 63 Inicializador de la clase VoltageReader para la lectura de datos de voltage.
64 64
65 65 Input:
66 66 dataOut : Objeto de la clase Voltage. Este objeto sera utilizado para
67 67 almacenar un perfil de datos cada vez que se haga un requerimiento
68 68 (getData). El perfil sera obtenido a partir del buffer de datos,
69 69 si el buffer esta vacio se hara un nuevo proceso de lectura de un
70 70 bloque de datos.
71 71 Si este parametro no es pasado se creara uno internamente.
72 72
73 73 Variables afectadas:
74 74 self.dataOut
75 75
76 76 Return:
77 77 None
78 78 """
79 79
80 80 ProcessingUnit.__init__(self)
81 81
82 82 self.isConfig = False
83 83
84 84 self.datablock = None
85 85
86 86 self.utc = 0
87 87
88 88 self.ext = ".r"
89 89
90 90 self.optchar = "D"
91 91
92 92 self.basicHeaderObj = BasicHeader(LOCALTIME)
93 93
94 94 self.systemHeaderObj = SystemHeader()
95 95
96 96 self.radarControllerHeaderObj = RadarControllerHeader()
97 97
98 98 self.processingHeaderObj = ProcessingHeader()
99 99
100 100 self.online = 0
101 101
102 102 self.fp = None
103 103
104 104 self.idFile = None
105 105
106 106 self.dtype = None
107 107
108 108 self.fileSizeByHeader = None
109 109
110 110 self.filenameList = []
111 111
112 112 self.filename = None
113 113
114 114 self.fileSize = None
115 115
116 116 self.firstHeaderSize = 0
117 117
118 118 self.basicHeaderSize = 24
119 119
120 120 self.pathList = []
121 121
122 122 self.filenameList = []
123 123
124 124 self.lastUTTime = 0
125 125
126 126 self.maxTimeStep = 30
127 127
128 128 self.flagNoMoreFiles = 0
129 129
130 130 self.set = 0
131 131
132 132 self.path = None
133 133
134 134 self.profileIndex = 2**32-1
135 135
136 136 self.delay = 3 #seconds
137 137
138 138 self.nTries = 3 #quantity tries
139 139
140 140 self.nFiles = 3 #number of files for searching
141 141
142 142 self.nReadBlocks = 0
143 143
144 144 self.flagIsNewFile = 1
145 145
146 146 self.__isFirstTimeOnline = 1
147 147
148 148 # self.ippSeconds = 0
149 149
150 150 self.flagDiscontinuousBlock = 0
151 151
152 152 self.flagIsNewBlock = 0
153 153
154 154 self.nTotalBlocks = 0
155 155
156 156 self.blocksize = 0
157 157
158 158 self.dataOut = self.createObjByDefault()
159 159
160 160 self.nTxs = 1
161 161
162 162 self.txIndex = 0
163 163
164 164 def createObjByDefault(self):
165 165
166 166 dataObj = Voltage()
167 167
168 168 return dataObj
169 169
170 170 def __hasNotDataInBuffer(self):
171 171
172 172 if self.profileIndex >= self.processingHeaderObj.profilesPerBlock*self.nTxs:
173 173 return 1
174 174
175 175 return 0
176 176
177 177
178 178 def getBlockDimension(self):
179 179 """
180 180 Obtiene la cantidad de puntos a leer por cada bloque de datos
181 181
182 182 Affected:
183 183 self.blocksize
184 184
185 185 Return:
186 186 None
187 187 """
188 188 pts2read = self.processingHeaderObj.profilesPerBlock * self.processingHeaderObj.nHeights * self.systemHeaderObj.nChannels
189 189 self.blocksize = pts2read
190 190
191 191
192 192 def readBlock(self):
193 193 """
194 194 readBlock lee el bloque de datos desde la posicion actual del puntero del archivo
195 195 (self.fp) y actualiza todos los parametros relacionados al bloque de datos
196 196 (metadata + data). La data leida es almacenada en el buffer y el contador del buffer
197 197 es seteado a 0
198 198
199 199 Inputs:
200 200 None
201 201
202 202 Return:
203 203 None
204 204
205 205 Affected:
206 206 self.profileIndex
207 207 self.datablock
208 208 self.flagIsNewFile
209 209 self.flagIsNewBlock
210 210 self.nTotalBlocks
211 211
212 212 Exceptions:
213 213 Si un bloque leido no es un bloque valido
214 214 """
215 215 current_pointer_location = self.fp.tell()
216 216 junk = numpy.fromfile( self.fp, self.dtype, self.blocksize )
217 217
218 218 try:
219 219 junk = junk.reshape( (self.processingHeaderObj.profilesPerBlock, self.processingHeaderObj.nHeights, self.systemHeaderObj.nChannels) )
220 220 except:
221 221 #print "The read block (%3d) has not enough data" %self.nReadBlocks
222 222
223 223 if self.waitDataBlock(pointer_location=current_pointer_location):
224 224 junk = numpy.fromfile( self.fp, self.dtype, self.blocksize )
225 225 junk = junk.reshape( (self.processingHeaderObj.profilesPerBlock, self.processingHeaderObj.nHeights, self.systemHeaderObj.nChannels) )
226 226 # return 0
227 227
228 228 #Dimensions : nChannels, nProfiles, nSamples
229 229
230 230 junk = numpy.transpose(junk, (2,0,1))
231 231 self.datablock = junk['real'] + junk['imag']*1j
232 232
233 233 self.profileIndex = 0
234 234
235 235 self.flagIsNewFile = 0
236 236 self.flagIsNewBlock = 1
237 237
238 238 self.nTotalBlocks += 1
239 239 self.nReadBlocks += 1
240 240
241 241 return 1
242 242
243 243 def getFirstHeader(self):
244 244
245 245 self.getBasicHeader()
246 246
247 247 self.dataOut.systemHeaderObj = self.systemHeaderObj.copy()
248 248
249 249 self.dataOut.radarControllerHeaderObj = self.radarControllerHeaderObj.copy()
250 250
251 251 if self.nTxs > 1:
252 252 self.dataOut.radarControllerHeaderObj.ippSeconds = self.radarControllerHeaderObj.ippSeconds/self.nTxs
253 253
254 254 #Time interval and code are propierties of dataOut. Its value depends of radarControllerHeaderObj.
255 255
256 256 # self.dataOut.timeInterval = self.radarControllerHeaderObj.ippSeconds * self.processingHeaderObj.nCohInt
257 257 #
258 258 # if self.radarControllerHeaderObj.code is not None:
259 259 #
260 260 # self.dataOut.nCode = self.radarControllerHeaderObj.nCode
261 261 #
262 262 # self.dataOut.nBaud = self.radarControllerHeaderObj.nBaud
263 263 #
264 264 # self.dataOut.code = self.radarControllerHeaderObj.code
265 265
266 266 self.dataOut.dtype = self.dtype
267 267
268 268 self.dataOut.nProfiles = self.processingHeaderObj.profilesPerBlock
269 269
270 270 self.dataOut.heightList = numpy.arange(self.processingHeaderObj.nHeights) *self.processingHeaderObj.deltaHeight + self.processingHeaderObj.firstHeight
271 271
272 272 self.dataOut.channelList = range(self.systemHeaderObj.nChannels)
273 273
274 274 self.dataOut.nCohInt = self.processingHeaderObj.nCohInt
275 275
276 276 self.dataOut.flagDecodeData = self.processingHeaderObj.flag_decode #asumo q la data no esta decodificada
277 277
278 278 self.dataOut.flagDeflipData = self.processingHeaderObj.flag_deflip #asumo q la data no esta sin flip
279 279
280 280 self.dataOut.flagShiftFFT = self.processingHeaderObj.shif_fft
281 281
282 282 def reshapeData(self):
283 283
284 284 if self.nTxs < 0:
285 285 return
286 286
287 287 if self.nTxs == 1:
288 288 return
289 289
290 290 if self.nTxs < 1 and self.processingHeaderObj.profilesPerBlock % (1./self.nTxs) != 0:
291 291 raise ValueError, "1./nTxs (=%f), should be a multiple of nProfiles (=%d)" %(1./self.nTxs, self.processingHeaderObj.profilesPerBlock)
292 292
293 293 if self.nTxs > 1 and self.processingHeaderObj.nHeights % self.nTxs != 0:
294 294 raise ValueError, "nTxs (=%d), should be a multiple of nHeights (=%d)" %(self.nTxs, self.processingHeaderObj.nHeights)
295 295
296 296 self.datablock = self.datablock.reshape((self.systemHeaderObj.nChannels, self.processingHeaderObj.profilesPerBlock*self.nTxs, self.processingHeaderObj.nHeights/self.nTxs))
297 297
298 298 self.dataOut.nProfiles = self.processingHeaderObj.profilesPerBlock*self.nTxs
299 299 self.dataOut.heightList = numpy.arange(self.processingHeaderObj.nHeights/self.nTxs) *self.processingHeaderObj.deltaHeight + self.processingHeaderObj.firstHeight
300 300 self.dataOut.radarControllerHeaderObj.ippSeconds = self.radarControllerHeaderObj.ippSeconds/self.nTxs
301 301
302 302 return
303 303
304 304 def getData(self):
305 305 """
306 306 getData obtiene una unidad de datos del buffer de lectura, un perfil, y la copia al objeto self.dataOut
307 307 del tipo "Voltage" con todos los parametros asociados a este (metadata). cuando no hay datos
308 308 en el buffer de lectura es necesario hacer una nueva lectura de los bloques de datos usando
309 309 "readNextBlock"
310 310
311 311 Ademas incrementa el contador del buffer "self.profileIndex" en 1.
312 312
313 313 Return:
314 314
315 315 Si el flag self.getByBlock ha sido seteado el bloque completo es copiado a self.dataOut y el self.profileIndex
316 316 es igual al total de perfiles leidos desde el archivo.
317 317
318 318 Si self.getByBlock == False:
319 319
320 320 self.dataOut.data = buffer[:, thisProfile, :]
321 321
322 322 shape = [nChannels, nHeis]
323 323
324 324 Si self.getByBlock == True:
325 325
326 326 self.dataOut.data = buffer[:, :, :]
327 327
328 328 shape = [nChannels, nProfiles, nHeis]
329 329
330 330 Variables afectadas:
331 331 self.dataOut
332 332 self.profileIndex
333 333
334 334 Affected:
335 335 self.dataOut
336 336 self.profileIndex
337 337 self.flagDiscontinuousBlock
338 338 self.flagIsNewBlock
339 339 """
340 340
341 341 if self.flagNoMoreFiles:
342 342 self.dataOut.flagNoData = True
343 343 print 'Process finished'
344 344 return 0
345 345
346 346 self.flagDiscontinuousBlock = 0
347 347 self.flagIsNewBlock = 0
348 348
349 349 if self.__hasNotDataInBuffer():
350 350
351 351 if not( self.readNextBlock() ):
352 352 return 0
353 353
354 354 self.getFirstHeader()
355 355
356 356 self.reshapeData()
357 357
358 358 if self.datablock is None:
359 359 self.dataOut.flagNoData = True
360 360 return 0
361 361
362 362 if not self.getByBlock:
363 363
364 364 """
365 365 Return profile by profile
366 366
367 367 If nTxs > 1 then one profile is divided by nTxs and number of total
368 368 blocks is increased by nTxs (nProfiles *= nTxs)
369 369 """
370 370 self.dataOut.flagDataAsBlock = False
371 371 self.dataOut.data = self.datablock[:,self.profileIndex,:]
372 372 self.dataOut.profileIndex = self.profileIndex
373 373
374 374 self.profileIndex += 1
375 375
376 376 # elif self.selBlocksize==None or self.selBlocksize==self.dataOut.nProfiles:
377 377 # """
378 378 # Return all block
379 379 # """
380 380 # self.dataOut.flagDataAsBlock = True
381 381 # self.dataOut.data = self.datablock
382 382 # self.dataOut.profileIndex = self.dataOut.nProfiles - 1
383 383 #
384 384 # self.profileIndex = self.dataOut.nProfiles
385 385
386 386 else:
387 387 """
388 388 Return a block
389 389 """
390 390 if self.selBlocksize == None: self.selBlocksize = self.dataOut.nProfiles
391 391 if self.selBlocktime != None: self.selBlocksize = int(self.dataOut.nProfiles*round(self.selBlocktime/(self.dataOut.ippSeconds*self.dataOut.nProfiles)))
392 392
393 393 self.dataOut.data = self.datablock[:,self.profileIndex:self.profileIndex+self.selBlocksize,:]
394 394 self.profileIndex += self.selBlocksize
395 395 datasize = self.dataOut.data.shape[1]
396 396
397 397 if datasize < self.selBlocksize:
398 398 buffer = numpy.zeros((self.dataOut.data.shape[0],self.selBlocksize,self.dataOut.data.shape[2]), dtype = 'complex')
399 399 buffer[:,:datasize,:] = self.dataOut.data
400 400
401 401 while datasize < self.selBlocksize: #Not enough profiles to fill the block
402 402 if not( self.readNextBlock() ):
403 403 return 0
404 404 self.getFirstHeader()
405 405 self.reshapeData()
406 406 if self.datablock is None:
407 407 self.dataOut.flagNoData = True
408 408 return 0
409 409 #stack data
410 410 blockIndex = self.selBlocksize - datasize
411 411 datablock1 = self.datablock[:,:blockIndex,:]
412 412
413 413 buffer[:,datasize:datasize+datablock1.shape[1],:] = datablock1
414 414 datasize += datablock1.shape[1]
415 415
416 416 self.dataOut.data = buffer
417 417 self.profileIndex = blockIndex
418 418
419 419 self.dataOut.flagDataAsBlock = True
420 self.dataOut.nProfiles = self.selBlocksize
420 self.dataOut.nProfiles = self.dataOut.data.shape[1]
421 421
422 422 self.dataOut.flagNoData = False
423 423
424 424 self.getBasicHeader()
425 425
426 426 self.dataOut.realtime = self.online
427 427
428 428 return self.dataOut.data
429 429
430 430 class VoltageWriter(JRODataWriter, Operation):
431 431 """
432 432 Esta clase permite escribir datos de voltajes a archivos procesados (.r). La escritura
433 433 de los datos siempre se realiza por bloques.
434 434 """
435 435
436 436 ext = ".r"
437 437
438 438 optchar = "D"
439 439
440 440 shapeBuffer = None
441 441
442 442
443 443 def __init__(self):
444 444 """
445 445 Inicializador de la clase VoltageWriter para la escritura de datos de espectros.
446 446
447 447 Affected:
448 448 self.dataOut
449 449
450 450 Return: None
451 451 """
452 452 Operation.__init__(self)
453 453
454 454 self.nTotalBlocks = 0
455 455
456 456 self.profileIndex = 0
457 457
458 458 self.isConfig = False
459 459
460 460 self.fp = None
461 461
462 462 self.flagIsNewFile = 1
463 463
464 464 self.blockIndex = 0
465 465
466 466 self.flagIsNewBlock = 0
467 467
468 468 self.setFile = None
469 469
470 470 self.dtype = None
471 471
472 472 self.path = None
473 473
474 474 self.filename = None
475 475
476 476 self.basicHeaderObj = BasicHeader(LOCALTIME)
477 477
478 478 self.systemHeaderObj = SystemHeader()
479 479
480 480 self.radarControllerHeaderObj = RadarControllerHeader()
481 481
482 482 self.processingHeaderObj = ProcessingHeader()
483 483
484 484 def hasAllDataInBuffer(self):
485 485 if self.profileIndex >= self.processingHeaderObj.profilesPerBlock:
486 486 return 1
487 487 return 0
488 488
489 489
490 490 def setBlockDimension(self):
491 491 """
492 492 Obtiene las formas dimensionales del los subbloques de datos que componen un bloque
493 493
494 494 Affected:
495 495 self.shape_spc_Buffer
496 496 self.shape_cspc_Buffer
497 497 self.shape_dc_Buffer
498 498
499 499 Return: None
500 500 """
501 501 self.shapeBuffer = (self.processingHeaderObj.profilesPerBlock,
502 502 self.processingHeaderObj.nHeights,
503 503 self.systemHeaderObj.nChannels)
504 504
505 505 self.datablock = numpy.zeros((self.systemHeaderObj.nChannels,
506 506 self.processingHeaderObj.profilesPerBlock,
507 507 self.processingHeaderObj.nHeights),
508 508 dtype=numpy.dtype('complex64'))
509 509
510 510 def writeBlock(self):
511 511 """
512 512 Escribe el buffer en el file designado
513 513
514 514 Affected:
515 515 self.profileIndex
516 516 self.flagIsNewFile
517 517 self.flagIsNewBlock
518 518 self.nTotalBlocks
519 519 self.blockIndex
520 520
521 521 Return: None
522 522 """
523 523 data = numpy.zeros( self.shapeBuffer, self.dtype )
524 524
525 525 junk = numpy.transpose(self.datablock, (1,2,0))
526 526
527 527 data['real'] = junk.real
528 528 data['imag'] = junk.imag
529 529
530 530 data = data.reshape( (-1) )
531 531
532 532 data.tofile( self.fp )
533 533
534 534 self.datablock.fill(0)
535 535
536 536 self.profileIndex = 0
537 537 self.flagIsNewFile = 0
538 538 self.flagIsNewBlock = 1
539 539
540 540 self.blockIndex += 1
541 541 self.nTotalBlocks += 1
542 542
543 543 # print "[Writing] Block = %04d" %self.blockIndex
544 544
545 545 def putData(self):
546 546 """
547 547 Setea un bloque de datos y luego los escribe en un file
548 548
549 549 Affected:
550 550 self.flagIsNewBlock
551 551 self.profileIndex
552 552
553 553 Return:
554 554 0 : Si no hay data o no hay mas files que puedan escribirse
555 555 1 : Si se escribio la data de un bloque en un file
556 556 """
557 557 if self.dataOut.flagNoData:
558 558 return 0
559 559
560 560 self.flagIsNewBlock = 0
561 561
562 562 if self.dataOut.flagDiscontinuousBlock:
563 563 self.datablock.fill(0)
564 564 self.profileIndex = 0
565 565 self.setNextFile()
566 566
567 567 if self.profileIndex == 0:
568 568 self.setBasicHeader()
569 569
570 570 self.datablock[:,self.profileIndex,:] = self.dataOut.data
571 571
572 572 self.profileIndex += 1
573 573
574 574 if self.hasAllDataInBuffer():
575 575 #if self.flagIsNewFile:
576 576 self.writeNextBlock()
577 577 # self.setFirstHeader()
578 578
579 579 return 1
580 580
581 581 def __getBlockSize(self):
582 582 '''
583 583 Este metodos determina el cantidad de bytes para un bloque de datos de tipo Voltage
584 584 '''
585 585
586 586 dtype_width = self.getDtypeWidth()
587 587
588 588 blocksize = int(self.dataOut.nHeights * self.dataOut.nChannels * self.profilesPerBlock * dtype_width * 2)
589 589
590 590 return blocksize
591 591
592 592 def setFirstHeader(self):
593 593
594 594 """
595 595 Obtiene una copia del First Header
596 596
597 597 Affected:
598 598 self.systemHeaderObj
599 599 self.radarControllerHeaderObj
600 600 self.dtype
601 601
602 602 Return:
603 603 None
604 604 """
605 605
606 606 self.systemHeaderObj = self.dataOut.systemHeaderObj.copy()
607 607 self.systemHeaderObj.nChannels = self.dataOut.nChannels
608 608 self.radarControllerHeaderObj = self.dataOut.radarControllerHeaderObj.copy()
609 609
610 610 self.processingHeaderObj.dtype = 0 # Voltage
611 611 self.processingHeaderObj.blockSize = self.__getBlockSize()
612 612 self.processingHeaderObj.profilesPerBlock = self.profilesPerBlock
613 613 self.processingHeaderObj.dataBlocksPerFile = self.blocksPerFile
614 614 self.processingHeaderObj.nWindows = 1 #podria ser 1 o self.dataOut.processingHeaderObj.nWindows
615 615 self.processingHeaderObj.nCohInt = self.dataOut.nCohInt
616 616 self.processingHeaderObj.nIncohInt = 1 # Cuando la data de origen es de tipo Voltage
617 617 self.processingHeaderObj.totalSpectra = 0 # Cuando la data de origen es de tipo Voltage
618 618
619 619 if self.dataOut.code is not None:
620 620 self.processingHeaderObj.code = self.dataOut.code
621 621 self.processingHeaderObj.nCode = self.dataOut.nCode
622 622 self.processingHeaderObj.nBaud = self.dataOut.nBaud
623 623
624 624 if self.processingHeaderObj.nWindows != 0:
625 625 self.processingHeaderObj.firstHeight = self.dataOut.heightList[0]
626 626 self.processingHeaderObj.deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
627 627 self.processingHeaderObj.nHeights = self.dataOut.nHeights
628 628 self.processingHeaderObj.samplesWin = self.dataOut.nHeights
629 629
630 630 self.processingHeaderObj.processFlags = self.getProcessFlags()
631 631
632 632 self.setBasicHeader() No newline at end of file
@@ -1,1262 +1,1276
1 1 import sys
2 2 import numpy
3 3
4 4 from jroproc_base import ProcessingUnit, Operation
5 5 from schainpy.model.data.jrodata import Voltage
6 6
7 7 class VoltageProc(ProcessingUnit):
8 8
9 9
10 10 def __init__(self):
11 11
12 12 ProcessingUnit.__init__(self)
13 13
14 14 # self.objectDict = {}
15 15 self.dataOut = Voltage()
16 16 self.flip = 1
17 17
18 18 def run(self):
19 19 if self.dataIn.type == 'AMISR':
20 20 self.__updateObjFromAmisrInput()
21 21
22 22 if self.dataIn.type == 'Voltage':
23 23 self.dataOut.copy(self.dataIn)
24 24
25 25 # self.dataOut.copy(self.dataIn)
26 26
27 27 def __updateObjFromAmisrInput(self):
28 28
29 29 self.dataOut.timeZone = self.dataIn.timeZone
30 30 self.dataOut.dstFlag = self.dataIn.dstFlag
31 31 self.dataOut.errorCount = self.dataIn.errorCount
32 32 self.dataOut.useLocalTime = self.dataIn.useLocalTime
33 33
34 34 self.dataOut.flagNoData = self.dataIn.flagNoData
35 35 self.dataOut.data = self.dataIn.data
36 36 self.dataOut.utctime = self.dataIn.utctime
37 37 self.dataOut.channelList = self.dataIn.channelList
38 38 # self.dataOut.timeInterval = self.dataIn.timeInterval
39 39 self.dataOut.heightList = self.dataIn.heightList
40 40 self.dataOut.nProfiles = self.dataIn.nProfiles
41 41
42 42 self.dataOut.nCohInt = self.dataIn.nCohInt
43 43 self.dataOut.ippSeconds = self.dataIn.ippSeconds
44 44 self.dataOut.frequency = self.dataIn.frequency
45 45
46 46 self.dataOut.azimuth = self.dataIn.azimuth
47 47 self.dataOut.zenith = self.dataIn.zenith
48 48
49 49 self.dataOut.beam.codeList = self.dataIn.beam.codeList
50 50 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
51 51 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
52 52 #
53 53 # pass#
54 54 #
55 55 # def init(self):
56 56 #
57 57 #
58 58 # if self.dataIn.type == 'AMISR':
59 59 # self.__updateObjFromAmisrInput()
60 60 #
61 61 # if self.dataIn.type == 'Voltage':
62 62 # self.dataOut.copy(self.dataIn)
63 63 # # No necesita copiar en cada init() los atributos de dataIn
64 64 # # la copia deberia hacerse por cada nuevo bloque de datos
65 65
66 66 def selectChannels(self, channelList):
67 67
68 68 channelIndexList = []
69 69
70 70 for channel in channelList:
71 71 if channel not in self.dataOut.channelList:
72 72 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
73 73
74 74 index = self.dataOut.channelList.index(channel)
75 75 channelIndexList.append(index)
76 76
77 77 self.selectChannelsByIndex(channelIndexList)
78 78
79 79 def selectChannelsByIndex(self, channelIndexList):
80 80 """
81 81 Selecciona un bloque de datos en base a canales segun el channelIndexList
82 82
83 83 Input:
84 84 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
85 85
86 86 Affected:
87 87 self.dataOut.data
88 88 self.dataOut.channelIndexList
89 89 self.dataOut.nChannels
90 90 self.dataOut.m_ProcessingHeader.totalSpectra
91 91 self.dataOut.systemHeaderObj.numChannels
92 92 self.dataOut.m_ProcessingHeader.blockSize
93 93
94 94 Return:
95 95 None
96 96 """
97 97
98 98 for channelIndex in channelIndexList:
99 99 if channelIndex not in self.dataOut.channelIndexList:
100 100 print channelIndexList
101 101 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
102 102
103 103 if self.dataOut.flagDataAsBlock:
104 104 """
105 105 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
106 106 """
107 107 data = self.dataOut.data[channelIndexList,:,:]
108 108 else:
109 109 data = self.dataOut.data[channelIndexList,:]
110 110
111 111 self.dataOut.data = data
112 112 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
113 113 # self.dataOut.nChannels = nChannels
114 114
115 115 return 1
116 116
117 117 def selectHeights(self, minHei=None, maxHei=None):
118 118 """
119 119 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
120 120 minHei <= height <= maxHei
121 121
122 122 Input:
123 123 minHei : valor minimo de altura a considerar
124 124 maxHei : valor maximo de altura a considerar
125 125
126 126 Affected:
127 127 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
128 128
129 129 Return:
130 130 1 si el metodo se ejecuto con exito caso contrario devuelve 0
131 131 """
132 132
133 133 if minHei == None:
134 134 minHei = self.dataOut.heightList[0]
135 135
136 136 if maxHei == None:
137 137 maxHei = self.dataOut.heightList[-1]
138 138
139 139 if (minHei < self.dataOut.heightList[0]):
140 140 minHei = self.dataOut.heightList[0]
141 141
142 142 if (maxHei > self.dataOut.heightList[-1]):
143 143 maxHei = self.dataOut.heightList[-1]
144 144
145 145 minIndex = 0
146 146 maxIndex = 0
147 147 heights = self.dataOut.heightList
148 148
149 149 inda = numpy.where(heights >= minHei)
150 150 indb = numpy.where(heights <= maxHei)
151 151
152 152 try:
153 153 minIndex = inda[0][0]
154 154 except:
155 155 minIndex = 0
156 156
157 157 try:
158 158 maxIndex = indb[0][-1]
159 159 except:
160 160 maxIndex = len(heights)
161 161
162 162 self.selectHeightsByIndex(minIndex, maxIndex)
163 163
164 164 return 1
165 165
166 166
167 167 def selectHeightsByIndex(self, minIndex, maxIndex):
168 168 """
169 169 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
170 170 minIndex <= index <= maxIndex
171 171
172 172 Input:
173 173 minIndex : valor de indice minimo de altura a considerar
174 174 maxIndex : valor de indice maximo de altura a considerar
175 175
176 176 Affected:
177 177 self.dataOut.data
178 178 self.dataOut.heightList
179 179
180 180 Return:
181 181 1 si el metodo se ejecuto con exito caso contrario devuelve 0
182 182 """
183 183
184 184 if (minIndex < 0) or (minIndex > maxIndex):
185 185 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
186 186
187 187 if (maxIndex >= self.dataOut.nHeights):
188 188 maxIndex = self.dataOut.nHeights
189 189
190 190 #voltage
191 191 if self.dataOut.flagDataAsBlock:
192 192 """
193 193 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
194 194 """
195 195 data = self.dataOut.data[:,:, minIndex:maxIndex]
196 196 else:
197 197 data = self.dataOut.data[:, minIndex:maxIndex]
198 198
199 199 # firstHeight = self.dataOut.heightList[minIndex]
200 200
201 201 self.dataOut.data = data
202 202 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
203 203
204 204 if self.dataOut.nHeights <= 1:
205 205 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
206 206
207 207 return 1
208 208
209 209
210 210 def filterByHeights(self, window):
211 211
212 212 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
213 213
214 214 if window == None:
215 215 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
216 216
217 217 newdelta = deltaHeight * window
218 218 r = self.dataOut.nHeights % window
219 219 newheights = (self.dataOut.nHeights-r)/window
220 220
221 221 if newheights <= 1:
222 222 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
223 223
224 224 if self.dataOut.flagDataAsBlock:
225 225 """
226 226 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
227 227 """
228 228 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
229 229 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
230 230 buffer = numpy.sum(buffer,3)
231 231
232 232 else:
233 233 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
234 234 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
235 235 buffer = numpy.sum(buffer,2)
236 236
237 237 self.dataOut.data = buffer
238 238 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
239 239 self.dataOut.windowOfFilter = window
240 240
241 241 def setH0(self, h0, deltaHeight = None):
242 242
243 243 if not deltaHeight:
244 244 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
245 245
246 246 nHeights = self.dataOut.nHeights
247 247
248 248 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
249 249
250 250 self.dataOut.heightList = newHeiRange
251 251
252 252 def deFlip(self, channelList = []):
253 253
254 254 data = self.dataOut.data.copy()
255 255
256 256 if self.dataOut.flagDataAsBlock:
257 257 flip = self.flip
258 258 profileList = range(self.dataOut.nProfiles)
259 259
260 260 if not channelList:
261 261 for thisProfile in profileList:
262 262 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
263 263 flip *= -1.0
264 264 else:
265 265 for thisChannel in channelList:
266 266 if thisChannel not in self.dataOut.channelList:
267 267 continue
268 268
269 269 for thisProfile in profileList:
270 270 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
271 271 flip *= -1.0
272 272
273 273 self.flip = flip
274 274
275 275 else:
276 276 if not channelList:
277 277 data[:,:] = data[:,:]*self.flip
278 278 else:
279 279 for thisChannel in channelList:
280 280 if thisChannel not in self.dataOut.channelList:
281 281 continue
282 282
283 283 data[thisChannel,:] = data[thisChannel,:]*self.flip
284 284
285 285 self.flip *= -1.
286 286
287 287 self.dataOut.data = data
288 288
289 289 def setRadarFrequency(self, frequency=None):
290 290
291 291 if frequency != None:
292 292 self.dataOut.frequency = frequency
293 293
294 294 return 1
295 295
296 def interpolateHeights(self, topLim, botLim):
297 #69 al 72 para julia
298 #82-84 para meteoros
299 if len(numpy.shape(self.dataOut.data))==2:
300 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
301 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
302 self.dataOut.data[:,botLim:limSup+1] = sampInterp
303 else:
304 sampInterp = (self.dataOut.data[:,:,botLim-1] + self.dataOut.data[:,:,topLim+1])/2
305 nInt = topLim - botLim + 1
306 for i in range(nInt):
307 self.dataOut.data[:,:,botLim+i] = sampInterp
308 # import collections
309
296 310 class CohInt(Operation):
297 311
298 312 isConfig = False
299 313
300 314 __profIndex = 0
301 315 __withOverapping = False
302 316
303 317 __byTime = False
304 318 __initime = None
305 319 __lastdatatime = None
306 320 __integrationtime = None
307 321
308 322 __buffer = None
309 323
310 324 __dataReady = False
311 325
312 326 n = None
313 327
314 328
315 329 def __init__(self):
316 330
317 331 Operation.__init__(self)
318 332
319 333 # self.isConfig = False
320 334
321 335 def setup(self, n=None, timeInterval=None, overlapping=False, byblock=False):
322 336 """
323 337 Set the parameters of the integration class.
324 338
325 339 Inputs:
326 340
327 341 n : Number of coherent integrations
328 342 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
329 343 overlapping :
330 344
331 345 """
332 346
333 347 self.__initime = None
334 348 self.__lastdatatime = 0
335 349 self.__buffer = None
336 350 self.__dataReady = False
337 351 self.byblock = byblock
338 352
339 353 if n == None and timeInterval == None:
340 354 raise ValueError, "n or timeInterval should be specified ..."
341 355
342 356 if n != None:
343 357 self.n = n
344 358 self.__byTime = False
345 359 else:
346 360 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
347 361 self.n = 9999
348 362 self.__byTime = True
349 363
350 364 if overlapping:
351 365 self.__withOverapping = True
352 366 self.__buffer = None
353 367 else:
354 368 self.__withOverapping = False
355 369 self.__buffer = 0
356 370
357 371 self.__profIndex = 0
358 372
359 373 def putData(self, data):
360 374
361 375 """
362 376 Add a profile to the __buffer and increase in one the __profileIndex
363 377
364 378 """
365 379
366 380 if not self.__withOverapping:
367 381 self.__buffer += data.copy()
368 382 self.__profIndex += 1
369 383 return
370 384
371 385 #Overlapping data
372 386 nChannels, nHeis = data.shape
373 387 data = numpy.reshape(data, (1, nChannels, nHeis))
374 388
375 389 #If the buffer is empty then it takes the data value
376 390 if self.__buffer is None:
377 391 self.__buffer = data
378 392 self.__profIndex += 1
379 393 return
380 394
381 395 #If the buffer length is lower than n then stakcing the data value
382 396 if self.__profIndex < self.n:
383 397 self.__buffer = numpy.vstack((self.__buffer, data))
384 398 self.__profIndex += 1
385 399 return
386 400
387 401 #If the buffer length is equal to n then replacing the last buffer value with the data value
388 402 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
389 403 self.__buffer[self.n-1] = data
390 404 self.__profIndex = self.n
391 405 return
392 406
393 407
394 408 def pushData(self):
395 409 """
396 410 Return the sum of the last profiles and the profiles used in the sum.
397 411
398 412 Affected:
399 413
400 414 self.__profileIndex
401 415
402 416 """
403 417
404 418 if not self.__withOverapping:
405 419 data = self.__buffer
406 420 n = self.__profIndex
407 421
408 422 self.__buffer = 0
409 423 self.__profIndex = 0
410 424
411 425 return data, n
412 426
413 427 #Integration with Overlapping
414 428 data = numpy.sum(self.__buffer, axis=0)
415 429 n = self.__profIndex
416 430
417 431 return data, n
418 432
419 433 def byProfiles(self, data):
420 434
421 435 self.__dataReady = False
422 436 avgdata = None
423 437 # n = None
424 438
425 439 self.putData(data)
426 440
427 441 if self.__profIndex == self.n:
428 442
429 443 avgdata, n = self.pushData()
430 444 self.__dataReady = True
431 445
432 446 return avgdata
433 447
434 448 def byTime(self, data, datatime):
435 449
436 450 self.__dataReady = False
437 451 avgdata = None
438 452 n = None
439 453
440 454 self.putData(data)
441 455
442 456 if (datatime - self.__initime) >= self.__integrationtime:
443 457 avgdata, n = self.pushData()
444 458 self.n = n
445 459 self.__dataReady = True
446 460
447 461 return avgdata
448 462
449 463 def integrate(self, data, datatime=None):
450 464
451 465 if self.__initime == None:
452 466 self.__initime = datatime
453 467
454 468 if self.__byTime:
455 469 avgdata = self.byTime(data, datatime)
456 470 else:
457 471 avgdata = self.byProfiles(data)
458 472
459 473
460 474 self.__lastdatatime = datatime
461 475
462 476 if avgdata is None:
463 477 return None, None
464 478
465 479 avgdatatime = self.__initime
466 480
467 481 deltatime = datatime -self.__lastdatatime
468 482
469 483 if not self.__withOverapping:
470 484 self.__initime = datatime
471 485 else:
472 486 self.__initime += deltatime
473 487
474 488 return avgdata, avgdatatime
475 489
476 490 def integrateByBlock(self, dataOut):
477 491
478 492 times = int(dataOut.data.shape[1]/self.n)
479 493 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
480 494
481 495 id_min = 0
482 496 id_max = self.n
483 497
484 498 for i in range(times):
485 499 junk = dataOut.data[:,id_min:id_max,:]
486 500 avgdata[:,i,:] = junk.sum(axis=1)
487 501 id_min += self.n
488 502 id_max += self.n
489 503
490 504 timeInterval = dataOut.ippSeconds*self.n
491 505 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
492 506 self.__dataReady = True
493 507 return avgdata, avgdatatime
494 508
495 509 def run(self, dataOut, **kwargs):
496 510
497 511 if not self.isConfig:
498 512 self.setup(**kwargs)
499 513 self.isConfig = True
500 514
501 515 if dataOut.flagDataAsBlock:
502 516 """
503 517 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
504 518 """
505 519 avgdata, avgdatatime = self.integrateByBlock(dataOut)
506 520 dataOut.nProfiles /= self.n
507 521 else:
508 522 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
509 523
510 524 # dataOut.timeInterval *= n
511 525 dataOut.flagNoData = True
512 526
513 527 if self.__dataReady:
514 528 dataOut.data = avgdata
515 529 dataOut.nCohInt *= self.n
516 530 dataOut.utctime = avgdatatime
517 531 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
518 532 dataOut.flagNoData = False
519 533
520 534 class Decoder(Operation):
521 535
522 536 isConfig = False
523 537 __profIndex = 0
524 538
525 539 code = None
526 540
527 541 nCode = None
528 542 nBaud = None
529 543
530 544
531 545 def __init__(self):
532 546
533 547 Operation.__init__(self)
534 548
535 549 self.times = None
536 550 self.osamp = None
537 551 # self.__setValues = False
538 552 self.isConfig = False
539 553
540 554 def setup(self, code, osamp, dataOut):
541 555
542 556 self.__profIndex = 0
543 557
544 558 self.code = code
545 559
546 560 self.nCode = len(code)
547 561 self.nBaud = len(code[0])
548 562
549 563 if (osamp != None) and (osamp >1):
550 564 self.osamp = osamp
551 565 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
552 566 self.nBaud = self.nBaud*self.osamp
553 567
554 568 self.__nChannels = dataOut.nChannels
555 569 self.__nProfiles = dataOut.nProfiles
556 570 self.__nHeis = dataOut.nHeights
557 571
558 572 if self.__nHeis < self.nBaud:
559 573 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
560 574
561 575 #Frequency
562 576 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
563 577
564 578 __codeBuffer[:,0:self.nBaud] = self.code
565 579
566 580 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
567 581
568 582 if dataOut.flagDataAsBlock:
569 583
570 584 self.ndatadec = self.__nHeis #- self.nBaud + 1
571 585
572 586 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
573 587
574 588 else:
575 589
576 590 #Time
577 591 self.ndatadec = self.__nHeis #- self.nBaud + 1
578 592
579 593 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
580 594
581 595 def __convolutionInFreq(self, data):
582 596
583 597 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
584 598
585 599 fft_data = numpy.fft.fft(data, axis=1)
586 600
587 601 conv = fft_data*fft_code
588 602
589 603 data = numpy.fft.ifft(conv,axis=1)
590 604
591 605 return data
592 606
593 607 def __convolutionInFreqOpt(self, data):
594 608
595 609 raise NotImplementedError
596 610
597 611 def __convolutionInTime(self, data):
598 612
599 613 code = self.code[self.__profIndex]
600 614
601 615 for i in range(self.__nChannels):
602 616 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
603 617
604 618 return self.datadecTime
605 619
606 620 def __convolutionByBlockInTime(self, data):
607 621
608 622 repetitions = self.__nProfiles / self.nCode
609 623
610 624 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
611 625 junk = junk.flatten()
612 626 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
613 627
614 628 for i in range(self.__nChannels):
615 629 for j in range(self.__nProfiles):
616 630 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
617 631
618 632 return self.datadecTime
619 633
620 634 def __convolutionByBlockInFreq(self, data):
621 635
622 636 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
623 637
624 638
625 639 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
626 640
627 641 fft_data = numpy.fft.fft(data, axis=2)
628 642
629 643 conv = fft_data*fft_code
630 644
631 645 data = numpy.fft.ifft(conv,axis=2)
632 646
633 647 return data
634 648
635 649 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
636 650
637 651 if dataOut.flagDecodeData:
638 652 print "This data is already decoded, recoding again ..."
639 653
640 654 if not self.isConfig:
641 655
642 656 if code is None:
643 657 if dataOut.code is None:
644 658 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
645 659
646 660 code = dataOut.code
647 661 else:
648 662 code = numpy.array(code).reshape(nCode,nBaud)
649 663
650 664 self.setup(code, osamp, dataOut)
651 665
652 666 self.isConfig = True
653 667
654 668 if mode == 3:
655 669 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
656 670
657 671 if times != None:
658 672 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
659 673
660 674 if self.code is None:
661 675 print "Fail decoding: Code is not defined."
662 676 return
663 677
664 678 datadec = None
665 679 if mode == 3:
666 680 mode = 0
667 681
668 682 if dataOut.flagDataAsBlock:
669 683 """
670 684 Decoding when data have been read as block,
671 685 """
672 686
673 687 if mode == 0:
674 688 datadec = self.__convolutionByBlockInTime(dataOut.data)
675 689 if mode == 1:
676 690 datadec = self.__convolutionByBlockInFreq(dataOut.data)
677 691 else:
678 692 """
679 693 Decoding when data have been read profile by profile
680 694 """
681 695 if mode == 0:
682 696 datadec = self.__convolutionInTime(dataOut.data)
683 697
684 698 if mode == 1:
685 699 datadec = self.__convolutionInFreq(dataOut.data)
686 700
687 701 if mode == 2:
688 702 datadec = self.__convolutionInFreqOpt(dataOut.data)
689 703
690 704 if datadec is None:
691 705 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
692 706
693 707 dataOut.code = self.code
694 708 dataOut.nCode = self.nCode
695 709 dataOut.nBaud = self.nBaud
696 710
697 711 dataOut.data = datadec
698 712
699 713 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
700 714
701 715 dataOut.flagDecodeData = True #asumo q la data esta decodificada
702 716
703 717 if self.__profIndex == self.nCode-1:
704 718 self.__profIndex = 0
705 719 return 1
706 720
707 721 self.__profIndex += 1
708 722
709 723 return 1
710 724 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
711 725
712 726
713 727 class ProfileConcat(Operation):
714 728
715 729 isConfig = False
716 730 buffer = None
717 731
718 732 def __init__(self):
719 733
720 734 Operation.__init__(self)
721 735 self.profileIndex = 0
722 736
723 737 def reset(self):
724 738 self.buffer = numpy.zeros_like(self.buffer)
725 739 self.start_index = 0
726 740 self.times = 1
727 741
728 742 def setup(self, data, m, n=1):
729 743 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
730 744 self.nHeights = data.shape[1]#.nHeights
731 745 self.start_index = 0
732 746 self.times = 1
733 747
734 748 def concat(self, data):
735 749
736 750 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
737 751 self.start_index = self.start_index + self.nHeights
738 752
739 753 def run(self, dataOut, m):
740 754
741 755 dataOut.flagNoData = True
742 756
743 757 if not self.isConfig:
744 758 self.setup(dataOut.data, m, 1)
745 759 self.isConfig = True
746 760
747 761 if dataOut.flagDataAsBlock:
748 762 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
749 763
750 764 else:
751 765 self.concat(dataOut.data)
752 766 self.times += 1
753 767 if self.times > m:
754 768 dataOut.data = self.buffer
755 769 self.reset()
756 770 dataOut.flagNoData = False
757 771 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
758 772 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
759 773 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
760 774 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
761 775 dataOut.ippSeconds *= m
762 776
763 777 class ProfileSelector(Operation):
764 778
765 779 profileIndex = None
766 780 # Tamanho total de los perfiles
767 781 nProfiles = None
768 782
769 783 def __init__(self):
770 784
771 785 Operation.__init__(self)
772 786 self.profileIndex = 0
773 787
774 788 def incProfileIndex(self):
775 789
776 790 self.profileIndex += 1
777 791
778 792 if self.profileIndex >= self.nProfiles:
779 793 self.profileIndex = 0
780 794
781 795 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
782 796
783 797 if profileIndex < minIndex:
784 798 return False
785 799
786 800 if profileIndex > maxIndex:
787 801 return False
788 802
789 803 return True
790 804
791 805 def isThisProfileInList(self, profileIndex, profileList):
792 806
793 807 if profileIndex not in profileList:
794 808 return False
795 809
796 810 return True
797 811
798 812 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
799 813
800 814 """
801 815 ProfileSelector:
802 816
803 817 Inputs:
804 818 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
805 819
806 820 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
807 821
808 822 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
809 823
810 824 """
811 825
812 826 if rangeList is not None:
813 827 if type(rangeList[0]) not in (tuple, list):
814 828 rangeList = [rangeList]
815 829
816 830 dataOut.flagNoData = True
817 831
818 832 if dataOut.flagDataAsBlock:
819 833 """
820 834 data dimension = [nChannels, nProfiles, nHeis]
821 835 """
822 836 if profileList != None:
823 837 dataOut.data = dataOut.data[:,profileList,:]
824 838
825 839 if profileRangeList != None:
826 840 minIndex = profileRangeList[0]
827 841 maxIndex = profileRangeList[1]
828 842 profileList = range(minIndex, maxIndex+1)
829 843
830 844 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
831 845
832 846 if rangeList != None:
833 847
834 848 profileList = []
835 849
836 850 for thisRange in rangeList:
837 851 minIndex = thisRange[0]
838 852 maxIndex = thisRange[1]
839 853
840 854 profileList.extend(range(minIndex, maxIndex+1))
841 855
842 856 dataOut.data = dataOut.data[:,profileList,:]
843 857
844 858 dataOut.nProfiles = len(profileList)
845 859 dataOut.profileIndex = dataOut.nProfiles - 1
846 860 dataOut.flagNoData = False
847 861
848 862 return True
849 863
850 864 """
851 865 data dimension = [nChannels, nHeis]
852 866 """
853 867
854 868 if profileList != None:
855 869
856 870 if self.isThisProfileInList(dataOut.profileIndex, profileList):
857 871
858 872 self.nProfiles = len(profileList)
859 873 dataOut.nProfiles = self.nProfiles
860 874 dataOut.profileIndex = self.profileIndex
861 875 dataOut.flagNoData = False
862 876
863 877 self.incProfileIndex()
864 878 return True
865 879
866 880 if profileRangeList != None:
867 881
868 882 minIndex = profileRangeList[0]
869 883 maxIndex = profileRangeList[1]
870 884
871 885 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
872 886
873 887 self.nProfiles = maxIndex - minIndex + 1
874 888 dataOut.nProfiles = self.nProfiles
875 889 dataOut.profileIndex = self.profileIndex
876 890 dataOut.flagNoData = False
877 891
878 892 self.incProfileIndex()
879 893 return True
880 894
881 895 if rangeList != None:
882 896
883 897 nProfiles = 0
884 898
885 899 for thisRange in rangeList:
886 900 minIndex = thisRange[0]
887 901 maxIndex = thisRange[1]
888 902
889 903 nProfiles += maxIndex - minIndex + 1
890 904
891 905 for thisRange in rangeList:
892 906
893 907 minIndex = thisRange[0]
894 908 maxIndex = thisRange[1]
895 909
896 910 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
897 911
898 912 self.nProfiles = nProfiles
899 913 dataOut.nProfiles = self.nProfiles
900 914 dataOut.profileIndex = self.profileIndex
901 915 dataOut.flagNoData = False
902 916
903 917 self.incProfileIndex()
904 918
905 919 break
906 920
907 921 return True
908 922
909 923
910 924 if beam != None: #beam is only for AMISR data
911 925 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
912 926 dataOut.flagNoData = False
913 927 dataOut.profileIndex = self.profileIndex
914 928
915 929 self.incProfileIndex()
916 930
917 931 return True
918 932
919 933 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
920 934
921 935 return False
922 936
923 937 class Reshaper(Operation):
924 938
925 939 def __init__(self):
926 940
927 941 Operation.__init__(self)
928 942
929 943 self.__buffer = None
930 944 self.__nitems = 0
931 945
932 946 def __appendProfile(self, dataOut, nTxs):
933 947
934 948 if self.__buffer is None:
935 949 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
936 950 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
937 951
938 952 ini = dataOut.nHeights * self.__nitems
939 953 end = ini + dataOut.nHeights
940 954
941 955 self.__buffer[:, ini:end] = dataOut.data
942 956
943 957 self.__nitems += 1
944 958
945 959 return int(self.__nitems*nTxs)
946 960
947 961 def __getBuffer(self):
948 962
949 963 if self.__nitems == int(1./self.__nTxs):
950 964
951 965 self.__nitems = 0
952 966
953 967 return self.__buffer.copy()
954 968
955 969 return None
956 970
957 971 def __checkInputs(self, dataOut, shape, nTxs):
958 972
959 973 if shape is None and nTxs is None:
960 974 raise ValueError, "Reshaper: shape of factor should be defined"
961 975
962 976 if nTxs:
963 977 if nTxs < 0:
964 978 raise ValueError, "nTxs should be greater than 0"
965 979
966 980 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
967 981 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
968 982
969 983 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
970 984
971 985 return shape, nTxs
972 986
973 987 if len(shape) != 2 and len(shape) != 3:
974 988 raise ValueError, "shape dimension should be equal to 2 or 3. shape = (nProfiles, nHeis) or (nChannels, nProfiles, nHeis). Actually shape = (%d, %d, %d)" %(dataOut.nChannels, dataOut.nProfiles, dataOut.nHeights)
975 989
976 990 if len(shape) == 2:
977 991 shape_tuple = [dataOut.nChannels]
978 992 shape_tuple.extend(shape)
979 993 else:
980 994 shape_tuple = list(shape)
981 995
982 996 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
983 997
984 998 return shape_tuple, nTxs
985 999
986 1000 def run(self, dataOut, shape=None, nTxs=None):
987 1001
988 1002 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
989 1003
990 1004 dataOut.flagNoData = True
991 1005 profileIndex = None
992 1006
993 1007 if dataOut.flagDataAsBlock:
994 1008
995 1009 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
996 1010 dataOut.flagNoData = False
997 1011
998 1012 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
999 1013
1000 1014 else:
1001 1015
1002 1016 if self.__nTxs < 1:
1003 1017
1004 1018 self.__appendProfile(dataOut, self.__nTxs)
1005 1019 new_data = self.__getBuffer()
1006 1020
1007 1021 if new_data is not None:
1008 1022 dataOut.data = new_data
1009 1023 dataOut.flagNoData = False
1010 1024
1011 1025 profileIndex = dataOut.profileIndex*nTxs
1012 1026
1013 1027 else:
1014 1028 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1015 1029
1016 1030 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1017 1031
1018 1032 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1019 1033
1020 1034 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1021 1035
1022 1036 dataOut.profileIndex = profileIndex
1023 1037
1024 1038 dataOut.ippSeconds /= self.__nTxs
1025 1039
1026 1040 class SplitProfiles(Operation):
1027 1041
1028 1042 def __init__(self):
1029 1043
1030 1044 Operation.__init__(self)
1031 1045
1032 1046 def run(self, dataOut, n):
1033 1047
1034 1048 dataOut.flagNoData = True
1035 1049 profileIndex = None
1036 1050
1037 1051 if dataOut.flagDataAsBlock:
1038 1052
1039 1053 #nchannels, nprofiles, nsamples
1040 1054 shape = dataOut.data.shape
1041 1055
1042 1056 if shape[2] % n != 0:
1043 1057 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1044 1058
1045 1059 new_shape = shape[0], shape[1]*n, shape[2]/n
1046 1060
1047 1061 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1048 1062 dataOut.flagNoData = False
1049 1063
1050 1064 profileIndex = int(dataOut.nProfiles/n) - 1
1051 1065
1052 1066 else:
1053 1067
1054 1068 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1055 1069
1056 1070 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1057 1071
1058 1072 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1059 1073
1060 1074 dataOut.nProfiles = int(dataOut.nProfiles*n)
1061 1075
1062 1076 dataOut.profileIndex = profileIndex
1063 1077
1064 1078 dataOut.ippSeconds /= n
1065 1079
1066 1080 class CombineProfiles(Operation):
1067 1081
1068 1082 def __init__(self):
1069 1083
1070 1084 Operation.__init__(self)
1071 1085
1072 1086 self.__remData = None
1073 1087 self.__profileIndex = 0
1074 1088
1075 1089 def run(self, dataOut, n):
1076 1090
1077 1091 dataOut.flagNoData = True
1078 1092 profileIndex = None
1079 1093
1080 1094 if dataOut.flagDataAsBlock:
1081 1095
1082 1096 #nchannels, nprofiles, nsamples
1083 1097 shape = dataOut.data.shape
1084 1098 new_shape = shape[0], shape[1]/n, shape[2]*n
1085 1099
1086 1100 if shape[1] % n != 0:
1087 1101 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1088 1102
1089 1103 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1090 1104 dataOut.flagNoData = False
1091 1105
1092 1106 profileIndex = int(dataOut.nProfiles*n) - 1
1093 1107
1094 1108 else:
1095 1109
1096 1110 #nchannels, nsamples
1097 1111 if self.__remData is None:
1098 1112 newData = dataOut.data
1099 1113 else:
1100 1114 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1101 1115
1102 1116 self.__profileIndex += 1
1103 1117
1104 1118 if self.__profileIndex < n:
1105 1119 self.__remData = newData
1106 1120 #continue
1107 1121 return
1108 1122
1109 1123 self.__profileIndex = 0
1110 1124 self.__remData = None
1111 1125
1112 1126 dataOut.data = newData
1113 1127 dataOut.flagNoData = False
1114 1128
1115 1129 profileIndex = dataOut.profileIndex/n
1116 1130
1117 1131
1118 1132 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1119 1133
1120 1134 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1121 1135
1122 1136 dataOut.nProfiles = int(dataOut.nProfiles/n)
1123 1137
1124 1138 dataOut.profileIndex = profileIndex
1125 1139
1126 1140 dataOut.ippSeconds *= n
1127 1141
1128 1142 # import collections
1129 1143 # from scipy.stats import mode
1130 1144 #
1131 1145 # class Synchronize(Operation):
1132 1146 #
1133 1147 # isConfig = False
1134 1148 # __profIndex = 0
1135 1149 #
1136 1150 # def __init__(self):
1137 1151 #
1138 1152 # Operation.__init__(self)
1139 1153 # # self.isConfig = False
1140 1154 # self.__powBuffer = None
1141 1155 # self.__startIndex = 0
1142 1156 # self.__pulseFound = False
1143 1157 #
1144 1158 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1145 1159 #
1146 1160 # #Read data
1147 1161 #
1148 1162 # powerdB = dataOut.getPower(channel = channel)
1149 1163 # noisedB = dataOut.getNoise(channel = channel)[0]
1150 1164 #
1151 1165 # self.__powBuffer.extend(powerdB.flatten())
1152 1166 #
1153 1167 # dataArray = numpy.array(self.__powBuffer)
1154 1168 #
1155 1169 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1156 1170 #
1157 1171 # maxValue = numpy.nanmax(filteredPower)
1158 1172 #
1159 1173 # if maxValue < noisedB + 10:
1160 1174 # #No se encuentra ningun pulso de transmision
1161 1175 # return None
1162 1176 #
1163 1177 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1164 1178 #
1165 1179 # if len(maxValuesIndex) < 2:
1166 1180 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1167 1181 # return None
1168 1182 #
1169 1183 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1170 1184 #
1171 1185 # #Seleccionar solo valores con un espaciamiento de nSamples
1172 1186 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1173 1187 #
1174 1188 # if len(pulseIndex) < 2:
1175 1189 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1176 1190 # return None
1177 1191 #
1178 1192 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1179 1193 #
1180 1194 # #remover senales que se distancien menos de 10 unidades o muestras
1181 1195 # #(No deberian existir IPP menor a 10 unidades)
1182 1196 #
1183 1197 # realIndex = numpy.where(spacing > 10 )[0]
1184 1198 #
1185 1199 # if len(realIndex) < 2:
1186 1200 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1187 1201 # return None
1188 1202 #
1189 1203 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1190 1204 # realPulseIndex = pulseIndex[realIndex]
1191 1205 #
1192 1206 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1193 1207 #
1194 1208 # print "IPP = %d samples" %period
1195 1209 #
1196 1210 # self.__newNSamples = dataOut.nHeights #int(period)
1197 1211 # self.__startIndex = int(realPulseIndex[0])
1198 1212 #
1199 1213 # return 1
1200 1214 #
1201 1215 #
1202 1216 # def setup(self, nSamples, nChannels, buffer_size = 4):
1203 1217 #
1204 1218 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1205 1219 # maxlen = buffer_size*nSamples)
1206 1220 #
1207 1221 # bufferList = []
1208 1222 #
1209 1223 # for i in range(nChannels):
1210 1224 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1211 1225 # maxlen = buffer_size*nSamples)
1212 1226 #
1213 1227 # bufferList.append(bufferByChannel)
1214 1228 #
1215 1229 # self.__nSamples = nSamples
1216 1230 # self.__nChannels = nChannels
1217 1231 # self.__bufferList = bufferList
1218 1232 #
1219 1233 # def run(self, dataOut, channel = 0):
1220 1234 #
1221 1235 # if not self.isConfig:
1222 1236 # nSamples = dataOut.nHeights
1223 1237 # nChannels = dataOut.nChannels
1224 1238 # self.setup(nSamples, nChannels)
1225 1239 # self.isConfig = True
1226 1240 #
1227 1241 # #Append new data to internal buffer
1228 1242 # for thisChannel in range(self.__nChannels):
1229 1243 # bufferByChannel = self.__bufferList[thisChannel]
1230 1244 # bufferByChannel.extend(dataOut.data[thisChannel])
1231 1245 #
1232 1246 # if self.__pulseFound:
1233 1247 # self.__startIndex -= self.__nSamples
1234 1248 #
1235 1249 # #Finding Tx Pulse
1236 1250 # if not self.__pulseFound:
1237 1251 # indexFound = self.__findTxPulse(dataOut, channel)
1238 1252 #
1239 1253 # if indexFound == None:
1240 1254 # dataOut.flagNoData = True
1241 1255 # return
1242 1256 #
1243 1257 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1244 1258 # self.__pulseFound = True
1245 1259 # self.__startIndex = indexFound
1246 1260 #
1247 1261 # #If pulse was found ...
1248 1262 # for thisChannel in range(self.__nChannels):
1249 1263 # bufferByChannel = self.__bufferList[thisChannel]
1250 1264 # #print self.__startIndex
1251 1265 # x = numpy.array(bufferByChannel)
1252 1266 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1253 1267 #
1254 1268 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1255 1269 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1256 1270 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1257 1271 #
1258 1272 # dataOut.data = self.__arrayBuffer
1259 1273 #
1260 1274 # self.__startIndex += self.__newNSamples
1261 1275 #
1262 1276 # return No newline at end of file
General Comments 0
You need to be logged in to leave comments. Login now