##// END OF EJS Templates
schain
RSS
Bug fixed in jroIO_usrp: nProfiles and profileIndex values were defined. (nProfiles = number of profiles in one second)
Miguel Valdez -
r823:459a29f11c52
parent child
Show More
Show file before | Show file after | Hide comments
@@ -1,595 +1,601
1 1 '''
2 2 Created on Jul 3, 2014
3 3
4 4 @author: roj-idl71
5 5 '''
6 6 import os
7 7 import datetime
8 8 import numpy
9 9
10 10 try:
11 11 from gevent import sleep
12 12 except:
13 13 from time import sleep
14 14
15 15 from schainpy.model.data.jroheaderIO import RadarControllerHeader, SystemHeader
16 16 from schainpy.model.data.jrodata import Voltage
17 17 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
18 18
19 19 try:
20 20 import digital_rf_hdf5
21 21 except:
22 22 print 'You should install "digital_rf_hdf5" module if you want to read USRP data'
23 23
24 24 class USRPReader(ProcessingUnit):
25 25 '''
26 26 classdocs
27 27 '''
28 28
29 29 def __init__(self):
30 30 '''
31 31 Constructor
32 32 '''
33 33
34 34 ProcessingUnit.__init__(self)
35 35
36 36 self.dataOut = Voltage()
37 37 self.__printInfo = True
38 38 self.__flagDiscontinuousBlock = False
39 39 self.__bufferIndex = 9999999
40 40
41 41 self.__ippKm = None
42 42 self.__codeType = 0
43 43 self.__nCode = None
44 44 self.__nBaud = None
45 45 self.__code = None
46 46
47 47 def __getCurrentSecond(self):
48 48
49 49 return self.__thisUnixSample/self.__sample_rate
50 50
51 51 thisSecond = property(__getCurrentSecond, "I'm the 'thisSecond' property.")
52 52
53 53 def __setFileHeader(self):
54 54 '''
55 55 In this method will be initialized every parameter of dataOut object (header, no data)
56 56 '''
57 nProfiles = self.__sample_rate #Number of profiles by second
57 ippSeconds = 1.0*self.__nSamples/self.__sample_rate
58
59 nProfiles = 1.0/ippSeconds #Number of profiles in one second
58 60
59 61 self.dataOut.radarControllerHeaderObj = RadarControllerHeader(ippKm=self.__ippKm,
60 62 txA=0,
61 63 txB=0,
62 64 nWindows=1,
63 65 nHeights=self.__nSamples,
64 66 firstHeight=self.__firstHeigth,
65 67 deltaHeight=self.__deltaHeigth,
66 68 codeType=self.__codeType,
67 69 nCode=self.__nCode, nBaud=self.__nBaud,
68 70 code = self.__code)
69 71
70 72 self.dataOut.systemHeaderObj = SystemHeader(nSamples=self.__nSamples,
71 73 nProfiles=nProfiles,
72 74 nChannels=len(self.__channelList),
73 75 adcResolution=14)
74 76
75 77 self.dataOut.type = "Voltage"
76 78
77 79 self.dataOut.data = None
78 80
79 81 self.dataOut.dtype = numpy.dtype([('real','<i8'),('imag','<i8')])
80 82
81 83 # self.dataOut.nChannels = 0
82 84
83 85 # self.dataOut.nHeights = 0
84 86
85 87 self.dataOut.nProfiles = nProfiles
86 88
87 89 self.dataOut.heightList = self.__firstHeigth + numpy.arange(self.__nSamples, dtype = numpy.float)*self.__deltaHeigth
88 90
89 91 self.dataOut.channelList = self.__channelList
90 92
91 93 self.dataOut.blocksize = self.dataOut.getNChannels() * self.dataOut.getNHeights()
92 94
93 95 # self.dataOut.channelIndexList = None
94 96
95 97 self.dataOut.flagNoData = True
96 98
97 99 #Set to TRUE if the data is discontinuous
98 100 self.dataOut.flagDiscontinuousBlock = False
99 101
100 102 self.dataOut.utctime = None
101 103
102 104 self.dataOut.timeZone = self.__timezone/60 #timezone like jroheader, difference in minutes between UTC and localtime
103 105
104 106 self.dataOut.dstFlag = 0
105 107
106 108 self.dataOut.errorCount = 0
107 109
108 110 self.dataOut.nCohInt = 1
109 111
110 112 self.dataOut.flagDecodeData = False #asumo que la data esta decodificada
111 113
112 114 self.dataOut.flagDeflipData = False #asumo que la data esta sin flip
113 115
114 116 self.dataOut.flagShiftFFT = False
115 117
116 self.dataOut.ippSeconds = 1.0*self.__nSamples/self.__sample_rate
118 self.dataOut.ippSeconds = ippSeconds
117 119
118 120 #Time interval between profiles
119 121 #self.dataOut.timeInterval = self.dataOut.ippSeconds * self.dataOut.nCohInt
120 122
121 123 self.dataOut.frequency = self.__frequency
122 124
123 125 self.dataOut.realtime = self.__online
124 126
125 127 def findDatafiles(self, path, startDate=None, endDate=None):
126 128
127 129 if not os.path.isdir(path):
128 130 return []
129 131
130 132 try:
131 133 digitalReadObj = digital_rf_hdf5.read_hdf5(path, load_all_metadata=True)
132 134 except:
133 135 digitalReadObj = digital_rf_hdf5.read_hdf5(path)
134 136
135 137 channelNameList = digitalReadObj.get_channels()
136 138
137 139 if not channelNameList:
138 140 return []
139 141
140 142 metadata_dict = digitalReadObj.get_rf_file_metadata(channelNameList[0])
141 143
142 144 sample_rate = metadata_dict['sample_rate'][0]
143 145
144 146 this_metadata_file = digitalReadObj.get_metadata(channelNameList[0])
145 147
146 148 try:
147 149 timezone = this_metadata_file['timezone'].value
148 150 except:
149 151 timezone = 0
150 152
151 153 startUTCSecond, endUTCSecond = digitalReadObj.get_bounds(channelNameList[0])/sample_rate - timezone
152 154
153 155 startDatetime = datetime.datetime.utcfromtimestamp(startUTCSecond)
154 156 endDatatime = datetime.datetime.utcfromtimestamp(endUTCSecond)
155 157
156 158 if not startDate:
157 159 startDate = startDatetime.date()
158 160
159 161 if not endDate:
160 162 endDate = endDatatime.date()
161 163
162 164 dateList = []
163 165
164 166 thisDatetime = startDatetime
165 167
166 168 while(thisDatetime<=endDatatime):
167 169
168 170 thisDate = thisDatetime.date()
169 171
170 172 if thisDate < startDate:
171 173 continue
172 174
173 175 if thisDate > endDate:
174 176 break
175 177
176 178 dateList.append(thisDate)
177 179 thisDatetime += datetime.timedelta(1)
178 180
179 181 return dateList
180 182
181 183 def setup(self, path = None,
182 184 startDate = None,
183 185 endDate = None,
184 186 startTime = datetime.time(0,0,0),
185 187 endTime = datetime.time(23,59,59),
186 188 channelList = None,
187 189 nSamples = None,
188 190 ippKm = 60,
189 191 online = False,
190 192 delay = 60,
191 193 buffer_size = 1024,
192 194 **kwargs):
193 195 '''
194 196 In this method we should set all initial parameters.
195 197
196 198 Inputs:
197 199 path
198 200 startDate
199 201 endDate
200 202 startTime
201 203 endTime
202 204 set
203 205 expLabel
204 206 ext
205 207 online
206 208 delay
207 209 '''
208 210
209 211 if not os.path.isdir(path):
210 212 raise ValueError, "[Reading] Directory %s does not exist" %path
211 213
212 214 try:
213 215 self.digitalReadObj = digital_rf_hdf5.read_hdf5(path, load_all_metadata=True)
214 216 except:
215 217 self.digitalReadObj = digital_rf_hdf5.read_hdf5(path)
216 218
217 219 channelNameList = self.digitalReadObj.get_channels()
218 220
219 221 if not channelNameList:
220 222 raise ValueError, "[Reading] Directory %s does not have any files" %path
221 223
222 224 if not channelList:
223 225 channelList = range(len(channelNameList))
224 226
225 227 ########## Reading metadata ######################
226 228
227 229 metadata_dict = self.digitalReadObj.get_rf_file_metadata(channelNameList[channelList[0]])
228 230
229 231 self.__sample_rate = metadata_dict['sample_rate'][0]
230 232 # self.__samples_per_file = metadata_dict['samples_per_file'][0]
231 233 self.__deltaHeigth = 1e6*0.15/self.__sample_rate
232 234
233 235 this_metadata_file = self.digitalReadObj.get_metadata(channelNameList[channelList[0]])
234 236
235 237 self.__frequency = None
236 238 try:
237 239 self.__frequency = this_metadata_file['center_frequencies'].value
238 240 except:
239 241 self.__frequency = this_metadata_file['fc'].value
240 242
241 243 if not self.__frequency:
242 244 raise ValueError, "Center Frequency is not defined in metadata file"
243 245
244 246 try:
245 247 self.__timezone = this_metadata_file['timezone'].value
246 248 except:
247 249 self.__timezone = 0
248 250
249 251 self.__firstHeigth = 0
250 252
251 253 try:
252 254 codeType = this_metadata_file['codeType'].value
253 255 except:
254 256 codeType = 0
255 257
256 258 nCode = 1
257 259 nBaud = 1
258 260 code = numpy.ones((nCode, nBaud), dtype=numpy.int)
259 261
260 262 if codeType:
261 263 nCode = this_metadata_file['nCode'].value
262 264 nBaud = this_metadata_file['nBaud'].value
263 265 code = this_metadata_file['code'].value
264 266
265 267 if not ippKm:
266 268 try:
267 269 #seconds to km
268 270 ippKm = 1e6*0.15*this_metadata_file['ipp'].value
269 271 except:
270 272 ippKm = None
271 273
272 274 ####################################################
273 275 startUTCSecond = None
274 276 endUTCSecond = None
275 277
276 278 if startDate:
277 279 startDatetime = datetime.datetime.combine(startDate, startTime)
278 280 startUTCSecond = (startDatetime-datetime.datetime(1970,1,1)).total_seconds() + self.__timezone
279 281
280 282 if endDate:
281 283 endDatetime = datetime.datetime.combine(endDate, endTime)
282 284 endUTCSecond = (endDatetime-datetime.datetime(1970,1,1)).total_seconds() + self.__timezone
283 285
284 286 start_index, end_index = self.digitalReadObj.get_bounds(channelNameList[channelList[0]])
285 287
286 288 if not startUTCSecond:
287 289 startUTCSecond = start_index/self.__sample_rate
288 290
289 291 if start_index > startUTCSecond*self.__sample_rate:
290 292 startUTCSecond = start_index/self.__sample_rate
291 293
292 294 if not endUTCSecond:
293 295 endUTCSecond = end_index/self.__sample_rate
294 296
295 297 if end_index < endUTCSecond*self.__sample_rate:
296 298 endUTCSecond = end_index/self.__sample_rate
297 299
298 300 if not nSamples:
299 301 if not ippKm:
300 302 raise ValueError, "[Reading] nSamples or ippKm should be defined"
301 303
302 304 nSamples = int(ippKm / (1e6*0.15/self.__sample_rate))
303 305
304 306 channelBoundList = []
305 307 channelNameListFiltered = []
306 308
307 309 for thisIndexChannel in channelList:
308 310 thisChannelName = channelNameList[thisIndexChannel]
309 311 start_index, end_index = self.digitalReadObj.get_bounds(thisChannelName)
310 312 channelBoundList.append((start_index, end_index))
311 313 channelNameListFiltered.append(thisChannelName)
312 314
313 315 self.profileIndex = 0
314 316
315 317 self.__delay = delay
316 318 self.__ippKm = ippKm
317 319 self.__codeType = codeType
318 320 self.__nCode = nCode
319 321 self.__nBaud = nBaud
320 322 self.__code = code
321 323
322 324 self.__datapath = path
323 325 self.__online = online
324 326 self.__channelList = channelList
325 327 self.__channelNameList = channelNameListFiltered
326 328 self.__channelBoundList = channelBoundList
327 329 self.__nSamples = nSamples
328 330 self.__samples_to_read = int(buffer_size*nSamples)
329 331 self.__nChannels = len(self.__channelList)
330 332
331 333 self.__startUTCSecond = startUTCSecond
332 334 self.__endUTCSecond = endUTCSecond
333 335
334 336 self.__timeInterval = 1.0 * self.__samples_to_read/self.__sample_rate #Time interval
335 337
336 338 if online:
337 339 # self.__thisUnixSample = int(endUTCSecond*self.__sample_rate - 4*self.__samples_to_read)
338 340 startUTCSecond = numpy.floor(endUTCSecond)
339 341
340 342 self.__thisUnixSample = int(startUTCSecond*self.__sample_rate) - self.__samples_to_read
341 343
342 344 self.__data_buffer = numpy.zeros((self.__nChannels, self.__samples_to_read), dtype = numpy.complex)
343 345
344 346 self.__setFileHeader()
345 347 self.isConfig = True
346 348
347 349 print "[Reading] USRP Data was found from %s to %s " %(
348 350 datetime.datetime.utcfromtimestamp(self.__startUTCSecond - self.__timezone),
349 351 datetime.datetime.utcfromtimestamp(self.__endUTCSecond - self.__timezone)
350 352 )
351 353
352 354 print "[Reading] Starting process from %s to %s" %(datetime.datetime.utcfromtimestamp(startUTCSecond - self.__timezone),
353 355 datetime.datetime.utcfromtimestamp(endUTCSecond - self.__timezone)
354 356 )
355 357
356 358 def __reload(self):
357 359
358 360 if not self.__online:
359 361 return
360 362
361 363 # print
362 364 # print "%s not in range [%s, %s]" %(
363 365 # datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
364 366 # datetime.datetime.utcfromtimestamp(self.__startUTCSecond - self.__timezone),
365 367 # datetime.datetime.utcfromtimestamp(self.__endUTCSecond - self.__timezone)
366 368 # )
367 369 print "[Reading] reloading metadata ..."
368 370
369 371 try:
370 372 self.digitalReadObj.reload(complete_update=True)
371 373 except:
372 374 self.digitalReadObj.reload()
373 375
374 376 start_index, end_index = self.digitalReadObj.get_bounds(self.__channelNameList[self.__channelList[0]])
375 377
376 378 if start_index > self.__startUTCSecond*self.__sample_rate:
377 379 self.__startUTCSecond = 1.0*start_index/self.__sample_rate
378 380
379 381 if end_index > self.__endUTCSecond*self.__sample_rate:
380 382 self.__endUTCSecond = 1.0*end_index/self.__sample_rate
381 383 print
382 384 print "[Reading] New timerange found [%s, %s] " %(
383 385 datetime.datetime.utcfromtimestamp(self.__startUTCSecond - self.__timezone),
384 386 datetime.datetime.utcfromtimestamp(self.__endUTCSecond - self.__timezone)
385 387 )
386 388
387 389 return True
388 390
389 391 return False
390 392
391 393 def __readNextBlock(self, seconds=30, volt_scale = 218776):
392 394 '''
393 395 '''
394 396
395 397 #Set the next data
396 398 self.__flagDiscontinuousBlock = False
397 399 self.__thisUnixSample += self.__samples_to_read
398 400
399 401 if self.__thisUnixSample + 2*self.__samples_to_read > self.__endUTCSecond*self.__sample_rate:
400 402 print "[Reading] There are no more data into selected time-range"
401 403
402 404 self.__reload()
403 405
404 406 if self.__thisUnixSample + 2*self.__samples_to_read > self.__endUTCSecond*self.__sample_rate:
405 407 self.__thisUnixSample -= self.__samples_to_read
406 408 return False
407 409
408 410 indexChannel = 0
409 411
410 412 dataOk = False
411 413
412 414 for thisChannelName in self.__channelNameList:
413 415
414 416 try:
415 417 result = self.digitalReadObj.read_vector_c81d(self.__thisUnixSample,
416 418 self.__samples_to_read,
417 419 thisChannelName)
418 420
419 421 except IOError, e:
420 422 #read next profile
421 423 self.__flagDiscontinuousBlock = True
422 424 print "[Reading] %s" %datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone), e
423 425 break
424 426
425 427 if result.shape[0] != self.__samples_to_read:
426 428 self.__flagDiscontinuousBlock = True
427 429 print "[Reading] %s: Too few samples were found, just %d/%d samples" %(datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
428 430 result.shape[0],
429 431 self.__samples_to_read)
430 432 break
431 433
432 434 self.__data_buffer[indexChannel,:] = result*volt_scale
433 435
434 436 indexChannel += 1
435 437
436 438 dataOk = True
437 439
438 440 self.__utctime = self.__thisUnixSample/self.__sample_rate
439 441
440 442 if not dataOk:
441 443 return False
442 444
443 445 print "[Reading] %s: %d samples <> %f sec" %(datetime.datetime.utcfromtimestamp(self.thisSecond - self.__timezone),
444 446 self.__samples_to_read,
445 447 self.__timeInterval)
446 448
447 449 self.__bufferIndex = 0
448 450
449 451 return True
450 452
451 453 def __isBufferEmpty(self):
452 454
453 455 if self.__bufferIndex <= self.__samples_to_read - self.__nSamples:
454 456 return False
455 457
456 458 return True
457 459
458 460 def getData(self, seconds=30, nTries=5):
459 461
460 462 '''
461 463 This method gets the data from files and put the data into the dataOut object
462 464
463 465 In addition, increase el the buffer counter in one.
464 466
465 467 Return:
466 468 data : retorna un perfil de voltages (alturas * canales) copiados desde el
467 469 buffer. Si no hay mas archivos a leer retorna None.
468 470
469 471 Affected:
470 472 self.dataOut
471 473 self.profileIndex
472 474 self.flagDiscontinuousBlock
473 475 self.flagIsNewBlock
474 476 '''
475 477
476 478 err_counter = 0
477 479 self.dataOut.flagNoData = True
478 480
479 481 if self.__isBufferEmpty():
480 482
481 483 self.__flagDiscontinuousBlock = False
482 484
483 485 while True:
484 486 if self.__readNextBlock():
485 487 break
486 488
487 489 if self.__thisUnixSample > self.__endUTCSecond*self.__sample_rate:
488 490 return False
489 491
490 492 if self.__flagDiscontinuousBlock:
491 493 print '[Reading] discontinuous block found ... continue with the next block'
492 494 continue
493 495
494 496 if not self.__online:
495 497 return False
496 498
497 499 err_counter += 1
498 500 if err_counter > nTries:
499 501 return False
500 502
501 503 print '[Reading] waiting %d seconds to read a new block' %seconds
502 504 sleep(seconds)
503 505
504 506 self.dataOut.data = self.__data_buffer[:,self.__bufferIndex:self.__bufferIndex+self.__nSamples]
505 507 self.dataOut.utctime = (self.__thisUnixSample + self.__bufferIndex)/self.__sample_rate
506 508 self.dataOut.flagNoData = False
507 509 self.dataOut.flagDiscontinuousBlock = self.__flagDiscontinuousBlock
510 self.dataOut.profileIndex = self.profileIndex
508 511
509 512 self.__bufferIndex += self.__nSamples
510 513 self.profileIndex += 1
511 514
515 if self.profileIndex == self.dataOut.nProfiles:
516 self.profileIndex = 0
517
512 518 return True
513 519
514 520 def printInfo(self):
515 521 '''
516 522 '''
517 523 if self.__printInfo == False:
518 524 return
519 525
520 526 # self.systemHeaderObj.printInfo()
521 527 # self.radarControllerHeaderObj.printInfo()
522 528
523 529 self.__printInfo = False
524 530
525 531 def printNumberOfBlock(self):
526 532 '''
527 533 '''
528 534
529 535 print self.profileIndex
530 536
531 537 def run(self, **kwargs):
532 538 '''
533 539 This method will be called many times so here you should put all your code
534 540 '''
535 541
536 542 if not self.isConfig:
537 543 self.setup(**kwargs)
538 544
539 545 self.getData(seconds=self.__delay)
540 546
541 547 return
542 548
543 549 class USRPWriter(Operation):
544 550 '''
545 551 classdocs
546 552 '''
547 553
548 554 def __init__(self):
549 555 '''
550 556 Constructor
551 557 '''
552 558 self.dataOut = None
553 559
554 560 def setup(self, dataIn, path, blocksPerFile, set=0, ext=None):
555 561 '''
556 562 In this method we should set all initial parameters.
557 563
558 564 Input:
559 565 dataIn : Input data will also be outputa data
560 566
561 567 '''
562 568 self.dataOut = dataIn
563 569
564 570
565 571
566 572
567 573
568 574 self.isConfig = True
569 575
570 576 return
571 577
572 578 def run(self, dataIn, **kwargs):
573 579 '''
574 580 This method will be called many times so here you should put all your code
575 581
576 582 Inputs:
577 583
578 584 dataIn : object with the data
579 585
580 586 '''
581 587
582 588 if not self.isConfig:
583 589 self.setup(dataIn, **kwargs)
584 590
585 591
586 592 if __name__ == '__main__':
587 593
588 594 readObj = USRPReader()
589 595
590 596 while True:
591 597 readObj.run(path='/Volumes/DATA/haystack/passive_radar/')
592 598 # readObj.printInfo()
593 599 readObj.printNumberOfBlock()
594 600
595 601 No newline at end of file
Show file before | Show file after | Hide comments
@@ -1,1162 +1,1262
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 296 class CohInt(Operation):
297 297
298 298 isConfig = False
299 299
300 300 __profIndex = 0
301 301 __withOverapping = False
302 302
303 303 __byTime = False
304 304 __initime = None
305 305 __lastdatatime = None
306 306 __integrationtime = None
307 307
308 308 __buffer = None
309 309
310 310 __dataReady = False
311 311
312 312 n = None
313 313
314 314
315 315 def __init__(self):
316 316
317 317 Operation.__init__(self)
318 318
319 319 # self.isConfig = False
320 320
321 321 def setup(self, n=None, timeInterval=None, overlapping=False, byblock=False):
322 322 """
323 323 Set the parameters of the integration class.
324 324
325 325 Inputs:
326 326
327 327 n : Number of coherent integrations
328 328 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
329 329 overlapping :
330 330
331 331 """
332 332
333 333 self.__initime = None
334 334 self.__lastdatatime = 0
335 335 self.__buffer = None
336 336 self.__dataReady = False
337 337 self.byblock = byblock
338 338
339 339 if n == None and timeInterval == None:
340 340 raise ValueError, "n or timeInterval should be specified ..."
341 341
342 342 if n != None:
343 343 self.n = n
344 344 self.__byTime = False
345 345 else:
346 346 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
347 347 self.n = 9999
348 348 self.__byTime = True
349 349
350 350 if overlapping:
351 351 self.__withOverapping = True
352 352 self.__buffer = None
353 353 else:
354 354 self.__withOverapping = False
355 355 self.__buffer = 0
356 356
357 357 self.__profIndex = 0
358 358
359 359 def putData(self, data):
360 360
361 361 """
362 362 Add a profile to the __buffer and increase in one the __profileIndex
363 363
364 364 """
365 365
366 366 if not self.__withOverapping:
367 367 self.__buffer += data.copy()
368 368 self.__profIndex += 1
369 369 return
370 370
371 371 #Overlapping data
372 372 nChannels, nHeis = data.shape
373 373 data = numpy.reshape(data, (1, nChannels, nHeis))
374 374
375 375 #If the buffer is empty then it takes the data value
376 376 if self.__buffer is None:
377 377 self.__buffer = data
378 378 self.__profIndex += 1
379 379 return
380 380
381 381 #If the buffer length is lower than n then stakcing the data value
382 382 if self.__profIndex < self.n:
383 383 self.__buffer = numpy.vstack((self.__buffer, data))
384 384 self.__profIndex += 1
385 385 return
386 386
387 387 #If the buffer length is equal to n then replacing the last buffer value with the data value
388 388 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
389 389 self.__buffer[self.n-1] = data
390 390 self.__profIndex = self.n
391 391 return
392 392
393 393
394 394 def pushData(self):
395 395 """
396 396 Return the sum of the last profiles and the profiles used in the sum.
397 397
398 398 Affected:
399 399
400 400 self.__profileIndex
401 401
402 402 """
403 403
404 404 if not self.__withOverapping:
405 405 data = self.__buffer
406 406 n = self.__profIndex
407 407
408 408 self.__buffer = 0
409 409 self.__profIndex = 0
410 410
411 411 return data, n
412 412
413 413 #Integration with Overlapping
414 414 data = numpy.sum(self.__buffer, axis=0)
415 415 n = self.__profIndex
416 416
417 417 return data, n
418 418
419 419 def byProfiles(self, data):
420 420
421 421 self.__dataReady = False
422 422 avgdata = None
423 423 # n = None
424 424
425 425 self.putData(data)
426 426
427 427 if self.__profIndex == self.n:
428 428
429 429 avgdata, n = self.pushData()
430 430 self.__dataReady = True
431 431
432 432 return avgdata
433 433
434 434 def byTime(self, data, datatime):
435 435
436 436 self.__dataReady = False
437 437 avgdata = None
438 438 n = None
439 439
440 440 self.putData(data)
441 441
442 442 if (datatime - self.__initime) >= self.__integrationtime:
443 443 avgdata, n = self.pushData()
444 444 self.n = n
445 445 self.__dataReady = True
446 446
447 447 return avgdata
448 448
449 449 def integrate(self, data, datatime=None):
450 450
451 451 if self.__initime == None:
452 452 self.__initime = datatime
453 453
454 454 if self.__byTime:
455 455 avgdata = self.byTime(data, datatime)
456 456 else:
457 457 avgdata = self.byProfiles(data)
458 458
459 459
460 460 self.__lastdatatime = datatime
461 461
462 462 if avgdata is None:
463 463 return None, None
464 464
465 465 avgdatatime = self.__initime
466 466
467 467 deltatime = datatime -self.__lastdatatime
468 468
469 469 if not self.__withOverapping:
470 470 self.__initime = datatime
471 471 else:
472 472 self.__initime += deltatime
473 473
474 474 return avgdata, avgdatatime
475 475
476 476 def integrateByBlock(self, dataOut):
477 477
478 478 times = int(dataOut.data.shape[1]/self.n)
479 479 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
480 480
481 481 id_min = 0
482 482 id_max = self.n
483 483
484 484 for i in range(times):
485 485 junk = dataOut.data[:,id_min:id_max,:]
486 486 avgdata[:,i,:] = junk.sum(axis=1)
487 487 id_min += self.n
488 488 id_max += self.n
489 489
490 490 timeInterval = dataOut.ippSeconds*self.n
491 491 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
492 492 self.__dataReady = True
493 493 return avgdata, avgdatatime
494 494
495 495 def run(self, dataOut, **kwargs):
496 496
497 497 if not self.isConfig:
498 498 self.setup(**kwargs)
499 499 self.isConfig = True
500 500
501 501 if dataOut.flagDataAsBlock:
502 502 """
503 503 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
504 504 """
505 505 avgdata, avgdatatime = self.integrateByBlock(dataOut)
506 506 dataOut.nProfiles /= self.n
507 507 else:
508 508 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
509 509
510 510 # dataOut.timeInterval *= n
511 511 dataOut.flagNoData = True
512 512
513 513 if self.__dataReady:
514 514 dataOut.data = avgdata
515 515 dataOut.nCohInt *= self.n
516 516 dataOut.utctime = avgdatatime
517 517 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
518 518 dataOut.flagNoData = False
519 519
520 520 class Decoder(Operation):
521 521
522 522 isConfig = False
523 523 __profIndex = 0
524 524
525 525 code = None
526 526
527 527 nCode = None
528 528 nBaud = None
529 529
530 530
531 531 def __init__(self):
532 532
533 533 Operation.__init__(self)
534 534
535 535 self.times = None
536 536 self.osamp = None
537 537 # self.__setValues = False
538 538 self.isConfig = False
539 539
540 540 def setup(self, code, osamp, dataOut):
541 541
542 542 self.__profIndex = 0
543 543
544 544 self.code = code
545 545
546 546 self.nCode = len(code)
547 547 self.nBaud = len(code[0])
548 548
549 549 if (osamp != None) and (osamp >1):
550 550 self.osamp = osamp
551 551 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
552 552 self.nBaud = self.nBaud*self.osamp
553 553
554 554 self.__nChannels = dataOut.nChannels
555 555 self.__nProfiles = dataOut.nProfiles
556 556 self.__nHeis = dataOut.nHeights
557 557
558 558 if self.__nHeis < self.nBaud:
559 559 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
560 560
561 561 #Frequency
562 562 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
563 563
564 564 __codeBuffer[:,0:self.nBaud] = self.code
565 565
566 566 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
567 567
568 568 if dataOut.flagDataAsBlock:
569 569
570 570 self.ndatadec = self.__nHeis #- self.nBaud + 1
571 571
572 572 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
573 573
574 574 else:
575 575
576 576 #Time
577 577 self.ndatadec = self.__nHeis #- self.nBaud + 1
578 578
579 579 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
580 580
581 581 def __convolutionInFreq(self, data):
582 582
583 583 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
584 584
585 585 fft_data = numpy.fft.fft(data, axis=1)
586 586
587 587 conv = fft_data*fft_code
588 588
589 589 data = numpy.fft.ifft(conv,axis=1)
590 590
591 591 return data
592 592
593 593 def __convolutionInFreqOpt(self, data):
594 594
595 595 raise NotImplementedError
596 596
597 597 def __convolutionInTime(self, data):
598 598
599 599 code = self.code[self.__profIndex]
600 600
601 601 for i in range(self.__nChannels):
602 602 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
603 603
604 604 return self.datadecTime
605 605
606 606 def __convolutionByBlockInTime(self, data):
607 607
608 608 repetitions = self.__nProfiles / self.nCode
609 609
610 610 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
611 611 junk = junk.flatten()
612 612 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
613 613
614 614 for i in range(self.__nChannels):
615 615 for j in range(self.__nProfiles):
616 616 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
617 617
618 618 return self.datadecTime
619 619
620 620 def __convolutionByBlockInFreq(self, data):
621 621
622 622 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
623 623
624 624
625 625 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
626 626
627 627 fft_data = numpy.fft.fft(data, axis=2)
628 628
629 629 conv = fft_data*fft_code
630 630
631 631 data = numpy.fft.ifft(conv,axis=2)
632 632
633 633 return data
634 634
635 635 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
636 636
637 637 if dataOut.flagDecodeData:
638 638 print "This data is already decoded, recoding again ..."
639 639
640 640 if not self.isConfig:
641 641
642 642 if code is None:
643 643 if dataOut.code is None:
644 644 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
645 645
646 646 code = dataOut.code
647 647 else:
648 648 code = numpy.array(code).reshape(nCode,nBaud)
649 649
650 650 self.setup(code, osamp, dataOut)
651 651
652 652 self.isConfig = True
653 653
654 654 if mode == 3:
655 655 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
656 656
657 657 if times != None:
658 658 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
659 659
660 660 if self.code is None:
661 661 print "Fail decoding: Code is not defined."
662 662 return
663 663
664 664 datadec = None
665 665 if mode == 3:
666 666 mode = 0
667 667
668 668 if dataOut.flagDataAsBlock:
669 669 """
670 670 Decoding when data have been read as block,
671 671 """
672 672
673 673 if mode == 0:
674 674 datadec = self.__convolutionByBlockInTime(dataOut.data)
675 675 if mode == 1:
676 676 datadec = self.__convolutionByBlockInFreq(dataOut.data)
677 677 else:
678 678 """
679 679 Decoding when data have been read profile by profile
680 680 """
681 681 if mode == 0:
682 682 datadec = self.__convolutionInTime(dataOut.data)
683 683
684 684 if mode == 1:
685 685 datadec = self.__convolutionInFreq(dataOut.data)
686 686
687 687 if mode == 2:
688 688 datadec = self.__convolutionInFreqOpt(dataOut.data)
689 689
690 690 if datadec is None:
691 691 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
692 692
693 693 dataOut.code = self.code
694 694 dataOut.nCode = self.nCode
695 695 dataOut.nBaud = self.nBaud
696 696
697 697 dataOut.data = datadec
698 698
699 699 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
700 700
701 701 dataOut.flagDecodeData = True #asumo q la data esta decodificada
702 702
703 703 if self.__profIndex == self.nCode-1:
704 704 self.__profIndex = 0
705 705 return 1
706 706
707 707 self.__profIndex += 1
708 708
709 709 return 1
710 710 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
711 711
712 712
713 713 class ProfileConcat(Operation):
714 714
715 715 isConfig = False
716 716 buffer = None
717 717
718 718 def __init__(self):
719 719
720 720 Operation.__init__(self)
721 721 self.profileIndex = 0
722 722
723 723 def reset(self):
724 724 self.buffer = numpy.zeros_like(self.buffer)
725 725 self.start_index = 0
726 726 self.times = 1
727 727
728 728 def setup(self, data, m, n=1):
729 729 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
730 730 self.nHeights = data.shape[1]#.nHeights
731 731 self.start_index = 0
732 732 self.times = 1
733 733
734 734 def concat(self, data):
735 735
736 736 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
737 737 self.start_index = self.start_index + self.nHeights
738 738
739 739 def run(self, dataOut, m):
740 740
741 741 dataOut.flagNoData = True
742 742
743 743 if not self.isConfig:
744 744 self.setup(dataOut.data, m, 1)
745 745 self.isConfig = True
746 746
747 747 if dataOut.flagDataAsBlock:
748 748 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
749 749
750 750 else:
751 751 self.concat(dataOut.data)
752 752 self.times += 1
753 753 if self.times > m:
754 754 dataOut.data = self.buffer
755 755 self.reset()
756 756 dataOut.flagNoData = False
757 757 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
758 758 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
759 759 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
760 760 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
761 761 dataOut.ippSeconds *= m
762 762
763 763 class ProfileSelector(Operation):
764 764
765 765 profileIndex = None
766 766 # Tamanho total de los perfiles
767 767 nProfiles = None
768 768
769 769 def __init__(self):
770 770
771 771 Operation.__init__(self)
772 772 self.profileIndex = 0
773 773
774 774 def incProfileIndex(self):
775 775
776 776 self.profileIndex += 1
777 777
778 778 if self.profileIndex >= self.nProfiles:
779 779 self.profileIndex = 0
780 780
781 781 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
782 782
783 783 if profileIndex < minIndex:
784 784 return False
785 785
786 786 if profileIndex > maxIndex:
787 787 return False
788 788
789 789 return True
790 790
791 791 def isThisProfileInList(self, profileIndex, profileList):
792 792
793 793 if profileIndex not in profileList:
794 794 return False
795 795
796 796 return True
797 797
798 798 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
799 799
800 800 """
801 801 ProfileSelector:
802 802
803 803 Inputs:
804 804 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
805 805
806 806 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
807 807
808 808 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
809 809
810 810 """
811 811
812 812 if rangeList is not None:
813 813 if type(rangeList[0]) not in (tuple, list):
814 814 rangeList = [rangeList]
815 815
816 816 dataOut.flagNoData = True
817 817
818 818 if dataOut.flagDataAsBlock:
819 819 """
820 820 data dimension = [nChannels, nProfiles, nHeis]
821 821 """
822 822 if profileList != None:
823 823 dataOut.data = dataOut.data[:,profileList,:]
824 824
825 825 if profileRangeList != None:
826 826 minIndex = profileRangeList[0]
827 827 maxIndex = profileRangeList[1]
828 828 profileList = range(minIndex, maxIndex+1)
829 829
830 830 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
831 831
832 832 if rangeList != None:
833 833
834 834 profileList = []
835 835
836 836 for thisRange in rangeList:
837 837 minIndex = thisRange[0]
838 838 maxIndex = thisRange[1]
839 839
840 840 profileList.extend(range(minIndex, maxIndex+1))
841 841
842 842 dataOut.data = dataOut.data[:,profileList,:]
843 843
844 844 dataOut.nProfiles = len(profileList)
845 845 dataOut.profileIndex = dataOut.nProfiles - 1
846 846 dataOut.flagNoData = False
847 847
848 848 return True
849 849
850 850 """
851 851 data dimension = [nChannels, nHeis]
852 852 """
853 853
854 854 if profileList != None:
855 855
856 856 if self.isThisProfileInList(dataOut.profileIndex, profileList):
857 857
858 858 self.nProfiles = len(profileList)
859 859 dataOut.nProfiles = self.nProfiles
860 860 dataOut.profileIndex = self.profileIndex
861 861 dataOut.flagNoData = False
862 862
863 863 self.incProfileIndex()
864 864 return True
865 865
866 866 if profileRangeList != None:
867 867
868 868 minIndex = profileRangeList[0]
869 869 maxIndex = profileRangeList[1]
870 870
871 871 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
872 872
873 873 self.nProfiles = maxIndex - minIndex + 1
874 874 dataOut.nProfiles = self.nProfiles
875 875 dataOut.profileIndex = self.profileIndex
876 876 dataOut.flagNoData = False
877 877
878 878 self.incProfileIndex()
879 879 return True
880 880
881 881 if rangeList != None:
882 882
883 883 nProfiles = 0
884 884
885 885 for thisRange in rangeList:
886 886 minIndex = thisRange[0]
887 887 maxIndex = thisRange[1]
888 888
889 889 nProfiles += maxIndex - minIndex + 1
890 890
891 891 for thisRange in rangeList:
892 892
893 893 minIndex = thisRange[0]
894 894 maxIndex = thisRange[1]
895 895
896 896 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
897 897
898 898 self.nProfiles = nProfiles
899 899 dataOut.nProfiles = self.nProfiles
900 900 dataOut.profileIndex = self.profileIndex
901 901 dataOut.flagNoData = False
902 902
903 903 self.incProfileIndex()
904 904
905 905 break
906 906
907 907 return True
908 908
909 909
910 910 if beam != None: #beam is only for AMISR data
911 911 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
912 912 dataOut.flagNoData = False
913 913 dataOut.profileIndex = self.profileIndex
914 914
915 915 self.incProfileIndex()
916 916
917 917 return True
918 918
919 919 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
920 920
921 921 return False
922 922
923
924
925 923 class Reshaper(Operation):
926 924
927 925 def __init__(self):
928 926
929 927 Operation.__init__(self)
930 928
931 929 self.__buffer = None
932 930 self.__nitems = 0
933 931
934 932 def __appendProfile(self, dataOut, nTxs):
935 933
936 934 if self.__buffer is None:
937 935 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
938 936 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
939 937
940 938 ini = dataOut.nHeights * self.__nitems
941 939 end = ini + dataOut.nHeights
942 940
943 941 self.__buffer[:, ini:end] = dataOut.data
944 942
945 943 self.__nitems += 1
946 944
947 945 return int(self.__nitems*nTxs)
948 946
949 947 def __getBuffer(self):
950 948
951 949 if self.__nitems == int(1./self.__nTxs):
952 950
953 951 self.__nitems = 0
954 952
955 953 return self.__buffer.copy()
956 954
957 955 return None
958 956
959 957 def __checkInputs(self, dataOut, shape, nTxs):
960 958
961 959 if shape is None and nTxs is None:
962 960 raise ValueError, "Reshaper: shape of factor should be defined"
963 961
964 962 if nTxs:
965 963 if nTxs < 0:
966 964 raise ValueError, "nTxs should be greater than 0"
967 965
968 966 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
969 967 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
970 968
971 969 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
972 970
973 971 return shape, nTxs
974 972
975 973 if len(shape) != 2 and len(shape) != 3:
976 974 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)
977 975
978 976 if len(shape) == 2:
979 977 shape_tuple = [dataOut.nChannels]
980 978 shape_tuple.extend(shape)
981 979 else:
982 980 shape_tuple = list(shape)
983 981
984 982 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
985 983
986 984 return shape_tuple, nTxs
987 985
988 986 def run(self, dataOut, shape=None, nTxs=None):
989 987
990 988 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
991 989
992 990 dataOut.flagNoData = True
993 991 profileIndex = None
994 992
995 993 if dataOut.flagDataAsBlock:
996 994
997 995 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
998 996 dataOut.flagNoData = False
999 997
1000 998 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1001 999
1002 1000 else:
1003 1001
1004 1002 if self.__nTxs < 1:
1005 1003
1006 1004 self.__appendProfile(dataOut, self.__nTxs)
1007 1005 new_data = self.__getBuffer()
1008 1006
1009 1007 if new_data is not None:
1010 1008 dataOut.data = new_data
1011 1009 dataOut.flagNoData = False
1012 1010
1013 1011 profileIndex = dataOut.profileIndex*nTxs
1014 1012
1015 1013 else:
1016 1014 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1017 1015
1018 1016 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1019 1017
1020 1018 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1021 1019
1022 1020 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1023 1021
1024 1022 dataOut.profileIndex = profileIndex
1025 1023
1026 1024 dataOut.ippSeconds /= self.__nTxs
1027 #
1025
1026 class SplitProfiles(Operation):
1027
1028 def __init__(self):
1029
1030 Operation.__init__(self)
1031
1032 def run(self, dataOut, n):
1033
1034 dataOut.flagNoData = True
1035 profileIndex = None
1036
1037 if dataOut.flagDataAsBlock:
1038
1039 #nchannels, nprofiles, nsamples
1040 shape = dataOut.data.shape
1041
1042 if shape[2] % n != 0:
1043 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1044
1045 new_shape = shape[0], shape[1]*n, shape[2]/n
1046
1047 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1048 dataOut.flagNoData = False
1049
1050 profileIndex = int(dataOut.nProfiles/n) - 1
1051
1052 else:
1053
1054 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1055
1056 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1057
1058 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1059
1060 dataOut.nProfiles = int(dataOut.nProfiles*n)
1061
1062 dataOut.profileIndex = profileIndex
1063
1064 dataOut.ippSeconds /= n
1065
1066 class CombineProfiles(Operation):
1067
1068 def __init__(self):
1069
1070 Operation.__init__(self)
1071
1072 self.__remData = None
1073 self.__profileIndex = 0
1074
1075 def run(self, dataOut, n):
1076
1077 dataOut.flagNoData = True
1078 profileIndex = None
1079
1080 if dataOut.flagDataAsBlock:
1081
1082 #nchannels, nprofiles, nsamples
1083 shape = dataOut.data.shape
1084 new_shape = shape[0], shape[1]/n, shape[2]*n
1085
1086 if shape[1] % n != 0:
1087 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1088
1089 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1090 dataOut.flagNoData = False
1091
1092 profileIndex = int(dataOut.nProfiles*n) - 1
1093
1094 else:
1095
1096 #nchannels, nsamples
1097 if self.__remData is None:
1098 newData = dataOut.data
1099 else:
1100 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1101
1102 self.__profileIndex += 1
1103
1104 if self.__profileIndex < n:
1105 self.__remData = newData
1106 #continue
1107 return
1108
1109 self.__profileIndex = 0
1110 self.__remData = None
1111
1112 dataOut.data = newData
1113 dataOut.flagNoData = False
1114
1115 profileIndex = dataOut.profileIndex/n
1116
1117
1118 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1119
1120 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1121
1122 dataOut.nProfiles = int(dataOut.nProfiles/n)
1123
1124 dataOut.profileIndex = profileIndex
1125
1126 dataOut.ippSeconds *= n
1127
1028 1128 # import collections
1029 1129 # from scipy.stats import mode
1030 1130 #
1031 1131 # class Synchronize(Operation):
1032 1132 #
1033 1133 # isConfig = False
1034 1134 # __profIndex = 0
1035 1135 #
1036 1136 # def __init__(self):
1037 1137 #
1038 1138 # Operation.__init__(self)
1039 1139 # # self.isConfig = False
1040 1140 # self.__powBuffer = None
1041 1141 # self.__startIndex = 0
1042 1142 # self.__pulseFound = False
1043 1143 #
1044 1144 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1045 1145 #
1046 1146 # #Read data
1047 1147 #
1048 1148 # powerdB = dataOut.getPower(channel = channel)
1049 1149 # noisedB = dataOut.getNoise(channel = channel)[0]
1050 1150 #
1051 1151 # self.__powBuffer.extend(powerdB.flatten())
1052 1152 #
1053 1153 # dataArray = numpy.array(self.__powBuffer)
1054 1154 #
1055 1155 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1056 1156 #
1057 1157 # maxValue = numpy.nanmax(filteredPower)
1058 1158 #
1059 1159 # if maxValue < noisedB + 10:
1060 1160 # #No se encuentra ningun pulso de transmision
1061 1161 # return None
1062 1162 #
1063 1163 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1064 1164 #
1065 1165 # if len(maxValuesIndex) < 2:
1066 1166 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1067 1167 # return None
1068 1168 #
1069 1169 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1070 1170 #
1071 1171 # #Seleccionar solo valores con un espaciamiento de nSamples
1072 1172 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1073 1173 #
1074 1174 # if len(pulseIndex) < 2:
1075 1175 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1076 1176 # return None
1077 1177 #
1078 1178 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1079 1179 #
1080 1180 # #remover senales que se distancien menos de 10 unidades o muestras
1081 1181 # #(No deberian existir IPP menor a 10 unidades)
1082 1182 #
1083 1183 # realIndex = numpy.where(spacing > 10 )[0]
1084 1184 #
1085 1185 # if len(realIndex) < 2:
1086 1186 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1087 1187 # return None
1088 1188 #
1089 1189 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1090 1190 # realPulseIndex = pulseIndex[realIndex]
1091 1191 #
1092 1192 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1093 1193 #
1094 1194 # print "IPP = %d samples" %period
1095 1195 #
1096 1196 # self.__newNSamples = dataOut.nHeights #int(period)
1097 1197 # self.__startIndex = int(realPulseIndex[0])
1098 1198 #
1099 1199 # return 1
1100 1200 #
1101 1201 #
1102 1202 # def setup(self, nSamples, nChannels, buffer_size = 4):
1103 1203 #
1104 1204 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1105 1205 # maxlen = buffer_size*nSamples)
1106 1206 #
1107 1207 # bufferList = []
1108 1208 #
1109 1209 # for i in range(nChannels):
1110 1210 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1111 1211 # maxlen = buffer_size*nSamples)
1112 1212 #
1113 1213 # bufferList.append(bufferByChannel)
1114 1214 #
1115 1215 # self.__nSamples = nSamples
1116 1216 # self.__nChannels = nChannels
1117 1217 # self.__bufferList = bufferList
1118 1218 #
1119 1219 # def run(self, dataOut, channel = 0):
1120 1220 #
1121 1221 # if not self.isConfig:
1122 1222 # nSamples = dataOut.nHeights
1123 1223 # nChannels = dataOut.nChannels
1124 1224 # self.setup(nSamples, nChannels)
1125 1225 # self.isConfig = True
1126 1226 #
1127 1227 # #Append new data to internal buffer
1128 1228 # for thisChannel in range(self.__nChannels):
1129 1229 # bufferByChannel = self.__bufferList[thisChannel]
1130 1230 # bufferByChannel.extend(dataOut.data[thisChannel])
1131 1231 #
1132 1232 # if self.__pulseFound:
1133 1233 # self.__startIndex -= self.__nSamples
1134 1234 #
1135 1235 # #Finding Tx Pulse
1136 1236 # if not self.__pulseFound:
1137 1237 # indexFound = self.__findTxPulse(dataOut, channel)
1138 1238 #
1139 1239 # if indexFound == None:
1140 1240 # dataOut.flagNoData = True
1141 1241 # return
1142 1242 #
1143 1243 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1144 1244 # self.__pulseFound = True
1145 1245 # self.__startIndex = indexFound
1146 1246 #
1147 1247 # #If pulse was found ...
1148 1248 # for thisChannel in range(self.__nChannels):
1149 1249 # bufferByChannel = self.__bufferList[thisChannel]
1150 1250 # #print self.__startIndex
1151 1251 # x = numpy.array(bufferByChannel)
1152 1252 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1153 1253 #
1154 1254 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1155 1255 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1156 1256 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1157 1257 #
1158 1258 # dataOut.data = self.__arrayBuffer
1159 1259 #
1160 1260 # self.__startIndex += self.__newNSamples
1161 1261 #
1162 1262 # return No newline at end of file
General Comments 0
You need to be logged in to leave comments. Login now