##// END OF EJS Templates
Update plots for realtime app
Juan C. Espinoza -
r1221:94b941a5310a
parent child
Show More
@@ -1,1366 +1,1365
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 import json
10 import json
11
11
12 from schainpy.utils import log
12 from schainpy.utils import log
13 from .jroheaderIO import SystemHeader, RadarControllerHeader
13 from .jroheaderIO import SystemHeader, RadarControllerHeader
14
14
15
15
16 def getNumpyDtype(dataTypeCode):
16 def getNumpyDtype(dataTypeCode):
17
17
18 if dataTypeCode == 0:
18 if dataTypeCode == 0:
19 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
19 numpyDtype = numpy.dtype([('real', '<i1'), ('imag', '<i1')])
20 elif dataTypeCode == 1:
20 elif dataTypeCode == 1:
21 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
21 numpyDtype = numpy.dtype([('real', '<i2'), ('imag', '<i2')])
22 elif dataTypeCode == 2:
22 elif dataTypeCode == 2:
23 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
23 numpyDtype = numpy.dtype([('real', '<i4'), ('imag', '<i4')])
24 elif dataTypeCode == 3:
24 elif dataTypeCode == 3:
25 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
25 numpyDtype = numpy.dtype([('real', '<i8'), ('imag', '<i8')])
26 elif dataTypeCode == 4:
26 elif dataTypeCode == 4:
27 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
27 numpyDtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
28 elif dataTypeCode == 5:
28 elif dataTypeCode == 5:
29 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
29 numpyDtype = numpy.dtype([('real', '<f8'), ('imag', '<f8')])
30 else:
30 else:
31 raise ValueError('dataTypeCode was not defined')
31 raise ValueError('dataTypeCode was not defined')
32
32
33 return numpyDtype
33 return numpyDtype
34
34
35
35
36 def getDataTypeCode(numpyDtype):
36 def getDataTypeCode(numpyDtype):
37
37
38 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
38 if numpyDtype == numpy.dtype([('real', '<i1'), ('imag', '<i1')]):
39 datatype = 0
39 datatype = 0
40 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
40 elif numpyDtype == numpy.dtype([('real', '<i2'), ('imag', '<i2')]):
41 datatype = 1
41 datatype = 1
42 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
42 elif numpyDtype == numpy.dtype([('real', '<i4'), ('imag', '<i4')]):
43 datatype = 2
43 datatype = 2
44 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
44 elif numpyDtype == numpy.dtype([('real', '<i8'), ('imag', '<i8')]):
45 datatype = 3
45 datatype = 3
46 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
46 elif numpyDtype == numpy.dtype([('real', '<f4'), ('imag', '<f4')]):
47 datatype = 4
47 datatype = 4
48 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
48 elif numpyDtype == numpy.dtype([('real', '<f8'), ('imag', '<f8')]):
49 datatype = 5
49 datatype = 5
50 else:
50 else:
51 datatype = None
51 datatype = None
52
52
53 return datatype
53 return datatype
54
54
55
55
56 def hildebrand_sekhon(data, navg):
56 def hildebrand_sekhon(data, navg):
57 """
57 """
58 This method is for the objective determination of the noise level in Doppler spectra. This
58 This method is for the objective determination of the noise level in Doppler spectra. This
59 implementation technique is based on the fact that the standard deviation of the spectral
59 implementation technique is based on the fact that the standard deviation of the spectral
60 densities is equal to the mean spectral density for white Gaussian noise
60 densities is equal to the mean spectral density for white Gaussian noise
61
61
62 Inputs:
62 Inputs:
63 Data : heights
63 Data : heights
64 navg : numbers of averages
64 navg : numbers of averages
65
65
66 Return:
66 Return:
67 mean : noise's level
67 mean : 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
75
76 nums_min = 5
76 nums_min = 5
77
77
78 sump = 0.
78 sump = 0.
79 sumq = 0.
79 sumq = 0.
80
80
81 j = 0
81 j = 0
82 cont = 1
82 cont = 1
83
83
84 while((cont == 1)and(j < lenOfData)):
84 while((cont == 1)and(j < lenOfData)):
85
85
86 sump += sortdata[j]
86 sump += sortdata[j]
87 sumq += sortdata[j]**2
87 sumq += sortdata[j]**2
88
88
89 if j > nums_min:
89 if j > nums_min:
90 rtest = float(j)/(j-1) + 1.0/navg
90 rtest = float(j)/(j-1) + 1.0/navg
91 if ((sumq*j) > (rtest*sump**2)):
91 if ((sumq*j) > (rtest*sump**2)):
92 j = j - 1
92 j = j - 1
93 sump = sump - sortdata[j]
93 sump = sump - sortdata[j]
94 sumq = sumq - sortdata[j]**2
94 sumq = sumq - sortdata[j]**2
95 cont = 0
95 cont = 0
96
96
97 j += 1
97 j += 1
98
98
99 lnoise = sump / j
99 lnoise = sump / j
100
100
101 return lnoise
101 return lnoise
102
102
103
103
104 class Beam:
104 class Beam:
105
105
106 def __init__(self):
106 def __init__(self):
107 self.codeList = []
107 self.codeList = []
108 self.azimuthList = []
108 self.azimuthList = []
109 self.zenithList = []
109 self.zenithList = []
110
110
111
111
112 class GenericData(object):
112 class GenericData(object):
113
113
114 flagNoData = True
114 flagNoData = True
115
115
116 def copy(self, inputObj=None):
116 def copy(self, inputObj=None):
117
117
118 if inputObj == None:
118 if inputObj == None:
119 return copy.deepcopy(self)
119 return copy.deepcopy(self)
120
120
121 for key in list(inputObj.__dict__.keys()):
121 for key in list(inputObj.__dict__.keys()):
122
122
123 attribute = inputObj.__dict__[key]
123 attribute = inputObj.__dict__[key]
124
124
125 # If this attribute is a tuple or list
125 # If this attribute is a tuple or list
126 if type(inputObj.__dict__[key]) in (tuple, list):
126 if type(inputObj.__dict__[key]) in (tuple, list):
127 self.__dict__[key] = attribute[:]
127 self.__dict__[key] = attribute[:]
128 continue
128 continue
129
129
130 # If this attribute is another object or instance
130 # If this attribute is another object or instance
131 if hasattr(attribute, '__dict__'):
131 if hasattr(attribute, '__dict__'):
132 self.__dict__[key] = attribute.copy()
132 self.__dict__[key] = attribute.copy()
133 continue
133 continue
134
134
135 self.__dict__[key] = inputObj.__dict__[key]
135 self.__dict__[key] = inputObj.__dict__[key]
136
136
137 def deepcopy(self):
137 def deepcopy(self):
138
138
139 return copy.deepcopy(self)
139 return copy.deepcopy(self)
140
140
141 def isEmpty(self):
141 def isEmpty(self):
142
142
143 return self.flagNoData
143 return self.flagNoData
144
144
145
145
146 class JROData(GenericData):
146 class JROData(GenericData):
147
147
148 # m_BasicHeader = BasicHeader()
148 # m_BasicHeader = BasicHeader()
149 # m_ProcessingHeader = ProcessingHeader()
149 # m_ProcessingHeader = ProcessingHeader()
150
150
151 systemHeaderObj = SystemHeader()
151 systemHeaderObj = SystemHeader()
152 radarControllerHeaderObj = RadarControllerHeader()
152 radarControllerHeaderObj = RadarControllerHeader()
153 # data = None
153 # data = None
154 type = None
154 type = None
155 datatype = None # dtype but in string
155 datatype = None # dtype but in string
156 # dtype = None
156 # dtype = None
157 # nChannels = None
157 # nChannels = None
158 # nHeights = None
158 # nHeights = None
159 nProfiles = None
159 nProfiles = None
160 heightList = None
160 heightList = None
161 channelList = None
161 channelList = None
162 flagDiscontinuousBlock = False
162 flagDiscontinuousBlock = False
163 useLocalTime = False
163 useLocalTime = False
164 utctime = None
164 utctime = None
165 timeZone = None
165 timeZone = None
166 dstFlag = None
166 dstFlag = None
167 errorCount = None
167 errorCount = None
168 blocksize = None
168 blocksize = None
169 # nCode = None
169 # nCode = None
170 # nBaud = None
170 # nBaud = None
171 # code = None
171 # code = None
172 flagDecodeData = False # asumo q la data no esta decodificada
172 flagDecodeData = False # asumo q la data no esta decodificada
173 flagDeflipData = False # asumo q la data no esta sin flip
173 flagDeflipData = False # asumo q la data no esta sin flip
174 flagShiftFFT = False
174 flagShiftFFT = False
175 # ippSeconds = None
175 # ippSeconds = None
176 # timeInterval = None
176 # timeInterval = None
177 nCohInt = None
177 nCohInt = None
178 # noise = None
178 # noise = None
179 windowOfFilter = 1
179 windowOfFilter = 1
180 # Speed of ligth
180 # Speed of ligth
181 C = 3e8
181 C = 3e8
182 frequency = 49.92e6
182 frequency = 49.92e6
183 realtime = False
183 realtime = False
184 beacon_heiIndexList = None
184 beacon_heiIndexList = None
185 last_block = None
185 last_block = None
186 blocknow = None
186 blocknow = None
187 azimuth = None
187 azimuth = None
188 zenith = None
188 zenith = None
189 beam = Beam()
189 beam = Beam()
190 profileIndex = None
190 profileIndex = None
191 error = None
191 error = None
192 data = None
192 data = None
193 nmodes = None
193 nmodes = None
194
194
195 def __str__(self):
195 def __str__(self):
196
196
197 return '{} - {}'.format(self.type, self.getDatatime())
197 return '{} - {}'.format(self.type, self.getDatatime())
198
198
199 def getNoise(self):
199 def getNoise(self):
200
200
201 raise NotImplementedError
201 raise NotImplementedError
202
202
203 def getNChannels(self):
203 def getNChannels(self):
204
204
205 return len(self.channelList)
205 return len(self.channelList)
206
206
207 def getChannelIndexList(self):
207 def getChannelIndexList(self):
208
208
209 return list(range(self.nChannels))
209 return list(range(self.nChannels))
210
210
211 def getNHeights(self):
211 def getNHeights(self):
212
212
213 return len(self.heightList)
213 return len(self.heightList)
214
214
215 def getHeiRange(self, extrapoints=0):
215 def getHeiRange(self, extrapoints=0):
216
216
217 heis = self.heightList
217 heis = self.heightList
218 # deltah = self.heightList[1] - self.heightList[0]
218 # deltah = self.heightList[1] - self.heightList[0]
219 #
219 #
220 # heis.append(self.heightList[-1])
220 # heis.append(self.heightList[-1])
221
221
222 return heis
222 return heis
223
223
224 def getDeltaH(self):
224 def getDeltaH(self):
225
225
226 delta = self.heightList[1] - self.heightList[0]
226 delta = self.heightList[1] - self.heightList[0]
227
227
228 return delta
228 return delta
229
229
230 def getltctime(self):
230 def getltctime(self):
231
231
232 if self.useLocalTime:
232 if self.useLocalTime:
233 return self.utctime - self.timeZone * 60
233 return self.utctime - self.timeZone * 60
234
234
235 return self.utctime
235 return self.utctime
236
236
237 def getDatatime(self):
237 def getDatatime(self):
238
238
239 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
239 datatimeValue = datetime.datetime.utcfromtimestamp(self.ltctime)
240 return datatimeValue
240 return datatimeValue
241
241
242 def getTimeRange(self):
242 def getTimeRange(self):
243
243
244 datatime = []
244 datatime = []
245
245
246 datatime.append(self.ltctime)
246 datatime.append(self.ltctime)
247 datatime.append(self.ltctime + self.timeInterval + 1)
247 datatime.append(self.ltctime + self.timeInterval + 1)
248
248
249 datatime = numpy.array(datatime)
249 datatime = numpy.array(datatime)
250
250
251 return datatime
251 return datatime
252
252
253 def getFmaxTimeResponse(self):
253 def getFmaxTimeResponse(self):
254
254
255 period = (10**-6) * self.getDeltaH() / (0.15)
255 period = (10**-6) * self.getDeltaH() / (0.15)
256
256
257 PRF = 1. / (period * self.nCohInt)
257 PRF = 1. / (period * self.nCohInt)
258
258
259 fmax = PRF
259 fmax = PRF
260
260
261 return fmax
261 return fmax
262
262
263 def getFmax(self):
263 def getFmax(self):
264 PRF = 1. / (self.ippSeconds * self.nCohInt)
264 PRF = 1. / (self.ippSeconds * self.nCohInt)
265
265
266 fmax = PRF
266 fmax = PRF
267 return fmax
267 return fmax
268
268
269 def getVmax(self):
269 def getVmax(self):
270
270
271 _lambda = self.C / self.frequency
271 _lambda = self.C / self.frequency
272
272
273 vmax = self.getFmax() * _lambda / 2
273 vmax = self.getFmax() * _lambda / 2
274
274
275 return vmax
275 return vmax
276
276
277 def get_ippSeconds(self):
277 def get_ippSeconds(self):
278 '''
278 '''
279 '''
279 '''
280 return self.radarControllerHeaderObj.ippSeconds
280 return self.radarControllerHeaderObj.ippSeconds
281
281
282 def set_ippSeconds(self, ippSeconds):
282 def set_ippSeconds(self, ippSeconds):
283 '''
283 '''
284 '''
284 '''
285
285
286 self.radarControllerHeaderObj.ippSeconds = ippSeconds
286 self.radarControllerHeaderObj.ippSeconds = ippSeconds
287
287
288 return
288 return
289
289
290 def get_dtype(self):
290 def get_dtype(self):
291 '''
291 '''
292 '''
292 '''
293 return getNumpyDtype(self.datatype)
293 return getNumpyDtype(self.datatype)
294
294
295 def set_dtype(self, numpyDtype):
295 def set_dtype(self, numpyDtype):
296 '''
296 '''
297 '''
297 '''
298
298
299 self.datatype = getDataTypeCode(numpyDtype)
299 self.datatype = getDataTypeCode(numpyDtype)
300
300
301 def get_code(self):
301 def get_code(self):
302 '''
302 '''
303 '''
303 '''
304 return self.radarControllerHeaderObj.code
304 return self.radarControllerHeaderObj.code
305
305
306 def set_code(self, code):
306 def set_code(self, code):
307 '''
307 '''
308 '''
308 '''
309 self.radarControllerHeaderObj.code = code
309 self.radarControllerHeaderObj.code = code
310
310
311 return
311 return
312
312
313 def get_ncode(self):
313 def get_ncode(self):
314 '''
314 '''
315 '''
315 '''
316 return self.radarControllerHeaderObj.nCode
316 return self.radarControllerHeaderObj.nCode
317
317
318 def set_ncode(self, nCode):
318 def set_ncode(self, nCode):
319 '''
319 '''
320 '''
320 '''
321 self.radarControllerHeaderObj.nCode = nCode
321 self.radarControllerHeaderObj.nCode = nCode
322
322
323 return
323 return
324
324
325 def get_nbaud(self):
325 def get_nbaud(self):
326 '''
326 '''
327 '''
327 '''
328 return self.radarControllerHeaderObj.nBaud
328 return self.radarControllerHeaderObj.nBaud
329
329
330 def set_nbaud(self, nBaud):
330 def set_nbaud(self, nBaud):
331 '''
331 '''
332 '''
332 '''
333 self.radarControllerHeaderObj.nBaud = nBaud
333 self.radarControllerHeaderObj.nBaud = nBaud
334
334
335 return
335 return
336
336
337 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
337 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
338 channelIndexList = property(
338 channelIndexList = property(
339 getChannelIndexList, "I'm the 'channelIndexList' property.")
339 getChannelIndexList, "I'm the 'channelIndexList' property.")
340 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
340 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
341 #noise = property(getNoise, "I'm the 'nHeights' property.")
341 #noise = property(getNoise, "I'm the 'nHeights' property.")
342 datatime = property(getDatatime, "I'm the 'datatime' property")
342 datatime = property(getDatatime, "I'm the 'datatime' property")
343 ltctime = property(getltctime, "I'm the 'ltctime' property")
343 ltctime = property(getltctime, "I'm the 'ltctime' property")
344 ippSeconds = property(get_ippSeconds, set_ippSeconds)
344 ippSeconds = property(get_ippSeconds, set_ippSeconds)
345 dtype = property(get_dtype, set_dtype)
345 dtype = property(get_dtype, set_dtype)
346 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
346 # timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
347 code = property(get_code, set_code)
347 code = property(get_code, set_code)
348 nCode = property(get_ncode, set_ncode)
348 nCode = property(get_ncode, set_ncode)
349 nBaud = property(get_nbaud, set_nbaud)
349 nBaud = property(get_nbaud, set_nbaud)
350
350
351
351
352 class Voltage(JROData):
352 class Voltage(JROData):
353
353
354 # data es un numpy array de 2 dmensiones (canales, alturas)
354 # data es un numpy array de 2 dmensiones (canales, alturas)
355 data = None
355 data = None
356
356
357 def __init__(self):
357 def __init__(self):
358 '''
358 '''
359 Constructor
359 Constructor
360 '''
360 '''
361
361
362 self.useLocalTime = True
362 self.useLocalTime = True
363 self.radarControllerHeaderObj = RadarControllerHeader()
363 self.radarControllerHeaderObj = RadarControllerHeader()
364 self.systemHeaderObj = SystemHeader()
364 self.systemHeaderObj = SystemHeader()
365 self.type = "Voltage"
365 self.type = "Voltage"
366 self.data = None
366 self.data = None
367 # self.dtype = None
367 # self.dtype = None
368 # self.nChannels = 0
368 # self.nChannels = 0
369 # self.nHeights = 0
369 # self.nHeights = 0
370 self.nProfiles = None
370 self.nProfiles = None
371 self.heightList = None
371 self.heightList = None
372 self.channelList = None
372 self.channelList = None
373 # self.channelIndexList = None
373 # self.channelIndexList = None
374 self.flagNoData = True
374 self.flagNoData = True
375 self.flagDiscontinuousBlock = False
375 self.flagDiscontinuousBlock = False
376 self.utctime = None
376 self.utctime = None
377 self.timeZone = None
377 self.timeZone = None
378 self.dstFlag = None
378 self.dstFlag = None
379 self.errorCount = None
379 self.errorCount = None
380 self.nCohInt = None
380 self.nCohInt = None
381 self.blocksize = None
381 self.blocksize = None
382 self.flagDecodeData = False # asumo q la data no esta decodificada
382 self.flagDecodeData = False # asumo q la data no esta decodificada
383 self.flagDeflipData = False # asumo q la data no esta sin flip
383 self.flagDeflipData = False # asumo q la data no esta sin flip
384 self.flagShiftFFT = False
384 self.flagShiftFFT = False
385 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
385 self.flagDataAsBlock = False # Asumo que la data es leida perfil a perfil
386 self.profileIndex = 0
386 self.profileIndex = 0
387
387
388 def getNoisebyHildebrand(self, channel=None):
388 def getNoisebyHildebrand(self, channel=None):
389 """
389 """
390 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
390 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
391
391
392 Return:
392 Return:
393 noiselevel
393 noiselevel
394 """
394 """
395
395
396 if channel != None:
396 if channel != None:
397 data = self.data[channel]
397 data = self.data[channel]
398 nChannels = 1
398 nChannels = 1
399 else:
399 else:
400 data = self.data
400 data = self.data
401 nChannels = self.nChannels
401 nChannels = self.nChannels
402
402
403 noise = numpy.zeros(nChannels)
403 noise = numpy.zeros(nChannels)
404 power = data * numpy.conjugate(data)
404 power = data * numpy.conjugate(data)
405
405
406 for thisChannel in range(nChannels):
406 for thisChannel in range(nChannels):
407 if nChannels == 1:
407 if nChannels == 1:
408 daux = power[:].real
408 daux = power[:].real
409 else:
409 else:
410 daux = power[thisChannel, :].real
410 daux = power[thisChannel, :].real
411 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
411 noise[thisChannel] = hildebrand_sekhon(daux, self.nCohInt)
412
412
413 return noise
413 return noise
414
414
415 def getNoise(self, type=1, channel=None):
415 def getNoise(self, type=1, channel=None):
416
416
417 if type == 1:
417 if type == 1:
418 noise = self.getNoisebyHildebrand(channel)
418 noise = self.getNoisebyHildebrand(channel)
419
419
420 return noise
420 return noise
421
421
422 def getPower(self, channel=None):
422 def getPower(self, channel=None):
423
423
424 if channel != None:
424 if channel != None:
425 data = self.data[channel]
425 data = self.data[channel]
426 else:
426 else:
427 data = self.data
427 data = self.data
428
428
429 power = data * numpy.conjugate(data)
429 power = data * numpy.conjugate(data)
430 powerdB = 10 * numpy.log10(power.real)
430 powerdB = 10 * numpy.log10(power.real)
431 powerdB = numpy.squeeze(powerdB)
431 powerdB = numpy.squeeze(powerdB)
432
432
433 return powerdB
433 return powerdB
434
434
435 def getTimeInterval(self):
435 def getTimeInterval(self):
436
436
437 timeInterval = self.ippSeconds * self.nCohInt
437 timeInterval = self.ippSeconds * self.nCohInt
438
438
439 return timeInterval
439 return timeInterval
440
440
441 noise = property(getNoise, "I'm the 'nHeights' property.")
441 noise = property(getNoise, "I'm the 'nHeights' property.")
442 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
442 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
443
443
444
444
445 class Spectra(JROData):
445 class Spectra(JROData):
446
446
447 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
447 # data spc es un numpy array de 2 dmensiones (canales, perfiles, alturas)
448 data_spc = None
448 data_spc = None
449 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
449 # data cspc es un numpy array de 2 dmensiones (canales, pares, alturas)
450 data_cspc = None
450 data_cspc = None
451 # data dc es un numpy array de 2 dmensiones (canales, alturas)
451 # data dc es un numpy array de 2 dmensiones (canales, alturas)
452 data_dc = None
452 data_dc = None
453 # data power
453 # data power
454 data_pwr = None
454 data_pwr = None
455 nFFTPoints = None
455 nFFTPoints = None
456 # nPairs = None
456 # nPairs = None
457 pairsList = None
457 pairsList = None
458 nIncohInt = None
458 nIncohInt = None
459 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
459 wavelength = None # Necesario para cacular el rango de velocidad desde la frecuencia
460 nCohInt = None # se requiere para determinar el valor de timeInterval
460 nCohInt = None # se requiere para determinar el valor de timeInterval
461 ippFactor = None
461 ippFactor = None
462 profileIndex = 0
462 profileIndex = 0
463 plotting = "spectra"
463 plotting = "spectra"
464
464
465 def __init__(self):
465 def __init__(self):
466 '''
466 '''
467 Constructor
467 Constructor
468 '''
468 '''
469
469
470 self.useLocalTime = True
470 self.useLocalTime = True
471 self.radarControllerHeaderObj = RadarControllerHeader()
471 self.radarControllerHeaderObj = RadarControllerHeader()
472 self.systemHeaderObj = SystemHeader()
472 self.systemHeaderObj = SystemHeader()
473 self.type = "Spectra"
473 self.type = "Spectra"
474 # self.data = None
474 # self.data = None
475 # self.dtype = None
475 # self.dtype = None
476 # self.nChannels = 0
476 # self.nChannels = 0
477 # self.nHeights = 0
477 # self.nHeights = 0
478 self.nProfiles = None
478 self.nProfiles = None
479 self.heightList = None
479 self.heightList = None
480 self.channelList = None
480 self.channelList = None
481 # self.channelIndexList = None
481 # self.channelIndexList = None
482 self.pairsList = None
482 self.pairsList = None
483 self.flagNoData = True
483 self.flagNoData = True
484 self.flagDiscontinuousBlock = False
484 self.flagDiscontinuousBlock = False
485 self.utctime = None
485 self.utctime = None
486 self.nCohInt = None
486 self.nCohInt = None
487 self.nIncohInt = None
487 self.nIncohInt = None
488 self.blocksize = None
488 self.blocksize = None
489 self.nFFTPoints = None
489 self.nFFTPoints = None
490 self.wavelength = None
490 self.wavelength = None
491 self.flagDecodeData = False # asumo q la data no esta decodificada
491 self.flagDecodeData = False # asumo q la data no esta decodificada
492 self.flagDeflipData = False # asumo q la data no esta sin flip
492 self.flagDeflipData = False # asumo q la data no esta sin flip
493 self.flagShiftFFT = False
493 self.flagShiftFFT = False
494 self.ippFactor = 1
494 self.ippFactor = 1
495 #self.noise = None
495 #self.noise = None
496 self.beacon_heiIndexList = []
496 self.beacon_heiIndexList = []
497 self.noise_estimation = None
497 self.noise_estimation = None
498
498
499 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
499 def getNoisebyHildebrand(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
500 """
500 """
501 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
501 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
502
502
503 Return:
503 Return:
504 noiselevel
504 noiselevel
505 """
505 """
506
506
507 noise = numpy.zeros(self.nChannels)
507 noise = numpy.zeros(self.nChannels)
508
508
509 for channel in range(self.nChannels):
509 for channel in range(self.nChannels):
510 daux = self.data_spc[channel,
510 daux = self.data_spc[channel,
511 xmin_index:xmax_index, ymin_index:ymax_index]
511 xmin_index:xmax_index, ymin_index:ymax_index]
512 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
512 noise[channel] = hildebrand_sekhon(daux, self.nIncohInt)
513
513
514 return noise
514 return noise
515
515
516 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
516 def getNoise(self, xmin_index=None, xmax_index=None, ymin_index=None, ymax_index=None):
517
517
518 if self.noise_estimation is not None:
518 if self.noise_estimation is not None:
519 # this was estimated by getNoise Operation defined in jroproc_spectra.py
519 # this was estimated by getNoise Operation defined in jroproc_spectra.py
520 return self.noise_estimation
520 return self.noise_estimation
521 else:
521 else:
522 noise = self.getNoisebyHildebrand(
522 noise = self.getNoisebyHildebrand(
523 xmin_index, xmax_index, ymin_index, ymax_index)
523 xmin_index, xmax_index, ymin_index, ymax_index)
524 return noise
524 return noise
525
525
526 def getFreqRangeTimeResponse(self, extrapoints=0):
526 def getFreqRangeTimeResponse(self, extrapoints=0):
527
527
528 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
528 deltafreq = self.getFmaxTimeResponse() / (self.nFFTPoints * self.ippFactor)
529 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
529 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
530
530
531 return freqrange
531 return freqrange
532
532
533 def getAcfRange(self, extrapoints=0):
533 def getAcfRange(self, extrapoints=0):
534
534
535 deltafreq = 10. / (self.getFmax() / (self.nFFTPoints * self.ippFactor))
535 deltafreq = 10. / (self.getFmax() / (self.nFFTPoints * self.ippFactor))
536 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
536 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
537
537
538 return freqrange
538 return freqrange
539
539
540 def getFreqRange(self, extrapoints=0):
540 def getFreqRange(self, extrapoints=0):
541
541
542 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
542 deltafreq = self.getFmax() / (self.nFFTPoints * self.ippFactor)
543 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
543 freqrange = deltafreq * (numpy.arange(self.nFFTPoints + extrapoints) -self.nFFTPoints / 2.) - deltafreq / 2
544
544
545 return freqrange
545 return freqrange
546
546
547 def getVelRange(self, extrapoints=0):
547 def getVelRange(self, extrapoints=0):
548
548
549 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
549 deltav = self.getVmax() / (self.nFFTPoints * self.ippFactor)
550 velrange = deltav * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.)
550 velrange = deltav * (numpy.arange(self.nFFTPoints + extrapoints) - self.nFFTPoints / 2.)
551
551
552 if self.nmodes:
552 if self.nmodes:
553 return velrange/self.nmodes
553 return velrange/self.nmodes
554 else:
554 else:
555 return velrange
555 return velrange
556
556
557 def getNPairs(self):
557 def getNPairs(self):
558
558
559 return len(self.pairsList)
559 return len(self.pairsList)
560
560
561 def getPairsIndexList(self):
561 def getPairsIndexList(self):
562
562
563 return list(range(self.nPairs))
563 return list(range(self.nPairs))
564
564
565 def getNormFactor(self):
565 def getNormFactor(self):
566
566
567 pwcode = 1
567 pwcode = 1
568
568
569 if self.flagDecodeData:
569 if self.flagDecodeData:
570 pwcode = numpy.sum(self.code[0]**2)
570 pwcode = numpy.sum(self.code[0]**2)
571 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
571 #normFactor = min(self.nFFTPoints,self.nProfiles)*self.nIncohInt*self.nCohInt*pwcode*self.windowOfFilter
572 normFactor = self.nProfiles * self.nIncohInt * self.nCohInt * pwcode * self.windowOfFilter
572 normFactor = self.nProfiles * self.nIncohInt * self.nCohInt * pwcode * self.windowOfFilter
573
573
574 return normFactor
574 return normFactor
575
575
576 def getFlagCspc(self):
576 def getFlagCspc(self):
577
577
578 if self.data_cspc is None:
578 if self.data_cspc is None:
579 return True
579 return True
580
580
581 return False
581 return False
582
582
583 def getFlagDc(self):
583 def getFlagDc(self):
584
584
585 if self.data_dc is None:
585 if self.data_dc is None:
586 return True
586 return True
587
587
588 return False
588 return False
589
589
590 def getTimeInterval(self):
590 def getTimeInterval(self):
591
591
592 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
592 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt * self.nProfiles * self.ippFactor
593 if self.nmodes:
593 if self.nmodes:
594 return self.nmodes*timeInterval
594 return self.nmodes*timeInterval
595 else:
595 else:
596 return timeInterval
596 return timeInterval
597
597
598 def getPower(self):
598 def getPower(self):
599
599
600 factor = self.normFactor
600 factor = self.normFactor
601 z = self.data_spc / factor
601 z = self.data_spc / factor
602 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
602 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
603 avg = numpy.average(z, axis=1)
603 avg = numpy.average(z, axis=1)
604
604
605 return 10 * numpy.log10(avg)
605 return 10 * numpy.log10(avg)
606
606
607 def getCoherence(self, pairsList=None, phase=False):
607 def getCoherence(self, pairsList=None, phase=False):
608
608
609 z = []
609 z = []
610 if pairsList is None:
610 if pairsList is None:
611 pairsIndexList = self.pairsIndexList
611 pairsIndexList = self.pairsIndexList
612 else:
612 else:
613 pairsIndexList = []
613 pairsIndexList = []
614 for pair in pairsList:
614 for pair in pairsList:
615 if pair not in self.pairsList:
615 if pair not in self.pairsList:
616 raise ValueError("Pair %s is not in dataOut.pairsList" % (
616 raise ValueError("Pair %s is not in dataOut.pairsList" % (
617 pair))
617 pair))
618 pairsIndexList.append(self.pairsList.index(pair))
618 pairsIndexList.append(self.pairsList.index(pair))
619 for i in range(len(pairsIndexList)):
619 for i in range(len(pairsIndexList)):
620 pair = self.pairsList[pairsIndexList[i]]
620 pair = self.pairsList[pairsIndexList[i]]
621 ccf = numpy.average(self.data_cspc[pairsIndexList[i], :, :], axis=0)
621 ccf = numpy.average(self.data_cspc[pairsIndexList[i], :, :], axis=0)
622 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
622 powa = numpy.average(self.data_spc[pair[0], :, :], axis=0)
623 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
623 powb = numpy.average(self.data_spc[pair[1], :, :], axis=0)
624 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
624 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
625 if phase:
625 if phase:
626 data = numpy.arctan2(avgcoherenceComplex.imag,
626 data = numpy.arctan2(avgcoherenceComplex.imag,
627 avgcoherenceComplex.real) * 180 / numpy.pi
627 avgcoherenceComplex.real) * 180 / numpy.pi
628 else:
628 else:
629 data = numpy.abs(avgcoherenceComplex)
629 data = numpy.abs(avgcoherenceComplex)
630
630
631 z.append(data)
631 z.append(data)
632
632
633 return numpy.array(z)
633 return numpy.array(z)
634
634
635 def setValue(self, value):
635 def setValue(self, value):
636
636
637 print("This property should not be initialized")
637 print("This property should not be initialized")
638
638
639 return
639 return
640
640
641 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
641 nPairs = property(getNPairs, setValue, "I'm the 'nPairs' property.")
642 pairsIndexList = property(
642 pairsIndexList = property(
643 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
643 getPairsIndexList, setValue, "I'm the 'pairsIndexList' property.")
644 normFactor = property(getNormFactor, setValue,
644 normFactor = property(getNormFactor, setValue,
645 "I'm the 'getNormFactor' property.")
645 "I'm the 'getNormFactor' property.")
646 flag_cspc = property(getFlagCspc, setValue)
646 flag_cspc = property(getFlagCspc, setValue)
647 flag_dc = property(getFlagDc, setValue)
647 flag_dc = property(getFlagDc, setValue)
648 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
648 noise = property(getNoise, setValue, "I'm the 'nHeights' property.")
649 timeInterval = property(getTimeInterval, setValue,
649 timeInterval = property(getTimeInterval, setValue,
650 "I'm the 'timeInterval' property")
650 "I'm the 'timeInterval' property")
651
651
652
652
653 class SpectraHeis(Spectra):
653 class SpectraHeis(Spectra):
654
654
655 data_spc = None
655 data_spc = None
656 data_cspc = None
656 data_cspc = None
657 data_dc = None
657 data_dc = None
658 nFFTPoints = None
658 nFFTPoints = None
659 # nPairs = None
659 # nPairs = None
660 pairsList = None
660 pairsList = None
661 nCohInt = None
661 nCohInt = None
662 nIncohInt = None
662 nIncohInt = None
663
663
664 def __init__(self):
664 def __init__(self):
665
665
666 self.radarControllerHeaderObj = RadarControllerHeader()
666 self.radarControllerHeaderObj = RadarControllerHeader()
667
667
668 self.systemHeaderObj = SystemHeader()
668 self.systemHeaderObj = SystemHeader()
669
669
670 self.type = "SpectraHeis"
670 self.type = "SpectraHeis"
671
671
672 # self.dtype = None
672 # self.dtype = None
673
673
674 # self.nChannels = 0
674 # self.nChannels = 0
675
675
676 # self.nHeights = 0
676 # self.nHeights = 0
677
677
678 self.nProfiles = None
678 self.nProfiles = None
679
679
680 self.heightList = None
680 self.heightList = None
681
681
682 self.channelList = None
682 self.channelList = None
683
683
684 # self.channelIndexList = None
684 # self.channelIndexList = None
685
685
686 self.flagNoData = True
686 self.flagNoData = True
687
687
688 self.flagDiscontinuousBlock = False
688 self.flagDiscontinuousBlock = False
689
689
690 # self.nPairs = 0
690 # self.nPairs = 0
691
691
692 self.utctime = None
692 self.utctime = None
693
693
694 self.blocksize = None
694 self.blocksize = None
695
695
696 self.profileIndex = 0
696 self.profileIndex = 0
697
697
698 self.nCohInt = 1
698 self.nCohInt = 1
699
699
700 self.nIncohInt = 1
700 self.nIncohInt = 1
701
701
702 def getNormFactor(self):
702 def getNormFactor(self):
703 pwcode = 1
703 pwcode = 1
704 if self.flagDecodeData:
704 if self.flagDecodeData:
705 pwcode = numpy.sum(self.code[0]**2)
705 pwcode = numpy.sum(self.code[0]**2)
706
706
707 normFactor = self.nIncohInt * self.nCohInt * pwcode
707 normFactor = self.nIncohInt * self.nCohInt * pwcode
708
708
709 return normFactor
709 return normFactor
710
710
711 def getTimeInterval(self):
711 def getTimeInterval(self):
712
712
713 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
713 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
714
714
715 return timeInterval
715 return timeInterval
716
716
717 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
717 normFactor = property(getNormFactor, "I'm the 'getNormFactor' property.")
718 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
718 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
719
719
720
720
721 class Fits(JROData):
721 class Fits(JROData):
722
722
723 heightList = None
723 heightList = None
724 channelList = None
724 channelList = None
725 flagNoData = True
725 flagNoData = True
726 flagDiscontinuousBlock = False
726 flagDiscontinuousBlock = False
727 useLocalTime = False
727 useLocalTime = False
728 utctime = None
728 utctime = None
729 timeZone = None
729 timeZone = None
730 # ippSeconds = None
730 # ippSeconds = None
731 # timeInterval = None
731 # timeInterval = None
732 nCohInt = None
732 nCohInt = None
733 nIncohInt = None
733 nIncohInt = None
734 noise = None
734 noise = None
735 windowOfFilter = 1
735 windowOfFilter = 1
736 # Speed of ligth
736 # Speed of ligth
737 C = 3e8
737 C = 3e8
738 frequency = 49.92e6
738 frequency = 49.92e6
739 realtime = False
739 realtime = False
740
740
741 def __init__(self):
741 def __init__(self):
742
742
743 self.type = "Fits"
743 self.type = "Fits"
744
744
745 self.nProfiles = None
745 self.nProfiles = None
746
746
747 self.heightList = None
747 self.heightList = None
748
748
749 self.channelList = None
749 self.channelList = None
750
750
751 # self.channelIndexList = None
751 # self.channelIndexList = None
752
752
753 self.flagNoData = True
753 self.flagNoData = True
754
754
755 self.utctime = None
755 self.utctime = None
756
756
757 self.nCohInt = 1
757 self.nCohInt = 1
758
758
759 self.nIncohInt = 1
759 self.nIncohInt = 1
760
760
761 self.useLocalTime = True
761 self.useLocalTime = True
762
762
763 self.profileIndex = 0
763 self.profileIndex = 0
764
764
765 # self.utctime = None
765 # self.utctime = None
766 # self.timeZone = None
766 # self.timeZone = None
767 # self.ltctime = None
767 # self.ltctime = None
768 # self.timeInterval = None
768 # self.timeInterval = None
769 # self.header = None
769 # self.header = None
770 # self.data_header = None
770 # self.data_header = None
771 # self.data = None
771 # self.data = None
772 # self.datatime = None
772 # self.datatime = None
773 # self.flagNoData = False
773 # self.flagNoData = False
774 # self.expName = ''
774 # self.expName = ''
775 # self.nChannels = None
775 # self.nChannels = None
776 # self.nSamples = None
776 # self.nSamples = None
777 # self.dataBlocksPerFile = None
777 # self.dataBlocksPerFile = None
778 # self.comments = ''
778 # self.comments = ''
779 #
779 #
780
780
781 def getltctime(self):
781 def getltctime(self):
782
782
783 if self.useLocalTime:
783 if self.useLocalTime:
784 return self.utctime - self.timeZone * 60
784 return self.utctime - self.timeZone * 60
785
785
786 return self.utctime
786 return self.utctime
787
787
788 def getDatatime(self):
788 def getDatatime(self):
789
789
790 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
790 datatime = datetime.datetime.utcfromtimestamp(self.ltctime)
791 return datatime
791 return datatime
792
792
793 def getTimeRange(self):
793 def getTimeRange(self):
794
794
795 datatime = []
795 datatime = []
796
796
797 datatime.append(self.ltctime)
797 datatime.append(self.ltctime)
798 datatime.append(self.ltctime + self.timeInterval)
798 datatime.append(self.ltctime + self.timeInterval)
799
799
800 datatime = numpy.array(datatime)
800 datatime = numpy.array(datatime)
801
801
802 return datatime
802 return datatime
803
803
804 def getHeiRange(self):
804 def getHeiRange(self):
805
805
806 heis = self.heightList
806 heis = self.heightList
807
807
808 return heis
808 return heis
809
809
810 def getNHeights(self):
810 def getNHeights(self):
811
811
812 return len(self.heightList)
812 return len(self.heightList)
813
813
814 def getNChannels(self):
814 def getNChannels(self):
815
815
816 return len(self.channelList)
816 return len(self.channelList)
817
817
818 def getChannelIndexList(self):
818 def getChannelIndexList(self):
819
819
820 return list(range(self.nChannels))
820 return list(range(self.nChannels))
821
821
822 def getNoise(self, type=1):
822 def getNoise(self, type=1):
823
823
824 #noise = numpy.zeros(self.nChannels)
824 #noise = numpy.zeros(self.nChannels)
825
825
826 if type == 1:
826 if type == 1:
827 noise = self.getNoisebyHildebrand()
827 noise = self.getNoisebyHildebrand()
828
828
829 if type == 2:
829 if type == 2:
830 noise = self.getNoisebySort()
830 noise = self.getNoisebySort()
831
831
832 if type == 3:
832 if type == 3:
833 noise = self.getNoisebyWindow()
833 noise = self.getNoisebyWindow()
834
834
835 return noise
835 return noise
836
836
837 def getTimeInterval(self):
837 def getTimeInterval(self):
838
838
839 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
839 timeInterval = self.ippSeconds * self.nCohInt * self.nIncohInt
840
840
841 return timeInterval
841 return timeInterval
842
842
843 def get_ippSeconds(self):
843 def get_ippSeconds(self):
844 '''
844 '''
845 '''
845 '''
846 return self.ipp_sec
846 return self.ipp_sec
847
847
848
848
849 datatime = property(getDatatime, "I'm the 'datatime' property")
849 datatime = property(getDatatime, "I'm the 'datatime' property")
850 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
850 nHeights = property(getNHeights, "I'm the 'nHeights' property.")
851 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
851 nChannels = property(getNChannels, "I'm the 'nChannel' property.")
852 channelIndexList = property(
852 channelIndexList = property(
853 getChannelIndexList, "I'm the 'channelIndexList' property.")
853 getChannelIndexList, "I'm the 'channelIndexList' property.")
854 noise = property(getNoise, "I'm the 'nHeights' property.")
854 noise = property(getNoise, "I'm the 'nHeights' property.")
855
855
856 ltctime = property(getltctime, "I'm the 'ltctime' property")
856 ltctime = property(getltctime, "I'm the 'ltctime' property")
857 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
857 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
858 ippSeconds = property(get_ippSeconds, '')
858 ippSeconds = property(get_ippSeconds, '')
859
859
860 class Correlation(JROData):
860 class Correlation(JROData):
861
861
862 noise = None
862 noise = None
863 SNR = None
863 SNR = None
864 #--------------------------------------------------
864 #--------------------------------------------------
865 mode = None
865 mode = None
866 split = False
866 split = False
867 data_cf = None
867 data_cf = None
868 lags = None
868 lags = None
869 lagRange = None
869 lagRange = None
870 pairsList = None
870 pairsList = None
871 normFactor = None
871 normFactor = None
872 #--------------------------------------------------
872 #--------------------------------------------------
873 # calculateVelocity = None
873 # calculateVelocity = None
874 nLags = None
874 nLags = None
875 nPairs = None
875 nPairs = None
876 nAvg = None
876 nAvg = None
877
877
878 def __init__(self):
878 def __init__(self):
879 '''
879 '''
880 Constructor
880 Constructor
881 '''
881 '''
882 self.radarControllerHeaderObj = RadarControllerHeader()
882 self.radarControllerHeaderObj = RadarControllerHeader()
883
883
884 self.systemHeaderObj = SystemHeader()
884 self.systemHeaderObj = SystemHeader()
885
885
886 self.type = "Correlation"
886 self.type = "Correlation"
887
887
888 self.data = None
888 self.data = None
889
889
890 self.dtype = None
890 self.dtype = None
891
891
892 self.nProfiles = None
892 self.nProfiles = None
893
893
894 self.heightList = None
894 self.heightList = None
895
895
896 self.channelList = None
896 self.channelList = None
897
897
898 self.flagNoData = True
898 self.flagNoData = True
899
899
900 self.flagDiscontinuousBlock = False
900 self.flagDiscontinuousBlock = False
901
901
902 self.utctime = None
902 self.utctime = None
903
903
904 self.timeZone = None
904 self.timeZone = None
905
905
906 self.dstFlag = None
906 self.dstFlag = None
907
907
908 self.errorCount = None
908 self.errorCount = None
909
909
910 self.blocksize = None
910 self.blocksize = None
911
911
912 self.flagDecodeData = False # asumo q la data no esta decodificada
912 self.flagDecodeData = False # asumo q la data no esta decodificada
913
913
914 self.flagDeflipData = False # asumo q la data no esta sin flip
914 self.flagDeflipData = False # asumo q la data no esta sin flip
915
915
916 self.pairsList = None
916 self.pairsList = None
917
917
918 self.nPoints = None
918 self.nPoints = None
919
919
920 def getPairsList(self):
920 def getPairsList(self):
921
921
922 return self.pairsList
922 return self.pairsList
923
923
924 def getNoise(self, mode=2):
924 def getNoise(self, mode=2):
925
925
926 indR = numpy.where(self.lagR == 0)[0][0]
926 indR = numpy.where(self.lagR == 0)[0][0]
927 indT = numpy.where(self.lagT == 0)[0][0]
927 indT = numpy.where(self.lagT == 0)[0][0]
928
928
929 jspectra0 = self.data_corr[:, :, indR, :]
929 jspectra0 = self.data_corr[:, :, indR, :]
930 jspectra = copy.copy(jspectra0)
930 jspectra = copy.copy(jspectra0)
931
931
932 num_chan = jspectra.shape[0]
932 num_chan = jspectra.shape[0]
933 num_hei = jspectra.shape[2]
933 num_hei = jspectra.shape[2]
934
934
935 freq_dc = jspectra.shape[1] / 2
935 freq_dc = jspectra.shape[1] / 2
936 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
936 ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
937
937
938 if ind_vel[0] < 0:
938 if ind_vel[0] < 0:
939 ind_vel[list(range(0, 1))] = ind_vel[list(
939 ind_vel[list(range(0, 1))] = ind_vel[list(
940 range(0, 1))] + self.num_prof
940 range(0, 1))] + self.num_prof
941
941
942 if mode == 1:
942 if mode == 1:
943 jspectra[:, freq_dc, :] = (
943 jspectra[:, freq_dc, :] = (
944 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
944 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
945
945
946 if mode == 2:
946 if mode == 2:
947
947
948 vel = numpy.array([-2, -1, 1, 2])
948 vel = numpy.array([-2, -1, 1, 2])
949 xx = numpy.zeros([4, 4])
949 xx = numpy.zeros([4, 4])
950
950
951 for fil in range(4):
951 for fil in range(4):
952 xx[fil, :] = vel[fil]**numpy.asarray(list(range(4)))
952 xx[fil, :] = vel[fil]**numpy.asarray(list(range(4)))
953
953
954 xx_inv = numpy.linalg.inv(xx)
954 xx_inv = numpy.linalg.inv(xx)
955 xx_aux = xx_inv[0, :]
955 xx_aux = xx_inv[0, :]
956
956
957 for ich in range(num_chan):
957 for ich in range(num_chan):
958 yy = jspectra[ich, ind_vel, :]
958 yy = jspectra[ich, ind_vel, :]
959 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
959 jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
960
960
961 junkid = jspectra[ich, freq_dc, :] <= 0
961 junkid = jspectra[ich, freq_dc, :] <= 0
962 cjunkid = sum(junkid)
962 cjunkid = sum(junkid)
963
963
964 if cjunkid.any():
964 if cjunkid.any():
965 jspectra[ich, freq_dc, junkid.nonzero()] = (
965 jspectra[ich, freq_dc, junkid.nonzero()] = (
966 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
966 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
967
967
968 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
968 noise = jspectra0[:, freq_dc, :] - jspectra[:, freq_dc, :]
969
969
970 return noise
970 return noise
971
971
972 def getTimeInterval(self):
972 def getTimeInterval(self):
973
973
974 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
974 timeInterval = self.ippSeconds * self.nCohInt * self.nProfiles
975
975
976 return timeInterval
976 return timeInterval
977
977
978 def splitFunctions(self):
978 def splitFunctions(self):
979
979
980 pairsList = self.pairsList
980 pairsList = self.pairsList
981 ccf_pairs = []
981 ccf_pairs = []
982 acf_pairs = []
982 acf_pairs = []
983 ccf_ind = []
983 ccf_ind = []
984 acf_ind = []
984 acf_ind = []
985 for l in range(len(pairsList)):
985 for l in range(len(pairsList)):
986 chan0 = pairsList[l][0]
986 chan0 = pairsList[l][0]
987 chan1 = pairsList[l][1]
987 chan1 = pairsList[l][1]
988
988
989 # Obteniendo pares de Autocorrelacion
989 # Obteniendo pares de Autocorrelacion
990 if chan0 == chan1:
990 if chan0 == chan1:
991 acf_pairs.append(chan0)
991 acf_pairs.append(chan0)
992 acf_ind.append(l)
992 acf_ind.append(l)
993 else:
993 else:
994 ccf_pairs.append(pairsList[l])
994 ccf_pairs.append(pairsList[l])
995 ccf_ind.append(l)
995 ccf_ind.append(l)
996
996
997 data_acf = self.data_cf[acf_ind]
997 data_acf = self.data_cf[acf_ind]
998 data_ccf = self.data_cf[ccf_ind]
998 data_ccf = self.data_cf[ccf_ind]
999
999
1000 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1000 return acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf
1001
1001
1002 def getNormFactor(self):
1002 def getNormFactor(self):
1003 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1003 acf_ind, ccf_ind, acf_pairs, ccf_pairs, data_acf, data_ccf = self.splitFunctions()
1004 acf_pairs = numpy.array(acf_pairs)
1004 acf_pairs = numpy.array(acf_pairs)
1005 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1005 normFactor = numpy.zeros((self.nPairs, self.nHeights))
1006
1006
1007 for p in range(self.nPairs):
1007 for p in range(self.nPairs):
1008 pair = self.pairsList[p]
1008 pair = self.pairsList[p]
1009
1009
1010 ch0 = pair[0]
1010 ch0 = pair[0]
1011 ch1 = pair[1]
1011 ch1 = pair[1]
1012
1012
1013 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1013 ch0_max = numpy.max(data_acf[acf_pairs == ch0, :, :], axis=1)
1014 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1014 ch1_max = numpy.max(data_acf[acf_pairs == ch1, :, :], axis=1)
1015 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1015 normFactor[p, :] = numpy.sqrt(ch0_max * ch1_max)
1016
1016
1017 return normFactor
1017 return normFactor
1018
1018
1019 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1019 timeInterval = property(getTimeInterval, "I'm the 'timeInterval' property")
1020 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1020 normFactor = property(getNormFactor, "I'm the 'normFactor property'")
1021
1021
1022
1022
1023 class Parameters(Spectra):
1023 class Parameters(Spectra):
1024
1024
1025 experimentInfo = None # Information about the experiment
1025 experimentInfo = None # Information about the experiment
1026 # Information from previous data
1026 # Information from previous data
1027 inputUnit = None # Type of data to be processed
1027 inputUnit = None # Type of data to be processed
1028 operation = None # Type of operation to parametrize
1028 operation = None # Type of operation to parametrize
1029 # normFactor = None #Normalization Factor
1029 # normFactor = None #Normalization Factor
1030 groupList = None # List of Pairs, Groups, etc
1030 groupList = None # List of Pairs, Groups, etc
1031 # Parameters
1031 # Parameters
1032 data_param = None # Parameters obtained
1032 data_param = None # Parameters obtained
1033 data_pre = None # Data Pre Parametrization
1033 data_pre = None # Data Pre Parametrization
1034 data_SNR = None # Signal to Noise Ratio
1034 data_SNR = None # Signal to Noise Ratio
1035 # heightRange = None #Heights
1035 # heightRange = None #Heights
1036 abscissaList = None # Abscissa, can be velocities, lags or time
1036 abscissaList = None # Abscissa, can be velocities, lags or time
1037 # noise = None #Noise Potency
1037 # noise = None #Noise Potency
1038 utctimeInit = None # Initial UTC time
1038 utctimeInit = None # Initial UTC time
1039 paramInterval = None # Time interval to calculate Parameters in seconds
1039 paramInterval = None # Time interval to calculate Parameters in seconds
1040 useLocalTime = True
1040 useLocalTime = True
1041 # Fitting
1041 # Fitting
1042 data_error = None # Error of the estimation
1042 data_error = None # Error of the estimation
1043 constants = None
1043 constants = None
1044 library = None
1044 library = None
1045 # Output signal
1045 # Output signal
1046 outputInterval = None # Time interval to calculate output signal in seconds
1046 outputInterval = None # Time interval to calculate output signal in seconds
1047 data_output = None # Out signal
1047 data_output = None # Out signal
1048 nAvg = None
1048 nAvg = None
1049 noise_estimation = None
1049 noise_estimation = None
1050 GauSPC = None # Fit gaussian SPC
1050 GauSPC = None # Fit gaussian SPC
1051
1051
1052 def __init__(self):
1052 def __init__(self):
1053 '''
1053 '''
1054 Constructor
1054 Constructor
1055 '''
1055 '''
1056 self.radarControllerHeaderObj = RadarControllerHeader()
1056 self.radarControllerHeaderObj = RadarControllerHeader()
1057
1057
1058 self.systemHeaderObj = SystemHeader()
1058 self.systemHeaderObj = SystemHeader()
1059
1059
1060 self.type = "Parameters"
1060 self.type = "Parameters"
1061
1061
1062 def getTimeRange1(self, interval):
1062 def getTimeRange1(self, interval):
1063
1063
1064 datatime = []
1064 datatime = []
1065
1065
1066 if self.useLocalTime:
1066 if self.useLocalTime:
1067 time1 = self.utctimeInit - self.timeZone * 60
1067 time1 = self.utctimeInit - self.timeZone * 60
1068 else:
1068 else:
1069 time1 = self.utctimeInit
1069 time1 = self.utctimeInit
1070
1070
1071 datatime.append(time1)
1071 datatime.append(time1)
1072 datatime.append(time1 + interval)
1072 datatime.append(time1 + interval)
1073 datatime = numpy.array(datatime)
1073 datatime = numpy.array(datatime)
1074
1074
1075 return datatime
1075 return datatime
1076
1076
1077 def getTimeInterval(self):
1077 def getTimeInterval(self):
1078
1078
1079 if hasattr(self, 'timeInterval1'):
1079 if hasattr(self, 'timeInterval1'):
1080 return self.timeInterval1
1080 return self.timeInterval1
1081 else:
1081 else:
1082 return self.paramInterval
1082 return self.paramInterval
1083
1083
1084 def setValue(self, value):
1084 def setValue(self, value):
1085
1085
1086 print("This property should not be initialized")
1086 print("This property should not be initialized")
1087
1087
1088 return
1088 return
1089
1089
1090 def getNoise(self):
1090 def getNoise(self):
1091
1091
1092 return self.spc_noise
1092 return self.spc_noise
1093
1093
1094 timeInterval = property(getTimeInterval)
1094 timeInterval = property(getTimeInterval)
1095 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
1095 noise = property(getNoise, setValue, "I'm the 'Noise' property.")
1096
1096
1097
1097
1098 class PlotterData(object):
1098 class PlotterData(object):
1099 '''
1099 '''
1100 Object to hold data to be plotted
1100 Object to hold data to be plotted
1101 '''
1101 '''
1102
1102
1103 MAXNUMX = 100
1103 MAXNUMX = 100
1104 MAXNUMY = 100
1104 MAXNUMY = 100
1105
1105
1106 def __init__(self, code, throttle_value, exp_code, buffering=True, snr=False):
1106 def __init__(self, code, throttle_value, exp_code, buffering=True, snr=False):
1107
1107
1108 self.key = code
1108 self.throttle = throttle_value
1109 self.throttle = throttle_value
1109 self.exp_code = exp_code
1110 self.exp_code = exp_code
1110 self.buffering = buffering
1111 self.buffering = buffering
1111 self.ready = False
1112 self.ready = False
1112 self.localtime = False
1113 self.localtime = False
1113 self.data = {}
1114 self.data = {}
1114 self.meta = {}
1115 self.meta = {}
1115 self.__times = []
1116 self.__times = []
1116 self.__heights = []
1117 self.__heights = []
1117
1118
1118 if 'snr' in code:
1119 if 'snr' in code:
1119 self.plottypes = ['snr']
1120 self.plottypes = ['snr']
1120 elif code == 'spc':
1121 elif code == 'spc':
1121 self.plottypes = ['spc', 'noise', 'rti']
1122 self.plottypes = ['spc', 'noise', 'rti']
1122 elif code == 'rti':
1123 elif code == 'rti':
1123 self.plottypes = ['noise', 'rti']
1124 self.plottypes = ['noise', 'rti']
1124 else:
1125 else:
1125 self.plottypes = [code]
1126 self.plottypes = [code]
1126
1127
1127 if 'snr' not in self.plottypes and snr:
1128 if 'snr' not in self.plottypes and snr:
1128 self.plottypes.append('snr')
1129 self.plottypes.append('snr')
1129
1130
1130 for plot in self.plottypes:
1131 for plot in self.plottypes:
1131 self.data[plot] = {}
1132 self.data[plot] = {}
1132
1133
1133 def __str__(self):
1134 def __str__(self):
1134 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1135 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1135 return 'Data[{}][{}]'.format(';'.join(dum), len(self.__times))
1136 return 'Data[{}][{}]'.format(';'.join(dum), len(self.__times))
1136
1137
1137 def __len__(self):
1138 def __len__(self):
1138 return len(self.__times)
1139 return len(self.__times)
1139
1140
1140 def __getitem__(self, key):
1141 def __getitem__(self, key):
1141
1142
1142 if key not in self.data:
1143 if key not in self.data:
1143 raise KeyError(log.error('Missing key: {}'.format(key)))
1144 raise KeyError(log.error('Missing key: {}'.format(key)))
1144 if 'spc' in key or not self.buffering:
1145 if 'spc' in key or not self.buffering:
1145 ret = self.data[key]
1146 ret = self.data[key]
1146 elif 'scope' in key:
1147 elif 'scope' in key:
1147 ret = numpy.array(self.data[key][float(self.tm)])
1148 ret = numpy.array(self.data[key][float(self.tm)])
1148 else:
1149 else:
1149 ret = numpy.array([self.data[key][x] for x in self.times])
1150 ret = numpy.array([self.data[key][x] for x in self.times])
1150 if ret.ndim > 1:
1151 if ret.ndim > 1:
1151 ret = numpy.swapaxes(ret, 0, 1)
1152 ret = numpy.swapaxes(ret, 0, 1)
1152 return ret
1153 return ret
1153
1154
1154 def __contains__(self, key):
1155 def __contains__(self, key):
1155 return key in self.data
1156 return key in self.data
1156
1157
1157 def setup(self):
1158 def setup(self):
1158 '''
1159 '''
1159 Configure object
1160 Configure object
1160 '''
1161 '''
1161
1162
1162 self.type = ''
1163 self.type = ''
1163 self.ready = False
1164 self.ready = False
1164 self.data = {}
1165 self.data = {}
1165 self.__times = []
1166 self.__times = []
1166 self.__heights = []
1167 self.__heights = []
1167 self.__all_heights = set()
1168 self.__all_heights = set()
1168 for plot in self.plottypes:
1169 for plot in self.plottypes:
1169 if 'snr' in plot:
1170 if 'snr' in plot:
1170 plot = 'snr'
1171 plot = 'snr'
1171 elif 'spc_moments' == plot:
1172 elif 'spc_moments' == plot:
1172 plot = 'moments'
1173 plot = 'moments'
1173 self.data[plot] = {}
1174 self.data[plot] = {}
1174
1175
1175 if 'spc' in self.data or 'rti' in self.data or 'cspc' in self.data or 'moments' in self.data:
1176 if 'spc' in self.data or 'rti' in self.data or 'cspc' in self.data or 'moments' in self.data:
1176 self.data['noise'] = {}
1177 self.data['noise'] = {}
1177 self.data['rti'] = {}
1178 self.data['rti'] = {}
1178 if 'noise' not in self.plottypes:
1179 if 'noise' not in self.plottypes:
1179 self.plottypes.append('noise')
1180 self.plottypes.append('noise')
1180 if 'rti' not in self.plottypes:
1181 if 'rti' not in self.plottypes:
1181 self.plottypes.append('rti')
1182 self.plottypes.append('rti')
1182
1183
1183 def shape(self, key):
1184 def shape(self, key):
1184 '''
1185 '''
1185 Get the shape of the one-element data for the given key
1186 Get the shape of the one-element data for the given key
1186 '''
1187 '''
1187
1188
1188 if len(self.data[key]):
1189 if len(self.data[key]):
1189 if 'spc' in key or not self.buffering:
1190 if 'spc' in key or not self.buffering:
1190 return self.data[key].shape
1191 return self.data[key].shape
1191 return self.data[key][self.__times[0]].shape
1192 return self.data[key][self.__times[0]].shape
1192 return (0,)
1193 return (0,)
1193
1194
1194 def update(self, dataOut, tm):
1195 def update(self, dataOut, tm):
1195 '''
1196 '''
1196 Update data object with new dataOut
1197 Update data object with new dataOut
1197 '''
1198 '''
1198
1199
1199 if tm in self.__times:
1200 if tm in self.__times:
1200 return
1201 return
1201 self.profileIndex = dataOut.profileIndex
1202 self.profileIndex = dataOut.profileIndex
1202 self.tm = tm
1203 self.tm = tm
1203 self.type = dataOut.type
1204 self.type = dataOut.type
1204 self.parameters = getattr(dataOut, 'parameters', [])
1205 self.parameters = getattr(dataOut, 'parameters', [])
1205 if hasattr(dataOut, 'pairsList'):
1206 self.pairs = dataOut.pairsList
1207 if hasattr(dataOut, 'meta'):
1206 if hasattr(dataOut, 'meta'):
1208 self.meta = dataOut.meta
1207 self.meta.update(dataOut.meta)
1209 self.channels = dataOut.channelList
1208 self.channels = dataOut.channelList
1210 self.interval = dataOut.getTimeInterval()
1209 self.interval = dataOut.getTimeInterval()
1211 self.localtime = dataOut.useLocalTime
1210 self.localtime = dataOut.useLocalTime
1212 if 'spc' in self.plottypes or 'cspc' in self.plottypes or 'spc_moments' in self.plottypes:
1211 if 'spc' in self.plottypes or 'cspc' in self.plottypes or 'spc_moments' in self.plottypes:
1213 self.xrange = (dataOut.getFreqRange(1)/1000.,
1212 self.xrange = (dataOut.getFreqRange(1)/1000.,
1214 dataOut.getAcfRange(1), dataOut.getVelRange(1))
1213 dataOut.getAcfRange(1), dataOut.getVelRange(1))
1215 self.factor = dataOut.normFactor
1214 self.factor = dataOut.normFactor
1216 self.__heights.append(dataOut.heightList)
1215 self.__heights.append(dataOut.heightList)
1217 self.__all_heights.update(dataOut.heightList)
1216 self.__all_heights.update(dataOut.heightList)
1218 self.__times.append(tm)
1217 self.__times.append(tm)
1219
1218
1220 for plot in self.plottypes:
1219 for plot in self.plottypes:
1221 if plot in ('spc', 'spc_moments'):
1220 if plot in ('spc', 'spc_moments'):
1222 z = dataOut.data_spc/dataOut.normFactor
1221 z = dataOut.data_spc/dataOut.normFactor
1223 buffer = 10*numpy.log10(z)
1222 buffer = 10*numpy.log10(z)
1224 if plot == 'cspc':
1223 if plot == 'cspc':
1225 z = dataOut.data_spc/dataOut.normFactor
1224 z = dataOut.data_spc/dataOut.normFactor
1226 buffer = (dataOut.data_spc, dataOut.data_cspc)
1225 buffer = (dataOut.data_spc, dataOut.data_cspc)
1227 if plot == 'noise':
1226 if plot == 'noise':
1228 buffer = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
1227 buffer = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
1229 if plot == 'rti':
1228 if plot == 'rti':
1230 buffer = dataOut.getPower()
1229 buffer = dataOut.getPower()
1231 if plot == 'snr_db':
1230 if plot == 'snr_db':
1232 buffer = dataOut.data_SNR
1231 buffer = dataOut.data_SNR
1233 if plot == 'snr':
1232 if plot == 'snr':
1234 buffer = 10*numpy.log10(dataOut.data_SNR)
1233 buffer = 10*numpy.log10(dataOut.data_SNR)
1235 if plot == 'dop':
1234 if plot == 'dop':
1236 buffer = 10*numpy.log10(dataOut.data_DOP)
1235 buffer = 10*numpy.log10(dataOut.data_DOP)
1237 if plot == 'mean':
1236 if plot == 'mean':
1238 buffer = dataOut.data_MEAN
1237 buffer = dataOut.data_MEAN
1239 if plot == 'std':
1238 if plot == 'std':
1240 buffer = dataOut.data_STD
1239 buffer = dataOut.data_STD
1241 if plot == 'coh':
1240 if plot == 'coh':
1242 buffer = dataOut.getCoherence()
1241 buffer = dataOut.getCoherence()
1243 if plot == 'phase':
1242 if plot == 'phase':
1244 buffer = dataOut.getCoherence(phase=True)
1243 buffer = dataOut.getCoherence(phase=True)
1245 if plot == 'output':
1244 if plot == 'output':
1246 buffer = dataOut.data_output
1245 buffer = dataOut.data_output
1247 if plot == 'param':
1246 if plot == 'param':
1248 buffer = dataOut.data_param
1247 buffer = dataOut.data_param
1249 if plot == 'scope':
1248 if plot == 'scope':
1250 buffer = dataOut.data
1249 buffer = dataOut.data
1251 self.flagDataAsBlock = dataOut.flagDataAsBlock
1250 self.flagDataAsBlock = dataOut.flagDataAsBlock
1252 self.nProfiles = dataOut.nProfiles
1251 self.nProfiles = dataOut.nProfiles
1253
1252
1254 if plot == 'spc':
1253 if plot == 'spc':
1255 self.data['spc'] = buffer
1254 self.data['spc'] = buffer
1256 elif plot == 'cspc':
1255 elif plot == 'cspc':
1257 self.data['spc'] = buffer[0]
1256 self.data['spc'] = buffer[0]
1258 self.data['cspc'] = buffer[1]
1257 self.data['cspc'] = buffer[1]
1259 elif plot == 'spc_moments':
1258 elif plot == 'spc_moments':
1260 self.data['spc'] = buffer
1259 self.data['spc'] = buffer
1261 self.data['moments'][tm] = dataOut.moments
1260 self.data['moments'][tm] = dataOut.moments
1262 else:
1261 else:
1263 if self.buffering:
1262 if self.buffering:
1264 self.data[plot][tm] = buffer
1263 self.data[plot][tm] = buffer
1265 else:
1264 else:
1266 self.data[plot] = buffer
1265 self.data[plot] = buffer
1267
1266
1268 def normalize_heights(self):
1267 def normalize_heights(self):
1269 '''
1268 '''
1270 Ensure same-dimension of the data for different heighList
1269 Ensure same-dimension of the data for different heighList
1271 '''
1270 '''
1272
1271
1273 H = numpy.array(list(self.__all_heights))
1272 H = numpy.array(list(self.__all_heights))
1274 H.sort()
1273 H.sort()
1275 for key in self.data:
1274 for key in self.data:
1276 shape = self.shape(key)[:-1] + H.shape
1275 shape = self.shape(key)[:-1] + H.shape
1277 for tm, obj in list(self.data[key].items()):
1276 for tm, obj in list(self.data[key].items()):
1278 h = self.__heights[self.__times.index(tm)]
1277 h = self.__heights[self.__times.index(tm)]
1279 if H.size == h.size:
1278 if H.size == h.size:
1280 continue
1279 continue
1281 index = numpy.where(numpy.in1d(H, h))[0]
1280 index = numpy.where(numpy.in1d(H, h))[0]
1282 dummy = numpy.zeros(shape) + numpy.nan
1281 dummy = numpy.zeros(shape) + numpy.nan
1283 if len(shape) == 2:
1282 if len(shape) == 2:
1284 dummy[:, index] = obj
1283 dummy[:, index] = obj
1285 else:
1284 else:
1286 dummy[index] = obj
1285 dummy[index] = obj
1287 self.data[key][tm] = dummy
1286 self.data[key][tm] = dummy
1288
1287
1289 self.__heights = [H for tm in self.__times]
1288 self.__heights = [H for tm in self.__times]
1290
1289
1291 def jsonify(self, decimate=False):
1290 def jsonify(self, plot_name, plot_type, decimate=False):
1292 '''
1291 '''
1293 Convert data to json
1292 Convert data to json
1294 '''
1293 '''
1295
1294
1296 data = {}
1297 tm = self.times[-1]
1295 tm = self.times[-1]
1298 dy = int(self.heights.size/self.MAXNUMY) + 1
1296 dy = int(self.heights.size/self.MAXNUMY) + 1
1299 for key in self.data:
1297 if self.key in ('spc', 'cspc') or not self.buffering:
1300 if key in ('spc', 'cspc') or not self.buffering:
1298 dx = int(self.data[self.key].shape[1]/self.MAXNUMX) + 1
1301 dx = int(self.data[key].shape[1]/self.MAXNUMX) + 1
1299 data = self.roundFloats(
1302 data[key] = self.roundFloats(
1300 self.data[self.key][::, ::dx, ::dy].tolist())
1303 self.data[key][::, ::dx, ::dy].tolist())
1304 else:
1301 else:
1305 data[key] = self.roundFloats(self.data[key][tm].tolist())
1302 data = self.roundFloats(self.data[self.key][tm].tolist())
1306
1303 if self.key is 'noise':
1307 ret = {'data': data}
1304 data = [[x] for x in data]
1308 ret['exp_code'] = self.exp_code
1305
1309 ret['time'] = float(tm)
1306 meta = {}
1310 ret['interval'] = float(self.interval)
1307 ret = {
1311 ret['localtime'] = self.localtime
1308 'plot': plot_name,
1312 ret['yrange'] = self.roundFloats(self.heights[::dy].tolist())
1309 'code': self.exp_code,
1310 'time': float(tm),
1311 'data': data,
1312 }
1313 meta['type'] = plot_type
1314 meta['interval'] = float(self.interval)
1315 meta['localtime'] = self.localtime
1316 meta['yrange'] = self.roundFloats(self.heights[::dy].tolist())
1313 if 'spc' in self.data or 'cspc' in self.data:
1317 if 'spc' in self.data or 'cspc' in self.data:
1314 ret['xrange'] = self.roundFloats(self.xrange[2][::dx].tolist())
1318 meta['xrange'] = self.roundFloats(self.xrange[2][::dx].tolist())
1315 else:
1316 ret['xrange'] = []
1317 if hasattr(self, 'pairs'):
1318 ret['pairs'] = [(int(p[0]), int(p[1])) for p in self.pairs]
1319 else:
1319 else:
1320 ret['pairs'] = []
1320 meta['xrange'] = []
1321
1322 for key, value in list(self.meta.items()):
1323 ret[key] = value
1324
1321
1322 meta.update(self.meta)
1323 ret['metadata'] = meta
1325 return json.dumps(ret)
1324 return json.dumps(ret)
1326
1325
1327 @property
1326 @property
1328 def times(self):
1327 def times(self):
1329 '''
1328 '''
1330 Return the list of times of the current data
1329 Return the list of times of the current data
1331 '''
1330 '''
1332
1331
1333 ret = numpy.array(self.__times)
1332 ret = numpy.array(self.__times)
1334 ret.sort()
1333 ret.sort()
1335 return ret
1334 return ret
1336
1335
1337 @property
1336 @property
1338 def min_time(self):
1337 def min_time(self):
1339 '''
1338 '''
1340 Return the minimun time value
1339 Return the minimun time value
1341 '''
1340 '''
1342
1341
1343 return self.times[0]
1342 return self.times[0]
1344
1343
1345 @property
1344 @property
1346 def max_time(self):
1345 def max_time(self):
1347 '''
1346 '''
1348 Return the maximun time value
1347 Return the maximun time value
1349 '''
1348 '''
1350
1349
1351 return self.times[-1]
1350 return self.times[-1]
1352
1351
1353 @property
1352 @property
1354 def heights(self):
1353 def heights(self):
1355 '''
1354 '''
1356 Return the list of heights of the current data
1355 Return the list of heights of the current data
1357 '''
1356 '''
1358
1357
1359 return numpy.array(self.__heights[-1])
1358 return numpy.array(self.__heights[-1])
1360
1359
1361 @staticmethod
1360 @staticmethod
1362 def roundFloats(obj):
1361 def roundFloats(obj):
1363 if isinstance(obj, list):
1362 if isinstance(obj, list):
1364 return list(map(PlotterData.roundFloats, obj))
1363 return list(map(PlotterData.roundFloats, obj))
1365 elif isinstance(obj, float):
1364 elif isinstance(obj, float):
1366 return round(obj, 2)
1365 return round(obj, 2)
@@ -1,806 +1,808
1
1
2 import os
2 import os
3 import sys
3 import sys
4 import zmq
4 import zmq
5 import time
5 import time
6 import numpy
6 import numpy
7 import datetime
7 import datetime
8 from functools import wraps
8 from functools import wraps
9 from threading import Thread
9 from threading import Thread
10 import matplotlib
10 import matplotlib
11
11
12 if 'BACKEND' in os.environ:
12 if 'BACKEND' in os.environ:
13 matplotlib.use(os.environ['BACKEND'])
13 matplotlib.use(os.environ['BACKEND'])
14 elif 'linux' in sys.platform:
14 elif 'linux' in sys.platform:
15 matplotlib.use("TkAgg")
15 matplotlib.use("TkAgg")
16 elif 'darwin' in sys.platform:
16 elif 'darwin' in sys.platform:
17 matplotlib.use('WxAgg')
17 matplotlib.use('WxAgg')
18 else:
18 else:
19 from schainpy.utils import log
19 from schainpy.utils import log
20 log.warning('Using default Backend="Agg"', 'INFO')
20 log.warning('Using default Backend="Agg"', 'INFO')
21 matplotlib.use('Agg')
21 matplotlib.use('Agg')
22
22
23 import matplotlib.pyplot as plt
23 import matplotlib.pyplot as plt
24 from matplotlib.patches import Polygon
24 from matplotlib.patches import Polygon
25 from mpl_toolkits.axes_grid1 import make_axes_locatable
25 from mpl_toolkits.axes_grid1 import make_axes_locatable
26 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
26 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
27
27
28 from schainpy.model.data.jrodata import PlotterData
28 from schainpy.model.data.jrodata import PlotterData
29 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
29 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
30 from schainpy.utils import log
30 from schainpy.utils import log
31
31
32 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
32 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
33 blu_values = matplotlib.pyplot.get_cmap(
33 blu_values = matplotlib.pyplot.get_cmap(
34 'seismic_r', 20)(numpy.arange(20))[10:15]
34 'seismic_r', 20)(numpy.arange(20))[10:15]
35 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
35 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
36 'jro', numpy.vstack((blu_values, jet_values)))
36 'jro', numpy.vstack((blu_values, jet_values)))
37 matplotlib.pyplot.register_cmap(cmap=ncmap)
37 matplotlib.pyplot.register_cmap(cmap=ncmap)
38
38
39 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis',
39 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis',
40 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm')]
40 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm')]
41
41
42 EARTH_RADIUS = 6.3710e3
42 EARTH_RADIUS = 6.3710e3
43
43
44 def ll2xy(lat1, lon1, lat2, lon2):
44 def ll2xy(lat1, lon1, lat2, lon2):
45
45
46 p = 0.017453292519943295
46 p = 0.017453292519943295
47 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
47 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
48 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
48 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
49 r = 12742 * numpy.arcsin(numpy.sqrt(a))
49 r = 12742 * numpy.arcsin(numpy.sqrt(a))
50 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
50 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
51 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
51 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
52 theta = -theta + numpy.pi/2
52 theta = -theta + numpy.pi/2
53 return r*numpy.cos(theta), r*numpy.sin(theta)
53 return r*numpy.cos(theta), r*numpy.sin(theta)
54
54
55
55
56 def km2deg(km):
56 def km2deg(km):
57 '''
57 '''
58 Convert distance in km to degrees
58 Convert distance in km to degrees
59 '''
59 '''
60
60
61 return numpy.rad2deg(km/EARTH_RADIUS)
61 return numpy.rad2deg(km/EARTH_RADIUS)
62
62
63
63
64 def figpause(interval):
64 def figpause(interval):
65 backend = plt.rcParams['backend']
65 backend = plt.rcParams['backend']
66 if backend in matplotlib.rcsetup.interactive_bk:
66 if backend in matplotlib.rcsetup.interactive_bk:
67 figManager = matplotlib._pylab_helpers.Gcf.get_active()
67 figManager = matplotlib._pylab_helpers.Gcf.get_active()
68 if figManager is not None:
68 if figManager is not None:
69 canvas = figManager.canvas
69 canvas = figManager.canvas
70 if canvas.figure.stale:
70 if canvas.figure.stale:
71 canvas.draw()
71 canvas.draw()
72 try:
72 try:
73 canvas.start_event_loop(interval)
73 canvas.start_event_loop(interval)
74 except:
74 except:
75 pass
75 pass
76 return
76 return
77
77
78
78
79 def popup(message):
79 def popup(message):
80 '''
80 '''
81 '''
81 '''
82
82
83 fig = plt.figure(figsize=(12, 8), facecolor='r')
83 fig = plt.figure(figsize=(12, 8), facecolor='r')
84 text = '\n'.join([s.strip() for s in message.split(':')])
84 text = '\n'.join([s.strip() for s in message.split(':')])
85 fig.text(0.01, 0.5, text, ha='left', va='center',
85 fig.text(0.01, 0.5, text, ha='left', va='center',
86 size='20', weight='heavy', color='w')
86 size='20', weight='heavy', color='w')
87 fig.show()
87 fig.show()
88 figpause(1000)
88 figpause(1000)
89
89
90
90
91 class Throttle(object):
91 class Throttle(object):
92 '''
92 '''
93 Decorator that prevents a function from being called more than once every
93 Decorator that prevents a function from being called more than once every
94 time period.
94 time period.
95 To create a function that cannot be called more than once a minute, but
95 To create a function that cannot be called more than once a minute, but
96 will sleep until it can be called:
96 will sleep until it can be called:
97 @Throttle(minutes=1)
97 @Throttle(minutes=1)
98 def foo():
98 def foo():
99 pass
99 pass
100
100
101 for i in range(10):
101 for i in range(10):
102 foo()
102 foo()
103 print "This function has run %s times." % i
103 print "This function has run %s times." % i
104 '''
104 '''
105
105
106 def __init__(self, seconds=0, minutes=0, hours=0):
106 def __init__(self, seconds=0, minutes=0, hours=0):
107 self.throttle_period = datetime.timedelta(
107 self.throttle_period = datetime.timedelta(
108 seconds=seconds, minutes=minutes, hours=hours
108 seconds=seconds, minutes=minutes, hours=hours
109 )
109 )
110
110
111 self.time_of_last_call = datetime.datetime.min
111 self.time_of_last_call = datetime.datetime.min
112
112
113 def __call__(self, fn):
113 def __call__(self, fn):
114 @wraps(fn)
114 @wraps(fn)
115 def wrapper(*args, **kwargs):
115 def wrapper(*args, **kwargs):
116 coerce = kwargs.pop('coerce', None)
116 coerce = kwargs.pop('coerce', None)
117 if coerce:
117 if coerce:
118 self.time_of_last_call = datetime.datetime.now()
118 self.time_of_last_call = datetime.datetime.now()
119 return fn(*args, **kwargs)
119 return fn(*args, **kwargs)
120 else:
120 else:
121 now = datetime.datetime.now()
121 now = datetime.datetime.now()
122 time_since_last_call = now - self.time_of_last_call
122 time_since_last_call = now - self.time_of_last_call
123 time_left = self.throttle_period - time_since_last_call
123 time_left = self.throttle_period - time_since_last_call
124
124
125 if time_left > datetime.timedelta(seconds=0):
125 if time_left > datetime.timedelta(seconds=0):
126 return
126 return
127
127
128 self.time_of_last_call = datetime.datetime.now()
128 self.time_of_last_call = datetime.datetime.now()
129 return fn(*args, **kwargs)
129 return fn(*args, **kwargs)
130
130
131 return wrapper
131 return wrapper
132
132
133 def apply_throttle(value):
133 def apply_throttle(value):
134
134
135 @Throttle(seconds=value)
135 @Throttle(seconds=value)
136 def fnThrottled(fn):
136 def fnThrottled(fn):
137 fn()
137 fn()
138
138
139 return fnThrottled
139 return fnThrottled
140
140
141
141
142 @MPDecorator
142 @MPDecorator
143 class Plot(Operation):
143 class Plot(Operation):
144 '''
144 '''
145 Base class for Schain plotting operations
145 Base class for Schain plotting operations
146 '''
146 '''
147
147
148 CODE = 'Figure'
148 CODE = 'Figure'
149 colormap = 'jet'
149 colormap = 'jet'
150 bgcolor = 'white'
150 bgcolor = 'white'
151 __missing = 1E30
151 __missing = 1E30
152
152
153 __attrs__ = ['show', 'save', 'xmin', 'xmax', 'ymin', 'ymax', 'zmin', 'zmax',
153 __attrs__ = ['show', 'save', 'xmin', 'xmax', 'ymin', 'ymax', 'zmin', 'zmax',
154 'zlimits', 'xlabel', 'ylabel', 'xaxis', 'cb_label', 'title',
154 'zlimits', 'xlabel', 'ylabel', 'xaxis', 'cb_label', 'title',
155 'colorbar', 'bgcolor', 'width', 'height', 'localtime', 'oneFigure',
155 'colorbar', 'bgcolor', 'width', 'height', 'localtime', 'oneFigure',
156 'showprofile', 'decimation', 'pause']
156 'showprofile', 'decimation', 'pause']
157
157
158 def __init__(self):
158 def __init__(self):
159
159
160 Operation.__init__(self)
160 Operation.__init__(self)
161 self.isConfig = False
161 self.isConfig = False
162 self.isPlotConfig = False
162 self.isPlotConfig = False
163 self.save_counter = 1
163 self.save_counter = 1
164 self.sender_counter = 1
164 self.sender_counter = 1
165 self.data = None
165 self.data = None
166
166
167 def __fmtTime(self, x, pos):
167 def __fmtTime(self, x, pos):
168 '''
168 '''
169 '''
169 '''
170
170
171 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
171 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
172
172
173 def __setup(self, **kwargs):
173 def __setup(self, **kwargs):
174 '''
174 '''
175 Initialize variables
175 Initialize variables
176 '''
176 '''
177
177
178 self.figures = []
178 self.figures = []
179 self.axes = []
179 self.axes = []
180 self.cb_axes = []
180 self.cb_axes = []
181 self.localtime = kwargs.pop('localtime', True)
181 self.localtime = kwargs.pop('localtime', True)
182 self.show = kwargs.get('show', True)
182 self.show = kwargs.get('show', True)
183 self.save = kwargs.get('save', False)
183 self.save = kwargs.get('save', False)
184 self.save_period = kwargs.get('save_period', 2)
184 self.save_period = kwargs.get('save_period', 2)
185 self.ftp = kwargs.get('ftp', False)
185 self.ftp = kwargs.get('ftp', False)
186 self.colormap = kwargs.get('colormap', self.colormap)
186 self.colormap = kwargs.get('colormap', self.colormap)
187 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
187 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
188 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
188 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
189 self.colormaps = kwargs.get('colormaps', None)
189 self.colormaps = kwargs.get('colormaps', None)
190 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
190 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
191 self.showprofile = kwargs.get('showprofile', False)
191 self.showprofile = kwargs.get('showprofile', False)
192 self.title = kwargs.get('wintitle', self.CODE.upper())
192 self.title = kwargs.get('wintitle', self.CODE.upper())
193 self.cb_label = kwargs.get('cb_label', None)
193 self.cb_label = kwargs.get('cb_label', None)
194 self.cb_labels = kwargs.get('cb_labels', None)
194 self.cb_labels = kwargs.get('cb_labels', None)
195 self.labels = kwargs.get('labels', None)
195 self.labels = kwargs.get('labels', None)
196 self.xaxis = kwargs.get('xaxis', 'frequency')
196 self.xaxis = kwargs.get('xaxis', 'frequency')
197 self.zmin = kwargs.get('zmin', None)
197 self.zmin = kwargs.get('zmin', None)
198 self.zmax = kwargs.get('zmax', None)
198 self.zmax = kwargs.get('zmax', None)
199 self.zlimits = kwargs.get('zlimits', None)
199 self.zlimits = kwargs.get('zlimits', None)
200 self.xmin = kwargs.get('xmin', None)
200 self.xmin = kwargs.get('xmin', None)
201 self.xmax = kwargs.get('xmax', None)
201 self.xmax = kwargs.get('xmax', None)
202 self.xrange = kwargs.get('xrange', 24)
202 self.xrange = kwargs.get('xrange', 24)
203 self.xscale = kwargs.get('xscale', None)
203 self.xscale = kwargs.get('xscale', None)
204 self.ymin = kwargs.get('ymin', None)
204 self.ymin = kwargs.get('ymin', None)
205 self.ymax = kwargs.get('ymax', None)
205 self.ymax = kwargs.get('ymax', None)
206 self.yscale = kwargs.get('yscale', None)
206 self.yscale = kwargs.get('yscale', None)
207 self.xlabel = kwargs.get('xlabel', None)
207 self.xlabel = kwargs.get('xlabel', None)
208 self.decimation = kwargs.get('decimation', None)
208 self.decimation = kwargs.get('decimation', None)
209 self.showSNR = kwargs.get('showSNR', False)
209 self.showSNR = kwargs.get('showSNR', False)
210 self.oneFigure = kwargs.get('oneFigure', True)
210 self.oneFigure = kwargs.get('oneFigure', True)
211 self.width = kwargs.get('width', None)
211 self.width = kwargs.get('width', None)
212 self.height = kwargs.get('height', None)
212 self.height = kwargs.get('height', None)
213 self.colorbar = kwargs.get('colorbar', True)
213 self.colorbar = kwargs.get('colorbar', True)
214 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
214 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
215 self.channels = kwargs.get('channels', None)
215 self.channels = kwargs.get('channels', None)
216 self.titles = kwargs.get('titles', [])
216 self.titles = kwargs.get('titles', [])
217 self.polar = False
217 self.polar = False
218 self.type = kwargs.get('type', 'iq')
218 self.type = kwargs.get('type', 'iq')
219 self.grid = kwargs.get('grid', False)
219 self.grid = kwargs.get('grid', False)
220 self.pause = kwargs.get('pause', False)
220 self.pause = kwargs.get('pause', False)
221 self.save_labels = kwargs.get('save_labels', None)
221 self.save_labels = kwargs.get('save_labels', None)
222 self.realtime = kwargs.get('realtime', True)
222 self.realtime = kwargs.get('realtime', True)
223 self.buffering = kwargs.get('buffering', True)
223 self.buffering = kwargs.get('buffering', True)
224 self.throttle = kwargs.get('throttle', 2)
224 self.throttle = kwargs.get('throttle', 2)
225 self.exp_code = kwargs.get('exp_code', None)
225 self.exp_code = kwargs.get('exp_code', None)
226 self.plot_server = kwargs.get('plot_server', False)
226 self.plot_server = kwargs.get('plot_server', False)
227 self.sender_period = kwargs.get('sender_period', 2)
227 self.sender_period = kwargs.get('sender_period', 2)
228 self.__throttle_plot = apply_throttle(self.throttle)
228 self.__throttle_plot = apply_throttle(self.throttle)
229 self.data = PlotterData(
229 self.data = PlotterData(
230 self.CODE, self.throttle, self.exp_code, self.buffering, snr=self.showSNR)
230 self.CODE, self.throttle, self.exp_code, self.buffering, snr=self.showSNR)
231
231
232 if self.plot_server:
232 if self.plot_server:
233 if not self.plot_server.startswith('tcp://'):
233 if not self.plot_server.startswith('tcp://'):
234 self.plot_server = 'tcp://{}'.format(self.plot_server)
234 self.plot_server = 'tcp://{}'.format(self.plot_server)
235 log.success(
235 log.success(
236 'Sending to server: {}'.format(self.plot_server),
236 'Sending to server: {}'.format(self.plot_server),
237 self.name
237 self.name
238 )
238 )
239 if 'plot_name' in kwargs:
240 self.plot_name = kwargs['plot_name']
239
241
240 def __setup_plot(self):
242 def __setup_plot(self):
241 '''
243 '''
242 Common setup for all figures, here figures and axes are created
244 Common setup for all figures, here figures and axes are created
243 '''
245 '''
244
246
245 self.setup()
247 self.setup()
246
248
247 self.time_label = 'LT' if self.localtime else 'UTC'
249 self.time_label = 'LT' if self.localtime else 'UTC'
248
250
249 if self.width is None:
251 if self.width is None:
250 self.width = 8
252 self.width = 8
251
253
252 self.figures = []
254 self.figures = []
253 self.axes = []
255 self.axes = []
254 self.cb_axes = []
256 self.cb_axes = []
255 self.pf_axes = []
257 self.pf_axes = []
256 self.cmaps = []
258 self.cmaps = []
257
259
258 size = '15%' if self.ncols == 1 else '30%'
260 size = '15%' if self.ncols == 1 else '30%'
259 pad = '4%' if self.ncols == 1 else '8%'
261 pad = '4%' if self.ncols == 1 else '8%'
260
262
261 if self.oneFigure:
263 if self.oneFigure:
262 if self.height is None:
264 if self.height is None:
263 self.height = 1.4 * self.nrows + 1
265 self.height = 1.4 * self.nrows + 1
264 fig = plt.figure(figsize=(self.width, self.height),
266 fig = plt.figure(figsize=(self.width, self.height),
265 edgecolor='k',
267 edgecolor='k',
266 facecolor='w')
268 facecolor='w')
267 self.figures.append(fig)
269 self.figures.append(fig)
268 for n in range(self.nplots):
270 for n in range(self.nplots):
269 ax = fig.add_subplot(self.nrows, self.ncols,
271 ax = fig.add_subplot(self.nrows, self.ncols,
270 n + 1, polar=self.polar)
272 n + 1, polar=self.polar)
271 ax.tick_params(labelsize=8)
273 ax.tick_params(labelsize=8)
272 ax.firsttime = True
274 ax.firsttime = True
273 ax.index = 0
275 ax.index = 0
274 ax.press = None
276 ax.press = None
275 self.axes.append(ax)
277 self.axes.append(ax)
276 if self.showprofile:
278 if self.showprofile:
277 cax = self.__add_axes(ax, size=size, pad=pad)
279 cax = self.__add_axes(ax, size=size, pad=pad)
278 cax.tick_params(labelsize=8)
280 cax.tick_params(labelsize=8)
279 self.pf_axes.append(cax)
281 self.pf_axes.append(cax)
280 else:
282 else:
281 if self.height is None:
283 if self.height is None:
282 self.height = 3
284 self.height = 3
283 for n in range(self.nplots):
285 for n in range(self.nplots):
284 fig = plt.figure(figsize=(self.width, self.height),
286 fig = plt.figure(figsize=(self.width, self.height),
285 edgecolor='k',
287 edgecolor='k',
286 facecolor='w')
288 facecolor='w')
287 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
289 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
288 ax.tick_params(labelsize=8)
290 ax.tick_params(labelsize=8)
289 ax.firsttime = True
291 ax.firsttime = True
290 ax.index = 0
292 ax.index = 0
291 ax.press = None
293 ax.press = None
292 self.figures.append(fig)
294 self.figures.append(fig)
293 self.axes.append(ax)
295 self.axes.append(ax)
294 if self.showprofile:
296 if self.showprofile:
295 cax = self.__add_axes(ax, size=size, pad=pad)
297 cax = self.__add_axes(ax, size=size, pad=pad)
296 cax.tick_params(labelsize=8)
298 cax.tick_params(labelsize=8)
297 self.pf_axes.append(cax)
299 self.pf_axes.append(cax)
298
300
299 for n in range(self.nrows):
301 for n in range(self.nrows):
300 if self.colormaps is not None:
302 if self.colormaps is not None:
301 cmap = plt.get_cmap(self.colormaps[n])
303 cmap = plt.get_cmap(self.colormaps[n])
302 else:
304 else:
303 cmap = plt.get_cmap(self.colormap)
305 cmap = plt.get_cmap(self.colormap)
304 cmap.set_bad(self.bgcolor, 1.)
306 cmap.set_bad(self.bgcolor, 1.)
305 self.cmaps.append(cmap)
307 self.cmaps.append(cmap)
306
308
307 for fig in self.figures:
309 for fig in self.figures:
308 fig.canvas.mpl_connect('key_press_event', self.OnKeyPress)
310 fig.canvas.mpl_connect('key_press_event', self.OnKeyPress)
309 fig.canvas.mpl_connect('scroll_event', self.OnBtnScroll)
311 fig.canvas.mpl_connect('scroll_event', self.OnBtnScroll)
310 fig.canvas.mpl_connect('button_press_event', self.onBtnPress)
312 fig.canvas.mpl_connect('button_press_event', self.onBtnPress)
311 fig.canvas.mpl_connect('motion_notify_event', self.onMotion)
313 fig.canvas.mpl_connect('motion_notify_event', self.onMotion)
312 fig.canvas.mpl_connect('button_release_event', self.onBtnRelease)
314 fig.canvas.mpl_connect('button_release_event', self.onBtnRelease)
313
315
314 def OnKeyPress(self, event):
316 def OnKeyPress(self, event):
315 '''
317 '''
316 Event for pressing keys (up, down) change colormap
318 Event for pressing keys (up, down) change colormap
317 '''
319 '''
318 ax = event.inaxes
320 ax = event.inaxes
319 if ax in self.axes:
321 if ax in self.axes:
320 if event.key == 'down':
322 if event.key == 'down':
321 ax.index += 1
323 ax.index += 1
322 elif event.key == 'up':
324 elif event.key == 'up':
323 ax.index -= 1
325 ax.index -= 1
324 if ax.index < 0:
326 if ax.index < 0:
325 ax.index = len(CMAPS) - 1
327 ax.index = len(CMAPS) - 1
326 elif ax.index == len(CMAPS):
328 elif ax.index == len(CMAPS):
327 ax.index = 0
329 ax.index = 0
328 cmap = CMAPS[ax.index]
330 cmap = CMAPS[ax.index]
329 ax.cbar.set_cmap(cmap)
331 ax.cbar.set_cmap(cmap)
330 ax.cbar.draw_all()
332 ax.cbar.draw_all()
331 ax.plt.set_cmap(cmap)
333 ax.plt.set_cmap(cmap)
332 ax.cbar.patch.figure.canvas.draw()
334 ax.cbar.patch.figure.canvas.draw()
333 self.colormap = cmap.name
335 self.colormap = cmap.name
334
336
335 def OnBtnScroll(self, event):
337 def OnBtnScroll(self, event):
336 '''
338 '''
337 Event for scrolling, scale figure
339 Event for scrolling, scale figure
338 '''
340 '''
339 cb_ax = event.inaxes
341 cb_ax = event.inaxes
340 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
342 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
341 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
343 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
342 pt = ax.cbar.ax.bbox.get_points()[:, 1]
344 pt = ax.cbar.ax.bbox.get_points()[:, 1]
343 nrm = ax.cbar.norm
345 nrm = ax.cbar.norm
344 vmin, vmax, p0, p1, pS = (
346 vmin, vmax, p0, p1, pS = (
345 nrm.vmin, nrm.vmax, pt[0], pt[1], event.y)
347 nrm.vmin, nrm.vmax, pt[0], pt[1], event.y)
346 scale = 2 if event.step == 1 else 0.5
348 scale = 2 if event.step == 1 else 0.5
347 point = vmin + (vmax - vmin) / (p1 - p0) * (pS - p0)
349 point = vmin + (vmax - vmin) / (p1 - p0) * (pS - p0)
348 ax.cbar.norm.vmin = point - scale * (point - vmin)
350 ax.cbar.norm.vmin = point - scale * (point - vmin)
349 ax.cbar.norm.vmax = point - scale * (point - vmax)
351 ax.cbar.norm.vmax = point - scale * (point - vmax)
350 ax.plt.set_norm(ax.cbar.norm)
352 ax.plt.set_norm(ax.cbar.norm)
351 ax.cbar.draw_all()
353 ax.cbar.draw_all()
352 ax.cbar.patch.figure.canvas.draw()
354 ax.cbar.patch.figure.canvas.draw()
353
355
354 def onBtnPress(self, event):
356 def onBtnPress(self, event):
355 '''
357 '''
356 Event for mouse button press
358 Event for mouse button press
357 '''
359 '''
358 cb_ax = event.inaxes
360 cb_ax = event.inaxes
359 if cb_ax is None:
361 if cb_ax is None:
360 return
362 return
361
363
362 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
364 if cb_ax in [ax.cbar.ax for ax in self.axes if ax.cbar]:
363 cb_ax.press = event.x, event.y
365 cb_ax.press = event.x, event.y
364 else:
366 else:
365 cb_ax.press = None
367 cb_ax.press = None
366
368
367 def onMotion(self, event):
369 def onMotion(self, event):
368 '''
370 '''
369 Event for move inside colorbar
371 Event for move inside colorbar
370 '''
372 '''
371 cb_ax = event.inaxes
373 cb_ax = event.inaxes
372 if cb_ax is None:
374 if cb_ax is None:
373 return
375 return
374 if cb_ax not in [ax.cbar.ax for ax in self.axes if ax.cbar]:
376 if cb_ax not in [ax.cbar.ax for ax in self.axes if ax.cbar]:
375 return
377 return
376 if cb_ax.press is None:
378 if cb_ax.press is None:
377 return
379 return
378
380
379 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
381 ax = [ax for ax in self.axes if cb_ax == ax.cbar.ax][0]
380 xprev, yprev = cb_ax.press
382 xprev, yprev = cb_ax.press
381 dx = event.x - xprev
383 dx = event.x - xprev
382 dy = event.y - yprev
384 dy = event.y - yprev
383 cb_ax.press = event.x, event.y
385 cb_ax.press = event.x, event.y
384 scale = ax.cbar.norm.vmax - ax.cbar.norm.vmin
386 scale = ax.cbar.norm.vmax - ax.cbar.norm.vmin
385 perc = 0.03
387 perc = 0.03
386
388
387 if event.button == 1:
389 if event.button == 1:
388 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
390 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
389 ax.cbar.norm.vmax -= (perc * scale) * numpy.sign(dy)
391 ax.cbar.norm.vmax -= (perc * scale) * numpy.sign(dy)
390 elif event.button == 3:
392 elif event.button == 3:
391 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
393 ax.cbar.norm.vmin -= (perc * scale) * numpy.sign(dy)
392 ax.cbar.norm.vmax += (perc * scale) * numpy.sign(dy)
394 ax.cbar.norm.vmax += (perc * scale) * numpy.sign(dy)
393
395
394 ax.cbar.draw_all()
396 ax.cbar.draw_all()
395 ax.plt.set_norm(ax.cbar.norm)
397 ax.plt.set_norm(ax.cbar.norm)
396 ax.cbar.patch.figure.canvas.draw()
398 ax.cbar.patch.figure.canvas.draw()
397
399
398 def onBtnRelease(self, event):
400 def onBtnRelease(self, event):
399 '''
401 '''
400 Event for mouse button release
402 Event for mouse button release
401 '''
403 '''
402 cb_ax = event.inaxes
404 cb_ax = event.inaxes
403 if cb_ax is not None:
405 if cb_ax is not None:
404 cb_ax.press = None
406 cb_ax.press = None
405
407
406 def __add_axes(self, ax, size='30%', pad='8%'):
408 def __add_axes(self, ax, size='30%', pad='8%'):
407 '''
409 '''
408 Add new axes to the given figure
410 Add new axes to the given figure
409 '''
411 '''
410 divider = make_axes_locatable(ax)
412 divider = make_axes_locatable(ax)
411 nax = divider.new_horizontal(size=size, pad=pad)
413 nax = divider.new_horizontal(size=size, pad=pad)
412 ax.figure.add_axes(nax)
414 ax.figure.add_axes(nax)
413 return nax
415 return nax
414
416
415 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
417 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
416 '''
418 '''
417 Create a masked array for missing data
419 Create a masked array for missing data
418 '''
420 '''
419 if x_buffer.shape[0] < 2:
421 if x_buffer.shape[0] < 2:
420 return x_buffer, y_buffer, z_buffer
422 return x_buffer, y_buffer, z_buffer
421
423
422 deltas = x_buffer[1:] - x_buffer[0:-1]
424 deltas = x_buffer[1:] - x_buffer[0:-1]
423 x_median = numpy.median(deltas)
425 x_median = numpy.median(deltas)
424
426
425 index = numpy.where(deltas > 5 * x_median)
427 index = numpy.where(deltas > 5 * x_median)
426
428
427 if len(index[0]) != 0:
429 if len(index[0]) != 0:
428 z_buffer[::, index[0], ::] = self.__missing
430 z_buffer[::, index[0], ::] = self.__missing
429 z_buffer = numpy.ma.masked_inside(z_buffer,
431 z_buffer = numpy.ma.masked_inside(z_buffer,
430 0.99 * self.__missing,
432 0.99 * self.__missing,
431 1.01 * self.__missing)
433 1.01 * self.__missing)
432
434
433 return x_buffer, y_buffer, z_buffer
435 return x_buffer, y_buffer, z_buffer
434
436
435 def decimate(self):
437 def decimate(self):
436
438
437 # dx = int(len(self.x)/self.__MAXNUMX) + 1
439 # dx = int(len(self.x)/self.__MAXNUMX) + 1
438 dy = int(len(self.y) / self.decimation) + 1
440 dy = int(len(self.y) / self.decimation) + 1
439
441
440 # x = self.x[::dx]
442 # x = self.x[::dx]
441 x = self.x
443 x = self.x
442 y = self.y[::dy]
444 y = self.y[::dy]
443 z = self.z[::, ::, ::dy]
445 z = self.z[::, ::, ::dy]
444
446
445 return x, y, z
447 return x, y, z
446
448
447 def format(self):
449 def format(self):
448 '''
450 '''
449 Set min and max values, labels, ticks and titles
451 Set min and max values, labels, ticks and titles
450 '''
452 '''
451
453
452 if self.xmin is None:
454 if self.xmin is None:
453 xmin = self.data.min_time
455 xmin = self.data.min_time
454 else:
456 else:
455 if self.xaxis is 'time':
457 if self.xaxis is 'time':
456 dt = self.getDateTime(self.data.min_time)
458 dt = self.getDateTime(self.data.min_time)
457 xmin = (dt.replace(hour=int(self.xmin), minute=0, second=0) -
459 xmin = (dt.replace(hour=int(self.xmin), minute=0, second=0) -
458 datetime.datetime(1970, 1, 1)).total_seconds()
460 datetime.datetime(1970, 1, 1)).total_seconds()
459 if self.data.localtime:
461 if self.data.localtime:
460 xmin += time.timezone
462 xmin += time.timezone
461 else:
463 else:
462 xmin = self.xmin
464 xmin = self.xmin
463
465
464 if self.xmax is None:
466 if self.xmax is None:
465 xmax = xmin + self.xrange * 60 * 60
467 xmax = xmin + self.xrange * 60 * 60
466 else:
468 else:
467 if self.xaxis is 'time':
469 if self.xaxis is 'time':
468 dt = self.getDateTime(self.data.max_time)
470 dt = self.getDateTime(self.data.max_time)
469 xmax = (dt.replace(hour=int(self.xmax), minute=59, second=59) -
471 xmax = (dt.replace(hour=int(self.xmax), minute=59, second=59) -
470 datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds()
472 datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds()
471 if self.data.localtime:
473 if self.data.localtime:
472 xmax += time.timezone
474 xmax += time.timezone
473 else:
475 else:
474 xmax = self.xmax
476 xmax = self.xmax
475
477
476 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
478 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
477 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
479 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
478 #Y = numpy.array([1, 2, 5, 10, 20, 50, 100, 200, 500, 1000, 2000, 5000, 10000, 20000, 50000])
480 #Y = numpy.array([1, 2, 5, 10, 20, 50, 100, 200, 500, 1000, 2000, 5000, 10000, 20000, 50000])
479
481
480 #i = 1 if numpy.where(
482 #i = 1 if numpy.where(
481 # abs(ymax-ymin) <= Y)[0][0] < 0 else numpy.where(abs(ymax-ymin) <= Y)[0][0]
483 # abs(ymax-ymin) <= Y)[0][0] < 0 else numpy.where(abs(ymax-ymin) <= Y)[0][0]
482 #ystep = Y[i] / 10.
484 #ystep = Y[i] / 10.
483 dig = int(numpy.log10(ymax))
485 dig = int(numpy.log10(ymax))
484 if dig == 0:
486 if dig == 0:
485 digD = len(str(ymax)) - 2
487 digD = len(str(ymax)) - 2
486 ydec = ymax*(10**digD)
488 ydec = ymax*(10**digD)
487
489
488 dig = int(numpy.log10(ydec))
490 dig = int(numpy.log10(ydec))
489 ystep = ((ydec + (10**(dig)))//10**(dig))*(10**(dig))
491 ystep = ((ydec + (10**(dig)))//10**(dig))*(10**(dig))
490 ystep = ystep/5
492 ystep = ystep/5
491 ystep = ystep/(10**digD)
493 ystep = ystep/(10**digD)
492
494
493 else:
495 else:
494 ystep = ((ymax + (10**(dig)))//10**(dig))*(10**(dig))
496 ystep = ((ymax + (10**(dig)))//10**(dig))*(10**(dig))
495 ystep = ystep/5
497 ystep = ystep/5
496
498
497 if self.xaxis is not 'time':
499 if self.xaxis is not 'time':
498
500
499 dig = int(numpy.log10(xmax))
501 dig = int(numpy.log10(xmax))
500
502
501 if dig <= 0:
503 if dig <= 0:
502 digD = len(str(xmax)) - 2
504 digD = len(str(xmax)) - 2
503 xdec = xmax*(10**digD)
505 xdec = xmax*(10**digD)
504
506
505 dig = int(numpy.log10(xdec))
507 dig = int(numpy.log10(xdec))
506 xstep = ((xdec + (10**(dig)))//10**(dig))*(10**(dig))
508 xstep = ((xdec + (10**(dig)))//10**(dig))*(10**(dig))
507 xstep = xstep*0.5
509 xstep = xstep*0.5
508 xstep = xstep/(10**digD)
510 xstep = xstep/(10**digD)
509
511
510 else:
512 else:
511 xstep = ((xmax + (10**(dig)))//10**(dig))*(10**(dig))
513 xstep = ((xmax + (10**(dig)))//10**(dig))*(10**(dig))
512 xstep = xstep/5
514 xstep = xstep/5
513
515
514 for n, ax in enumerate(self.axes):
516 for n, ax in enumerate(self.axes):
515 if ax.firsttime:
517 if ax.firsttime:
516 ax.set_facecolor(self.bgcolor)
518 ax.set_facecolor(self.bgcolor)
517 ax.yaxis.set_major_locator(MultipleLocator(ystep))
519 ax.yaxis.set_major_locator(MultipleLocator(ystep))
518 if self.xscale:
520 if self.xscale:
519 ax.xaxis.set_major_formatter(FuncFormatter(
521 ax.xaxis.set_major_formatter(FuncFormatter(
520 lambda x, pos: '{0:g}'.format(x*self.xscale)))
522 lambda x, pos: '{0:g}'.format(x*self.xscale)))
521 if self.xscale:
523 if self.xscale:
522 ax.yaxis.set_major_formatter(FuncFormatter(
524 ax.yaxis.set_major_formatter(FuncFormatter(
523 lambda x, pos: '{0:g}'.format(x*self.yscale)))
525 lambda x, pos: '{0:g}'.format(x*self.yscale)))
524 if self.xaxis is 'time':
526 if self.xaxis is 'time':
525 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
527 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
526 ax.xaxis.set_major_locator(LinearLocator(9))
528 ax.xaxis.set_major_locator(LinearLocator(9))
527 else:
529 else:
528 ax.xaxis.set_major_locator(MultipleLocator(xstep))
530 ax.xaxis.set_major_locator(MultipleLocator(xstep))
529 if self.xlabel is not None:
531 if self.xlabel is not None:
530 ax.set_xlabel(self.xlabel)
532 ax.set_xlabel(self.xlabel)
531 ax.set_ylabel(self.ylabel)
533 ax.set_ylabel(self.ylabel)
532 ax.firsttime = False
534 ax.firsttime = False
533 if self.showprofile:
535 if self.showprofile:
534 self.pf_axes[n].set_ylim(ymin, ymax)
536 self.pf_axes[n].set_ylim(ymin, ymax)
535 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
537 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
536 self.pf_axes[n].set_xlabel('dB')
538 self.pf_axes[n].set_xlabel('dB')
537 self.pf_axes[n].grid(b=True, axis='x')
539 self.pf_axes[n].grid(b=True, axis='x')
538 [tick.set_visible(False)
540 [tick.set_visible(False)
539 for tick in self.pf_axes[n].get_yticklabels()]
541 for tick in self.pf_axes[n].get_yticklabels()]
540 if self.colorbar:
542 if self.colorbar:
541 ax.cbar = plt.colorbar(
543 ax.cbar = plt.colorbar(
542 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
544 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
543 ax.cbar.ax.tick_params(labelsize=8)
545 ax.cbar.ax.tick_params(labelsize=8)
544 ax.cbar.ax.press = None
546 ax.cbar.ax.press = None
545 if self.cb_label:
547 if self.cb_label:
546 ax.cbar.set_label(self.cb_label, size=8)
548 ax.cbar.set_label(self.cb_label, size=8)
547 elif self.cb_labels:
549 elif self.cb_labels:
548 ax.cbar.set_label(self.cb_labels[n], size=8)
550 ax.cbar.set_label(self.cb_labels[n], size=8)
549 else:
551 else:
550 ax.cbar = None
552 ax.cbar = None
551 if self.grid:
553 if self.grid:
552 ax.grid(True)
554 ax.grid(True)
553
555
554 if not self.polar:
556 if not self.polar:
555 ax.set_xlim(xmin, xmax)
557 ax.set_xlim(xmin, xmax)
556 ax.set_ylim(ymin, ymax)
558 ax.set_ylim(ymin, ymax)
557 ax.set_title('{} {} {}'.format(
559 ax.set_title('{} {} {}'.format(
558 self.titles[n],
560 self.titles[n],
559 self.getDateTime(self.data.max_time).strftime(
561 self.getDateTime(self.data.max_time).strftime(
560 '%Y-%m-%d %H:%M:%S'),
562 '%Y-%m-%d %H:%M:%S'),
561 self.time_label),
563 self.time_label),
562 size=8)
564 size=8)
563 else:
565 else:
564 ax.set_title('{}'.format(self.titles[n]), size=8)
566 ax.set_title('{}'.format(self.titles[n]), size=8)
565 ax.set_ylim(0, 90)
567 ax.set_ylim(0, 90)
566 ax.set_yticks(numpy.arange(0, 90, 20))
568 ax.set_yticks(numpy.arange(0, 90, 20))
567 ax.yaxis.labelpad = 40
569 ax.yaxis.labelpad = 40
568
570
569 def clear_figures(self):
571 def clear_figures(self):
570 '''
572 '''
571 Reset axes for redraw plots
573 Reset axes for redraw plots
572 '''
574 '''
573
575
574 for ax in self.axes:
576 for ax in self.axes:
575 ax.clear()
577 ax.clear()
576 ax.firsttime = True
578 ax.firsttime = True
577 if ax.cbar:
579 if ax.cbar:
578 ax.cbar.remove()
580 ax.cbar.remove()
579
581
580 def __plot(self):
582 def __plot(self):
581 '''
583 '''
582 Main function to plot, format and save figures
584 Main function to plot, format and save figures
583 '''
585 '''
584
586
585 try:
587 try:
586 self.plot()
588 self.plot()
587 self.format()
589 self.format()
588 except Exception as e:
590 except Exception as e:
589 log.warning('{} Plot could not be updated... check data'.format(
591 log.warning('{} Plot could not be updated... check data'.format(
590 self.CODE), self.name)
592 self.CODE), self.name)
591 log.error(str(e), '')
593 log.error(str(e), '')
592 return
594 return
593
595
594 for n, fig in enumerate(self.figures):
596 for n, fig in enumerate(self.figures):
595 if self.nrows == 0 or self.nplots == 0:
597 if self.nrows == 0 or self.nplots == 0:
596 log.warning('No data', self.name)
598 log.warning('No data', self.name)
597 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
599 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
598 fig.canvas.manager.set_window_title(self.CODE)
600 fig.canvas.manager.set_window_title(self.CODE)
599 continue
601 continue
600
602
601 fig.tight_layout()
603 fig.tight_layout()
602 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
604 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
603 self.getDateTime(self.data.max_time).strftime('%Y/%m/%d')))
605 self.getDateTime(self.data.max_time).strftime('%Y/%m/%d')))
604 fig.canvas.draw()
606 fig.canvas.draw()
605 if self.show:
607 if self.show:
606 fig.show()
608 fig.show()
607 figpause(0.1)
609 figpause(0.1)
608
610
609 if self.save:
611 if self.save:
610 self.save_figure(n)
612 self.save_figure(n)
611
613
612 if self.plot_server:
614 if self.plot_server:
613 self.send_to_server()
615 self.send_to_server()
614 # t = Thread(target=self.send_to_server)
616 # t = Thread(target=self.send_to_server)
615 # t.start()
617 # t.start()
616
618
617 def save_figure(self, n):
619 def save_figure(self, n):
618 '''
620 '''
619 '''
621 '''
620
622
621 if self.save_counter < self.save_period:
623 if self.save_counter < self.save_period:
622 self.save_counter += 1
624 self.save_counter += 1
623 return
625 return
624
626
625 self.save_counter = 1
627 self.save_counter = 1
626
628
627 fig = self.figures[n]
629 fig = self.figures[n]
628
630
629 if self.save_labels:
631 if self.save_labels:
630 labels = self.save_labels
632 labels = self.save_labels
631 else:
633 else:
632 labels = list(range(self.nrows))
634 labels = list(range(self.nrows))
633
635
634 if self.oneFigure:
636 if self.oneFigure:
635 label = ''
637 label = ''
636 else:
638 else:
637 label = '-{}'.format(labels[n])
639 label = '-{}'.format(labels[n])
638 figname = os.path.join(
640 figname = os.path.join(
639 self.save,
641 self.save,
640 self.CODE,
642 self.CODE,
641 '{}{}_{}.png'.format(
643 '{}{}_{}.png'.format(
642 self.CODE,
644 self.CODE,
643 label,
645 label,
644 self.getDateTime(self.data.max_time).strftime(
646 self.getDateTime(self.data.max_time).strftime(
645 '%Y%m%d_%H%M%S'
647 '%Y%m%d_%H%M%S'
646 ),
648 ),
647 )
649 )
648 )
650 )
649 log.log('Saving figure: {}'.format(figname), self.name)
651 log.log('Saving figure: {}'.format(figname), self.name)
650 if not os.path.isdir(os.path.dirname(figname)):
652 if not os.path.isdir(os.path.dirname(figname)):
651 os.makedirs(os.path.dirname(figname))
653 os.makedirs(os.path.dirname(figname))
652 fig.savefig(figname)
654 fig.savefig(figname)
653
655
654 if self.realtime:
656 if self.realtime:
655 figname = os.path.join(
657 figname = os.path.join(
656 self.save,
658 self.save,
657 '{}{}_{}.png'.format(
659 '{}{}_{}.png'.format(
658 self.CODE,
660 self.CODE,
659 label,
661 label,
660 self.getDateTime(self.data.min_time).strftime(
662 self.getDateTime(self.data.min_time).strftime(
661 '%Y%m%d'
663 '%Y%m%d'
662 ),
664 ),
663 )
665 )
664 )
666 )
665 fig.savefig(figname)
667 fig.savefig(figname)
666
668
667 def send_to_server(self):
669 def send_to_server(self):
668 '''
670 '''
669 '''
671 '''
670
672
671 if self.sender_counter < self.sender_period:
673 if self.sender_counter < self.sender_period:
672 self.sender_counter += 1
674 self.sender_counter += 1
673
675
674 self.sender_counter = 1
676 self.sender_counter = 1
675
677 self.data.meta['titles'] = self.titles
676 retries = 2
678 retries = 2
677 while True:
679 while True:
678 self.socket.send_string(self.data.jsonify())
680 self.socket.send_string(self.data.jsonify(self.plot_name, self.plot_type))
679 socks = dict(self.poll.poll(5000))
681 socks = dict(self.poll.poll(5000))
680 if socks.get(self.socket) == zmq.POLLIN:
682 if socks.get(self.socket) == zmq.POLLIN:
681 reply = self.socket.recv_string()
683 reply = self.socket.recv_string()
682 if reply == 'ok':
684 if reply == 'ok':
683 log.log("Response from server ok", self.name)
685 log.log("Response from server ok", self.name)
684 break
686 break
685 else:
687 else:
686 log.warning(
688 log.warning(
687 "Malformed reply from server: {}".format(reply), self.name)
689 "Malformed reply from server: {}".format(reply), self.name)
688
690
689 else:
691 else:
690 log.warning(
692 log.warning(
691 "No response from server, retrying...", self.name)
693 "No response from server, retrying...", self.name)
692 self.socket.setsockopt(zmq.LINGER, 0)
694 self.socket.setsockopt(zmq.LINGER, 0)
693 self.socket.close()
695 self.socket.close()
694 self.poll.unregister(self.socket)
696 self.poll.unregister(self.socket)
695 retries -= 1
697 retries -= 1
696 if retries == 0:
698 if retries == 0:
697 log.error(
699 log.error(
698 "Server seems to be offline, abandoning", self.name)
700 "Server seems to be offline, abandoning", self.name)
699 self.socket = self.context.socket(zmq.REQ)
701 self.socket = self.context.socket(zmq.REQ)
700 self.socket.connect(self.plot_server)
702 self.socket.connect(self.plot_server)
701 self.poll.register(self.socket, zmq.POLLIN)
703 self.poll.register(self.socket, zmq.POLLIN)
702 time.sleep(1)
704 time.sleep(1)
703 break
705 break
704 self.socket = self.context.socket(zmq.REQ)
706 self.socket = self.context.socket(zmq.REQ)
705 self.socket.connect(self.plot_server)
707 self.socket.connect(self.plot_server)
706 self.poll.register(self.socket, zmq.POLLIN)
708 self.poll.register(self.socket, zmq.POLLIN)
707 time.sleep(0.5)
709 time.sleep(0.5)
708
710
709 def setup(self):
711 def setup(self):
710 '''
712 '''
711 This method should be implemented in the child class, the following
713 This method should be implemented in the child class, the following
712 attributes should be set:
714 attributes should be set:
713
715
714 self.nrows: number of rows
716 self.nrows: number of rows
715 self.ncols: number of cols
717 self.ncols: number of cols
716 self.nplots: number of plots (channels or pairs)
718 self.nplots: number of plots (channels or pairs)
717 self.ylabel: label for Y axes
719 self.ylabel: label for Y axes
718 self.titles: list of axes title
720 self.titles: list of axes title
719
721
720 '''
722 '''
721 raise NotImplementedError
723 raise NotImplementedError
722
724
723 def plot(self):
725 def plot(self):
724 '''
726 '''
725 Must be defined in the child class
727 Must be defined in the child class
726 '''
728 '''
727 raise NotImplementedError
729 raise NotImplementedError
728
730
729 def run(self, dataOut, **kwargs):
731 def run(self, dataOut, **kwargs):
730 '''
732 '''
731 Main plotting routine
733 Main plotting routine
732 '''
734 '''
733
735
734 if self.isConfig is False:
736 if self.isConfig is False:
735 self.__setup(**kwargs)
737 self.__setup(**kwargs)
736 if dataOut.type == 'Parameters':
738 if dataOut.type == 'Parameters':
737 t = dataOut.utctimeInit
739 t = dataOut.utctimeInit
738 else:
740 else:
739 t = dataOut.utctime
741 t = dataOut.utctime
740
742
741 if dataOut.useLocalTime:
743 if dataOut.useLocalTime:
742 self.getDateTime = datetime.datetime.fromtimestamp
744 self.getDateTime = datetime.datetime.fromtimestamp
743 if not self.localtime:
745 if not self.localtime:
744 t += time.timezone
746 t += time.timezone
745 else:
747 else:
746 self.getDateTime = datetime.datetime.utcfromtimestamp
748 self.getDateTime = datetime.datetime.utcfromtimestamp
747 if self.localtime:
749 if self.localtime:
748 t -= time.timezone
750 t -= time.timezone
749
751
750 if self.xmin is None:
752 if self.xmin is None:
751 self.tmin = t
753 self.tmin = t
752 else:
754 else:
753 self.tmin = (
755 self.tmin = (
754 self.getDateTime(t).replace(
756 self.getDateTime(t).replace(
755 hour=self.xmin,
757 hour=self.xmin,
756 minute=0,
758 minute=0,
757 second=0) - self.getDateTime(0)).total_seconds()
759 second=0) - self.getDateTime(0)).total_seconds()
758
760
759 self.data.setup()
761 self.data.setup()
760 self.isConfig = True
762 self.isConfig = True
761 if self.plot_server:
763 if self.plot_server:
762 self.context = zmq.Context()
764 self.context = zmq.Context()
763 self.socket = self.context.socket(zmq.REQ)
765 self.socket = self.context.socket(zmq.REQ)
764 self.socket.connect(self.plot_server)
766 self.socket.connect(self.plot_server)
765 self.poll = zmq.Poller()
767 self.poll = zmq.Poller()
766 self.poll.register(self.socket, zmq.POLLIN)
768 self.poll.register(self.socket, zmq.POLLIN)
767
769
768 if dataOut.type == 'Parameters':
770 if dataOut.type == 'Parameters':
769 tm = dataOut.utctimeInit
771 tm = dataOut.utctimeInit
770 else:
772 else:
771 tm = dataOut.utctime
773 tm = dataOut.utctime
772
774
773 if not dataOut.useLocalTime and self.localtime:
775 if not dataOut.useLocalTime and self.localtime:
774 tm -= time.timezone
776 tm -= time.timezone
775 if dataOut.useLocalTime and not self.localtime:
777 if dataOut.useLocalTime and not self.localtime:
776 tm += time.timezone
778 tm += time.timezone
777
779
778 if self.xaxis is 'time' and self.data and (tm - self.tmin) >= self.xrange*60*60:
780 if self.xaxis is 'time' and self.data and (tm - self.tmin) >= self.xrange*60*60:
779 self.save_counter = self.save_period
781 self.save_counter = self.save_period
780 self.__plot()
782 self.__plot()
781 self.xmin += self.xrange
783 self.xmin += self.xrange
782 if self.xmin >= 24:
784 if self.xmin >= 24:
783 self.xmin -= 24
785 self.xmin -= 24
784 self.tmin += self.xrange*60*60
786 self.tmin += self.xrange*60*60
785 self.data.setup()
787 self.data.setup()
786 self.clear_figures()
788 self.clear_figures()
787
789
788 self.data.update(dataOut, tm)
790 self.data.update(dataOut, tm)
789
791
790 if self.isPlotConfig is False:
792 if self.isPlotConfig is False:
791 self.__setup_plot()
793 self.__setup_plot()
792 self.isPlotConfig = True
794 self.isPlotConfig = True
793
795
794 if self.realtime:
796 if self.realtime:
795 self.__plot()
797 self.__plot()
796 else:
798 else:
797 self.__throttle_plot(self.__plot)#, coerce=coerce)
799 self.__throttle_plot(self.__plot)#, coerce=coerce)
798
800
799 def close(self):
801 def close(self):
800
802
801 if self.data:
803 if self.data:
802 self.save_counter = self.save_period
804 self.save_counter = self.save_period
803 self.__plot()
805 self.__plot()
804 if self.data and self.pause:
806 if self.data and self.pause:
805 figpause(10)
807 figpause(10)
806
808
@@ -1,748 +1,759
1 '''
1 '''
2 New Plots Operations
2 New Plots Operations
3
3
4 @author: juan.espinoza@jro.igp.gob.pe
4 @author: juan.espinoza@jro.igp.gob.pe
5 '''
5 '''
6
6
7
7
8 import time
8 import time
9 import datetime
9 import datetime
10 import numpy
10 import numpy
11
11
12 from schainpy.model.graphics.jroplot_base import Plot, plt
12 from schainpy.model.graphics.jroplot_base import Plot, plt
13 from schainpy.utils import log
13 from schainpy.utils import log
14
14
15 EARTH_RADIUS = 6.3710e3
15 EARTH_RADIUS = 6.3710e3
16
16
17
17
18 def ll2xy(lat1, lon1, lat2, lon2):
18 def ll2xy(lat1, lon1, lat2, lon2):
19
19
20 p = 0.017453292519943295
20 p = 0.017453292519943295
21 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
21 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
22 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
22 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
23 r = 12742 * numpy.arcsin(numpy.sqrt(a))
23 r = 12742 * numpy.arcsin(numpy.sqrt(a))
24 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
24 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
25 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
25 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
26 theta = -theta + numpy.pi/2
26 theta = -theta + numpy.pi/2
27 return r*numpy.cos(theta), r*numpy.sin(theta)
27 return r*numpy.cos(theta), r*numpy.sin(theta)
28
28
29
29
30 def km2deg(km):
30 def km2deg(km):
31 '''
31 '''
32 Convert distance in km to degrees
32 Convert distance in km to degrees
33 '''
33 '''
34
34
35 return numpy.rad2deg(km/EARTH_RADIUS)
35 return numpy.rad2deg(km/EARTH_RADIUS)
36
36
37
37
38 class SpectraPlot(Plot):
38 class SpectraPlot(Plot):
39 '''
39 '''
40 Plot for Spectra data
40 Plot for Spectra data
41 '''
41 '''
42
42
43 CODE = 'spc'
43 CODE = 'spc'
44 colormap = 'jro'
44 colormap = 'jro'
45 plot_name = 'Spectra'
46 plot_type = 'pcolor'
45
47
46 def setup(self):
48 def setup(self):
47 self.nplots = len(self.data.channels)
49 self.nplots = len(self.data.channels)
48 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
50 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
49 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
51 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
50 self.width = 3.4 * self.ncols
52 self.width = 3.4 * self.ncols
51 self.height = 3 * self.nrows
53 self.height = 3 * self.nrows
52 self.cb_label = 'dB'
54 self.cb_label = 'dB'
53 if self.showprofile:
55 if self.showprofile:
54 self.width += 0.8 * self.ncols
56 self.width += 0.8 * self.ncols
55
57
56 self.ylabel = 'Range [km]'
58 self.ylabel = 'Range [km]'
57
59
58 def plot(self):
60 def plot(self):
59 if self.xaxis == "frequency":
61 if self.xaxis == "frequency":
60 x = self.data.xrange[0]
62 x = self.data.xrange[0]
61 self.xlabel = "Frequency (kHz)"
63 self.xlabel = "Frequency (kHz)"
62 elif self.xaxis == "time":
64 elif self.xaxis == "time":
63 x = self.data.xrange[1]
65 x = self.data.xrange[1]
64 self.xlabel = "Time (ms)"
66 self.xlabel = "Time (ms)"
65 else:
67 else:
66 x = self.data.xrange[2]
68 x = self.data.xrange[2]
67 self.xlabel = "Velocity (m/s)"
69 self.xlabel = "Velocity (m/s)"
68
70
69 if self.CODE == 'spc_moments':
71 if self.CODE == 'spc_moments':
70 x = self.data.xrange[2]
72 x = self.data.xrange[2]
71 self.xlabel = "Velocity (m/s)"
73 self.xlabel = "Velocity (m/s)"
72
74
73 self.titles = []
75 self.titles = []
74
76
75 y = self.data.heights
77 y = self.data.heights
76 self.y = y
78 self.y = y
77 z = self.data['spc']
79 z = self.data['spc']
78
80
79 for n, ax in enumerate(self.axes):
81 for n, ax in enumerate(self.axes):
80 noise = self.data['noise'][n][-1]
82 noise = self.data['noise'][n][-1]
81 if self.CODE == 'spc_moments':
83 if self.CODE == 'spc_moments':
82 mean = self.data['moments'][n, :, 1, :][-1]
84 mean = self.data['moments'][n, :, 1, :][-1]
83 if ax.firsttime:
85 if ax.firsttime:
84 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
86 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
85 self.xmin = self.xmin if self.xmin else -self.xmax
87 self.xmin = self.xmin if self.xmin else -self.xmax
86 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
88 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
87 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
89 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
88 ax.plt = ax.pcolormesh(x, y, z[n].T,
90 ax.plt = ax.pcolormesh(x, y, z[n].T,
89 vmin=self.zmin,
91 vmin=self.zmin,
90 vmax=self.zmax,
92 vmax=self.zmax,
91 cmap=plt.get_cmap(self.colormap)
93 cmap=plt.get_cmap(self.colormap)
92 )
94 )
93
95
94 if self.showprofile:
96 if self.showprofile:
95 ax.plt_profile = self.pf_axes[n].plot(
97 ax.plt_profile = self.pf_axes[n].plot(
96 self.data['rti'][n][-1], y)[0]
98 self.data['rti'][n][-1], y)[0]
97 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
99 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
98 color="k", linestyle="dashed", lw=1)[0]
100 color="k", linestyle="dashed", lw=1)[0]
99 if self.CODE == 'spc_moments':
101 if self.CODE == 'spc_moments':
100 ax.plt_mean = ax.plot(mean, y, color='k')[0]
102 ax.plt_mean = ax.plot(mean, y, color='k')[0]
101 else:
103 else:
102 ax.plt.set_array(z[n].T.ravel())
104 ax.plt.set_array(z[n].T.ravel())
103 if self.showprofile:
105 if self.showprofile:
104 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
106 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
105 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
107 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
106 if self.CODE == 'spc_moments':
108 if self.CODE == 'spc_moments':
107 ax.plt_mean.set_data(mean, y)
109 ax.plt_mean.set_data(mean, y)
108 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
110 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
109
111
110
112
111 class CrossSpectraPlot(Plot):
113 class CrossSpectraPlot(Plot):
112
114
113 CODE = 'cspc'
115 CODE = 'cspc'
114 colormap = 'jet'
116 colormap = 'jet'
117 plot_name = 'CrossSpectra'
118 plot_type = 'pcolor'
115 zmin_coh = None
119 zmin_coh = None
116 zmax_coh = None
120 zmax_coh = None
117 zmin_phase = None
121 zmin_phase = None
118 zmax_phase = None
122 zmax_phase = None
119
123
120 def setup(self):
124 def setup(self):
121
125
122 self.ncols = 4
126 self.ncols = 4
123 self.nrows = len(self.data.pairs)
127 self.nrows = len(self.data.pairs)
124 self.nplots = self.nrows * 4
128 self.nplots = self.nrows * 4
125 self.width = 3.4 * self.ncols
129 self.width = 3.4 * self.ncols
126 self.height = 3 * self.nrows
130 self.height = 3 * self.nrows
127 self.ylabel = 'Range [km]'
131 self.ylabel = 'Range [km]'
128 self.showprofile = False
132 self.showprofile = False
129
133
130 def plot(self):
134 def plot(self):
131
135
132 if self.xaxis == "frequency":
136 if self.xaxis == "frequency":
133 x = self.data.xrange[0]
137 x = self.data.xrange[0]
134 self.xlabel = "Frequency (kHz)"
138 self.xlabel = "Frequency (kHz)"
135 elif self.xaxis == "time":
139 elif self.xaxis == "time":
136 x = self.data.xrange[1]
140 x = self.data.xrange[1]
137 self.xlabel = "Time (ms)"
141 self.xlabel = "Time (ms)"
138 else:
142 else:
139 x = self.data.xrange[2]
143 x = self.data.xrange[2]
140 self.xlabel = "Velocity (m/s)"
144 self.xlabel = "Velocity (m/s)"
141
145
142 self.titles = []
146 self.titles = []
143
147
144 y = self.data.heights
148 y = self.data.heights
145 self.y = y
149 self.y = y
146 spc = self.data['spc']
150 spc = self.data['spc']
147 cspc = self.data['cspc']
151 cspc = self.data['cspc']
148
152
149 for n in range(self.nrows):
153 for n in range(self.nrows):
150 noise = self.data['noise'][n][-1]
154 noise = self.data['noise'][n][-1]
151 pair = self.data.pairs[n]
155 pair = self.data.pairs[n]
152 ax = self.axes[4 * n]
156 ax = self.axes[4 * n]
153 spc0 = 10.*numpy.log10(spc[pair[0]]/self.data.factor)
157 spc0 = 10.*numpy.log10(spc[pair[0]]/self.data.factor)
154 if ax.firsttime:
158 if ax.firsttime:
155 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
159 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
156 self.xmin = self.xmin if self.xmin else -self.xmax
160 self.xmin = self.xmin if self.xmin else -self.xmax
157 self.zmin = self.zmin if self.zmin else numpy.nanmin(spc)
161 self.zmin = self.zmin if self.zmin else numpy.nanmin(spc)
158 self.zmax = self.zmax if self.zmax else numpy.nanmax(spc)
162 self.zmax = self.zmax if self.zmax else numpy.nanmax(spc)
159 ax.plt = ax.pcolormesh(x , y , spc0.T,
163 ax.plt = ax.pcolormesh(x , y , spc0.T,
160 vmin=self.zmin,
164 vmin=self.zmin,
161 vmax=self.zmax,
165 vmax=self.zmax,
162 cmap=plt.get_cmap(self.colormap)
166 cmap=plt.get_cmap(self.colormap)
163 )
167 )
164 else:
168 else:
165 ax.plt.set_array(spc0.T.ravel())
169 ax.plt.set_array(spc0.T.ravel())
166 self.titles.append('CH {}: {:3.2f}dB'.format(pair[0], noise))
170 self.titles.append('CH {}: {:3.2f}dB'.format(pair[0], noise))
167
171
168 ax = self.axes[4 * n + 1]
172 ax = self.axes[4 * n + 1]
169 spc1 = 10.*numpy.log10(spc[pair[1]]/self.data.factor)
173 spc1 = 10.*numpy.log10(spc[pair[1]]/self.data.factor)
170 if ax.firsttime:
174 if ax.firsttime:
171 ax.plt = ax.pcolormesh(x , y, spc1.T,
175 ax.plt = ax.pcolormesh(x , y, spc1.T,
172 vmin=self.zmin,
176 vmin=self.zmin,
173 vmax=self.zmax,
177 vmax=self.zmax,
174 cmap=plt.get_cmap(self.colormap)
178 cmap=plt.get_cmap(self.colormap)
175 )
179 )
176 else:
180 else:
177 ax.plt.set_array(spc1.T.ravel())
181 ax.plt.set_array(spc1.T.ravel())
178 self.titles.append('CH {}: {:3.2f}dB'.format(pair[1], noise))
182 self.titles.append('CH {}: {:3.2f}dB'.format(pair[1], noise))
179
183
180 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
184 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
181 coh = numpy.abs(out)
185 coh = numpy.abs(out)
182 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
186 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
183
187
184 ax = self.axes[4 * n + 2]
188 ax = self.axes[4 * n + 2]
185 if ax.firsttime:
189 if ax.firsttime:
186 ax.plt = ax.pcolormesh(x, y, coh.T,
190 ax.plt = ax.pcolormesh(x, y, coh.T,
187 vmin=0,
191 vmin=0,
188 vmax=1,
192 vmax=1,
189 cmap=plt.get_cmap(self.colormap_coh)
193 cmap=plt.get_cmap(self.colormap_coh)
190 )
194 )
191 else:
195 else:
192 ax.plt.set_array(coh.T.ravel())
196 ax.plt.set_array(coh.T.ravel())
193 self.titles.append(
197 self.titles.append(
194 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
198 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
195
199
196 ax = self.axes[4 * n + 3]
200 ax = self.axes[4 * n + 3]
197 if ax.firsttime:
201 if ax.firsttime:
198 ax.plt = ax.pcolormesh(x, y, phase.T,
202 ax.plt = ax.pcolormesh(x, y, phase.T,
199 vmin=-180,
203 vmin=-180,
200 vmax=180,
204 vmax=180,
201 cmap=plt.get_cmap(self.colormap_phase)
205 cmap=plt.get_cmap(self.colormap_phase)
202 )
206 )
203 else:
207 else:
204 ax.plt.set_array(phase.T.ravel())
208 ax.plt.set_array(phase.T.ravel())
205 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
209 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
206
210
207
211
208 class SpectralMomentsPlot(SpectraPlot):
212 class SpectralMomentsPlot(SpectraPlot):
209 '''
213 '''
210 Plot for Spectral Moments
214 Plot for Spectral Moments
211 '''
215 '''
212 CODE = 'spc_moments'
216 CODE = 'spc_moments'
213 colormap = 'jro'
217 colormap = 'jro'
218 plot_name = 'SpectralMoments'
219 plot_type = 'pcolor'
214
220
215
221
216 class RTIPlot(Plot):
222 class RTIPlot(Plot):
217 '''
223 '''
218 Plot for RTI data
224 Plot for RTI data
219 '''
225 '''
220
226
221 CODE = 'rti'
227 CODE = 'rti'
222 colormap = 'jro'
228 colormap = 'jro'
229 plot_name = 'RTI'
230 plot_type = 'pcolorbuffer'
223
231
224 def setup(self):
232 def setup(self):
225 self.xaxis = 'time'
233 self.xaxis = 'time'
226 self.ncols = 1
234 self.ncols = 1
227 self.nrows = len(self.data.channels)
235 self.nrows = len(self.data.channels)
228 self.nplots = len(self.data.channels)
236 self.nplots = len(self.data.channels)
229 self.ylabel = 'Range [km]'
237 self.ylabel = 'Range [km]'
230 self.cb_label = 'dB'
238 self.cb_label = 'dB'
231 self.titles = ['{} Channel {}'.format(
239 self.titles = ['{} Channel {}'.format(
232 self.CODE.upper(), x) for x in range(self.nrows)]
240 self.CODE.upper(), x) for x in range(self.nrows)]
233
241
234 def plot(self):
242 def plot(self):
235 self.x = self.data.times
243 self.x = self.data.times
236 self.y = self.data.heights
244 self.y = self.data.heights
237 self.z = self.data[self.CODE]
245 self.z = self.data[self.CODE]
238 self.z = numpy.ma.masked_invalid(self.z)
246 self.z = numpy.ma.masked_invalid(self.z)
239
247
240 if self.decimation is None:
248 if self.decimation is None:
241 x, y, z = self.fill_gaps(self.x, self.y, self.z)
249 x, y, z = self.fill_gaps(self.x, self.y, self.z)
242 else:
250 else:
243 x, y, z = self.fill_gaps(*self.decimate())
251 x, y, z = self.fill_gaps(*self.decimate())
244
252
245 for n, ax in enumerate(self.axes):
253 for n, ax in enumerate(self.axes):
246 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
254 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
247 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
255 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
248 if ax.firsttime:
256 if ax.firsttime:
249 ax.plt = ax.pcolormesh(x, y, z[n].T,
257 ax.plt = ax.pcolormesh(x, y, z[n].T,
250 vmin=self.zmin,
258 vmin=self.zmin,
251 vmax=self.zmax,
259 vmax=self.zmax,
252 cmap=plt.get_cmap(self.colormap)
260 cmap=plt.get_cmap(self.colormap)
253 )
261 )
254 if self.showprofile:
262 if self.showprofile:
255 ax.plot_profile = self.pf_axes[n].plot(
263 ax.plot_profile = self.pf_axes[n].plot(
256 self.data['rti'][n][-1], self.y)[0]
264 self.data['rti'][n][-1], self.y)[0]
257 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
265 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
258 color="k", linestyle="dashed", lw=1)[0]
266 color="k", linestyle="dashed", lw=1)[0]
259 else:
267 else:
260 ax.collections.remove(ax.collections[0])
268 ax.collections.remove(ax.collections[0])
261 ax.plt = ax.pcolormesh(x, y, z[n].T,
269 ax.plt = ax.pcolormesh(x, y, z[n].T,
262 vmin=self.zmin,
270 vmin=self.zmin,
263 vmax=self.zmax,
271 vmax=self.zmax,
264 cmap=plt.get_cmap(self.colormap)
272 cmap=plt.get_cmap(self.colormap)
265 )
273 )
266 if self.showprofile:
274 if self.showprofile:
267 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
275 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
268 ax.plot_noise.set_data(numpy.repeat(
276 ax.plot_noise.set_data(numpy.repeat(
269 self.data['noise'][n][-1], len(self.y)), self.y)
277 self.data['noise'][n][-1], len(self.y)), self.y)
270
278
271
279
272 class CoherencePlot(RTIPlot):
280 class CoherencePlot(RTIPlot):
273 '''
281 '''
274 Plot for Coherence data
282 Plot for Coherence data
275 '''
283 '''
276
284
277 CODE = 'coh'
285 CODE = 'coh'
286 plot_name = 'Coherence'
278
287
279 def setup(self):
288 def setup(self):
280 self.xaxis = 'time'
289 self.xaxis = 'time'
281 self.ncols = 1
290 self.ncols = 1
282 self.nrows = len(self.data.pairs)
291 self.nrows = len(self.data.pairs)
283 self.nplots = len(self.data.pairs)
292 self.nplots = len(self.data.pairs)
284 self.ylabel = 'Range [km]'
293 self.ylabel = 'Range [km]'
285 if self.CODE == 'coh':
294 if self.CODE == 'coh':
286 self.cb_label = ''
295 self.cb_label = ''
287 self.titles = [
296 self.titles = [
288 'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
297 'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
289 else:
298 else:
290 self.cb_label = 'Degrees'
299 self.cb_label = 'Degrees'
291 self.titles = [
300 self.titles = [
292 'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
301 'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
293
302
294
303
295 class PhasePlot(CoherencePlot):
304 class PhasePlot(CoherencePlot):
296 '''
305 '''
297 Plot for Phase map data
306 Plot for Phase map data
298 '''
307 '''
299
308
300 CODE = 'phase'
309 CODE = 'phase'
301 colormap = 'seismic'
310 colormap = 'seismic'
311 plot_name = 'Phase'
302
312
303
313
304 class NoisePlot(Plot):
314 class NoisePlot(Plot):
305 '''
315 '''
306 Plot for noise
316 Plot for noise
307 '''
317 '''
308
318
309 CODE = 'noise'
319 CODE = 'noise'
320 plot_name = 'Noise'
321 plot_type = 'scatterbuffer'
322
310
323
311 def setup(self):
324 def setup(self):
312 self.xaxis = 'time'
325 self.xaxis = 'time'
313 self.ncols = 1
326 self.ncols = 1
314 self.nrows = 1
327 self.nrows = 1
315 self.nplots = 1
328 self.nplots = 1
316 self.ylabel = 'Intensity [dB]'
329 self.ylabel = 'Intensity [dB]'
317 self.titles = ['Noise']
330 self.titles = ['Noise']
318 self.colorbar = False
331 self.colorbar = False
319
332
320 def plot(self):
333 def plot(self):
321
334
322 x = self.data.times
335 x = self.data.times
323 xmin = self.data.min_time
336 xmin = self.data.min_time
324 xmax = xmin + self.xrange * 60 * 60
337 xmax = xmin + self.xrange * 60 * 60
325 Y = self.data[self.CODE]
338 Y = self.data[self.CODE]
326
339
327 if self.axes[0].firsttime:
340 if self.axes[0].firsttime:
328 for ch in self.data.channels:
341 for ch in self.data.channels:
329 y = Y[ch]
342 y = Y[ch]
330 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
343 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
331 plt.legend()
344 plt.legend()
332 else:
345 else:
333 for ch in self.data.channels:
346 for ch in self.data.channels:
334 y = Y[ch]
347 y = Y[ch]
335 self.axes[0].lines[ch].set_data(x, y)
348 self.axes[0].lines[ch].set_data(x, y)
336
349
337 self.ymin = numpy.nanmin(Y) - 5
350 self.ymin = numpy.nanmin(Y) - 5
338 self.ymax = numpy.nanmax(Y) + 5
351 self.ymax = numpy.nanmax(Y) + 5
339
352
340
353
341 class SnrPlot(RTIPlot):
354 class SnrPlot(RTIPlot):
342 '''
355 '''
343 Plot for SNR Data
356 Plot for SNR Data
344 '''
357 '''
345
358
346 CODE = 'snr'
359 CODE = 'snr'
347 colormap = 'jet'
360 colormap = 'jet'
361 plot_name = 'SNR'
348
362
349
363
350 class DopplerPlot(RTIPlot):
364 class DopplerPlot(RTIPlot):
351 '''
365 '''
352 Plot for DOPPLER Data
366 Plot for DOPPLER Data
353 '''
367 '''
354
368
355 CODE = 'dop'
369 CODE = 'dop'
356 colormap = 'jet'
370 colormap = 'jet'
371 plot_name = 'Doppler'
357
372
358
373
359 class SkyMapPlot(Plot):
374 class SkyMapPlot(Plot):
360 '''
375 '''
361 Plot for meteors detection data
376 Plot for meteors detection data
362 '''
377 '''
363
378
364 CODE = 'param'
379 CODE = 'param'
365
380
366 def setup(self):
381 def setup(self):
367
382
368 self.ncols = 1
383 self.ncols = 1
369 self.nrows = 1
384 self.nrows = 1
370 self.width = 7.2
385 self.width = 7.2
371 self.height = 7.2
386 self.height = 7.2
372 self.nplots = 1
387 self.nplots = 1
373 self.xlabel = 'Zonal Zenith Angle (deg)'
388 self.xlabel = 'Zonal Zenith Angle (deg)'
374 self.ylabel = 'Meridional Zenith Angle (deg)'
389 self.ylabel = 'Meridional Zenith Angle (deg)'
375 self.polar = True
390 self.polar = True
376 self.ymin = -180
391 self.ymin = -180
377 self.ymax = 180
392 self.ymax = 180
378 self.colorbar = False
393 self.colorbar = False
379
394
380 def plot(self):
395 def plot(self):
381
396
382 arrayParameters = numpy.concatenate(self.data['param'])
397 arrayParameters = numpy.concatenate(self.data['param'])
383 error = arrayParameters[:, -1]
398 error = arrayParameters[:, -1]
384 indValid = numpy.where(error == 0)[0]
399 indValid = numpy.where(error == 0)[0]
385 finalMeteor = arrayParameters[indValid, :]
400 finalMeteor = arrayParameters[indValid, :]
386 finalAzimuth = finalMeteor[:, 3]
401 finalAzimuth = finalMeteor[:, 3]
387 finalZenith = finalMeteor[:, 4]
402 finalZenith = finalMeteor[:, 4]
388
403
389 x = finalAzimuth * numpy.pi / 180
404 x = finalAzimuth * numpy.pi / 180
390 y = finalZenith
405 y = finalZenith
391
406
392 ax = self.axes[0]
407 ax = self.axes[0]
393
408
394 if ax.firsttime:
409 if ax.firsttime:
395 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
410 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
396 else:
411 else:
397 ax.plot.set_data(x, y)
412 ax.plot.set_data(x, y)
398
413
399 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
414 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
400 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
415 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
401 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
416 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
402 dt2,
417 dt2,
403 len(x))
418 len(x))
404 self.titles[0] = title
419 self.titles[0] = title
405
420
406
421
407 class ParametersPlot(RTIPlot):
422 class ParametersPlot(RTIPlot):
408 '''
423 '''
409 Plot for data_param object
424 Plot for data_param object
410 '''
425 '''
411
426
412 CODE = 'param'
427 CODE = 'param'
413 colormap = 'seismic'
428 colormap = 'seismic'
429 plot_name = 'Parameters'
414
430
415 def setup(self):
431 def setup(self):
416 self.xaxis = 'time'
432 self.xaxis = 'time'
417 self.ncols = 1
433 self.ncols = 1
418 self.nrows = self.data.shape(self.CODE)[0]
434 self.nrows = self.data.shape(self.CODE)[0]
419 self.nplots = self.nrows
435 self.nplots = self.nrows
420 if self.showSNR:
436 if self.showSNR:
421 self.nrows += 1
437 self.nrows += 1
422 self.nplots += 1
438 self.nplots += 1
423
439
424 self.ylabel = 'Height [km]'
440 self.ylabel = 'Height [km]'
425 if not self.titles:
441 if not self.titles:
426 self.titles = self.data.parameters \
442 self.titles = self.data.parameters \
427 if self.data.parameters else ['Param {}'.format(x) for x in range(self.nrows)]
443 if self.data.parameters else ['Param {}'.format(x) for x in range(self.nrows)]
428 if self.showSNR:
444 if self.showSNR:
429 self.titles.append('SNR')
445 self.titles.append('SNR')
430
446
431 def plot(self):
447 def plot(self):
432 self.data.normalize_heights()
448 self.data.normalize_heights()
433 self.x = self.data.times
449 self.x = self.data.times
434 self.y = self.data.heights
450 self.y = self.data.heights
435 if self.showSNR:
451 if self.showSNR:
436 self.z = numpy.concatenate(
452 self.z = numpy.concatenate(
437 (self.data[self.CODE], self.data['snr'])
453 (self.data[self.CODE], self.data['snr'])
438 )
454 )
439 else:
455 else:
440 self.z = self.data[self.CODE]
456 self.z = self.data[self.CODE]
441
457
442 self.z = numpy.ma.masked_invalid(self.z)
458 self.z = numpy.ma.masked_invalid(self.z)
443
459
444 if self.decimation is None:
460 if self.decimation is None:
445 x, y, z = self.fill_gaps(self.x, self.y, self.z)
461 x, y, z = self.fill_gaps(self.x, self.y, self.z)
446 else:
462 else:
447 x, y, z = self.fill_gaps(*self.decimate())
463 x, y, z = self.fill_gaps(*self.decimate())
448
464
449 for n, ax in enumerate(self.axes):
465 for n, ax in enumerate(self.axes):
450
466
451 self.zmax = self.zmax if self.zmax is not None else numpy.max(
467 self.zmax = self.zmax if self.zmax is not None else numpy.max(
452 self.z[n])
468 self.z[n])
453 self.zmin = self.zmin if self.zmin is not None else numpy.min(
469 self.zmin = self.zmin if self.zmin is not None else numpy.min(
454 self.z[n])
470 self.z[n])
455
471
456 if ax.firsttime:
472 if ax.firsttime:
457 if self.zlimits is not None:
473 if self.zlimits is not None:
458 self.zmin, self.zmax = self.zlimits[n]
474 self.zmin, self.zmax = self.zlimits[n]
459
475
460 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
476 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
461 vmin=self.zmin,
477 vmin=self.zmin,
462 vmax=self.zmax,
478 vmax=self.zmax,
463 cmap=self.cmaps[n]
479 cmap=self.cmaps[n]
464 )
480 )
465 else:
481 else:
466 if self.zlimits is not None:
482 if self.zlimits is not None:
467 self.zmin, self.zmax = self.zlimits[n]
483 self.zmin, self.zmax = self.zlimits[n]
468 ax.collections.remove(ax.collections[0])
484 ax.collections.remove(ax.collections[0])
469 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
485 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
470 vmin=self.zmin,
486 vmin=self.zmin,
471 vmax=self.zmax,
487 vmax=self.zmax,
472 cmap=self.cmaps[n]
488 cmap=self.cmaps[n]
473 )
489 )
474
490
475
491
476 class OutputPlot(ParametersPlot):
492 class OutputPlot(ParametersPlot):
477 '''
493 '''
478 Plot data_output object
494 Plot data_output object
479 '''
495 '''
480
496
481 CODE = 'output'
497 CODE = 'output'
482 colormap = 'seismic'
498 colormap = 'seismic'
499 plot_name = 'Output'
483
500
484
501
485 class PolarMapPlot(Plot):
502 class PolarMapPlot(Plot):
486 '''
503 '''
487 Plot for weather radar
504 Plot for weather radar
488 '''
505 '''
489
506
490 CODE = 'param'
507 CODE = 'param'
491 colormap = 'seismic'
508 colormap = 'seismic'
492
509
493 def setup(self):
510 def setup(self):
494 self.ncols = 1
511 self.ncols = 1
495 self.nrows = 1
512 self.nrows = 1
496 self.width = 9
513 self.width = 9
497 self.height = 8
514 self.height = 8
498 self.mode = self.data.meta['mode']
515 self.mode = self.data.meta['mode']
499 if self.channels is not None:
516 if self.channels is not None:
500 self.nplots = len(self.channels)
517 self.nplots = len(self.channels)
501 self.nrows = len(self.channels)
518 self.nrows = len(self.channels)
502 else:
519 else:
503 self.nplots = self.data.shape(self.CODE)[0]
520 self.nplots = self.data.shape(self.CODE)[0]
504 self.nrows = self.nplots
521 self.nrows = self.nplots
505 self.channels = list(range(self.nplots))
522 self.channels = list(range(self.nplots))
506 if self.mode == 'E':
523 if self.mode == 'E':
507 self.xlabel = 'Longitude'
524 self.xlabel = 'Longitude'
508 self.ylabel = 'Latitude'
525 self.ylabel = 'Latitude'
509 else:
526 else:
510 self.xlabel = 'Range (km)'
527 self.xlabel = 'Range (km)'
511 self.ylabel = 'Height (km)'
528 self.ylabel = 'Height (km)'
512 self.bgcolor = 'white'
529 self.bgcolor = 'white'
513 self.cb_labels = self.data.meta['units']
530 self.cb_labels = self.data.meta['units']
514 self.lat = self.data.meta['latitude']
531 self.lat = self.data.meta['latitude']
515 self.lon = self.data.meta['longitude']
532 self.lon = self.data.meta['longitude']
516 self.xmin, self.xmax = float(
533 self.xmin, self.xmax = float(
517 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
534 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
518 self.ymin, self.ymax = float(
535 self.ymin, self.ymax = float(
519 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
536 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
520 # self.polar = True
537 # self.polar = True
521
538
522 def plot(self):
539 def plot(self):
523
540
524 for n, ax in enumerate(self.axes):
541 for n, ax in enumerate(self.axes):
525 data = self.data['param'][self.channels[n]]
542 data = self.data['param'][self.channels[n]]
526
543
527 zeniths = numpy.linspace(
544 zeniths = numpy.linspace(
528 0, self.data.meta['max_range'], data.shape[1])
545 0, self.data.meta['max_range'], data.shape[1])
529 if self.mode == 'E':
546 if self.mode == 'E':
530 azimuths = -numpy.radians(self.data.heights)+numpy.pi/2
547 azimuths = -numpy.radians(self.data.heights)+numpy.pi/2
531 r, theta = numpy.meshgrid(zeniths, azimuths)
548 r, theta = numpy.meshgrid(zeniths, azimuths)
532 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
549 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
533 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
550 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
534 x = km2deg(x) + self.lon
551 x = km2deg(x) + self.lon
535 y = km2deg(y) + self.lat
552 y = km2deg(y) + self.lat
536 else:
553 else:
537 azimuths = numpy.radians(self.data.heights)
554 azimuths = numpy.radians(self.data.heights)
538 r, theta = numpy.meshgrid(zeniths, azimuths)
555 r, theta = numpy.meshgrid(zeniths, azimuths)
539 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
556 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
540 self.y = zeniths
557 self.y = zeniths
541
558
542 if ax.firsttime:
559 if ax.firsttime:
543 if self.zlimits is not None:
560 if self.zlimits is not None:
544 self.zmin, self.zmax = self.zlimits[n]
561 self.zmin, self.zmax = self.zlimits[n]
545 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
562 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
546 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
563 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
547 vmin=self.zmin,
564 vmin=self.zmin,
548 vmax=self.zmax,
565 vmax=self.zmax,
549 cmap=self.cmaps[n])
566 cmap=self.cmaps[n])
550 else:
567 else:
551 if self.zlimits is not None:
568 if self.zlimits is not None:
552 self.zmin, self.zmax = self.zlimits[n]
569 self.zmin, self.zmax = self.zlimits[n]
553 ax.collections.remove(ax.collections[0])
570 ax.collections.remove(ax.collections[0])
554 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
571 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
555 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
572 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
556 vmin=self.zmin,
573 vmin=self.zmin,
557 vmax=self.zmax,
574 vmax=self.zmax,
558 cmap=self.cmaps[n])
575 cmap=self.cmaps[n])
559
576
560 if self.mode == 'A':
577 if self.mode == 'A':
561 continue
578 continue
562
579
563 # plot district names
580 # plot district names
564 f = open('/data/workspace/schain_scripts/distrito.csv')
581 f = open('/data/workspace/schain_scripts/distrito.csv')
565 for line in f:
582 for line in f:
566 label, lon, lat = [s.strip() for s in line.split(',') if s]
583 label, lon, lat = [s.strip() for s in line.split(',') if s]
567 lat = float(lat)
584 lat = float(lat)
568 lon = float(lon)
585 lon = float(lon)
569 # ax.plot(lon, lat, '.b', ms=2)
586 # ax.plot(lon, lat, '.b', ms=2)
570 ax.text(lon, lat, label.decode('utf8'), ha='center',
587 ax.text(lon, lat, label.decode('utf8'), ha='center',
571 va='bottom', size='8', color='black')
588 va='bottom', size='8', color='black')
572
589
573 # plot limites
590 # plot limites
574 limites = []
591 limites = []
575 tmp = []
592 tmp = []
576 for line in open('/data/workspace/schain_scripts/lima.csv'):
593 for line in open('/data/workspace/schain_scripts/lima.csv'):
577 if '#' in line:
594 if '#' in line:
578 if tmp:
595 if tmp:
579 limites.append(tmp)
596 limites.append(tmp)
580 tmp = []
597 tmp = []
581 continue
598 continue
582 values = line.strip().split(',')
599 values = line.strip().split(',')
583 tmp.append((float(values[0]), float(values[1])))
600 tmp.append((float(values[0]), float(values[1])))
584 for points in limites:
601 for points in limites:
585 ax.add_patch(
602 ax.add_patch(
586 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
603 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
587
604
588 # plot Cuencas
605 # plot Cuencas
589 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
606 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
590 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
607 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
591 values = [line.strip().split(',') for line in f]
608 values = [line.strip().split(',') for line in f]
592 points = [(float(s[0]), float(s[1])) for s in values]
609 points = [(float(s[0]), float(s[1])) for s in values]
593 ax.add_patch(Polygon(points, ec='b', fc='none'))
610 ax.add_patch(Polygon(points, ec='b', fc='none'))
594
611
595 # plot grid
612 # plot grid
596 for r in (15, 30, 45, 60):
613 for r in (15, 30, 45, 60):
597 ax.add_artist(plt.Circle((self.lon, self.lat),
614 ax.add_artist(plt.Circle((self.lon, self.lat),
598 km2deg(r), color='0.6', fill=False, lw=0.2))
615 km2deg(r), color='0.6', fill=False, lw=0.2))
599 ax.text(
616 ax.text(
600 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
617 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
601 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
618 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
602 '{}km'.format(r),
619 '{}km'.format(r),
603 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
620 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
604
621
605 if self.mode == 'E':
622 if self.mode == 'E':
606 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
623 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
607 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
624 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
608 else:
625 else:
609 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
626 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
610 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
627 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
611
628
612 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
629 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
613 self.titles = ['{} {}'.format(
630 self.titles = ['{} {}'.format(
614 self.data.parameters[x], title) for x in self.channels]
631 self.data.parameters[x], title) for x in self.channels]
615
632
616
633
617 class ScopePlot(Plot):
634 class ScopePlot(Plot):
618
635
619 '''
636 '''
620 Plot for Scope
637 Plot for Scope
621 '''
638 '''
622
639
623 CODE = 'scope'
640 CODE = 'scope'
641 plot_name = 'Scope'
642 plot_type = 'scatter'
624
643
625 def setup(self):
644 def setup(self):
626
645
627 self.xaxis = 'Range (Km)'
646 self.xaxis = 'Range (Km)'
628 self.ncols = 1
647 self.ncols = 1
629 self.nrows = 1
648 self.nrows = 1
630 self.nplots = 1
649 self.nplots = 1
631 self.ylabel = 'Intensity [dB]'
650 self.ylabel = 'Intensity [dB]'
632 self.titles = ['Scope']
651 self.titles = ['Scope']
633 self.colorbar = False
652 self.colorbar = False
634 colspan = 3
653 colspan = 3
635 rowspan = 1
654 rowspan = 1
636
655
637 def plot_iq(self, x, y, channelIndexList, thisDatetime, wintitle):
656 def plot_iq(self, x, y, channelIndexList, thisDatetime, wintitle):
638
657
639 yreal = y[channelIndexList,:].real
658 yreal = y[channelIndexList,:].real
640 yimag = y[channelIndexList,:].imag
659 yimag = y[channelIndexList,:].imag
641 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
660 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
642 self.xlabel = "Range (Km)"
661 self.xlabel = "Range (Km)"
643 self.ylabel = "Intensity - IQ"
662 self.ylabel = "Intensity - IQ"
644
663
645 self.y = yreal
664 self.y = yreal
646 self.x = x
665 self.x = x
647 self.xmin = min(x)
666 self.xmin = min(x)
648 self.xmax = max(x)
667 self.xmax = max(x)
649
668
650
669
651 self.titles[0] = title
670 self.titles[0] = title
652
671
653 for i,ax in enumerate(self.axes):
672 for i,ax in enumerate(self.axes):
654 title = "Channel %d" %(i)
673 title = "Channel %d" %(i)
655 if ax.firsttime:
674 if ax.firsttime:
656 ax.plt_r = ax.plot(x, yreal[i,:], color='b')[0]
675 ax.plt_r = ax.plot(x, yreal[i,:], color='b')[0]
657 ax.plt_i = ax.plot(x, yimag[i,:], color='r')[0]
676 ax.plt_i = ax.plot(x, yimag[i,:], color='r')[0]
658 else:
677 else:
659 #pass
678 #pass
660 ax.plt_r.set_data(x, yreal[i,:])
679 ax.plt_r.set_data(x, yreal[i,:])
661 ax.plt_i.set_data(x, yimag[i,:])
680 ax.plt_i.set_data(x, yimag[i,:])
662
681
663 def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle):
682 def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle):
664 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
683 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
665 yreal = y.real
684 yreal = y.real
666 self.y = yreal
685 self.y = yreal
667 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
686 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
668 self.xlabel = "Range (Km)"
687 self.xlabel = "Range (Km)"
669 self.ylabel = "Intensity"
688 self.ylabel = "Intensity"
670 self.xmin = min(x)
689 self.xmin = min(x)
671 self.xmax = max(x)
690 self.xmax = max(x)
672
691
673
692
674 self.titles[0] = title
693 self.titles[0] = title
675
694
676 for i,ax in enumerate(self.axes):
695 for i,ax in enumerate(self.axes):
677 title = "Channel %d" %(i)
696 title = "Channel %d" %(i)
678
697
679 ychannel = yreal[i,:]
698 ychannel = yreal[i,:]
680
699
681 if ax.firsttime:
700 if ax.firsttime:
682 ax.plt_r = ax.plot(x, ychannel)[0]
701 ax.plt_r = ax.plot(x, ychannel)[0]
683 else:
702 else:
684 #pass
703 #pass
685 ax.plt_r.set_data(x, ychannel)
704 ax.plt_r.set_data(x, ychannel)
686
705
687
706
688 def plot(self):
707 def plot(self):
689
708
690 if self.channels:
709 if self.channels:
691 channels = self.channels
710 channels = self.channels
692 else:
711 else:
693 channels = self.data.channels
712 channels = self.data.channels
694
713
695
714
696
715
697 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1])
716 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1])
698
717
699 scope = self.data['scope']
718 scope = self.data['scope']
700
719
701
720
702 if self.data.flagDataAsBlock:
721 if self.data.flagDataAsBlock:
703
722
704 for i in range(self.data.nProfiles):
723 for i in range(self.data.nProfiles):
705
724
706 wintitle1 = " [Profile = %d] " %i
725 wintitle1 = " [Profile = %d] " %i
707
726
708 if self.type == "power":
727 if self.type == "power":
709 self.plot_power(self.data.heights,
728 self.plot_power(self.data.heights,
710 scope[:,i,:],
729 scope[:,i,:],
711 channels,
730 channels,
712 thisDatetime,
731 thisDatetime,
713 wintitle1
732 wintitle1
714 )
733 )
715
734
716 if self.type == "iq":
735 if self.type == "iq":
717 self.plot_iq(self.data.heights,
736 self.plot_iq(self.data.heights,
718 scope[:,i,:],
737 scope[:,i,:],
719 channels,
738 channels,
720 thisDatetime,
739 thisDatetime,
721 wintitle1
740 wintitle1
722 )
741 )
723
724
725
726
727
728 else:
742 else:
729 wintitle = " [Profile = %d] " %self.data.profileIndex
743 wintitle = " [Profile = %d] " %self.data.profileIndex
730
744
731 if self.type == "power":
745 if self.type == "power":
732 self.plot_power(self.data.heights,
746 self.plot_power(self.data.heights,
733 scope,
747 scope,
734 channels,
748 channels,
735 thisDatetime,
749 thisDatetime,
736 wintitle
750 wintitle
737 )
751 )
738
752
739 if self.type == "iq":
753 if self.type == "iq":
740 self.plot_iq(self.data.heights,
754 self.plot_iq(self.data.heights,
741 scope,
755 scope,
742 channels,
756 channels,
743 thisDatetime,
757 thisDatetime,
744 wintitle
758 wintitle
745 )
759 )
746
747
748 No newline at end of file
General Comments 0
You need to be logged in to leave comments. Login now