##// END OF EJS Templates
profile concat parameter
avaldez -
r1247:9ada11315bac
parent child
Show More
@@ -1,1251 +1,1261
1 '''
1 '''
2
2
3 $Author: murco $
3 $Author: murco $
4 $Id: JROData.py 173 2012-11-20 15:06:21Z murco $
4 $Id: JROData.py 173 2012-11-20 15:06:21Z murco $
5 '''
5 '''
6
6
7 import copy
7 import copy
8 import numpy
8 import numpy
9 import datetime
9 import datetime
10
10
11 from jroheaderIO import SystemHeader, RadarControllerHeader
11 from jroheaderIO import SystemHeader, RadarControllerHeader
12 from schainpy import cSchain
12 from schainpy import cSchain
13
13
14
14
15 def getNumpyDtype(dataTypeCode):
15 def getNumpyDtype(dataTypeCode):
16
16
17 if dataTypeCode == 0:
17 if dataTypeCode == 0:
18 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
18 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
19 elif dataTypeCode == 1:
19 elif dataTypeCode == 1:
20 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
20 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
21 elif dataTypeCode == 2:
21 elif dataTypeCode == 2:
22 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
22 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
23 elif dataTypeCode == 3:
23 elif dataTypeCode == 3:
24 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
24 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
25 elif dataTypeCode == 4:
25 elif dataTypeCode == 4:
26 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
26 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
27 elif dataTypeCode == 5:
27 elif dataTypeCode == 5:
28 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
28 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
29 else:
29 else:
30 raise ValueError, 'dataTypeCode was not defined'
30 raise ValueError, 'dataTypeCode was not defined'
31
31
32 return numpyDtype
32 return numpyDtype
33
33
34
34
35 def getDataTypeCode(numpyDtype):
35 def getDataTypeCode(numpyDtype):
36
36
37 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
37 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
38 datatype = 0
38 datatype = 0
39 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
39 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
40 datatype = 1
40 datatype = 1
41 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
41 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
42 datatype = 2
42 datatype = 2
43 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
43 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
44 datatype = 3
44 datatype = 3
45 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
45 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
46 datatype = 4
46 datatype = 4
47 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
47 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
48 datatype = 5
48 datatype = 5
49 else:
49 else:
50 datatype = None
50 datatype = None
51
51
52 return datatype
52 return datatype
53
53
54
54
55 def hildebrand_sekhon(data, navg):
55 def hildebrand_sekhon(data, navg):
56 """
56 """
57 This method is for the objective determination of the noise level in Doppler spectra. This
57 This method is for the objective determination of the noise level in Doppler spectra. This
58 implementation technique is based on the fact that the standard deviation of the spectral
58 implementation technique is based on the fact that the standard deviation of the spectral
59 densities is equal to the mean spectral density for white Gaussian noise
59 densities is equal to the mean spectral density for white Gaussian noise
60
60
61 Inputs:
61 Inputs:
62 Data : heights
62 Data : heights
63 navg : numbers of averages
63 navg : numbers of averages
64
64
65 Return:
65 Return:
66 -1 : any error
66 -1 : any error
67 anoise : noise's level
67 anoise : noise's level
68 """
68 """
69
69
70 sortdata = numpy.sort(data, axis=None)
70 sortdata = numpy.sort(data, axis=None)
71 # lenOfData = len(sortdata)
71 # lenOfData = len(sortdata)
72 # nums_min = lenOfData*0.2
72 # nums_min = lenOfData*0.2
73 #
73 #
74 # if nums_min <= 5:
74 # if nums_min <= 5:
75 # nums_min = 5
75 # nums_min = 5
76 #
76 #
77 # sump = 0.
77 # sump = 0.
78 #
78 #
79 # sumq = 0.
79 # sumq = 0.
80 #
80 #
81 # j = 0
81 # j = 0
82 #
82 #
83 # cont = 1
83 # cont = 1
84 #
84 #
85 # while((cont==1)and(j<lenOfData)):
85 # while((cont==1)and(j<lenOfData)):
86 #
86 #
87 # sump += sortdata[j]
87 # sump += sortdata[j]
88 #
88 #
89 # sumq += sortdata[j]**2
89 # sumq += sortdata[j]**2
90 #
90 #
91 # if j > nums_min:
91 # if j > nums_min:
92 # rtest = float(j)/(j-1) + 1.0/navg
92 # rtest = float(j)/(j-1) + 1.0/navg
93 # if ((sumq*j) > (rtest*sump**2)):
93 # if ((sumq*j) > (rtest*sump**2)):
94 # j = j - 1
94 # j = j - 1
95 # sump = sump - sortdata[j]
95 # sump = sump - sortdata[j]
96 # sumq = sumq - sortdata[j]**2
96 # sumq = sumq - sortdata[j]**2
97 # cont = 0
97 # cont = 0
98 #
98 #
99 # j += 1
99 # j += 1
100 #
100 #
101 # lnoise = sump /j
101 # lnoise = sump /j
102 #
102 #
103 # return lnoise
103 # return lnoise
104
104
105 return cSchain.hildebrand_sekhon(sortdata, navg)
105 return cSchain.hildebrand_sekhon(sortdata, navg)
106
106
107
107
108 class Beam:
108 class Beam:
109
109
110 def __init__(self):
110 def __init__(self):
111 self.codeList = []
111 self.codeList = []
112 self.azimuthList = []
112 self.azimuthList = []
113 self.zenithList = []
113 self.zenithList = []
114
114
115
115
116 class GenericData(object):
116 class GenericData(object):
117
117
118 flagNoData = True
118 flagNoData = True
119
119
120 def copy(self, inputObj=None):
120 def copy(self, inputObj=None):
121
121
122 if inputObj == None:
122 if inputObj == None:
123 return copy.deepcopy(self)
123 return copy.deepcopy(self)
124
124
125 for key in inputObj.__dict__.keys():
125 for key in inputObj.__dict__.keys():
126
126
127 attribute = inputObj.__dict__[key]
127 attribute = inputObj.__dict__[key]
128
128
129 # If this attribute is a tuple or list
129 # If this attribute is a tuple or list
130 if type(inputObj.__dict__[key]) in (tuple, list):
130 if type(inputObj.__dict__[key]) in (tuple, list):
131 self.__dict__[key] = attribute[:]
131 self.__dict__[key] = attribute[:]
132 continue
132 continue
133
133
134 # If this attribute is another object or instance
134 # If this attribute is another object or instance
135 if hasattr(attribute, '__dict__'):
135 if hasattr(attribute, '__dict__'):
136 self.__dict__[key] = attribute.copy()
136 self.__dict__[key] = attribute.copy()
137 continue
137 continue
138
138
139 self.__dict__[key] = inputObj.__dict__[key]
139 self.__dict__[key] = inputObj.__dict__[key]
140
140
141 def deepcopy(self):
141 def deepcopy(self):
142
142
143 return copy.deepcopy(self)
143 return copy.deepcopy(self)
144
144
145 def isEmpty(self):
145 def isEmpty(self):
146
146
147 return self.flagNoData
147 return self.flagNoData
148
148
149
149
150 class JROData(GenericData):
150 class JROData(GenericData):
151
151
152 # m_BasicHeader = BasicHeader()
152 # m_BasicHeader = BasicHeader()
153 # m_ProcessingHeader = ProcessingHeader()
153 # m_ProcessingHeader = ProcessingHeader()
154
154
155 systemHeaderObj = SystemHeader()
155 systemHeaderObj = SystemHeader()
156
156
157 radarControllerHeaderObj = RadarControllerHeader()
157 radarControllerHeaderObj = RadarControllerHeader()
158
158
159 # data = None
159 # data = None
160
160
161 type = None
161 type = None
162
162
163 datatype = None # dtype but in string
163 datatype = None # dtype but in string
164
164
165 # dtype = None
165 # dtype = None
166
166
167 # nChannels = None
167 # nChannels = None
168
168
169 # nHeights = None
169 # nHeights = None
170
170
171 nProfiles = None
171 nProfiles = None
172
172
173 heightList = None
173 heightList = None
174
174
175 channelList = None
175 channelList = None
176
176
177 flagDiscontinuousBlock = False
177 flagDiscontinuousBlock = False
178
178
179 useLocalTime = False
179 useLocalTime = False
180
180
181 utctime = None
181 utctime = None
182
182
183 timeZone = None
183 timeZone = None
184
184
185 dstFlag = None
185 dstFlag = None
186
186
187 errorCount = None
187 errorCount = None
188
188
189 blocksize = None
189 blocksize = None
190
190
191 # nCode = None
191 # nCode = None
192 #
192 #
193 # nBaud = None
193 # nBaud = None
194 #
194 #
195 # code = None
195 # code = None
196
196
197 flagDecodeData = False # asumo q la data no esta decodificada
197 flagDecodeData = False # asumo q la data no esta decodificada
198
198
199 flagDeflipData = False # asumo q la data no esta sin flip
199 flagDeflipData = False # asumo q la data no esta sin flip
200
200
201 flagShiftFFT = False
201 flagShiftFFT = False
202
202
203 # ippSeconds = None
203 # ippSeconds = None
204
204
205 # timeInterval = None
205 # timeInterval = None
206
206
207 nCohInt = None
207 nCohInt = None
208
208
209 # noise = None
209 # noise = None
210
210
211 windowOfFilter = 1
211 windowOfFilter = 1
212
212
213 # Speed of ligth
213 # Speed of ligth
214 C = 3e8
214 C = 3e8
215
215
216 frequency = 49.92e6
216 frequency = 49.92e6
217
217
218 realtime = False
218 realtime = False
219
219
220 beacon_heiIndexList = None
220 beacon_heiIndexList = None
221
221
222 last_block = None
222 last_block = None
223
223
224 blocknow = None
224 blocknow = None
225
225
226 azimuth = None
226 azimuth = None
227
227
228 zenith = None
228 zenith = None
229
229
230 beam = Beam()
230 beam = Beam()
231
231
232 profileIndex = None
232 profileIndex = None
233
233
234 def getNoise(self):
234 def getNoise(self):
235
235
236 raise NotImplementedError
236 raise NotImplementedError
237
237
238 def getNChannels(self):
238 def getNChannels(self):
239
239
240 return len(self.channelList)
240 return len(self.channelList)
241
241
242 def getChannelIndexList(self):
242 def getChannelIndexList(self):
243
243
244 return range(self.nChannels)
244 return range(self.nChannels)
245
245
246 def getNHeights(self):
246 def getNHeights(self):
247
247
248 return len(self.heightList)
248 return len(self.heightList)
249
249
250 def getHeiRange(self, extrapoints=0):
250 def getHeiRange(self, extrapoints=0):
251
251
252 heis = self.heightList
252 heis = self.heightList
253 # deltah = self.heightList[1] - self.heightList[0]
253 # deltah = self.heightList[1] - self.heightList[0]
254 #
254 #
255 # heis.append(self.heightList[-1])
255 # heis.append(self.heightList[-1])
256
256
257 return heis
257 return heis
258
258
259 def getDeltaH(self):
259 def getDeltaH(self):
260
260
261 delta = self.heightList[1] - self.heightList[0]
261 delta = self.heightList[1] - self.heightList[0]
262
262
263 return delta
263 return delta
264
264
265 def getltctime(self):
265 def getltctime(self):
266
266
267 if self.useLocalTime:
267 if self.useLocalTime:
268 return self.utctime - self.timeZone * 60
268 return self.utctime - self.timeZone * 60
269
269
270 return self.utctime
270 return self.utctime
271
271
272 def getDatatime(self):
272 def getDatatime(self):
273
273
274 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
274 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
275 return datatimeValue
275 return datatimeValue
276
276
277 def getTimeRange(self):
277 def getTimeRange(self):
278
278
279 datatime = []
279 datatime = []
280
280
281 datatime.append(self.ltctime)
281 datatime.append(self.ltctime)
282 datatime.append(self.ltctime + self.timeInterval + 1)
282 datatime.append(self.ltctime + self.timeInterval + 1)
283
283
284 datatime = numpy.array(datatime)
284 datatime = numpy.array(datatime)
285
285
286 return datatime
286 return datatime
287
287
288 def getFmaxTimeResponse(self):
288 def getFmaxTimeResponse(self):
289
289
290 period = (10**-6) * self.getDeltaH() / (0.15)
290 period = (10**-6) * self.getDeltaH() / (0.15)
291
291
292 PRF = 1. / (period * self.nCohInt)
292 PRF = 1. / (period * self.nCohInt)
293
293
294 fmax = PRF
294 fmax = PRF
295
295
296 return fmax
296 return fmax
297
297
298 def getFmax(self):
298 def getFmax(self):
299 PRF = 1. / (self.ippSeconds * self.nCohInt)
299 PRF = 1. / (self.ippSeconds * self.nCohInt)
300 #print "ippSeconds",self.ippSeconds
301 #print "nCohInt",self.nIncohInt
302 #print PRF
303 #import time
304 #time.sleep(30)
300 fmax = PRF
305 fmax = PRF
301 return fmax
306 return fmax
302
307
303 def getVmax(self):
308 def getVmax(self):
304
309
305 _lambda = self.C / self.frequency
310 _lambda = self.C / self.frequency
306
311
307 vmax = self.getFmax() * _lambda / 2
312 vmax = self.getFmax() * _lambda / 2
308
313
309 return vmax
314 return vmax
310
315
311 def get_ippSeconds(self):
316 def get_ippSeconds(self):
312 '''
317 '''
313 '''
318 '''
314 return self.radarControllerHeaderObj.ippSeconds
319 return self.radarControllerHeaderObj.ippSeconds
315
320
316 def set_ippSeconds(self, ippSeconds):
321 def set_ippSeconds(self, ippSeconds):
317 '''
322 '''
318 '''
323 '''
319
324
320 self.radarControllerHeaderObj.ippSeconds = ippSeconds
325 self.radarControllerHeaderObj.ippSeconds = ippSeconds
321
326
322 return
327 return
323
328
324 def get_dtype(self):
329 def get_dtype(self):
325 '''
330 '''
326 '''
331 '''
327 return getNumpyDtype(self.datatype)
332 return getNumpyDtype(self.datatype)
328
333
329 def set_dtype(self, numpyDtype):
334 def set_dtype(self, numpyDtype):
330 '''
335 '''
331 '''
336 '''
332
337
333 self.datatype = getDataTypeCode(numpyDtype)
338 self.datatype = getDataTypeCode(numpyDtype)
334
339
335 def get_code(self):
340 def get_code(self):
336 '''
341 '''
337 '''
342 '''
338 return self.radarControllerHeaderObj.code
343 return self.radarControllerHeaderObj.code
339
344
340 def set_code(self, code):
345 def set_code(self, code):
341 '''
346 '''
342 '''
347 '''
343 self.radarControllerHeaderObj.code = code
348 self.radarControllerHeaderObj.code = code
344
349
345 return
350 return
346
351
347 def get_ncode(self):
352 def get_ncode(self):
348 '''
353 '''
349 '''
354 '''
350 return self.radarControllerHeaderObj.nCode
355 return self.radarControllerHeaderObj.nCode
351
356
352 def set_ncode(self, nCode):
357 def set_ncode(self, nCode):
353 '''
358 '''
354 '''
359 '''
355 self.radarControllerHeaderObj.nCode = nCode
360 self.radarControllerHeaderObj.nCode = nCode
356
361
357 return
362 return
358
363
359 def get_nbaud(self):
364 def get_nbaud(self):
360 '''
365 '''
361 '''
366 '''
362 return self.radarControllerHeaderObj.nBaud
367 return self.radarControllerHeaderObj.nBaud
363
368
364 def set_nbaud(self, nBaud):
369 def set_nbaud(self, nBaud):
365 '''
370 '''
366 '''
371 '''
367 self.radarControllerHeaderObj.nBaud = nBaud
372 self.radarControllerHeaderObj.nBaud = nBaud
368
373
369 return
374 return
370
375
371 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
376 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
372 channelIndexList = property(
377 channelIndexList = property(
373 getChannelIndexList, "I'm the 'channelIndexList' property.")
378 getChannelIndexList, "I'm the 'channelIndexList' property.")
374 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
379 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
375 #noise = property(getNoise, "I'm the 'nHeights' property.")
380 #noise = property(getNoise, "I'm the 'nHeights' property.")
376 datatime = property(getDatatime, "I'm the 'datatime' property")
381 datatime = property(getDatatime, "I'm the 'datatime' property")
377 ltctime = property(getltctime, "I'm the 'ltctime' property")
382 ltctime = property(getltctime, "I'm the 'ltctime' property")
378 ippSeconds = property(get_ippSeconds, set_ippSeconds)
383 ippSeconds = property(get_ippSeconds, set_ippSeconds)
379 dtype = property(get_dtype, set_dtype)
384 dtype = property(get_dtype, set_dtype)
380 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
385 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
381 code = property(get_code, set_code)
386 code = property(get_code, set_code)
382 nCode = property(get_ncode, set_ncode)
387 nCode = property(get_ncode, set_ncode)
383 nBaud = property(get_nbaud, set_nbaud)
388 nBaud = property(get_nbaud, set_nbaud)
384
389
385
390
386 class Voltage(JROData):
391 class Voltage(JROData):
387
392
388 # data es un numpy array de 2 dmensiones (canales, alturas)
393 # data es un numpy array de 2 dmensiones (canales, alturas)
389 data = None
394 data = None
390
395
391 def __init__(self):
396 def __init__(self):
392 '''
397 '''
393 Constructor
398 Constructor
394 '''
399 '''
395
400
396 self.useLocalTime = True
401 self.useLocalTime = True
397
402
398 self.radarControllerHeaderObj = RadarControllerHeader()
403 self.radarControllerHeaderObj = RadarControllerHeader()
399
404
400 self.systemHeaderObj = SystemHeader()
405 self.systemHeaderObj = SystemHeader()
401
406
402 self.type = "Voltage"
407 self.type = "Voltage"
403
408
404 self.data = None
409 self.data = None
405
410
406 # self.dtype = None
411 # self.dtype = None
407
412
408 # self.nChannels = 0
413 # self.nChannels = 0
409
414
410 # self.nHeights = 0
415 # self.nHeights = 0
411
416
412 self.nProfiles = None
417 self.nProfiles = None
413
418
414 self.heightList = None
419 self.heightList = None
415
420
416 self.channelList = None
421 self.channelList = None
417
422
418 # self.channelIndexList = None
423 # self.channelIndexList = None
419
424
420 self.flagNoData = True
425 self.flagNoData = True
421
426
422 self.flagDiscontinuousBlock = False
427 self.flagDiscontinuousBlock = False
423
428
424 self.utctime = None
429 self.utctime = None
425
430
426 self.timeZone = None
431 self.timeZone = None
427
432
428 self.dstFlag = None
433 self.dstFlag = None
429
434
430 self.errorCount = None
435 self.errorCount = None
431
436
432 self.nCohInt = None
437 self.nCohInt = None
433
438
434 self.blocksize = None
439 self.blocksize = None
435
440
436 self.flagDecodeData = False # asumo q la data no esta decodificada
441 self.flagDecodeData = False # asumo q la data no esta decodificada
437
442
438 self.flagDeflipData = False # asumo q la data no esta sin flip
443 self.flagDeflipData = False # asumo q la data no esta sin flip
439
444
440 self.flagShiftFFT = False
445 self.flagShiftFFT = False
441
446
442 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
447 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
443
448
444 self.profileIndex = 0
449 self.profileIndex = 0
445
450
446 def getNoisebyHildebrand(self, channel=None):
451 def getNoisebyHildebrand(self, channel=None):
447 """
452 """
448 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
453 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
449
454
450 Return:
455 Return:
451 noiselevel
456 noiselevel
452 """
457 """
453
458
454 if channel != None:
459 if channel != None:
455 data = self.data[channel]
460 data = self.data[channel]
456 nChannels = 1
461 nChannels = 1
457 else:
462 else:
458 data = self.data
463 data = self.data
459 nChannels = self.nChannels
464 nChannels = self.nChannels
460
465
461 noise = numpy.zeros(nChannels)
466 noise = numpy.zeros(nChannels)
462 power = data * numpy.conjugate(data)
467 power = data * numpy.conjugate(data)
463
468
464 for thisChannel in range(nChannels):
469 for thisChannel in range(nChannels):
465 if nChannels == 1:
470 if nChannels == 1:
466 daux = power[:].real
471 daux = power[:].real
467 else:
472 else:
468 daux = power[thisChannel, :].real
473 daux = power[thisChannel, :].real
469 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
474 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
470
475
471 return noise
476 return noise
472
477
473 def getNoise(self, type=1, channel=None):
478 def getNoise(self, type=1, channel=None):
474
479
475 if type == 1:
480 if type == 1:
476 noise = self.getNoisebyHildebrand(channel)
481 noise = self.getNoisebyHildebrand(channel)
477
482
478 return noise
483 return noise
479
484
480 def getPower(self, channel=None):
485 def getPower(self, channel=None):
481
486
482 if channel != None:
487 if channel != None:
483 data = self.data[channel]
488 data = self.data[channel]
484 else:
489 else:
485 data = self.data
490 data = self.data
486
491
487 power = data * numpy.conjugate(data)
492 power = data * numpy.conjugate(data)
488 powerdB = 10 * numpy.log10(power.real)
493 powerdB = 10 * numpy.log10(power.real)
489 powerdB = numpy.squeeze(powerdB)
494 powerdB = numpy.squeeze(powerdB)
490
495
491 return powerdB
496 return powerdB
492
497
493 def getTimeInterval(self):
498 def getTimeInterval(self):
494
499
495 timeInterval = self.ippSeconds * self.nCohInt
500 timeInterval = self.ippSeconds * self.nCohInt
496
501
497 return timeInterval
502 return timeInterval
498
503
499 noise = property(getNoise, "I'm the 'nHeights' property.")
504 noise = property(getNoise, "I'm the 'nHeights' property.")
500 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
505 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
501
506
502
507
503 class Spectra(JROData):
508 class Spectra(JROData):
504
509
505 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
510 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
506 data_spc = None
511 data_spc = None
507
512
508 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
513 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
509 data_cspc = None
514 data_cspc = None
510
515
511 # data dc es un numpy array de 2 dmensiones (canales, alturas)
516 # data dc es un numpy array de 2 dmensiones (canales, alturas)
512 data_dc = None
517 data_dc = None
513
518
514 # data power
519 # data power
515 data_pwr = None
520 data_pwr = None
516
521
517 nFFTPoints = None
522 nFFTPoints = None
518
523
519 # nPairs = None
524 # nPairs = None
520
525
521 pairsList = None
526 pairsList = None
522
527
523 nIncohInt = None
528 nIncohInt = None
524
529
525 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
530 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
526
531
527 nCohInt = None # se requiere para determinar el valor de timeInterval
532 nCohInt = None # se requiere para determinar el valor de timeInterval
528
533
529 ippFactor = None
534 ippFactor = None
530
535
531 profileIndex = 0
536 profileIndex = 0
532
537
533 plotting = "spectra"
538 plotting = "spectra"
534
539
535 def __init__(self):
540 def __init__(self):
536 '''
541 '''
537 Constructor
542 Constructor
538 '''
543 '''
539
544
540 self.useLocalTime = True
545 self.useLocalTime = True
541
546
542 self.radarControllerHeaderObj = RadarControllerHeader()
547 self.radarControllerHeaderObj = RadarControllerHeader()
543
548
544 self.systemHeaderObj = SystemHeader()
549 self.systemHeaderObj = SystemHeader()
545
550
546 self.type = "Spectra"
551 self.type = "Spectra"
547
552
548 # self.data = None
553 # self.data = None
549
554
550 # self.dtype = None
555 # self.dtype = None
551
556
552 # self.nChannels = 0
557 # self.nChannels = 0
553
558
554 # self.nHeights = 0
559 # self.nHeights = 0
555
560
556 self.nProfiles = None
561 self.nProfiles = None
557
562
558 self.heightList = None
563 self.heightList = None
559
564
560 self.channelList = None
565 self.channelList = None
561
566
562 # self.channelIndexList = None
567 # self.channelIndexList = None
563
568
564 self.pairsList = None
569 self.pairsList = None
565
570
566 self.flagNoData = True
571 self.flagNoData = True
567
572
568 self.flagDiscontinuousBlock = False
573 self.flagDiscontinuousBlock = False
569
574
570 self.utctime = None
575 self.utctime = None
571
576
572 self.nCohInt = None
577 self.nCohInt = None
573
578
574 self.nIncohInt = None
579 self.nIncohInt = None
575
580
576 self.blocksize = None
581 self.blocksize = None
577
582
578 self.nFFTPoints = None
583 self.nFFTPoints = None
579
584
580 self.wavelength = None
585 self.wavelength = None
581
586
582 self.flagDecodeData = False # asumo q la data no esta decodificada
587 self.flagDecodeData = False # asumo q la data no esta decodificada
583
588
584 self.flagDeflipData = False # asumo q la data no esta sin flip
589 self.flagDeflipData = False # asumo q la data no esta sin flip
585
590
586 self.flagShiftFFT = False
591 self.flagShiftFFT = False
587
592
588 self.ippFactor = 1
593 self.ippFactor = 1
589
594
590 #self.noise = None
595 #self.noise = None
591
596
592 self.beacon_heiIndexList = []
597 self.beacon_heiIndexList = []
593
598
594 self.noise_estimation = None
599 self.noise_estimation = None
595
600
596 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
601 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
597 """
602 """
598 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
603 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
599
604
600 Return:
605 Return:
601 noiselevel
606 noiselevel
602 """
607 """
603
608
604 noise = numpy.zeros(self.nChannels)
609 noise = numpy.zeros(self.nChannels)
605
610
606 for channel in range(self.nChannels):
611 for channel in range(self.nChannels):
607 daux = self.data_spc[channel,
612 daux = self.data_spc[channel,
608 xmin_index:xmax_index, ymin_index:ymax_index]
613 xmin_index:xmax_index, ymin_index:ymax_index]
609
614
610 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
615 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
611
616
612 return noise
617 return noise
613
618
614 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
619 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
615
620
616 if self.noise_estimation is not None:
621 if self.noise_estimation is not None:
617 # this was estimated by getNoise Operation defined in jroproc_spectra.py
622 # this was estimated by getNoise Operation defined in jroproc_spectra.py
618 return self.noise_estimation
623 return self.noise_estimation
619 else:
624 else:
620 noise = self.getNoisebyHildebrand(
625 noise = self.getNoisebyHildebrand(
621 xmin_index, xmax_index, ymin_index, ymax_index)
626 xmin_index, xmax_index, ymin_index, ymax_index)
622 return noise
627 return noise
623
628
624 def getFreqRangeTimeResponse(self, extrapoints=0):
629 def getFreqRangeTimeResponse(self, extrapoints=0):
625
630
626 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
631 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
627 freqrange = deltafreq * \
632 freqrange = deltafreq * \
628 (numpy.arange(self.nFFTPoints + extrapoints) -
633 (numpy.arange(self.nFFTPoints + extrapoints) -
629 self.nFFTPoints / 2.) - deltafreq / 2
634 self.nFFTPoints / 2.) - deltafreq / 2
630
635
631 return freqrange
636 return freqrange
632
637
633 def getAcfRange(self, extrapoints=0):
638 def getAcfRange(self, extrapoints=0):
639 #print "NFFTPoints",self.nFFTPoints
640 #print "IPPFactor", self.ippFactor
634 deltafreq = 10. / ( self.getFmax() / (self.nFFTPoints * self.ippFactor) )
641 deltafreq = 10. / ( self.getFmax() / (self.nFFTPoints * self.ippFactor) )
642 #print "deltafreq",deltafreq
643 #import time
644 #time.sleep(30)
635 freqrange = deltafreq * \
645 freqrange = deltafreq * \
636 (numpy.arange(self.nFFTPoints + extrapoints) -
646 (numpy.arange(self.nFFTPoints + extrapoints) -
637 self.nFFTPoints / 2.) - deltafreq / 2
647 self.nFFTPoints / 2.) - deltafreq / 2
638
648
639 return freqrange
649 return freqrange
640
650
641 def getFreqRange(self, extrapoints=0):
651 def getFreqRange(self, extrapoints=0):
642
652
643 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
653 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
644 #print "deltafreq", deltafreq
654 #print "deltafreq", deltafreq
645 freqrange = deltafreq * \
655 freqrange = deltafreq * \
646 (numpy.arange(self.nFFTPoints + extrapoints) -
656 (numpy.arange(self.nFFTPoints + extrapoints) -
647 self.nFFTPoints / 2.) - deltafreq / 2
657 self.nFFTPoints / 2.) - deltafreq / 2
648 #print "freqrange",freqrange
658 #print "freqrange",freqrange
649 return freqrange
659 return freqrange
650
660
651 def getVelRange(self, extrapoints=0):
661 def getVelRange(self, extrapoints=0):
652
662
653 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
663 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
654 velrange = deltav * (numpy.arange(self.nFFTPoints +
664 velrange = deltav * (numpy.arange(self.nFFTPoints +
655 extrapoints) - self.nFFTPoints / 2.) # - deltav/2
665 extrapoints) - self.nFFTPoints / 2.) # - deltav/2
656
666
657 return velrange
667 return velrange
658
668
659 def getNPairs(self):
669 def getNPairs(self):
660
670
661 return len(self.pairsList)
671 return len(self.pairsList)
662
672
663 def getPairsIndexList(self):
673 def getPairsIndexList(self):
664
674
665 return range(self.nPairs)
675 return range(self.nPairs)
666
676
667 def getNormFactor(self):
677 def getNormFactor(self):
668
678
669 pwcode = 1
679 pwcode = 1
670
680
671 if self.flagDecodeData:
681 if self.flagDecodeData:
672 pwcode = numpy.sum(self.code[0]**2)
682 pwcode = numpy.sum(self.code[0]**2)
673 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
683 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
674 normFactor = self.nProfiles * self.nIncohInt * \
684 normFactor = self.nProfiles * self.nIncohInt * \
675 self.nCohInt * pwcode * self.windowOfFilter
685 self.nCohInt * pwcode * self.windowOfFilter
676
686
677 return normFactor
687 return normFactor
678
688
679 def getFlagCspc(self):
689 def getFlagCspc(self):
680
690
681 if self.data_cspc is None:
691 if self.data_cspc is None:
682 return True
692 return True
683
693
684 return False
694 return False
685
695
686 def getFlagDc(self):
696 def getFlagDc(self):
687
697
688 if self.data_dc is None:
698 if self.data_dc is None:
689 return True
699 return True
690
700
691 return False
701 return False
692
702
693 def getTimeInterval(self):
703 def getTimeInterval(self):
694
704
695 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
705 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
696
706
697 return timeInterval
707 return timeInterval
698
708
699 def getPower(self):
709 def getPower(self):
700
710
701 factor = self.normFactor
711 factor = self.normFactor
702 z = self.data_spc / factor
712 z = self.data_spc / factor
703 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
713 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
704 avg = numpy.average(z, axis=1)
714 avg = numpy.average(z, axis=1)
705
715
706 return 10 * numpy.log10(avg)
716 return 10 * numpy.log10(avg)
707
717
708 def getCoherence(self, pairsList=None, phase=False):
718 def getCoherence(self, pairsList=None, phase=False):
709
719
710 z = []
720 z = []
711 if pairsList is None:
721 if pairsList is None:
712 pairsIndexList = self.pairsIndexList
722 pairsIndexList = self.pairsIndexList
713 else:
723 else:
714 pairsIndexList = []
724 pairsIndexList = []
715 for pair in pairsList:
725 for pair in pairsList:
716 if pair not in self.pairsList:
726 if pair not in self.pairsList:
717 raise ValueError, "Pair %s is not in dataOut.pairsList" % (
727 raise ValueError, "Pair %s is not in dataOut.pairsList" % (
718 pair)
728 pair)
719 pairsIndexList.append(self.pairsList.index(pair))
729 pairsIndexList.append(self.pairsList.index(pair))
720 for i in range(len(pairsIndexList)):
730 for i in range(len(pairsIndexList)):
721 pair = self.pairsList[pairsIndexList[i]]
731 pair = self.pairsList[pairsIndexList[i]]
722 ccf = numpy.average(
732 ccf = numpy.average(
723 self.data_cspc[pairsIndexList[i], :, :], axis=0)
733 self.data_cspc[pairsIndexList[i], :, :], axis=0)
724 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
734 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
725 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
735 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
726 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
736 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
727 if phase:
737 if phase:
728 data = numpy.arctan2(avgcoherenceComplex.imag,
738 data = numpy.arctan2(avgcoherenceComplex.imag,
729 avgcoherenceComplex.real) * 180 / numpy.pi
739 avgcoherenceComplex.real) * 180 / numpy.pi
730 else:
740 else:
731 data = numpy.abs(avgcoherenceComplex)
741 data = numpy.abs(avgcoherenceComplex)
732
742
733 z.append(data)
743 z.append(data)
734
744
735 return numpy.array(z)
745 return numpy.array(z)
736
746
737 def setValue(self, value):
747 def setValue(self, value):
738
748
739 print "This property should not be initialized"
749 print "This property should not be initialized"
740
750
741 return
751 return
742
752
743 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
753 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
744 pairsIndexList = property(
754 pairsIndexList = property(
745 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
755 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
746 normFactor = property(getNormFactor, setValue,
756 normFactor = property(getNormFactor, setValue,
747 "I'm the 'getNormFactor' property.")
757 "I'm the 'getNormFactor' property.")
748 flag_cspc = property(getFlagCspc, setValue)
758 flag_cspc = property(getFlagCspc, setValue)
749 flag_dc = property(getFlagDc, setValue)
759 flag_dc = property(getFlagDc, setValue)
750 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
760 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
751 timeInterval = property(getTimeInterval, setValue,
761 timeInterval = property(getTimeInterval, setValue,
752 "I'm the 'timeInterval' property")
762 "I'm the 'timeInterval' property")
753
763
754
764
755 class SpectraHeis(Spectra):
765 class SpectraHeis(Spectra):
756
766
757 data_spc = None
767 data_spc = None
758
768
759 data_cspc = None
769 data_cspc = None
760
770
761 data_dc = None
771 data_dc = None
762
772
763 nFFTPoints = None
773 nFFTPoints = None
764
774
765 # nPairs = None
775 # nPairs = None
766
776
767 pairsList = None
777 pairsList = None
768
778
769 nCohInt = None
779 nCohInt = None
770
780
771 nIncohInt = None
781 nIncohInt = None
772
782
773 def __init__(self):
783 def __init__(self):
774
784
775 self.radarControllerHeaderObj = RadarControllerHeader()
785 self.radarControllerHeaderObj = RadarControllerHeader()
776
786
777 self.systemHeaderObj = SystemHeader()
787 self.systemHeaderObj = SystemHeader()
778
788
779 self.type = "SpectraHeis"
789 self.type = "SpectraHeis"
780
790
781 # self.dtype = None
791 # self.dtype = None
782
792
783 # self.nChannels = 0
793 # self.nChannels = 0
784
794
785 # self.nHeights = 0
795 # self.nHeights = 0
786
796
787 self.nProfiles = None
797 self.nProfiles = None
788
798
789 self.heightList = None
799 self.heightList = None
790
800
791 self.channelList = None
801 self.channelList = None
792
802
793 # self.channelIndexList = None
803 # self.channelIndexList = None
794
804
795 self.flagNoData = True
805 self.flagNoData = True
796
806
797 self.flagDiscontinuousBlock = False
807 self.flagDiscontinuousBlock = False
798
808
799 # self.nPairs = 0
809 # self.nPairs = 0
800
810
801 self.utctime = None
811 self.utctime = None
802
812
803 self.blocksize = None
813 self.blocksize = None
804
814
805 self.profileIndex = 0
815 self.profileIndex = 0
806
816
807 self.nCohInt = 1
817 self.nCohInt = 1
808
818
809 self.nIncohInt = 1
819 self.nIncohInt = 1
810
820
811 def getNormFactor(self):
821 def getNormFactor(self):
812 pwcode = 1
822 pwcode = 1
813 if self.flagDecodeData:
823 if self.flagDecodeData:
814 pwcode = numpy.sum(self.code[0]**2)
824 pwcode = numpy.sum(self.code[0]**2)
815
825
816 normFactor = self.nIncohInt * self.nCohInt * pwcode
826 normFactor = self.nIncohInt * self.nCohInt * pwcode
817
827
818 return normFactor
828 return normFactor
819
829
820 def getTimeInterval(self):
830 def getTimeInterval(self):
821
831
822 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
832 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
823
833
824 return timeInterval
834 return timeInterval
825
835
826 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
836 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
827 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
837 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
828
838
829
839
830 class Fits(JROData):
840 class Fits(JROData):
831
841
832 heightList = None
842 heightList = None
833
843
834 channelList = None
844 channelList = None
835
845
836 flagNoData = True
846 flagNoData = True
837
847
838 flagDiscontinuousBlock = False
848 flagDiscontinuousBlock = False
839
849
840 useLocalTime = False
850 useLocalTime = False
841
851
842 utctime = None
852 utctime = None
843
853
844 timeZone = None
854 timeZone = None
845
855
846 # ippSeconds = None
856 # ippSeconds = None
847
857
848 # timeInterval = None
858 # timeInterval = None
849
859
850 nCohInt = None
860 nCohInt = None
851
861
852 nIncohInt = None
862 nIncohInt = None
853
863
854 noise = None
864 noise = None
855
865
856 windowOfFilter = 1
866 windowOfFilter = 1
857
867
858 # Speed of ligth
868 # Speed of ligth
859 C = 3e8
869 C = 3e8
860
870
861 frequency = 49.92e6
871 frequency = 49.92e6
862
872
863 realtime = False
873 realtime = False
864
874
865 def __init__(self):
875 def __init__(self):
866
876
867 self.type = "Fits"
877 self.type = "Fits"
868
878
869 self.nProfiles = None
879 self.nProfiles = None
870
880
871 self.heightList = None
881 self.heightList = None
872
882
873 self.channelList = None
883 self.channelList = None
874
884
875 # self.channelIndexList = None
885 # self.channelIndexList = None
876
886
877 self.flagNoData = True
887 self.flagNoData = True
878
888
879 self.utctime = None
889 self.utctime = None
880
890
881 self.nCohInt = 1
891 self.nCohInt = 1
882
892
883 self.nIncohInt = 1
893 self.nIncohInt = 1
884
894
885 self.useLocalTime = True
895 self.useLocalTime = True
886
896
887 self.profileIndex = 0
897 self.profileIndex = 0
888
898
889 # self.utctime = None
899 # self.utctime = None
890 # self.timeZone = None
900 # self.timeZone = None
891 # self.ltctime = None
901 # self.ltctime = None
892 # self.timeInterval = None
902 # self.timeInterval = None
893 # self.header = None
903 # self.header = None
894 # self.data_header = None
904 # self.data_header = None
895 # self.data = None
905 # self.data = None
896 # self.datatime = None
906 # self.datatime = None
897 # self.flagNoData = False
907 # self.flagNoData = False
898 # self.expName = ''
908 # self.expName = ''
899 # self.nChannels = None
909 # self.nChannels = None
900 # self.nSamples = None
910 # self.nSamples = None
901 # self.dataBlocksPerFile = None
911 # self.dataBlocksPerFile = None
902 # self.comments = ''
912 # self.comments = ''
903 #
913 #
904
914
905 def getltctime(self):
915 def getltctime(self):
906
916
907 if self.useLocalTime:
917 if self.useLocalTime:
908 return self.utctime - self.timeZone * 60
918 return self.utctime - self.timeZone * 60
909
919
910 return self.utctime
920 return self.utctime
911
921
912 def getDatatime(self):
922 def getDatatime(self):
913
923
914 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
924 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
915 return datatime
925 return datatime
916
926
917 def getTimeRange(self):
927 def getTimeRange(self):
918
928
919 datatime = []
929 datatime = []
920
930
921 datatime.append(self.ltctime)
931 datatime.append(self.ltctime)
922 datatime.append(self.ltctime + self.timeInterval)
932 datatime.append(self.ltctime + self.timeInterval)
923
933
924 datatime = numpy.array(datatime)
934 datatime = numpy.array(datatime)
925
935
926 return datatime
936 return datatime
927
937
928 def getHeiRange(self):
938 def getHeiRange(self):
929
939
930 heis = self.heightList
940 heis = self.heightList
931
941
932 return heis
942 return heis
933
943
934 def getNHeights(self):
944 def getNHeights(self):
935
945
936 return len(self.heightList)
946 return len(self.heightList)
937
947
938 def getNChannels(self):
948 def getNChannels(self):
939
949
940 return len(self.channelList)
950 return len(self.channelList)
941
951
942 def getChannelIndexList(self):
952 def getChannelIndexList(self):
943
953
944 return range(self.nChannels)
954 return range(self.nChannels)
945
955
946 def getNoise(self, type=1):
956 def getNoise(self, type=1):
947
957
948 #noise = numpy.zeros(self.nChannels)
958 #noise = numpy.zeros(self.nChannels)
949
959
950 if type == 1:
960 if type == 1:
951 noise = self.getNoisebyHildebrand()
961 noise = self.getNoisebyHildebrand()
952
962
953 if type == 2:
963 if type == 2:
954 noise = self.getNoisebySort()
964 noise = self.getNoisebySort()
955
965
956 if type == 3:
966 if type == 3:
957 noise = self.getNoisebyWindow()
967 noise = self.getNoisebyWindow()
958
968
959 return noise
969 return noise
960
970
961 def getTimeInterval(self):
971 def getTimeInterval(self):
962
972
963 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
973 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
964
974
965 return timeInterval
975 return timeInterval
966
976
967 datatime = property(getDatatime, "I'm the 'datatime' property")
977 datatime = property(getDatatime, "I'm the 'datatime' property")
968 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
978 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
969 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
979 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
970 channelIndexList = property(
980 channelIndexList = property(
971 getChannelIndexList, "I'm the 'channelIndexList' property.")
981 getChannelIndexList, "I'm the 'channelIndexList' property.")
972 noise = property(getNoise, "I'm the 'nHeights' property.")
982 noise = property(getNoise, "I'm the 'nHeights' property.")
973
983
974 ltctime = property(getltctime, "I'm the 'ltctime' property")
984 ltctime = property(getltctime, "I'm the 'ltctime' property")
975 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
985 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
976
986
977
987
978 class Correlation(JROData):
988 class Correlation(JROData):
979
989
980 noise = None
990 noise = None
981
991
982 SNR = None
992 SNR = None
983
993
984 #--------------------------------------------------
994 #--------------------------------------------------
985
995
986 mode = None
996 mode = None
987
997
988 split = False
998 split = False
989
999
990 data_cf = None
1000 data_cf = None
991
1001
992 lags = None
1002 lags = None
993
1003
994 lagRange = None
1004 lagRange = None
995
1005
996 pairsList = None
1006 pairsList = None
997
1007
998 normFactor = None
1008 normFactor = None
999
1009
1000 #--------------------------------------------------
1010 #--------------------------------------------------
1001
1011
1002 # calculateVelocity = None
1012 # calculateVelocity = None
1003
1013
1004 nLags = None
1014 nLags = None
1005
1015
1006 nPairs = None
1016 nPairs = None
1007
1017
1008 nAvg = None
1018 nAvg = None
1009
1019
1010 def __init__(self):
1020 def __init__(self):
1011 '''
1021 '''
1012 Constructor
1022 Constructor
1013 '''
1023 '''
1014 self.radarControllerHeaderObj = RadarControllerHeader()
1024 self.radarControllerHeaderObj = RadarControllerHeader()
1015
1025
1016 self.systemHeaderObj = SystemHeader()
1026 self.systemHeaderObj = SystemHeader()
1017
1027
1018 self.type = "Correlation"
1028 self.type = "Correlation"
1019
1029
1020 self.data = None
1030 self.data = None
1021
1031
1022 self.dtype = None
1032 self.dtype = None
1023
1033
1024 self.nProfiles = None
1034 self.nProfiles = None
1025
1035
1026 self.heightList = None
1036 self.heightList = None
1027
1037
1028 self.channelList = None
1038 self.channelList = None
1029
1039
1030 self.flagNoData = True
1040 self.flagNoData = True
1031
1041
1032 self.flagDiscontinuousBlock = False
1042 self.flagDiscontinuousBlock = False
1033
1043
1034 self.utctime = None
1044 self.utctime = None
1035
1045
1036 self.timeZone = None
1046 self.timeZone = None
1037
1047
1038 self.dstFlag = None
1048 self.dstFlag = None
1039
1049
1040 self.errorCount = None
1050 self.errorCount = None
1041
1051
1042 self.blocksize = None
1052 self.blocksize = None
1043
1053
1044 self.flagDecodeData = False # asumo q la data no esta decodificada
1054 self.flagDecodeData = False # asumo q la data no esta decodificada
1045
1055
1046 self.flagDeflipData = False # asumo q la data no esta sin flip
1056 self.flagDeflipData = False # asumo q la data no esta sin flip
1047
1057
1048 self.pairsList = None
1058 self.pairsList = None
1049
1059
1050 self.nPoints = None
1060 self.nPoints = None
1051
1061
1052 def getPairsList(self):
1062 def getPairsList(self):
1053
1063
1054 return self.pairsList
1064 return self.pairsList
1055
1065
1056 def getNoise(self, mode=2):
1066 def getNoise(self, mode=2):
1057
1067
1058 indR = numpy.where(self.lagR == 0)[0][0]
1068 indR = numpy.where(self.lagR == 0)[0][0]
1059 indT = numpy.where(self.lagT == 0)[0][0]
1069 indT = numpy.where(self.lagT == 0)[0][0]
1060
1070
1061 jspectra0 = self.data_corr[:, :, indR, :]
1071 jspectra0 = self.data_corr[:, :, indR, :]
1062 jspectra = copy.copy(jspectra0)
1072 jspectra = copy.copy(jspectra0)
1063
1073
1064 num_chan = jspectra.shape[0]
1074 num_chan = jspectra.shape[0]
1065 num_hei = jspectra.shape[2]
1075 num_hei = jspectra.shape[2]
1066
1076
1067 freq_dc = jspectra.shape[1] / 2
1077 freq_dc = jspectra.shape[1] / 2
1068 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
1078 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
1069
1079
1070 if ind_vel[0] < 0:
1080 if ind_vel[0] < 0:
1071 ind_vel[range(0, 1)] = ind_vel[range(0, 1)] + self.num_prof
1081 ind_vel[range(0, 1)] = ind_vel[range(0, 1)] + self.num_prof
1072
1082
1073 if mode == 1:
1083 if mode == 1:
1074 jspectra[:, freq_dc, :] = (
1084 jspectra[:, freq_dc, :] = (
1075 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
1085 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
1076
1086
1077 if mode == 2:
1087 if mode == 2:
1078
1088
1079 vel = numpy.array([-2, -1, 1, 2])
1089 vel = numpy.array([-2, -1, 1, 2])
1080 xx = numpy.zeros([4, 4])
1090 xx = numpy.zeros([4, 4])
1081
1091
1082 for fil in range(4):
1092 for fil in range(4):
1083 xx[fil, :] = vel[fil]**numpy.asarray(range(4))
1093 xx[fil, :] = vel[fil]**numpy.asarray(range(4))
1084
1094
1085 xx_inv = numpy.linalg.inv(xx)
1095 xx_inv = numpy.linalg.inv(xx)
1086 xx_aux = xx_inv[0, :]
1096 xx_aux = xx_inv[0, :]
1087
1097
1088 for ich in range(num_chan):
1098 for ich in range(num_chan):
1089 yy = jspectra[ich, ind_vel, :]
1099 yy = jspectra[ich, ind_vel, :]
1090 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
1100 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
1091
1101
1092 junkid = jspectra[ich, freq_dc, :] <= 0
1102 junkid = jspectra[ich, freq_dc, :] <= 0
1093 cjunkid = sum(junkid)
1103 cjunkid = sum(junkid)
1094
1104
1095 if cjunkid.any():
1105 if cjunkid.any():
1096 jspectra[ich, freq_dc, junkid.nonzero()] = (
1106 jspectra[ich, freq_dc, junkid.nonzero()] = (
1097 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
1107 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
1098
1108
1099 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
1109 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
1100
1110
1101 return noise
1111 return noise
1102
1112
1103 def getTimeInterval(self):
1113 def getTimeInterval(self):
1104
1114
1105 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
1115 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
1106
1116
1107 return timeInterval
1117 return timeInterval
1108
1118
1109 def splitFunctions(self):
1119 def splitFunctions(self):
1110
1120
1111 pairsList = self.pairsList
1121 pairsList = self.pairsList
1112 ccf_pairs = []
1122 ccf_pairs = []
1113 acf_pairs = []
1123 acf_pairs = []
1114 ccf_ind = []
1124 ccf_ind = []
1115 acf_ind = []
1125 acf_ind = []
1116 for l in range(len(pairsList)):
1126 for l in range(len(pairsList)):
1117 chan0 = pairsList[l][0]
1127 chan0 = pairsList[l][0]
1118 chan1 = pairsList[l][1]
1128 chan1 = pairsList[l][1]
1119
1129
1120 # Obteniendo pares de Autocorrelacion
1130 # Obteniendo pares de Autocorrelacion
1121 if chan0 == chan1:
1131 if chan0 == chan1:
1122 acf_pairs.append(chan0)
1132 acf_pairs.append(chan0)
1123 acf_ind.append(l)
1133 acf_ind.append(l)
1124 else:
1134 else:
1125 ccf_pairs.append(pairsList[l])
1135 ccf_pairs.append(pairsList[l])
1126 ccf_ind.append(l)
1136 ccf_ind.append(l)
1127
1137
1128 data_acf = self.data_cf[acf_ind]
1138 data_acf = self.data_cf[acf_ind]
1129 data_ccf = self.data_cf[ccf_ind]
1139 data_ccf = self.data_cf[ccf_ind]
1130
1140
1131 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1141 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1132
1142
1133 def getNormFactor(self):
1143 def getNormFactor(self):
1134 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1144 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1135 acf_pairs = numpy.array(acf_pairs)
1145 acf_pairs = numpy.array(acf_pairs)
1136 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1146 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1137
1147
1138 for p in range(self.nPairs):
1148 for p in range(self.nPairs):
1139 pair = self.pairsList[p]
1149 pair = self.pairsList[p]
1140
1150
1141 ch0 = pair[0]
1151 ch0 = pair[0]
1142 ch1 = pair[1]
1152 ch1 = pair[1]
1143
1153
1144 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1154 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1145 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1155 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1146 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1156 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1147
1157
1148 return normFactor
1158 return normFactor
1149
1159
1150 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1160 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1151 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1161 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1152
1162
1153
1163
1154 class Parameters(Spectra):
1164 class Parameters(Spectra):
1155
1165
1156 experimentInfo = None # Information about the experiment
1166 experimentInfo = None # Information about the experiment
1157
1167
1158 # Information from previous data
1168 # Information from previous data
1159
1169
1160 inputUnit = None # Type of data to be processed
1170 inputUnit = None # Type of data to be processed
1161
1171
1162 operation = None # Type of operation to parametrize
1172 operation = None # Type of operation to parametrize
1163
1173
1164 # normFactor = None #Normalization Factor
1174 # normFactor = None #Normalization Factor
1165
1175
1166 groupList = None # List of Pairs, Groups, etc
1176 groupList = None # List of Pairs, Groups, etc
1167
1177
1168 # Parameters
1178 # Parameters
1169
1179
1170 data_param = None # Parameters obtained
1180 data_param = None # Parameters obtained
1171
1181
1172 data_pre = None # Data Pre Parametrization
1182 data_pre = None # Data Pre Parametrization
1173
1183
1174 data_SNR = None # Signal to Noise Ratio
1184 data_SNR = None # Signal to Noise Ratio
1175
1185
1176 # heightRange = None #Heights
1186 # heightRange = None #Heights
1177
1187
1178 abscissaList = None # Abscissa, can be velocities, lags or time
1188 abscissaList = None # Abscissa, can be velocities, lags or time
1179
1189
1180 # noise = None #Noise Potency
1190 # noise = None #Noise Potency
1181
1191
1182 utctimeInit = None # Initial UTC time
1192 utctimeInit = None # Initial UTC time
1183
1193
1184 paramInterval = None # Time interval to calculate Parameters in seconds
1194 paramInterval = None # Time interval to calculate Parameters in seconds
1185
1195
1186 useLocalTime = True
1196 useLocalTime = True
1187
1197
1188 # Fitting
1198 # Fitting
1189
1199
1190 data_error = None # Error of the estimation
1200 data_error = None # Error of the estimation
1191
1201
1192 constants = None
1202 constants = None
1193
1203
1194 library = None
1204 library = None
1195
1205
1196 # Output signal
1206 # Output signal
1197
1207
1198 outputInterval = None # Time interval to calculate output signal in seconds
1208 outputInterval = None # Time interval to calculate output signal in seconds
1199
1209
1200 data_output = None # Out signal
1210 data_output = None # Out signal
1201
1211
1202 nAvg = None
1212 nAvg = None
1203
1213
1204 noise_estimation = None
1214 noise_estimation = None
1205
1215
1206 GauSPC = None # Fit gaussian SPC
1216 GauSPC = None # Fit gaussian SPC
1207
1217
1208 def __init__(self):
1218 def __init__(self):
1209 '''
1219 '''
1210 Constructor
1220 Constructor
1211 '''
1221 '''
1212 self.radarControllerHeaderObj = RadarControllerHeader()
1222 self.radarControllerHeaderObj = RadarControllerHeader()
1213
1223
1214 self.systemHeaderObj = SystemHeader()
1224 self.systemHeaderObj = SystemHeader()
1215
1225
1216 self.type = "Parameters"
1226 self.type = "Parameters"
1217
1227
1218 def getTimeRange1(self, interval):
1228 def getTimeRange1(self, interval):
1219
1229
1220 datatime = []
1230 datatime = []
1221
1231
1222 if self.useLocalTime:
1232 if self.useLocalTime:
1223 time1 = self.utctimeInit - self.timeZone * 60
1233 time1 = self.utctimeInit - self.timeZone * 60
1224 else:
1234 else:
1225 time1 = self.utctimeInit
1235 time1 = self.utctimeInit
1226
1236
1227 datatime.append(time1)
1237 datatime.append(time1)
1228 datatime.append(time1 + interval)
1238 datatime.append(time1 + interval)
1229 datatime = numpy.array(datatime)
1239 datatime = numpy.array(datatime)
1230
1240
1231 return datatime
1241 return datatime
1232
1242
1233 def getTimeInterval(self):
1243 def getTimeInterval(self):
1234
1244
1235 if hasattr(self, 'timeInterval1'):
1245 if hasattr(self, 'timeInterval1'):
1236 return self.timeInterval1
1246 return self.timeInterval1
1237 else:
1247 else:
1238 return self.paramInterval
1248 return self.paramInterval
1239
1249
1240 def setValue(self, value):
1250 def setValue(self, value):
1241
1251
1242 print "This property should not be initialized"
1252 print "This property should not be initialized"
1243
1253
1244 return
1254 return
1245
1255
1246 def getNoise(self):
1256 def getNoise(self):
1247
1257
1248 return self.spc_noise
1258 return self.spc_noise
1249
1259
1250 timeInterval = property(getTimeInterval)
1260 timeInterval = property(getTimeInterval)
1251 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
1261 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
@@ -1,1520 +1,1530
1 import sys
1 import sys
2 import numpy
2 import numpy
3 from scipy import interpolate
3 from scipy import interpolate
4 from schainpy import cSchain
4 from schainpy import cSchain
5 from jroproc_base import ProcessingUnit, Operation
5 from jroproc_base import ProcessingUnit, Operation
6 from schainpy.model.data.jrodata import Voltage
6 from schainpy.model.data.jrodata import Voltage
7 from time import time
7 from time import time
8
8
9 import math
9 import math
10
10
11 def rep_seq(x, rep=10):
11 def rep_seq(x, rep=10):
12 L = len(x) * rep
12 L = len(x) * rep
13 res = numpy.zeros(L, dtype=x.dtype)
13 res = numpy.zeros(L, dtype=x.dtype)
14 idx = numpy.arange(len(x)) * rep
14 idx = numpy.arange(len(x)) * rep
15 for i in numpy.arange(rep):
15 for i in numpy.arange(rep):
16 res[idx + i] = x
16 res[idx + i] = x
17 return(res)
17 return(res)
18
18
19
19
20 def create_pseudo_random_code(clen=10000, seed=0):
20 def create_pseudo_random_code(clen=10000, seed=0):
21 """
21 """
22 seed is a way of reproducing the random code without
22 seed is a way of reproducing the random code without
23 having to store all actual codes. the seed can then
23 having to store all actual codes. the seed can then
24 act as a sort of station_id.
24 act as a sort of station_id.
25
25
26 """
26 """
27 numpy.random.seed(seed)
27 numpy.random.seed(seed)
28 phases = numpy.array(
28 phases = numpy.array(
29 numpy.exp(1.0j * 2.0 * math.pi * numpy.random.random(clen)),
29 numpy.exp(1.0j * 2.0 * math.pi * numpy.random.random(clen)),
30 dtype=numpy.complex64,
30 dtype=numpy.complex64,
31 )
31 )
32 return(phases)
32 return(phases)
33
33
34
34
35 def periodic_convolution_matrix(envelope, rmin=0, rmax=100):
35 def periodic_convolution_matrix(envelope, rmin=0, rmax=100):
36 """
36 """
37 we imply that the number of measurements is equal to the number of elements
37 we imply that the number of measurements is equal to the number of elements
38 in code
38 in code
39
39
40 """
40 """
41 L = len(envelope)
41 L = len(envelope)
42 ridx = numpy.arange(rmin, rmax)
42 ridx = numpy.arange(rmin, rmax)
43 A = numpy.zeros([L, rmax-rmin], dtype=numpy.complex64)
43 A = numpy.zeros([L, rmax-rmin], dtype=numpy.complex64)
44 for i in numpy.arange(L):
44 for i in numpy.arange(L):
45 A[i, :] = envelope[(i-ridx) % L]
45 A[i, :] = envelope[(i-ridx) % L]
46 result = {}
46 result = {}
47 result['A'] = A
47 result['A'] = A
48 result['ridx'] = ridx
48 result['ridx'] = ridx
49 return(result)
49 return(result)
50
50
51
51
52 B_cache = 0
52 B_cache = 0
53 r_cache = 0
53 r_cache = 0
54 B_cached = False
54 B_cached = False
55 def create_estimation_matrix(code, rmin=0, rmax=1000, cache=True):
55 def create_estimation_matrix(code, rmin=0, rmax=1000, cache=True):
56 global B_cache
56 global B_cache
57 global r_cache
57 global r_cache
58 global B_cached
58 global B_cached
59
59
60 if not cache or not B_cached:
60 if not cache or not B_cached:
61 r_cache = periodic_convolution_matrix(
61 r_cache = periodic_convolution_matrix(
62 envelope=code, rmin=rmin, rmax=rmax,
62 envelope=code, rmin=rmin, rmax=rmax,
63 )
63 )
64 A = r_cache['A']
64 A = r_cache['A']
65 Ah = numpy.transpose(numpy.conjugate(A))
65 Ah = numpy.transpose(numpy.conjugate(A))
66 B_cache = numpy.dot(numpy.linalg.inv(numpy.dot(Ah, A)), Ah)
66 B_cache = numpy.dot(numpy.linalg.inv(numpy.dot(Ah, A)), Ah)
67 r_cache['B'] = B_cache
67 r_cache['B'] = B_cache
68 B_cached = True
68 B_cached = True
69 return(r_cache)
69 return(r_cache)
70 else:
70 else:
71 return(r_cache)
71 return(r_cache)
72
72
73 class VoltageProc(ProcessingUnit):
73 class VoltageProc(ProcessingUnit):
74
74
75
75
76 def __init__(self, **kwargs):
76 def __init__(self, **kwargs):
77
77
78 ProcessingUnit.__init__(self, **kwargs)
78 ProcessingUnit.__init__(self, **kwargs)
79
79
80 # self.objectDict = {}
80 # self.objectDict = {}
81 self.dataOut = Voltage()
81 self.dataOut = Voltage()
82 self.flip = 1
82 self.flip = 1
83
83
84 def run(self):
84 def run(self):
85 if self.dataIn.type == 'AMISR':
85 if self.dataIn.type == 'AMISR':
86 self.__updateObjFromAmisrInput()
86 self.__updateObjFromAmisrInput()
87
87
88 if self.dataIn.type == 'Voltage':
88 if self.dataIn.type == 'Voltage':
89 self.dataOut.copy(self.dataIn)
89 self.dataOut.copy(self.dataIn)
90
90
91 # self.dataOut.copy(self.dataIn)
91 # self.dataOut.copy(self.dataIn)
92
92
93 def __updateObjFromAmisrInput(self):
93 def __updateObjFromAmisrInput(self):
94
94
95 self.dataOut.timeZone = self.dataIn.timeZone
95 self.dataOut.timeZone = self.dataIn.timeZone
96 self.dataOut.dstFlag = self.dataIn.dstFlag
96 self.dataOut.dstFlag = self.dataIn.dstFlag
97 self.dataOut.errorCount = self.dataIn.errorCount
97 self.dataOut.errorCount = self.dataIn.errorCount
98 self.dataOut.useLocalTime = self.dataIn.useLocalTime
98 self.dataOut.useLocalTime = self.dataIn.useLocalTime
99
99
100 self.dataOut.flagNoData = self.dataIn.flagNoData
100 self.dataOut.flagNoData = self.dataIn.flagNoData
101 self.dataOut.data = self.dataIn.data
101 self.dataOut.data = self.dataIn.data
102 self.dataOut.utctime = self.dataIn.utctime
102 self.dataOut.utctime = self.dataIn.utctime
103 self.dataOut.channelList = self.dataIn.channelList
103 self.dataOut.channelList = self.dataIn.channelList
104 # self.dataOut.timeInterval = self.dataIn.timeInterval
104 # self.dataOut.timeInterval = self.dataIn.timeInterval
105 self.dataOut.heightList = self.dataIn.heightList
105 self.dataOut.heightList = self.dataIn.heightList
106 self.dataOut.nProfiles = self.dataIn.nProfiles
106 self.dataOut.nProfiles = self.dataIn.nProfiles
107
107
108 self.dataOut.nCohInt = self.dataIn.nCohInt
108 self.dataOut.nCohInt = self.dataIn.nCohInt
109 self.dataOut.ippSeconds = self.dataIn.ippSeconds
109 self.dataOut.ippSeconds = self.dataIn.ippSeconds
110 self.dataOut.frequency = self.dataIn.frequency
110 self.dataOut.frequency = self.dataIn.frequency
111
111
112 self.dataOut.azimuth = self.dataIn.azimuth
112 self.dataOut.azimuth = self.dataIn.azimuth
113 self.dataOut.zenith = self.dataIn.zenith
113 self.dataOut.zenith = self.dataIn.zenith
114
114
115 self.dataOut.beam.codeList = self.dataIn.beam.codeList
115 self.dataOut.beam.codeList = self.dataIn.beam.codeList
116 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
116 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
117 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
117 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
118 #
118 #
119 # pass#
119 # pass#
120 #
120 #
121 # def init(self):
121 # def init(self):
122 #
122 #
123 #
123 #
124 # if self.dataIn.type == 'AMISR':
124 # if self.dataIn.type == 'AMISR':
125 # self.__updateObjFromAmisrInput()
125 # self.__updateObjFromAmisrInput()
126 #
126 #
127 # if self.dataIn.type == 'Voltage':
127 # if self.dataIn.type == 'Voltage':
128 # self.dataOut.copy(self.dataIn)
128 # self.dataOut.copy(self.dataIn)
129 # # No necesita copiar en cada init() los atributos de dataIn
129 # # No necesita copiar en cada init() los atributos de dataIn
130 # # la copia deberia hacerse por cada nuevo bloque de datos
130 # # la copia deberia hacerse por cada nuevo bloque de datos
131
131
132 def selectChannels(self, channelList):
132 def selectChannels(self, channelList):
133
133
134 channelIndexList = []
134 channelIndexList = []
135
135
136 for channel in channelList:
136 for channel in channelList:
137 if channel not in self.dataOut.channelList:
137 if channel not in self.dataOut.channelList:
138 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
138 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
139
139
140 index = self.dataOut.channelList.index(channel)
140 index = self.dataOut.channelList.index(channel)
141 channelIndexList.append(index)
141 channelIndexList.append(index)
142
142
143 self.selectChannelsByIndex(channelIndexList)
143 self.selectChannelsByIndex(channelIndexList)
144
144
145 def selectChannelsByIndex(self, channelIndexList):
145 def selectChannelsByIndex(self, channelIndexList):
146 """
146 """
147 Selecciona un bloque de datos en base a canales segun el channelIndexList
147 Selecciona un bloque de datos en base a canales segun el channelIndexList
148
148
149 Input:
149 Input:
150 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
150 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
151
151
152 Affected:
152 Affected:
153 self.dataOut.data
153 self.dataOut.data
154 self.dataOut.channelIndexList
154 self.dataOut.channelIndexList
155 self.dataOut.nChannels
155 self.dataOut.nChannels
156 self.dataOut.m_ProcessingHeader.totalSpectra
156 self.dataOut.m_ProcessingHeader.totalSpectra
157 self.dataOut.systemHeaderObj.numChannels
157 self.dataOut.systemHeaderObj.numChannels
158 self.dataOut.m_ProcessingHeader.blockSize
158 self.dataOut.m_ProcessingHeader.blockSize
159
159
160 Return:
160 Return:
161 None
161 None
162 """
162 """
163
163
164 for channelIndex in channelIndexList:
164 for channelIndex in channelIndexList:
165 if channelIndex not in self.dataOut.channelIndexList:
165 if channelIndex not in self.dataOut.channelIndexList:
166 print channelIndexList
166 print channelIndexList
167 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
167 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
168
168
169 if self.dataOut.flagDataAsBlock:
169 if self.dataOut.flagDataAsBlock:
170 """
170 """
171 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
171 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
172 """
172 """
173 data = self.dataOut.data[channelIndexList,:,:]
173 data = self.dataOut.data[channelIndexList,:,:]
174 else:
174 else:
175 data = self.dataOut.data[channelIndexList,:]
175 data = self.dataOut.data[channelIndexList,:]
176
176
177 self.dataOut.data = data
177 self.dataOut.data = data
178 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
178 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
179 # self.dataOut.nChannels = nChannels
179 # self.dataOut.nChannels = nChannels
180
180
181 return 1
181 return 1
182
182
183 def selectHeights(self, minHei=None, maxHei=None):
183 def selectHeights(self, minHei=None, maxHei=None):
184 """
184 """
185 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
185 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
186 minHei <= height <= maxHei
186 minHei <= height <= maxHei
187
187
188 Input:
188 Input:
189 minHei : valor minimo de altura a considerar
189 minHei : valor minimo de altura a considerar
190 maxHei : valor maximo de altura a considerar
190 maxHei : valor maximo de altura a considerar
191
191
192 Affected:
192 Affected:
193 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
193 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
194
194
195 Return:
195 Return:
196 1 si el metodo se ejecuto con exito caso contrario devuelve 0
196 1 si el metodo se ejecuto con exito caso contrario devuelve 0
197 """
197 """
198
198
199 if minHei == None:
199 if minHei == None:
200 minHei = self.dataOut.heightList[0]
200 minHei = self.dataOut.heightList[0]
201
201
202 if maxHei == None:
202 if maxHei == None:
203 maxHei = self.dataOut.heightList[-1]
203 maxHei = self.dataOut.heightList[-1]
204
204
205 if (minHei < self.dataOut.heightList[0]):
205 if (minHei < self.dataOut.heightList[0]):
206 minHei = self.dataOut.heightList[0]
206 minHei = self.dataOut.heightList[0]
207
207
208 if (maxHei > self.dataOut.heightList[-1]):
208 if (maxHei > self.dataOut.heightList[-1]):
209 maxHei = self.dataOut.heightList[-1]
209 maxHei = self.dataOut.heightList[-1]
210
210
211 minIndex = 0
211 minIndex = 0
212 maxIndex = 0
212 maxIndex = 0
213 heights = self.dataOut.heightList
213 heights = self.dataOut.heightList
214
214
215 inda = numpy.where(heights >= minHei)
215 inda = numpy.where(heights >= minHei)
216 indb = numpy.where(heights <= maxHei)
216 indb = numpy.where(heights <= maxHei)
217
217
218 try:
218 try:
219 minIndex = inda[0][0]
219 minIndex = inda[0][0]
220 except:
220 except:
221 minIndex = 0
221 minIndex = 0
222
222
223 try:
223 try:
224 maxIndex = indb[0][-1]
224 maxIndex = indb[0][-1]
225 except:
225 except:
226 maxIndex = len(heights)
226 maxIndex = len(heights)
227
227
228 self.selectHeightsByIndex(minIndex, maxIndex)
228 self.selectHeightsByIndex(minIndex, maxIndex)
229
229
230 return 1
230 return 1
231
231
232
232
233 def selectHeightsByIndex(self, minIndex, maxIndex):
233 def selectHeightsByIndex(self, minIndex, maxIndex):
234 """
234 """
235 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
235 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
236 minIndex <= index <= maxIndex
236 minIndex <= index <= maxIndex
237
237
238 Input:
238 Input:
239 minIndex : valor de indice minimo de altura a considerar
239 minIndex : valor de indice minimo de altura a considerar
240 maxIndex : valor de indice maximo de altura a considerar
240 maxIndex : valor de indice maximo de altura a considerar
241
241
242 Affected:
242 Affected:
243 self.dataOut.data
243 self.dataOut.data
244 self.dataOut.heightList
244 self.dataOut.heightList
245
245
246 Return:
246 Return:
247 1 si el metodo se ejecuto con exito caso contrario devuelve 0
247 1 si el metodo se ejecuto con exito caso contrario devuelve 0
248 """
248 """
249
249
250 if (minIndex < 0) or (minIndex > maxIndex):
250 if (minIndex < 0) or (minIndex > maxIndex):
251 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
251 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
252
252
253 if (maxIndex >= self.dataOut.nHeights):
253 if (maxIndex >= self.dataOut.nHeights):
254 maxIndex = self.dataOut.nHeights
254 maxIndex = self.dataOut.nHeights
255
255
256 #voltage
256 #voltage
257 if self.dataOut.flagDataAsBlock:
257 if self.dataOut.flagDataAsBlock:
258 """
258 """
259 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
259 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
260 """
260 """
261 data = self.dataOut.data[:,:, minIndex:maxIndex]
261 data = self.dataOut.data[:,:, minIndex:maxIndex]
262 else:
262 else:
263 data = self.dataOut.data[:, minIndex:maxIndex]
263 data = self.dataOut.data[:, minIndex:maxIndex]
264
264
265 # firstHeight = self.dataOut.heightList[minIndex]
265 # firstHeight = self.dataOut.heightList[minIndex]
266
266
267 self.dataOut.data = data
267 self.dataOut.data = data
268 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
268 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
269
269
270 if self.dataOut.nHeights <= 1:
270 if self.dataOut.nHeights <= 1:
271 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
271 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
272
272
273 return 1
273 return 1
274
274
275
275
276 def filterByHeights(self, window):
276 def filterByHeights(self, window):
277
277
278 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
278 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
279
279
280 if window == None:
280 if window == None:
281 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
281 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
282
282
283 newdelta = deltaHeight * window
283 newdelta = deltaHeight * window
284 r = self.dataOut.nHeights % window
284 r = self.dataOut.nHeights % window
285 newheights = (self.dataOut.nHeights-r)/window
285 newheights = (self.dataOut.nHeights-r)/window
286
286
287 if newheights <= 1:
287 if newheights <= 1:
288 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
288 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
289
289
290 if self.dataOut.flagDataAsBlock:
290 if self.dataOut.flagDataAsBlock:
291 """
291 """
292 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
292 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
293 """
293 """
294 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
294 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
295 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
295 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
296 buffer = numpy.sum(buffer,3)
296 buffer = numpy.sum(buffer,3)
297
297
298 else:
298 else:
299 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
299 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
300 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
300 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
301 buffer = numpy.sum(buffer,2)
301 buffer = numpy.sum(buffer,2)
302
302
303 self.dataOut.data = buffer
303 self.dataOut.data = buffer
304 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
304 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
305 self.dataOut.windowOfFilter = window
305 self.dataOut.windowOfFilter = window
306
306
307 def setH0(self, h0, deltaHeight = None):
307 def setH0(self, h0, deltaHeight = None):
308
308
309 if not deltaHeight:
309 if not deltaHeight:
310 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
310 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
311
311
312 nHeights = self.dataOut.nHeights
312 nHeights = self.dataOut.nHeights
313
313
314 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
314 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
315
315
316 self.dataOut.heightList = newHeiRange
316 self.dataOut.heightList = newHeiRange
317
317
318 def deFlip(self, channelList = []):
318 def deFlip(self, channelList = []):
319
319
320 data = self.dataOut.data.copy()
320 data = self.dataOut.data.copy()
321
321
322 if self.dataOut.flagDataAsBlock:
322 if self.dataOut.flagDataAsBlock:
323 flip = self.flip
323 flip = self.flip
324 profileList = range(self.dataOut.nProfiles)
324 profileList = range(self.dataOut.nProfiles)
325
325
326 if not channelList:
326 if not channelList:
327 for thisProfile in profileList:
327 for thisProfile in profileList:
328 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
328 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
329 flip *= -1.0
329 flip *= -1.0
330 else:
330 else:
331 for thisChannel in channelList:
331 for thisChannel in channelList:
332 if thisChannel not in self.dataOut.channelList:
332 if thisChannel not in self.dataOut.channelList:
333 continue
333 continue
334
334
335 for thisProfile in profileList:
335 for thisProfile in profileList:
336 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
336 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
337 flip *= -1.0
337 flip *= -1.0
338
338
339 self.flip = flip
339 self.flip = flip
340
340
341 else:
341 else:
342 if not channelList:
342 if not channelList:
343 data[:,:] = data[:,:]*self.flip
343 data[:,:] = data[:,:]*self.flip
344 else:
344 else:
345 for thisChannel in channelList:
345 for thisChannel in channelList:
346 if thisChannel not in self.dataOut.channelList:
346 if thisChannel not in self.dataOut.channelList:
347 continue
347 continue
348
348
349 data[thisChannel,:] = data[thisChannel,:]*self.flip
349 data[thisChannel,:] = data[thisChannel,:]*self.flip
350
350
351 self.flip *= -1.
351 self.flip *= -1.
352
352
353 self.dataOut.data = data
353 self.dataOut.data = data
354
354
355 def setRadarFrequency(self, frequency=None):
355 def setRadarFrequency(self, frequency=None):
356
356
357 if frequency != None:
357 if frequency != None:
358 self.dataOut.frequency = frequency
358 self.dataOut.frequency = frequency
359
359
360 return 1
360 return 1
361
361
362 def interpolateHeights(self, topLim, botLim):
362 def interpolateHeights(self, topLim, botLim):
363 #69 al 72 para julia
363 #69 al 72 para julia
364 #82-84 para meteoros
364 #82-84 para meteoros
365 if len(numpy.shape(self.dataOut.data))==2:
365 if len(numpy.shape(self.dataOut.data))==2:
366 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
366 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
367 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
367 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
368 #self.dataOut.data[:,botLim:limSup+1] = sampInterp
368 #self.dataOut.data[:,botLim:limSup+1] = sampInterp
369 self.dataOut.data[:,botLim:topLim+1] = sampInterp
369 self.dataOut.data[:,botLim:topLim+1] = sampInterp
370 else:
370 else:
371 nHeights = self.dataOut.data.shape[2]
371 nHeights = self.dataOut.data.shape[2]
372 x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
372 x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
373 y = self.dataOut.data[:,:,range(botLim)+range(topLim+1,nHeights)]
373 y = self.dataOut.data[:,:,range(botLim)+range(topLim+1,nHeights)]
374 f = interpolate.interp1d(x, y, axis = 2)
374 f = interpolate.interp1d(x, y, axis = 2)
375 xnew = numpy.arange(botLim,topLim+1)
375 xnew = numpy.arange(botLim,topLim+1)
376 ynew = f(xnew)
376 ynew = f(xnew)
377
377
378 self.dataOut.data[:,:,botLim:topLim+1] = ynew
378 self.dataOut.data[:,:,botLim:topLim+1] = ynew
379
379
380 # import collections
380 # import collections
381
381
382 class CohInt(Operation):
382 class CohInt(Operation):
383
383
384 isConfig = False
384 isConfig = False
385 __profIndex = 0
385 __profIndex = 0
386 __byTime = False
386 __byTime = False
387 __initime = None
387 __initime = None
388 __lastdatatime = None
388 __lastdatatime = None
389 __integrationtime = None
389 __integrationtime = None
390 __buffer = None
390 __buffer = None
391 __bufferStride = []
391 __bufferStride = []
392 __dataReady = False
392 __dataReady = False
393 __profIndexStride = 0
393 __profIndexStride = 0
394 __dataToPutStride = False
394 __dataToPutStride = False
395 n = None
395 n = None
396
396
397 def __init__(self, **kwargs):
397 def __init__(self, **kwargs):
398
398
399 Operation.__init__(self, **kwargs)
399 Operation.__init__(self, **kwargs)
400
400
401 # self.isConfig = False
401 # self.isConfig = False
402
402
403 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False):
403 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False):
404 """
404 """
405 Set the parameters of the integration class.
405 Set the parameters of the integration class.
406
406
407 Inputs:
407 Inputs:
408
408
409 n : Number of coherent integrations
409 n : Number of coherent integrations
410 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
410 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
411 overlapping :
411 overlapping :
412 """
412 """
413
413
414 self.__initime = None
414 self.__initime = None
415 self.__lastdatatime = 0
415 self.__lastdatatime = 0
416 self.__buffer = None
416 self.__buffer = None
417 self.__dataReady = False
417 self.__dataReady = False
418 self.byblock = byblock
418 self.byblock = byblock
419 self.stride = stride
419 self.stride = stride
420
420
421 if n == None and timeInterval == None:
421 if n == None and timeInterval == None:
422 raise ValueError, "n or timeInterval should be specified ..."
422 raise ValueError, "n or timeInterval should be specified ..."
423
423
424 if n != None:
424 if n != None:
425 self.n = n
425 self.n = n
426 self.__byTime = False
426 self.__byTime = False
427 else:
427 else:
428 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
428 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
429 self.n = 9999
429 self.n = 9999
430 self.__byTime = True
430 self.__byTime = True
431
431
432 if overlapping:
432 if overlapping:
433 self.__withOverlapping = True
433 self.__withOverlapping = True
434 self.__buffer = None
434 self.__buffer = None
435 else:
435 else:
436 self.__withOverlapping = False
436 self.__withOverlapping = False
437 self.__buffer = 0
437 self.__buffer = 0
438
438
439 self.__profIndex = 0
439 self.__profIndex = 0
440
440
441 def putData(self, data):
441 def putData(self, data):
442
442
443 """
443 """
444 Add a profile to the __buffer and increase in one the __profileIndex
444 Add a profile to the __buffer and increase in one the __profileIndex
445
445
446 """
446 """
447
447
448 if not self.__withOverlapping:
448 if not self.__withOverlapping:
449 self.__buffer += data.copy()
449 self.__buffer += data.copy()
450 self.__profIndex += 1
450 self.__profIndex += 1
451 return
451 return
452
452
453 #Overlapping data
453 #Overlapping data
454 nChannels, nHeis = data.shape
454 nChannels, nHeis = data.shape
455 data = numpy.reshape(data, (1, nChannels, nHeis))
455 data = numpy.reshape(data, (1, nChannels, nHeis))
456
456
457 #If the buffer is empty then it takes the data value
457 #If the buffer is empty then it takes the data value
458 if self.__buffer is None:
458 if self.__buffer is None:
459 self.__buffer = data
459 self.__buffer = data
460 self.__profIndex += 1
460 self.__profIndex += 1
461 return
461 return
462
462
463 #If the buffer length is lower than n then stakcing the data value
463 #If the buffer length is lower than n then stakcing the data value
464 if self.__profIndex < self.n:
464 if self.__profIndex < self.n:
465 self.__buffer = numpy.vstack((self.__buffer, data))
465 self.__buffer = numpy.vstack((self.__buffer, data))
466 self.__profIndex += 1
466 self.__profIndex += 1
467 return
467 return
468
468
469 #If the buffer length is equal to n then replacing the last buffer value with the data value
469 #If the buffer length is equal to n then replacing the last buffer value with the data value
470 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
470 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
471 self.__buffer[self.n-1] = data
471 self.__buffer[self.n-1] = data
472 self.__profIndex = self.n
472 self.__profIndex = self.n
473 return
473 return
474
474
475
475
476 def pushData(self):
476 def pushData(self):
477 """
477 """
478 Return the sum of the last profiles and the profiles used in the sum.
478 Return the sum of the last profiles and the profiles used in the sum.
479
479
480 Affected:
480 Affected:
481
481
482 self.__profileIndex
482 self.__profileIndex
483
483
484 """
484 """
485
485
486 if not self.__withOverlapping:
486 if not self.__withOverlapping:
487 data = self.__buffer
487 data = self.__buffer
488 n = self.__profIndex
488 n = self.__profIndex
489
489
490 self.__buffer = 0
490 self.__buffer = 0
491 self.__profIndex = 0
491 self.__profIndex = 0
492
492
493 return data, n
493 return data, n
494
494
495 #Integration with Overlapping
495 #Integration with Overlapping
496 data = numpy.sum(self.__buffer, axis=0)
496 data = numpy.sum(self.__buffer, axis=0)
497 # print data
497 # print data
498 # raise
498 # raise
499 n = self.__profIndex
499 n = self.__profIndex
500
500
501 return data, n
501 return data, n
502
502
503 def byProfiles(self, data):
503 def byProfiles(self, data):
504
504
505 self.__dataReady = False
505 self.__dataReady = False
506 avgdata = None
506 avgdata = None
507 # n = None
507 # n = None
508 # print data
508 # print data
509 # raise
509 # raise
510 self.putData(data)
510 self.putData(data)
511
511
512 if self.__profIndex == self.n:
512 if self.__profIndex == self.n:
513 avgdata, n = self.pushData()
513 avgdata, n = self.pushData()
514 self.__dataReady = True
514 self.__dataReady = True
515
515
516 return avgdata
516 return avgdata
517
517
518 def byTime(self, data, datatime):
518 def byTime(self, data, datatime):
519
519
520 self.__dataReady = False
520 self.__dataReady = False
521 avgdata = None
521 avgdata = None
522 n = None
522 n = None
523
523
524 self.putData(data)
524 self.putData(data)
525
525
526 if (datatime - self.__initime) >= self.__integrationtime:
526 if (datatime - self.__initime) >= self.__integrationtime:
527 avgdata, n = self.pushData()
527 avgdata, n = self.pushData()
528 self.n = n
528 self.n = n
529 self.__dataReady = True
529 self.__dataReady = True
530
530
531 return avgdata
531 return avgdata
532
532
533 def integrateByStride(self, data, datatime):
533 def integrateByStride(self, data, datatime):
534 # print data
534 # print data
535 if self.__profIndex == 0:
535 if self.__profIndex == 0:
536 self.__buffer = [[data.copy(), datatime]]
536 self.__buffer = [[data.copy(), datatime]]
537 else:
537 else:
538 self.__buffer.append([data.copy(),datatime])
538 self.__buffer.append([data.copy(),datatime])
539 self.__profIndex += 1
539 self.__profIndex += 1
540 self.__dataReady = False
540 self.__dataReady = False
541
541
542 if self.__profIndex == self.n * self.stride :
542 if self.__profIndex == self.n * self.stride :
543 self.__dataToPutStride = True
543 self.__dataToPutStride = True
544 self.__profIndexStride = 0
544 self.__profIndexStride = 0
545 self.__profIndex = 0
545 self.__profIndex = 0
546 self.__bufferStride = []
546 self.__bufferStride = []
547 for i in range(self.stride):
547 for i in range(self.stride):
548 current = self.__buffer[i::self.stride]
548 current = self.__buffer[i::self.stride]
549 data = numpy.sum([t[0] for t in current], axis=0)
549 data = numpy.sum([t[0] for t in current], axis=0)
550 avgdatatime = numpy.average([t[1] for t in current])
550 avgdatatime = numpy.average([t[1] for t in current])
551 # print data
551 # print data
552 self.__bufferStride.append((data, avgdatatime))
552 self.__bufferStride.append((data, avgdatatime))
553
553
554 if self.__dataToPutStride:
554 if self.__dataToPutStride:
555 self.__dataReady = True
555 self.__dataReady = True
556 self.__profIndexStride += 1
556 self.__profIndexStride += 1
557 if self.__profIndexStride == self.stride:
557 if self.__profIndexStride == self.stride:
558 self.__dataToPutStride = False
558 self.__dataToPutStride = False
559 # print self.__bufferStride[self.__profIndexStride - 1]
559 # print self.__bufferStride[self.__profIndexStride - 1]
560 # raise
560 # raise
561 return self.__bufferStride[self.__profIndexStride - 1]
561 return self.__bufferStride[self.__profIndexStride - 1]
562
562
563
563
564 return None, None
564 return None, None
565
565
566 def integrate(self, data, datatime=None):
566 def integrate(self, data, datatime=None):
567
567
568 if self.__initime == None:
568 if self.__initime == None:
569 self.__initime = datatime
569 self.__initime = datatime
570
570
571 if self.__byTime:
571 if self.__byTime:
572 avgdata = self.byTime(data, datatime)
572 avgdata = self.byTime(data, datatime)
573 else:
573 else:
574 avgdata = self.byProfiles(data)
574 avgdata = self.byProfiles(data)
575
575
576
576
577 self.__lastdatatime = datatime
577 self.__lastdatatime = datatime
578
578
579 if avgdata is None:
579 if avgdata is None:
580 return None, None
580 return None, None
581
581
582 avgdatatime = self.__initime
582 avgdatatime = self.__initime
583
583
584 deltatime = datatime - self.__lastdatatime
584 deltatime = datatime - self.__lastdatatime
585
585
586 if not self.__withOverlapping:
586 if not self.__withOverlapping:
587 self.__initime = datatime
587 self.__initime = datatime
588 else:
588 else:
589 self.__initime += deltatime
589 self.__initime += deltatime
590
590
591 return avgdata, avgdatatime
591 return avgdata, avgdatatime
592
592
593 def integrateByBlock(self, dataOut):
593 def integrateByBlock(self, dataOut):
594
594
595 times = int(dataOut.data.shape[1]/self.n)
595 times = int(dataOut.data.shape[1]/self.n)
596 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
596 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
597
597
598 id_min = 0
598 id_min = 0
599 id_max = self.n
599 id_max = self.n
600
600
601 for i in range(times):
601 for i in range(times):
602 junk = dataOut.data[:,id_min:id_max,:]
602 junk = dataOut.data[:,id_min:id_max,:]
603 avgdata[:,i,:] = junk.sum(axis=1)
603 avgdata[:,i,:] = junk.sum(axis=1)
604 id_min += self.n
604 id_min += self.n
605 id_max += self.n
605 id_max += self.n
606
606
607 timeInterval = dataOut.ippSeconds*self.n
607 timeInterval = dataOut.ippSeconds*self.n
608 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
608 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
609 self.__dataReady = True
609 self.__dataReady = True
610 return avgdata, avgdatatime
610 return avgdata, avgdatatime
611
611
612 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
612 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
613 if not self.isConfig:
613 if not self.isConfig:
614 self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
614 self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
615 self.isConfig = True
615 self.isConfig = True
616
616
617 if dataOut.flagDataAsBlock:
617 if dataOut.flagDataAsBlock:
618 """
618 """
619 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
619 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
620 """
620 """
621 avgdata, avgdatatime = self.integrateByBlock(dataOut)
621 avgdata, avgdatatime = self.integrateByBlock(dataOut)
622 dataOut.nProfiles /= self.n
622 dataOut.nProfiles /= self.n
623 else:
623 else:
624 if stride is None:
624 if stride is None:
625 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
625 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
626 else:
626 else:
627 avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
627 avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
628
628
629
629
630 # dataOut.timeInterval *= n
630 # dataOut.timeInterval *= n
631 dataOut.flagNoData = True
631 dataOut.flagNoData = True
632
632
633 if self.__dataReady:
633 if self.__dataReady:
634 dataOut.data = avgdata
634 dataOut.data = avgdata
635 dataOut.nCohInt *= self.n
635 dataOut.nCohInt *= self.n
636 dataOut.utctime = avgdatatime
636 dataOut.utctime = avgdatatime
637 # print avgdata, avgdatatime
637 # print avgdata, avgdatatime
638 # raise
638 # raise
639 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
639 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
640 dataOut.flagNoData = False
640 dataOut.flagNoData = False
641
641
642 class Decoder(Operation):
642 class Decoder(Operation):
643
643
644 isConfig = False
644 isConfig = False
645 __profIndex = 0
645 __profIndex = 0
646
646
647 code = None
647 code = None
648
648
649 nCode = None
649 nCode = None
650 nBaud = None
650 nBaud = None
651
651
652 def __init__(self, **kwargs):
652 def __init__(self, **kwargs):
653
653
654 Operation.__init__(self, **kwargs)
654 Operation.__init__(self, **kwargs)
655
655
656 self.times = None
656 self.times = None
657 self.osamp = None
657 self.osamp = None
658 # self.__setValues = False
658 # self.__setValues = False
659 self.isConfig = False
659 self.isConfig = False
660
660
661 def setup(self, code, osamp, dataOut):
661 def setup(self, code, osamp, dataOut):
662
662
663 self.__profIndex = 0
663 self.__profIndex = 0
664
664
665 self.code = code
665 self.code = code
666
666
667 self.nCode = len(code)
667 self.nCode = len(code)
668 self.nBaud = len(code[0])
668 self.nBaud = len(code[0])
669
669
670 if (osamp != None) and (osamp >1):
670 if (osamp != None) and (osamp >1):
671 self.osamp = osamp
671 self.osamp = osamp
672 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
672 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
673 self.nBaud = self.nBaud*self.osamp
673 self.nBaud = self.nBaud*self.osamp
674
674
675 self.__nChannels = dataOut.nChannels
675 self.__nChannels = dataOut.nChannels
676 self.__nProfiles = dataOut.nProfiles
676 self.__nProfiles = dataOut.nProfiles
677 self.__nHeis = dataOut.nHeights
677 self.__nHeis = dataOut.nHeights
678
678
679 if self.__nHeis < self.nBaud:
679 if self.__nHeis < self.nBaud:
680 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
680 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
681
681
682 #Frequency
682 #Frequency
683 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
683 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
684
684
685 __codeBuffer[:,0:self.nBaud] = self.code
685 __codeBuffer[:,0:self.nBaud] = self.code
686
686
687 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
687 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
688
688
689 if dataOut.flagDataAsBlock:
689 if dataOut.flagDataAsBlock:
690
690
691 self.ndatadec = self.__nHeis #- self.nBaud + 1
691 self.ndatadec = self.__nHeis #- self.nBaud + 1
692
692
693 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
693 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
694
694
695 else:
695 else:
696
696
697 #Time
697 #Time
698 self.ndatadec = self.__nHeis #- self.nBaud + 1
698 self.ndatadec = self.__nHeis #- self.nBaud + 1
699
699
700 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
700 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
701
701
702 def __convolutionInFreq(self, data):
702 def __convolutionInFreq(self, data):
703
703
704 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
704 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
705
705
706 fft_data = numpy.fft.fft(data, axis=1)
706 fft_data = numpy.fft.fft(data, axis=1)
707
707
708 conv = fft_data*fft_code
708 conv = fft_data*fft_code
709
709
710 data = numpy.fft.ifft(conv,axis=1)
710 data = numpy.fft.ifft(conv,axis=1)
711
711
712 return data
712 return data
713
713
714 def __convolutionInFreqOpt(self, data):
714 def __convolutionInFreqOpt(self, data):
715
715
716 raise NotImplementedError
716 raise NotImplementedError
717
717
718 def __convolutionInTime(self, data):
718 def __convolutionInTime(self, data):
719
719
720 code = self.code[self.__profIndex]
720 code = self.code[self.__profIndex]
721 for i in range(self.__nChannels):
721 for i in range(self.__nChannels):
722 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
722 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
723
723
724 return self.datadecTime
724 return self.datadecTime
725
725
726 def __convolutionByBlockInTime(self, data):
726 def __convolutionByBlockInTime(self, data):
727
727
728 repetitions = self.__nProfiles / self.nCode
728 repetitions = self.__nProfiles / self.nCode
729
729
730 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
730 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
731 junk = junk.flatten()
731 junk = junk.flatten()
732 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
732 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
733 profilesList = xrange(self.__nProfiles)
733 profilesList = xrange(self.__nProfiles)
734
734
735 for i in range(self.__nChannels):
735 for i in range(self.__nChannels):
736 for j in profilesList:
736 for j in profilesList:
737 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
737 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
738 return self.datadecTime
738 return self.datadecTime
739
739
740 def __convolutionByBlockInFreq(self, data):
740 def __convolutionByBlockInFreq(self, data):
741
741
742 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
742 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
743
743
744
744
745 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
745 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
746
746
747 fft_data = numpy.fft.fft(data, axis=2)
747 fft_data = numpy.fft.fft(data, axis=2)
748
748
749 conv = fft_data*fft_code
749 conv = fft_data*fft_code
750
750
751 data = numpy.fft.ifft(conv,axis=2)
751 data = numpy.fft.ifft(conv,axis=2)
752
752
753 return data
753 return data
754
754
755
755
756 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
756 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
757
757
758 if dataOut.flagDecodeData:
758 if dataOut.flagDecodeData:
759 print "This data is already decoded, recoding again ..."
759 print "This data is already decoded, recoding again ..."
760
760
761 if not self.isConfig:
761 if not self.isConfig:
762
762
763 if code is None:
763 if code is None:
764 if dataOut.code is None:
764 if dataOut.code is None:
765 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
765 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
766
766
767 code = dataOut.code
767 code = dataOut.code
768 else:
768 else:
769 code = numpy.array(code).reshape(nCode,nBaud)
769 code = numpy.array(code).reshape(nCode,nBaud)
770 self.setup(code, osamp, dataOut)
770 self.setup(code, osamp, dataOut)
771
771
772 self.isConfig = True
772 self.isConfig = True
773
773
774 if mode == 3:
774 if mode == 3:
775 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
775 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
776
776
777 if times != None:
777 if times != None:
778 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
778 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
779
779
780 if self.code is None:
780 if self.code is None:
781 print "Fail decoding: Code is not defined."
781 print "Fail decoding: Code is not defined."
782 return
782 return
783
783
784 self.__nProfiles = dataOut.nProfiles
784 self.__nProfiles = dataOut.nProfiles
785 datadec = None
785 datadec = None
786
786
787 if mode == 3:
787 if mode == 3:
788 mode = 0
788 mode = 0
789
789
790 if dataOut.flagDataAsBlock:
790 if dataOut.flagDataAsBlock:
791 """
791 """
792 Decoding when data have been read as block,
792 Decoding when data have been read as block,
793 """
793 """
794
794
795 if mode == 0:
795 if mode == 0:
796 datadec = self.__convolutionByBlockInTime(dataOut.data)
796 datadec = self.__convolutionByBlockInTime(dataOut.data)
797 if mode == 1:
797 if mode == 1:
798 datadec = self.__convolutionByBlockInFreq(dataOut.data)
798 datadec = self.__convolutionByBlockInFreq(dataOut.data)
799 else:
799 else:
800 """
800 """
801 Decoding when data have been read profile by profile
801 Decoding when data have been read profile by profile
802 """
802 """
803 if mode == 0:
803 if mode == 0:
804 datadec = self.__convolutionInTime(dataOut.data)
804 datadec = self.__convolutionInTime(dataOut.data)
805
805
806 if mode == 1:
806 if mode == 1:
807 datadec = self.__convolutionInFreq(dataOut.data)
807 datadec = self.__convolutionInFreq(dataOut.data)
808
808
809 if mode == 2:
809 if mode == 2:
810 datadec = self.__convolutionInFreqOpt(dataOut.data)
810 datadec = self.__convolutionInFreqOpt(dataOut.data)
811
811
812 if datadec is None:
812 if datadec is None:
813 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
813 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
814
814
815 dataOut.code = self.code
815 dataOut.code = self.code
816 dataOut.nCode = self.nCode
816 dataOut.nCode = self.nCode
817 dataOut.nBaud = self.nBaud
817 dataOut.nBaud = self.nBaud
818
818
819 dataOut.data = datadec
819 dataOut.data = datadec
820
820
821 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
821 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
822
822
823 dataOut.flagDecodeData = True #asumo q la data esta decodificada
823 dataOut.flagDecodeData = True #asumo q la data esta decodificada
824
824
825 if self.__profIndex == self.nCode-1:
825 if self.__profIndex == self.nCode-1:
826 self.__profIndex = 0
826 self.__profIndex = 0
827 return 1
827 return 1
828
828
829 self.__profIndex += 1
829 self.__profIndex += 1
830
830
831 return 1
831 return 1
832 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
832 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
833
833
834
834
835 class ProfileConcat(Operation):
835 class ProfileConcat(Operation):
836
836
837 isConfig = False
837 isConfig = False
838 buffer = None
838 buffer = None
839 concat_m =None
839
840
840 def __init__(self, **kwargs):
841 def __init__(self, **kwargs):
841
842
842 Operation.__init__(self, **kwargs)
843 Operation.__init__(self, **kwargs)
843 self.profileIndex = 0
844 self.profileIndex = 0
844
845
845 def reset(self):
846 def reset(self):
846 self.buffer = numpy.zeros_like(self.buffer)
847 self.buffer = numpy.zeros_like(self.buffer)
847 self.start_index = 0
848 self.start_index = 0
848 self.times = 1
849 self.times = 1
849
850
850 def setup(self, data, m, n=1):
851 def setup(self, data, m, n=1):
851 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
852 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
852 self.nHeights = data.shape[1]#.nHeights
853 self.nHeights = data.shape[1]#.nHeights
853 self.start_index = 0
854 self.start_index = 0
854 self.times = 1
855 self.times = 1
855
856
856 def concat(self, data):
857 def concat(self, data):
857
858
858 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
859 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
859 self.start_index = self.start_index + self.nHeights
860 self.start_index = self.start_index + self.nHeights
860
861
861 def run(self, dataOut, m):
862 def run(self, dataOut, m):
862
863
864 self.concat_m= m
863 dataOut.flagNoData = True
865 dataOut.flagNoData = True
864
866
865 if not self.isConfig:
867 if not self.isConfig:
866 self.setup(dataOut.data, m, 1)
868 self.setup(dataOut.data, m, 1)
867 self.isConfig = True
869 self.isConfig = True
868
870
869 if dataOut.flagDataAsBlock:
871 if dataOut.flagDataAsBlock:
870 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
872 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
871
873
872 else:
874 else:
873 self.concat(dataOut.data)
875 self.concat(dataOut.data)
874 self.times += 1
876 self.times += 1
875 if self.times > m:
877 if self.times > m:
876 dataOut.data = self.buffer
878 dataOut.data = self.buffer
877 self.reset()
879 self.reset()
878 dataOut.flagNoData = False
880 dataOut.flagNoData = False
879 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
881 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
880 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
882 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
881 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
883 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
882 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
884 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
883 dataOut.ippSeconds *= m
885 dataOut.ippSeconds *= m
886 dataOut.concat_m = int(m)
884
887
885 class ProfileSelector(Operation):
888 class ProfileSelector(Operation):
886
889
887 profileIndex = None
890 profileIndex = None
888 # Tamanho total de los perfiles
891 # Tamanho total de los perfiles
889 nProfiles = None
892 nProfiles = None
890
893
891 def __init__(self, **kwargs):
894 def __init__(self, **kwargs):
892
895
893 Operation.__init__(self, **kwargs)
896 Operation.__init__(self, **kwargs)
894 self.profileIndex = 0
897 self.profileIndex = 0
895
898
896 def incProfileIndex(self):
899 def incProfileIndex(self):
897
900
898 self.profileIndex += 1
901 self.profileIndex += 1
899
902
900 if self.profileIndex >= self.nProfiles:
903 if self.profileIndex >= self.nProfiles:
901 self.profileIndex = 0
904 self.profileIndex = 0
902
905
903 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
906 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
904
907
905 if profileIndex < minIndex:
908 if profileIndex < minIndex:
906 return False
909 return False
907
910
908 if profileIndex > maxIndex:
911 if profileIndex > maxIndex:
909 return False
912 return False
910
913
911 return True
914 return True
912
915
913 def isThisProfileInList(self, profileIndex, profileList):
916 def isThisProfileInList(self, profileIndex, profileList):
914
917
915 if profileIndex not in profileList:
918 if profileIndex not in profileList:
916 return False
919 return False
917
920
918 return True
921 return True
919
922
920 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
923 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
921
924
922 """
925 """
923 ProfileSelector:
926 ProfileSelector:
924
927
925 Inputs:
928 Inputs:
926 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
929 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
927
930
928 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
931 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
929
932
930 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
933 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
931
934
932 """
935 """
933
936
934 if rangeList is not None:
937 if rangeList is not None:
935 if type(rangeList[0]) not in (tuple, list):
938 if type(rangeList[0]) not in (tuple, list):
936 rangeList = [rangeList]
939 rangeList = [rangeList]
937
940
938 dataOut.flagNoData = True
941 dataOut.flagNoData = True
939
942
940 if dataOut.flagDataAsBlock:
943 if dataOut.flagDataAsBlock:
941 """
944 """
942 data dimension = [nChannels, nProfiles, nHeis]
945 data dimension = [nChannels, nProfiles, nHeis]
943 """
946 """
944 if profileList != None:
947 if profileList != None:
945 dataOut.data = dataOut.data[:,profileList,:]
948 dataOut.data = dataOut.data[:,profileList,:]
946
949
947 if profileRangeList != None:
950 if profileRangeList != None:
948 minIndex = profileRangeList[0]
951 minIndex = profileRangeList[0]
949 maxIndex = profileRangeList[1]
952 maxIndex = profileRangeList[1]
950 profileList = range(minIndex, maxIndex+1)
953 profileList = range(minIndex, maxIndex+1)
951
954
952 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
955 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
953
956
954 if rangeList != None:
957 if rangeList != None:
955
958
956 profileList = []
959 profileList = []
957
960
958 for thisRange in rangeList:
961 for thisRange in rangeList:
959 minIndex = thisRange[0]
962 minIndex = thisRange[0]
960 maxIndex = thisRange[1]
963 maxIndex = thisRange[1]
961
964
962 profileList.extend(range(minIndex, maxIndex+1))
965 profileList.extend(range(minIndex, maxIndex+1))
963
966
964 dataOut.data = dataOut.data[:,profileList,:]
967 dataOut.data = dataOut.data[:,profileList,:]
965
968
966 dataOut.nProfiles = len(profileList)
969 dataOut.nProfiles = len(profileList)
967 dataOut.profileIndex = dataOut.nProfiles - 1
970 dataOut.profileIndex = dataOut.nProfiles - 1
968 dataOut.flagNoData = False
971 dataOut.flagNoData = False
969
972
970 return True
973 return True
971
974
972 """
975 """
973 data dimension = [nChannels, nHeis]
976 data dimension = [nChannels, nHeis]
974 """
977 """
975
978
976 if profileList != None:
979 if profileList != None:
977
980
978 if self.isThisProfileInList(dataOut.profileIndex, profileList):
981 if self.isThisProfileInList(dataOut.profileIndex, profileList):
979
982
980 self.nProfiles = len(profileList)
983 self.nProfiles = len(profileList)
981 dataOut.nProfiles = self.nProfiles
984 dataOut.nProfiles = self.nProfiles
982 dataOut.profileIndex = self.profileIndex
985 dataOut.profileIndex = self.profileIndex
983 dataOut.flagNoData = False
986 dataOut.flagNoData = False
984
987
985 self.incProfileIndex()
988 self.incProfileIndex()
986 return True
989 return True
987
990
988 if profileRangeList != None:
991 if profileRangeList != None:
989
992
990 minIndex = profileRangeList[0]
993 minIndex = profileRangeList[0]
991 maxIndex = profileRangeList[1]
994 maxIndex = profileRangeList[1]
992
995
993 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
996 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
994
997
995 self.nProfiles = maxIndex - minIndex + 1
998 self.nProfiles = maxIndex - minIndex + 1
996 dataOut.nProfiles = self.nProfiles
999 dataOut.nProfiles = self.nProfiles
997 dataOut.profileIndex = self.profileIndex
1000 dataOut.profileIndex = self.profileIndex
998 dataOut.flagNoData = False
1001 dataOut.flagNoData = False
999
1002
1000 self.incProfileIndex()
1003 self.incProfileIndex()
1001 return True
1004 return True
1002
1005
1003 if rangeList != None:
1006 if rangeList != None:
1004
1007
1005 nProfiles = 0
1008 nProfiles = 0
1006
1009
1007 for thisRange in rangeList:
1010 for thisRange in rangeList:
1008 minIndex = thisRange[0]
1011 minIndex = thisRange[0]
1009 maxIndex = thisRange[1]
1012 maxIndex = thisRange[1]
1010
1013
1011 nProfiles += maxIndex - minIndex + 1
1014 nProfiles += maxIndex - minIndex + 1
1012
1015
1013 for thisRange in rangeList:
1016 for thisRange in rangeList:
1014
1017
1015 minIndex = thisRange[0]
1018 minIndex = thisRange[0]
1016 maxIndex = thisRange[1]
1019 maxIndex = thisRange[1]
1017
1020
1018 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
1021 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
1019
1022
1020 self.nProfiles = nProfiles
1023 self.nProfiles = nProfiles
1021 dataOut.nProfiles = self.nProfiles
1024 dataOut.nProfiles = self.nProfiles
1022 dataOut.profileIndex = self.profileIndex
1025 dataOut.profileIndex = self.profileIndex
1023 dataOut.flagNoData = False
1026 dataOut.flagNoData = False
1024
1027
1025 self.incProfileIndex()
1028 self.incProfileIndex()
1026
1029
1027 break
1030 break
1028
1031
1029 return True
1032 return True
1030
1033
1031
1034
1032 if beam != None: #beam is only for AMISR data
1035 if beam != None: #beam is only for AMISR data
1033 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
1036 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
1034 dataOut.flagNoData = False
1037 dataOut.flagNoData = False
1035 dataOut.profileIndex = self.profileIndex
1038 dataOut.profileIndex = self.profileIndex
1036
1039
1037 self.incProfileIndex()
1040 self.incProfileIndex()
1038
1041
1039 return True
1042 return True
1040
1043
1041 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
1044 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
1042
1045
1043 return False
1046 return False
1044
1047
1045 class Reshaper(Operation):
1048 class Reshaper(Operation):
1046
1049
1047 def __init__(self, **kwargs):
1050 def __init__(self, **kwargs):
1048
1051
1049 Operation.__init__(self, **kwargs)
1052 Operation.__init__(self, **kwargs)
1050
1053
1051 self.__buffer = None
1054 self.__buffer = None
1052 self.__nitems = 0
1055 self.__nitems = 0
1053
1056
1054 def __appendProfile(self, dataOut, nTxs):
1057 def __appendProfile(self, dataOut, nTxs):
1055
1058
1056 if self.__buffer is None:
1059 if self.__buffer is None:
1057 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
1060 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
1058 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
1061 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
1059
1062
1060 ini = dataOut.nHeights * self.__nitems
1063 ini = dataOut.nHeights * self.__nitems
1061 end = ini + dataOut.nHeights
1064 end = ini + dataOut.nHeights
1062
1065
1063 self.__buffer[:, ini:end] = dataOut.data
1066 self.__buffer[:, ini:end] = dataOut.data
1064
1067
1065 self.__nitems += 1
1068 self.__nitems += 1
1066
1069
1067 return int(self.__nitems*nTxs)
1070 return int(self.__nitems*nTxs)
1068
1071
1069 def __getBuffer(self):
1072 def __getBuffer(self):
1070
1073
1071 if self.__nitems == int(1./self.__nTxs):
1074 if self.__nitems == int(1./self.__nTxs):
1072
1075
1073 self.__nitems = 0
1076 self.__nitems = 0
1074
1077
1075 return self.__buffer.copy()
1078 return self.__buffer.copy()
1076
1079
1077 return None
1080 return None
1078
1081
1079 def __checkInputs(self, dataOut, shape, nTxs):
1082 def __checkInputs(self, dataOut, shape, nTxs):
1080
1083
1081 if shape is None and nTxs is None:
1084 if shape is None and nTxs is None:
1082 raise ValueError, "Reshaper: shape of factor should be defined"
1085 raise ValueError, "Reshaper: shape of factor should be defined"
1083
1086
1084 if nTxs:
1087 if nTxs:
1085 if nTxs < 0:
1088 if nTxs < 0:
1086 raise ValueError, "nTxs should be greater than 0"
1089 raise ValueError, "nTxs should be greater than 0"
1087
1090
1088 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
1091 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
1089 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
1092 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
1090
1093
1091 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
1094 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
1092
1095
1093 return shape, nTxs
1096 return shape, nTxs
1094
1097
1095 if len(shape) != 2 and len(shape) != 3:
1098 if len(shape) != 2 and len(shape) != 3:
1096 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)
1099 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)
1097
1100
1098 if len(shape) == 2:
1101 if len(shape) == 2:
1099 shape_tuple = [dataOut.nChannels]
1102 shape_tuple = [dataOut.nChannels]
1100 shape_tuple.extend(shape)
1103 shape_tuple.extend(shape)
1101 else:
1104 else:
1102 shape_tuple = list(shape)
1105 shape_tuple = list(shape)
1103
1106
1104 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
1107 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
1105
1108
1106 return shape_tuple, nTxs
1109 return shape_tuple, nTxs
1107
1110
1108 def run(self, dataOut, shape=None, nTxs=None):
1111 def run(self, dataOut, shape=None, nTxs=None):
1109
1112
1110 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
1113 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
1111
1114
1112 dataOut.flagNoData = True
1115 dataOut.flagNoData = True
1113 profileIndex = None
1116 profileIndex = None
1114
1117
1115 if dataOut.flagDataAsBlock:
1118 if dataOut.flagDataAsBlock:
1116
1119
1117 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
1120 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
1118 dataOut.flagNoData = False
1121 dataOut.flagNoData = False
1119
1122
1120 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1123 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1121
1124
1122 else:
1125 else:
1123
1126
1124 if self.__nTxs < 1:
1127 if self.__nTxs < 1:
1125
1128
1126 self.__appendProfile(dataOut, self.__nTxs)
1129 self.__appendProfile(dataOut, self.__nTxs)
1127 new_data = self.__getBuffer()
1130 new_data = self.__getBuffer()
1128
1131
1129 if new_data is not None:
1132 if new_data is not None:
1130 dataOut.data = new_data
1133 dataOut.data = new_data
1131 dataOut.flagNoData = False
1134 dataOut.flagNoData = False
1132
1135
1133 profileIndex = dataOut.profileIndex*nTxs
1136 profileIndex = dataOut.profileIndex*nTxs
1134
1137
1135 else:
1138 else:
1136 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1139 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1137
1140
1138 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1141 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1139
1142
1140 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1143 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1141
1144
1142 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1145 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1143
1146
1144 dataOut.profileIndex = profileIndex
1147 dataOut.profileIndex = profileIndex
1145
1148
1146 dataOut.ippSeconds /= self.__nTxs
1149 dataOut.ippSeconds /= self.__nTxs
1147
1150
1148 class SplitProfiles(Operation):
1151 class SplitProfiles(Operation):
1149
1152
1150 def __init__(self, **kwargs):
1153 def __init__(self, **kwargs):
1151
1154
1152 Operation.__init__(self, **kwargs)
1155 Operation.__init__(self, **kwargs)
1153
1156
1154 def run(self, dataOut, n):
1157 def run(self, dataOut, n):
1155
1158
1156 dataOut.flagNoData = True
1159 dataOut.flagNoData = True
1157 profileIndex = None
1160 profileIndex = None
1158
1161
1159 if dataOut.flagDataAsBlock:
1162 if dataOut.flagDataAsBlock:
1160
1163
1161 #nchannels, nprofiles, nsamples
1164 #nchannels, nprofiles, nsamples
1162 shape = dataOut.data.shape
1165 shape = dataOut.data.shape
1163
1166
1164 if shape[2] % n != 0:
1167 if shape[2] % n != 0:
1165 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1168 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1166
1169
1167 new_shape = shape[0], shape[1]*n, shape[2]/n
1170 new_shape = shape[0], shape[1]*n, shape[2]/n
1168
1171
1169 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1172 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1170 dataOut.flagNoData = False
1173 dataOut.flagNoData = False
1171
1174
1172 profileIndex = int(dataOut.nProfiles/n) - 1
1175 profileIndex = int(dataOut.nProfiles/n) - 1
1173
1176
1174 else:
1177 else:
1175
1178
1176 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1179 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1177
1180
1178 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1181 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1179
1182
1180 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1183 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1181
1184
1182 dataOut.nProfiles = int(dataOut.nProfiles*n)
1185 dataOut.nProfiles = int(dataOut.nProfiles*n)
1183
1186
1184 dataOut.profileIndex = profileIndex
1187 dataOut.profileIndex = profileIndex
1185
1188
1186 dataOut.ippSeconds /= n
1189 dataOut.ippSeconds /= n
1187
1190
1188 class CombineProfiles(Operation):
1191 class CombineProfiles(Operation):
1189
1192
1190 def __init__(self, **kwargs):
1193 def __init__(self, **kwargs):
1191
1194
1192 Operation.__init__(self, **kwargs)
1195 Operation.__init__(self, **kwargs)
1193
1196
1194 self.__remData = None
1197 self.__remData = None
1195 self.__profileIndex = 0
1198 self.__profileIndex = 0
1196
1199
1197 def run(self, dataOut, n):
1200 def run(self, dataOut, n):
1198
1201
1199 dataOut.flagNoData = True
1202 dataOut.flagNoData = True
1200 profileIndex = None
1203 profileIndex = None
1201
1204
1202 if dataOut.flagDataAsBlock:
1205 if dataOut.flagDataAsBlock:
1203
1206
1204 #nchannels, nprofiles, nsamples
1207 #nchannels, nprofiles, nsamples
1205 shape = dataOut.data.shape
1208 shape = dataOut.data.shape
1206 new_shape = shape[0], shape[1]/n, shape[2]*n
1209 new_shape = shape[0], shape[1]/n, shape[2]*n
1207
1210
1208 if shape[1] % n != 0:
1211 if shape[1] % n != 0:
1209 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1212 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1210
1213
1211 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1214 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1212 dataOut.flagNoData = False
1215 dataOut.flagNoData = False
1213
1216
1214 profileIndex = int(dataOut.nProfiles*n) - 1
1217 profileIndex = int(dataOut.nProfiles*n) - 1
1215
1218
1216 else:
1219 else:
1217
1220
1218 #nchannels, nsamples
1221 #nchannels, nsamples
1219 if self.__remData is None:
1222 if self.__remData is None:
1220 newData = dataOut.data
1223 newData = dataOut.data
1221 else:
1224 else:
1222 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1225 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1223
1226
1224 self.__profileIndex += 1
1227 self.__profileIndex += 1
1225
1228
1226 if self.__profileIndex < n:
1229 if self.__profileIndex < n:
1227 self.__remData = newData
1230 self.__remData = newData
1228 #continue
1231 #continue
1229 return
1232 return
1230
1233
1231 self.__profileIndex = 0
1234 self.__profileIndex = 0
1232 self.__remData = None
1235 self.__remData = None
1233
1236
1234 dataOut.data = newData
1237 dataOut.data = newData
1235 dataOut.flagNoData = False
1238 dataOut.flagNoData = False
1236
1239
1237 profileIndex = dataOut.profileIndex/n
1240 profileIndex = dataOut.profileIndex/n
1238
1241
1239
1242
1240 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1243 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1241
1244
1242 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1245 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1243
1246
1244 dataOut.nProfiles = int(dataOut.nProfiles/n)
1247 dataOut.nProfiles = int(dataOut.nProfiles/n)
1245
1248
1246 dataOut.profileIndex = profileIndex
1249 dataOut.profileIndex = profileIndex
1247
1250
1248 dataOut.ippSeconds *= n
1251 dataOut.ippSeconds *= n
1249
1252
1250
1253
1251 class SSheightProfiles(Operation):
1254 class SSheightProfiles(Operation):
1252
1255
1253 step = None
1256 step = None
1254 nsamples = None
1257 nsamples = None
1255 bufferShape = None
1258 bufferShape = None
1256 profileShape = None
1259 profileShape = None
1257 sshProfiles = None
1260 sshProfiles = None
1258 profileIndex = None
1261 profileIndex = None
1259
1262
1260 def __init__(self, **kwargs):
1263 def __init__(self, **kwargs):
1261
1264
1262 Operation.__init__(self, **kwargs)
1265 Operation.__init__(self, **kwargs)
1263 self.isConfig = False
1266 self.isConfig = False
1264
1267
1265 def setup(self,dataOut ,step = None , nsamples = None):
1268 def setup(self,dataOut ,step = None , nsamples = None):
1266
1269
1267 if step == None and nsamples == None:
1270 if step == None and nsamples == None:
1268 raise ValueError, "step or nheights should be specified ..."
1271 raise ValueError, "step or nheights should be specified ..."
1269
1272
1270 self.step = step
1273 self.step = step
1271 self.nsamples = nsamples
1274 self.nsamples = nsamples
1272 self.__nChannels = dataOut.nChannels
1275 self.__nChannels = dataOut.nChannels
1273 self.__nProfiles = dataOut.nProfiles
1276 self.__nProfiles = dataOut.nProfiles
1274 self.__nHeis = dataOut.nHeights
1277 self.__nHeis = dataOut.nHeights
1275 shape = dataOut.data.shape #nchannels, nprofiles, nsamples
1278 shape = dataOut.data.shape #nchannels, nprofiles, nsamples
1276
1279
1277
1280
1278 residue = (shape[1] - self.nsamples) % self.step
1281 residue = (shape[1] - self.nsamples) % self.step
1279 if residue != 0:
1282 if residue != 0:
1280 print "The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue)
1283 print "The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue)
1281
1284
1282 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1285 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1283 numberProfile = self.nsamples
1286 numberProfile = self.nsamples
1284 numberSamples = (shape[1] - self.nsamples)/self.step
1287 numberSamples = (shape[1] - self.nsamples)/self.step
1285
1288
1286 print "New number of profile: %d, number of height: %d, Resolution %d Km"%(numberProfile,numberSamples,deltaHeight*self.step)
1289 print "New number of profile: %d, number of height: %d, Resolution %d Km"%(numberProfile,numberSamples,deltaHeight*self.step)
1287
1290
1288 self.bufferShape = shape[0], numberSamples, numberProfile # nchannels, nsamples , nprofiles
1291 self.bufferShape = shape[0], numberSamples, numberProfile # nchannels, nsamples , nprofiles
1289 self.profileShape = shape[0], numberProfile, numberSamples # nchannels, nprofiles, nsamples
1292 self.profileShape = shape[0], numberProfile, numberSamples # nchannels, nprofiles, nsamples
1290
1293
1291 self.buffer = numpy.zeros(self.bufferShape , dtype=numpy.complex)
1294 self.buffer = numpy.zeros(self.bufferShape , dtype=numpy.complex)
1292 self.sshProfiles = numpy.zeros(self.profileShape, dtype=numpy.complex)
1295 self.sshProfiles = numpy.zeros(self.profileShape, dtype=numpy.complex)
1293
1296
1294 def run(self, dataOut, step, nsamples):
1297 def run(self, dataOut, step, nsamples):
1295
1298
1296 dataOut.flagNoData = True
1299 dataOut.flagNoData = True
1297 dataOut.flagDataAsBlock = False
1300 dataOut.flagDataAsBlock = False
1298 profileIndex = None
1301 profileIndex = None
1299
1302
1303
1300 if not self.isConfig:
1304 if not self.isConfig:
1301 self.setup(dataOut, step=step , nsamples=nsamples)
1305 self.setup(dataOut, step=step , nsamples=nsamples)
1302 self.isConfig = True
1306 self.isConfig = True
1303
1307
1304 for i in range(self.buffer.shape[1]):
1308 for i in range(self.buffer.shape[1]):
1305 self.buffer[:,i] = numpy.flip(dataOut.data[:,i*self.step:i*self.step + self.nsamples])
1309 self.buffer[:,i] = numpy.flip(dataOut.data[:,i*self.step:i*self.step + self.nsamples])
1306 #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
1310 #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
1307
1311
1308 for j in range(self.buffer.shape[0]):
1312 for j in range(self.buffer.shape[0]):
1309 self.sshProfiles[j] = numpy.transpose(self.buffer[j])
1313 self.sshProfiles[j] = numpy.transpose(self.buffer[j])
1310
1314
1311 profileIndex = self.nsamples
1315 profileIndex = self.nsamples
1312 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1316 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1313 ippSeconds = (deltaHeight*1.0e-6)/(0.15)
1317 ippSeconds = (deltaHeight*1.0e-6)/(0.15)
1318 try:
1319 if dataOut.concat_m is not None:
1320 ippSeconds= ippSeconds/float(dataOut.concat_m)
1321 #print "Profile concat %d"%dataOut.concat_m
1322 except:
1323 pass
1314
1324
1315 dataOut.data = self.sshProfiles
1325 dataOut.data = self.sshProfiles
1316 dataOut.flagNoData = False
1326 dataOut.flagNoData = False
1317 dataOut.heightList = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
1327 dataOut.heightList = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
1318 dataOut.nProfiles = int(dataOut.nProfiles*self.nsamples)
1328 dataOut.nProfiles = int(dataOut.nProfiles*self.nsamples)
1319 dataOut.profileIndex = profileIndex
1329 dataOut.profileIndex = profileIndex
1320 dataOut.flagDataAsBlock = True
1330 dataOut.flagDataAsBlock = True
1321 dataOut.ippSeconds = ippSeconds
1331 dataOut.ippSeconds = ippSeconds
1322 dataOut.step = self.step
1332 dataOut.step = self.step
1323
1333
1324
1334
1325 import time
1335 import time
1326 #################################################
1336 #################################################
1327
1337
1328 class decoPseudorandom(Operation):
1338 class decoPseudorandom(Operation):
1329
1339
1330 nProfiles= 0
1340 nProfiles= 0
1331 buffer= None
1341 buffer= None
1332 isConfig = False
1342 isConfig = False
1333
1343
1334 def setup(self, clen= 10000,seed= 0,Nranges= 1000,oversample=1):
1344 def setup(self, clen= 10000,seed= 0,Nranges= 1000,oversample=1):
1335 #code = create_pseudo_random_code(clen=clen, seed=seed)
1345 #code = create_pseudo_random_code(clen=clen, seed=seed)
1336 code= rep_seq(create_pseudo_random_code(clen=clen, seed=seed),rep=oversample)
1346 code= rep_seq(create_pseudo_random_code(clen=clen, seed=seed),rep=oversample)
1337 #print ("code_rx", code.shape)
1347 #print ("code_rx", code.shape)
1338 #N = int(an_len/clen) # 100
1348 #N = int(an_len/clen) # 100
1339 B_cache = 0
1349 B_cache = 0
1340 r_cache = 0
1350 r_cache = 0
1341 B_cached = False
1351 B_cached = False
1342 r = create_estimation_matrix(code=code, cache=True, rmax=Nranges)
1352 r = create_estimation_matrix(code=code, cache=True, rmax=Nranges)
1343 #print ("code shape", code.shape)
1353 #print ("code shape", code.shape)
1344 #print ("seed",seed)
1354 #print ("seed",seed)
1345 #print ("Code", code[0:10])
1355 #print ("Code", code[0:10])
1346 self.B = r['B']
1356 self.B = r['B']
1347
1357
1348
1358
1349 def run (self,dataOut,length_code= 10000,seed= 0,Nranges= 1000,oversample=1):
1359 def run (self,dataOut,length_code= 10000,seed= 0,Nranges= 1000,oversample=1):
1350 #print((dataOut.data.shape))
1360 #print((dataOut.data.shape))
1351 if not self.isConfig:
1361 if not self.isConfig:
1352 self.setup(clen= length_code,seed= seed,Nranges= Nranges,oversample=oversample)
1362 self.setup(clen= length_code,seed= seed,Nranges= Nranges,oversample=oversample)
1353 self.isConfig = True
1363 self.isConfig = True
1354
1364
1355 dataOut.flagNoData = True
1365 dataOut.flagNoData = True
1356 data =dataOut.data
1366 data =dataOut.data
1357 #print "length_CODE",length_code
1367 #print "length_CODE",length_code
1358 data_shape = (data.shape[1])
1368 data_shape = (data.shape[1])
1359 #print "data_shape",data_shape
1369 #print "data_shape",data_shape
1360 n = (length_code /data_shape)
1370 n = (length_code /data_shape)
1361 #print "we need this number of sample",n
1371 #print "we need this number of sample",n
1362
1372
1363 if n>0 and self.buffer is None:
1373 if n>0 and self.buffer is None:
1364 self.buffer = numpy.zeros([1, length_code], dtype=numpy.complex64)
1374 self.buffer = numpy.zeros([1, length_code], dtype=numpy.complex64)
1365 self.buffer[0][0:data_shape] = data[0]
1375 self.buffer[0][0:data_shape] = data[0]
1366 #print "FIRST CREATION",self.buffer.shape
1376 #print "FIRST CREATION",self.buffer.shape
1367
1377
1368 else:
1378 else:
1369 self.buffer[0][self.nProfiles*data_shape:(self.nProfiles+1)*data_shape]=data[0]
1379 self.buffer[0][self.nProfiles*data_shape:(self.nProfiles+1)*data_shape]=data[0]
1370
1380
1371 #print "buffer_shape",(self.buffer.shape)
1381 #print "buffer_shape",(self.buffer.shape)
1372 self.nProfiles += 1
1382 self.nProfiles += 1
1373 #print "count",self.nProfiles
1383 #print "count",self.nProfiles
1374
1384
1375 if self.nProfiles== n:
1385 if self.nProfiles== n:
1376 temporal = numpy.dot(self.B, numpy.transpose(self.buffer))
1386 temporal = numpy.dot(self.B, numpy.transpose(self.buffer))
1377 #print temporal.shape
1387 #print temporal.shape
1378 #import time
1388 #import time
1379 #time.sleep(40)
1389 #time.sleep(40)
1380 dataOut.data=numpy.transpose(temporal)
1390 dataOut.data=numpy.transpose(temporal)
1381
1391
1382 dataOut.flagNoData = False
1392 dataOut.flagNoData = False
1383 self.buffer= None
1393 self.buffer= None
1384 self.nProfiles = 0
1394 self.nProfiles = 0
1385
1395
1386 # import collections
1396 # import collections
1387 # from scipy.stats import mode
1397 # from scipy.stats import mode
1388 #
1398 #
1389 # class Synchronize(Operation):
1399 # class Synchronize(Operation):
1390 #
1400 #
1391 # isConfig = False
1401 # isConfig = False
1392 # __profIndex = 0
1402 # __profIndex = 0
1393 #
1403 #
1394 # def __init__(self, **kwargs):
1404 # def __init__(self, **kwargs):
1395 #
1405 #
1396 # Operation.__init__(self, **kwargs)
1406 # Operation.__init__(self, **kwargs)
1397 # # self.isConfig = False
1407 # # self.isConfig = False
1398 # self.__powBuffer = None
1408 # self.__powBuffer = None
1399 # self.__startIndex = 0
1409 # self.__startIndex = 0
1400 # self.__pulseFound = False
1410 # self.__pulseFound = False
1401 #
1411 #
1402 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1412 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1403 #
1413 #
1404 # #Read data
1414 # #Read data
1405 #
1415 #
1406 # powerdB = dataOut.getPower(channel = channel)
1416 # powerdB = dataOut.getPower(channel = channel)
1407 # noisedB = dataOut.getNoise(channel = channel)[0]
1417 # noisedB = dataOut.getNoise(channel = channel)[0]
1408 #
1418 #
1409 # self.__powBuffer.extend(powerdB.flatten())
1419 # self.__powBuffer.extend(powerdB.flatten())
1410 #
1420 #
1411 # dataArray = numpy.array(self.__powBuffer)
1421 # dataArray = numpy.array(self.__powBuffer)
1412 #
1422 #
1413 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1423 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1414 #
1424 #
1415 # maxValue = numpy.nanmax(filteredPower)
1425 # maxValue = numpy.nanmax(filteredPower)
1416 #
1426 #
1417 # if maxValue < noisedB + 10:
1427 # if maxValue < noisedB + 10:
1418 # #No se encuentra ningun pulso de transmision
1428 # #No se encuentra ningun pulso de transmision
1419 # return None
1429 # return None
1420 #
1430 #
1421 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1431 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1422 #
1432 #
1423 # if len(maxValuesIndex) < 2:
1433 # if len(maxValuesIndex) < 2:
1424 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1434 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1425 # return None
1435 # return None
1426 #
1436 #
1427 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1437 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1428 #
1438 #
1429 # #Seleccionar solo valores con un espaciamiento de nSamples
1439 # #Seleccionar solo valores con un espaciamiento de nSamples
1430 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1440 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1431 #
1441 #
1432 # if len(pulseIndex) < 2:
1442 # if len(pulseIndex) < 2:
1433 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1443 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1434 # return None
1444 # return None
1435 #
1445 #
1436 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1446 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1437 #
1447 #
1438 # #remover senales que se distancien menos de 10 unidades o muestras
1448 # #remover senales que se distancien menos de 10 unidades o muestras
1439 # #(No deberian existir IPP menor a 10 unidades)
1449 # #(No deberian existir IPP menor a 10 unidades)
1440 #
1450 #
1441 # realIndex = numpy.where(spacing > 10 )[0]
1451 # realIndex = numpy.where(spacing > 10 )[0]
1442 #
1452 #
1443 # if len(realIndex) < 2:
1453 # if len(realIndex) < 2:
1444 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1454 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1445 # return None
1455 # return None
1446 #
1456 #
1447 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1457 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1448 # realPulseIndex = pulseIndex[realIndex]
1458 # realPulseIndex = pulseIndex[realIndex]
1449 #
1459 #
1450 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1460 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1451 #
1461 #
1452 # print "IPP = %d samples" %period
1462 # print "IPP = %d samples" %period
1453 #
1463 #
1454 # self.__newNSamples = dataOut.nHeights #int(period)
1464 # self.__newNSamples = dataOut.nHeights #int(period)
1455 # self.__startIndex = int(realPulseIndex[0])
1465 # self.__startIndex = int(realPulseIndex[0])
1456 #
1466 #
1457 # return 1
1467 # return 1
1458 #
1468 #
1459 #
1469 #
1460 # def setup(self, nSamples, nChannels, buffer_size = 4):
1470 # def setup(self, nSamples, nChannels, buffer_size = 4):
1461 #
1471 #
1462 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1472 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1463 # maxlen = buffer_size*nSamples)
1473 # maxlen = buffer_size*nSamples)
1464 #
1474 #
1465 # bufferList = []
1475 # bufferList = []
1466 #
1476 #
1467 # for i in range(nChannels):
1477 # for i in range(nChannels):
1468 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1478 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1469 # maxlen = buffer_size*nSamples)
1479 # maxlen = buffer_size*nSamples)
1470 #
1480 #
1471 # bufferList.append(bufferByChannel)
1481 # bufferList.append(bufferByChannel)
1472 #
1482 #
1473 # self.__nSamples = nSamples
1483 # self.__nSamples = nSamples
1474 # self.__nChannels = nChannels
1484 # self.__nChannels = nChannels
1475 # self.__bufferList = bufferList
1485 # self.__bufferList = bufferList
1476 #
1486 #
1477 # def run(self, dataOut, channel = 0):
1487 # def run(self, dataOut, channel = 0):
1478 #
1488 #
1479 # if not self.isConfig:
1489 # if not self.isConfig:
1480 # nSamples = dataOut.nHeights
1490 # nSamples = dataOut.nHeights
1481 # nChannels = dataOut.nChannels
1491 # nChannels = dataOut.nChannels
1482 # self.setup(nSamples, nChannels)
1492 # self.setup(nSamples, nChannels)
1483 # self.isConfig = True
1493 # self.isConfig = True
1484 #
1494 #
1485 # #Append new data to internal buffer
1495 # #Append new data to internal buffer
1486 # for thisChannel in range(self.__nChannels):
1496 # for thisChannel in range(self.__nChannels):
1487 # bufferByChannel = self.__bufferList[thisChannel]
1497 # bufferByChannel = self.__bufferList[thisChannel]
1488 # bufferByChannel.extend(dataOut.data[thisChannel])
1498 # bufferByChannel.extend(dataOut.data[thisChannel])
1489 #
1499 #
1490 # if self.__pulseFound:
1500 # if self.__pulseFound:
1491 # self.__startIndex -= self.__nSamples
1501 # self.__startIndex -= self.__nSamples
1492 #
1502 #
1493 # #Finding Tx Pulse
1503 # #Finding Tx Pulse
1494 # if not self.__pulseFound:
1504 # if not self.__pulseFound:
1495 # indexFound = self.__findTxPulse(dataOut, channel)
1505 # indexFound = self.__findTxPulse(dataOut, channel)
1496 #
1506 #
1497 # if indexFound == None:
1507 # if indexFound == None:
1498 # dataOut.flagNoData = True
1508 # dataOut.flagNoData = True
1499 # return
1509 # return
1500 #
1510 #
1501 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1511 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1502 # self.__pulseFound = True
1512 # self.__pulseFound = True
1503 # self.__startIndex = indexFound
1513 # self.__startIndex = indexFound
1504 #
1514 #
1505 # #If pulse was found ...
1515 # #If pulse was found ...
1506 # for thisChannel in range(self.__nChannels):
1516 # for thisChannel in range(self.__nChannels):
1507 # bufferByChannel = self.__bufferList[thisChannel]
1517 # bufferByChannel = self.__bufferList[thisChannel]
1508 # #print self.__startIndex
1518 # #print self.__startIndex
1509 # x = numpy.array(bufferByChannel)
1519 # x = numpy.array(bufferByChannel)
1510 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1520 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1511 #
1521 #
1512 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1522 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1513 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1523 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1514 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1524 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1515 #
1525 #
1516 # dataOut.data = self.__arrayBuffer
1526 # dataOut.data = self.__arrayBuffer
1517 #
1527 #
1518 # self.__startIndex += self.__newNSamples
1528 # self.__startIndex += self.__newNSamples
1519 #
1529 #
1520 # return
1530 # return
General Comments 0
You need to be logged in to leave comments. Login now