##// END OF EJS Templates
schain
RSS
Fix publish and plots operations issue #929
Juan C. Espinoza -
r1062:8048843f4edf
parent child
Show More
Side by Side Unified
Context file:
r1062:8048843f4edf Loading diff...:

The requested changes are too big and content was truncated. Show full diff

Show file before | Show file after | Hide comments
1 NO CONTENT: modified file
Show file before | Show file after | Hide comments
This diff has been collapsed as it changes many lines, (1122 lines changed) Show them Hide them
@@ -1,964 +1,782
1 1
2 2 import os
3 import zmq
4 3 import time
5 import numpy
4 import glob
6 5 import datetime
7 import numpy as np
6 from multiprocessing import Process
7
8 import zmq
9 import numpy
8 10 import matplotlib
9 import glob
10 matplotlib.use('TkAgg')
11 11 import matplotlib.pyplot as plt
12 12 from mpl_toolkits.axes_grid1 import make_axes_locatable
13 from matplotlib.ticker import FuncFormatter, LinearLocator
14 from multiprocessing import Process
13 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
15 14
16 15 from schainpy.model.proc.jroproc_base import Operation
17
18 plt.ion()
16 from schainpy.utils import log
19 17
20 18 func = lambda x, pos: ('%s') %(datetime.datetime.fromtimestamp(x).strftime('%H:%M'))
21 fromtimestamp = lambda x, mintime : (datetime.datetime.utcfromtimestamp(mintime).replace(hour=(x + 5), minute=0) - d1970).total_seconds()
22
23 19
24 20 d1970 = datetime.datetime(1970,1,1)
25 21
22
26 23 class PlotData(Operation, Process):
24 '''
25 Base class for Schain plotting operations
26 '''
27 27
28 28 CODE = 'Figure'
29 29 colormap = 'jro'
30 bgcolor = 'white'
30 31 CONFLATE = False
31 32 __MAXNUMX = 80
32 33 __missing = 1E30
33 34
34 35 def __init__(self, **kwargs):
35 36
36 37 Operation.__init__(self, plot=True, **kwargs)
37 38 Process.__init__(self)
38 39 self.kwargs['code'] = self.CODE
39 40 self.mp = False
40 self.dataOut = None
41 self.data = None
41 42 self.isConfig = False
42 self.figure = None
43 self.figures = []
43 44 self.axes = []
45 self.cb_axes = []
44 46 self.localtime = kwargs.pop('localtime', True)
45 47 self.show = kwargs.get('show', True)
46 48 self.save = kwargs.get('save', False)
47 49 self.colormap = kwargs.get('colormap', self.colormap)
48 50 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
49 51 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
50 self.showprofile = kwargs.get('showprofile', True)
51 self.title = kwargs.get('wintitle', '')
52 self.colormaps = kwargs.get('colormaps', None)
53 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
54 self.showprofile = kwargs.get('showprofile', False)
55 self.title = kwargs.get('wintitle', self.CODE.upper())
56 self.cb_label = kwargs.get('cb_label', None)
57 self.cb_labels = kwargs.get('cb_labels', None)
52 58 self.xaxis = kwargs.get('xaxis', 'frequency')
53 59 self.zmin = kwargs.get('zmin', None)
54 60 self.zmax = kwargs.get('zmax', None)
61 self.zlimits = kwargs.get('zlimits', None)
55 62 self.xmin = kwargs.get('xmin', None)
63 if self.xmin is not None:
64 self.xmin += 5
56 65 self.xmax = kwargs.get('xmax', None)
57 66 self.xrange = kwargs.get('xrange', 24)
58 67 self.ymin = kwargs.get('ymin', None)
59 68 self.ymax = kwargs.get('ymax', None)
60 self.__MAXNUMY = kwargs.get('decimation', 5000)
61 self.throttle_value = 5
62 self.times = []
63 #self.interactive = self.kwargs['parent']
69 self.xlabel = kwargs.get('xlabel', None)
70 self.__MAXNUMY = kwargs.get('decimation', 100)
71 self.showSNR = kwargs.get('showSNR', False)
72 self.oneFigure = kwargs.get('oneFigure', True)
73 self.width = kwargs.get('width', None)
74 self.height = kwargs.get('height', None)
75 self.colorbar = kwargs.get('colorbar', True)
76 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
77 self.titles = ['' for __ in range(16)]
78
79 def __setup(self):
80 '''
81 Common setup for all figures, here figures and axes are created
82 '''
83
84 self.setup()
85
86 if self.width is None:
87 self.width = 8
88
89 self.figures = []
90 self.axes = []
91 self.cb_axes = []
92 self.pf_axes = []
93 self.cmaps = []
64 94
95 size = '15%' if self.ncols==1 else '30%'
96 pad = '4%' if self.ncols==1 else '8%'
97
98 if self.oneFigure:
99 if self.height is None:
100 self.height = 1.4*self.nrows + 1
101 fig = plt.figure(figsize=(self.width, self.height),
102 edgecolor='k',
103 facecolor='w')
104 self.figures.append(fig)
105 for n in range(self.nplots):
106 ax = fig.add_subplot(self.nrows, self.ncols, n+1)
107 ax.tick_params(labelsize=8)
108 ax.firsttime = True
109 self.axes.append(ax)
110 if self.showprofile:
111 cax = self.__add_axes(ax, size=size, pad=pad)
112 cax.tick_params(labelsize=8)
113 self.pf_axes.append(cax)
114 else:
115 if self.height is None:
116 self.height = 3
117 for n in range(self.nplots):
118 fig = plt.figure(figsize=(self.width, self.height),
119 edgecolor='k',
120 facecolor='w')
121 ax = fig.add_subplot(1, 1, 1)
122 ax.tick_params(labelsize=8)
123 ax.firsttime = True
124 self.figures.append(fig)
125 self.axes.append(ax)
126 if self.showprofile:
127 cax = self.__add_axes(ax, size=size, pad=pad)
128 cax.tick_params(labelsize=8)
129 self.pf_axes.append(cax)
130
131 for n in range(self.nrows):
132 if self.colormaps is not None:
133 cmap = plt.get_cmap(self.colormaps[n])
134 else:
135 cmap = plt.get_cmap(self.colormap)
136 cmap.set_bad(self.bgcolor, 1.)
137 self.cmaps.append(cmap)
138
139 def __add_axes(self, ax, size='30%', pad='8%'):
65 140 '''
66 this new parameter is created to plot data from varius channels at different figures
67 1. crear una lista de figuras donde se puedan plotear las figuras,
68 2. dar las opciones de configuracion a cada figura, estas opciones son iguales para ambas figuras
69 3. probar?
141 Add new axes to the given figure
70 142 '''
71 self.ind_plt_ch = kwargs.get('ind_plt_ch', False)
72 self.figurelist = None
143 divider = make_axes_locatable(ax)
144 nax = divider.new_horizontal(size=size, pad=pad)
145 ax.figure.add_axes(nax)
146 return nax
73 147
74 148
75 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
149 def setup(self):
150 '''
151 This method should be implemented in the child class, the following
152 attributes should be set:
76 153
154 self.nrows: number of rows
155 self.ncols: number of cols
156 self.nplots: number of plots (channels or pairs)
157 self.ylabel: label for Y axes
158 self.titles: list of axes title
159
160 '''
161 raise(NotImplementedError, 'Implement this method in child class')
162
163 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
164 '''
165 Create a masked array for missing data
166 '''
77 167 if x_buffer.shape[0] < 2:
78 168 return x_buffer, y_buffer, z_buffer
79 169
80 170 deltas = x_buffer[1:] - x_buffer[0:-1]
81 x_median = np.median(deltas)
171 x_median = numpy.median(deltas)
82 172
83 index = np.where(deltas > 5*x_median)
173 index = numpy.where(deltas > 5*x_median)
84 174
85 175 if len(index[0]) != 0:
86 176 z_buffer[::, index[0], ::] = self.__missing
87 z_buffer = np.ma.masked_inside(z_buffer,
177 z_buffer = numpy.ma.masked_inside(z_buffer,
88 178 0.99*self.__missing,
89 179 1.01*self.__missing)
90 180
91 181 return x_buffer, y_buffer, z_buffer
92 182
93 183 def decimate(self):
94 184
95 185 # dx = int(len(self.x)/self.__MAXNUMX) + 1
96 186 dy = int(len(self.y)/self.__MAXNUMY) + 1
97 187
98 188 # x = self.x[::dx]
99 189 x = self.x
100 190 y = self.y[::dy]
101 191 z = self.z[::, ::, ::dy]
102 192
103 193 return x, y, z
104 194
195 def format(self):
105 196 '''
106 JM:
107 elimana las otras imagenes generadas debido a que lso workers no llegan en orden y le pueden
108 poner otro tiempo a la figura q no necesariamente es el ultimo.
109 Solo se realiza cuando termina la imagen.
110 Problemas:
111
112 File "/home/ci-81/workspace/schainv2.3/schainpy/model/graphics/jroplot_data.py", line 145, in __plot
113 for n, eachfigure in enumerate(self.figurelist):
114 TypeError: 'NoneType' object is not iterable
115
197 Set min and max values, labels, ticks and titles
116 198 '''
117 def deleteanotherfiles(self):
118 figurenames=[]
119 if self.figurelist != None:
120 for n, eachfigure in enumerate(self.figurelist):
121 #add specific name for each channel in channelList
122 ghostfigname = os.path.join(self.save, '{}_{}_{}'.format(self.titles[n].replace(' ',''),self.CODE,
123 datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d')))
124 figname = os.path.join(self.save, '{}_{}_{}.png'.format(self.titles[n].replace(' ',''),self.CODE,
125 datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
126 199
127 for ghostfigure in glob.glob(ghostfigname+'*'): #ghostfigure will adopt all posible names of figures
128 if ghostfigure != figname:
129 os.remove(ghostfigure)
130 print 'Removing GhostFigures:' , figname
200 if self.xmin is None:
201 xmin = self.min_time
202 else:
203 if self.xaxis is 'time':
204 dt = datetime.datetime.fromtimestamp(self.min_time)
205 xmin = (datetime.datetime.combine(dt.date(),
206 datetime.time(int(self.xmin), 0, 0))-d1970).total_seconds()
131 207 else :
132 '''Erasing ghost images for just on******************'''
133 ghostfigname = os.path.join(self.save, '{}_{}'.format(self.CODE,datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d')))
134 figname = os.path.join(self.save, '{}_{}.png'.format(self.CODE,datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
135 for ghostfigure in glob.glob(ghostfigname+'*'): #ghostfigure will adopt all posible names of figures
136 if ghostfigure != figname:
137 os.remove(ghostfigure)
138 print 'Removing GhostFigures:' , figname
208 xmin = self.xmin
209
210 if self.xmax is None:
211 xmax = xmin+self.xrange*60*60
212 else:
213 if self.xaxis is 'time':
214 dt = datetime.datetime.fromtimestamp(self.min_time)
215 xmax = (datetime.datetime.combine(dt.date(),
216 datetime.time(int(self.xmax), 0, 0))-d1970).total_seconds()
217 else:
218 xmax = self.xmax
219
220 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
221 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
222
223 ystep = 200 if ymax>= 800 else 100 if ymax>=400 else 50 if ymax>=200 else 20
224
225 for n, ax in enumerate(self.axes):
226 if ax.firsttime:
227 ax.set_facecolor(self.bgcolor)
228 ax.yaxis.set_major_locator(MultipleLocator(ystep))
229 if self.xaxis is 'time':
230 ax.xaxis.set_major_formatter(FuncFormatter(func))
231 ax.xaxis.set_major_locator(LinearLocator(9))
232 if self.xlabel is not None:
233 ax.set_xlabel(self.xlabel)
234 ax.set_ylabel(self.ylabel)
235 ax.firsttime = False
236 if self.showprofile:
237 self.pf_axes[n].set_ylim(ymin, ymax)
238 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
239 self.pf_axes[n].set_xlabel('dB')
240 self.pf_axes[n].grid(b=True, axis='x')
241 [tick.set_visible(False) for tick in self.pf_axes[n].get_yticklabels()]
242 if self.colorbar:
243 cb = plt.colorbar(ax.plt, ax=ax, pad=0.02)
244 cb.ax.tick_params(labelsize=8)
245 if self.cb_label:
246 cb.set_label(self.cb_label, size=8)
247 elif self.cb_labels:
248 cb.set_label(self.cb_labels[n], size=8)
249
250 ax.set_title('{} - {} UTC'.format(
251 self.titles[n],
252 datetime.datetime.fromtimestamp(self.max_time).strftime('%H:%M:%S')),
253 size=8)
254 ax.set_xlim(xmin, xmax)
255 ax.set_ylim(ymin, ymax)
256
139 257
140 258 def __plot(self):
259 '''
260 '''
261 log.success('Plotting', self.name)
141 262
142 print 'plotting...{}'.format(self.CODE)
143 if self.ind_plt_ch is False : #standard
144 if self.show:
145 self.figure.show()
146 263 self.plot()
147 plt.tight_layout()
148 self.figure.canvas.manager.set_window_title('{} {} - {}'.format(self.title, self.CODE.upper(),
149 datetime.datetime.fromtimestamp(self.max_time).strftime('%Y/%m/%d')))
150 else :
151 print 'len(self.figurelist): ',len(self.figurelist)
152 for n, eachfigure in enumerate(self.figurelist):
264 self.format()
265
266 for n, fig in enumerate(self.figures):
267 if self.nrows == 0 or self.nplots == 0:
268 log.warning('No data', self.name)
269 continue
153 270 if self.show:
154 eachfigure.show()
271 fig.show()
155 272
156 self.plot()
157 eachfigure.tight_layout() # ajuste de cada subplot
158 eachfigure.canvas.manager.set_window_title('{} {} - {}'.format(self.title[n], self.CODE.upper(),
273 fig.tight_layout()
274 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
159 275 datetime.datetime.fromtimestamp(self.max_time).strftime('%Y/%m/%d')))
276 # fig.canvas.draw()
160 277
161 # if self.save:
162 # if self.ind_plt_ch is False : #standard
163 # figname = os.path.join(self.save, '{}_{}.png'.format(self.CODE,
164 # datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
165 # print 'Saving figure: {}'.format(figname)
166 # self.figure.savefig(figname)
167 # else :
168 # for n, eachfigure in enumerate(self.figurelist):
169 # #add specific name for each channel in channelList
170 # figname = os.path.join(self.save, '{}_{}_{}.png'.format(self.titles[n],self.CODE,
171 # datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
172 #
173 # print 'Saving figure: {}'.format(figname)
174 # eachfigure.savefig(figname)
175
176 if self.ind_plt_ch is False :
177 self.figure.canvas.draw()
178 else :
179 for eachfigure in self.figurelist:
180 eachfigure.canvas.draw()
181
182 if self.save:
183 if self.ind_plt_ch is False : #standard
184 figname = os.path.join(self.save, '{}_{}.png'.format(self.CODE,
185 datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
186 print 'Saving figure: {}'.format(figname)
187 self.figure.savefig(figname)
278 if self.save and self.data.ended:
279 channels = range(self.nrows)
280 if self.oneFigure:
281 label = ''
188 282 else :
189 for n, eachfigure in enumerate(self.figurelist):
190 #add specific name for each channel in channelList
191 figname = os.path.join(self.save, '{}_{}_{}.png'.format(self.titles[n].replace(' ',''),self.CODE,
192 datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')))
193
283 label = '_{}'.format(channels[n])
284 figname = os.path.join(
285 self.save,
286 '{}{}_{}.png'.format(
287 self.CODE,
288 label,
289 datetime.datetime.fromtimestamp(self.saveTime).strftime('%y%m%d_%H%M%S')
290 )
291 )
194 292 print 'Saving figure: {}'.format(figname)
195 eachfigure.savefig(figname)
196
293 fig.savefig(figname)
197 294
198 295 def plot(self):
199
200 print 'plotting...{}'.format(self.CODE.upper())
201 return
296 '''
297 '''
298 raise(NotImplementedError, 'Implement this method in child class')
202 299
203 300 def run(self):
204 301
205 print '[Starting] {}'.format(self.name)
302 log.success('Starting', self.name)
206 303
207 304 context = zmq.Context()
208 305 receiver = context.socket(zmq.SUB)
209 306 receiver.setsockopt(zmq.SUBSCRIBE, '')
210 307 receiver.setsockopt(zmq.CONFLATE, self.CONFLATE)
211 308
212 309 if 'server' in self.kwargs['parent']:
213 310 receiver.connect('ipc:///tmp/{}.plots'.format(self.kwargs['parent']['server']))
214 311 else:
215 312 receiver.connect("ipc:///tmp/zmq.plots")
216 313
217 seconds_passed = 0
218
219 314 while True:
220 315 try:
221 self.data = receiver.recv_pyobj(flags=zmq.NOBLOCK)#flags=zmq.NOBLOCK
222 self.started = self.data['STARTED']
223 self.dataOut = self.data['dataOut']
316 self.data = receiver.recv_pyobj(flags=zmq.NOBLOCK)
224 317
225 if (len(self.times) < len(self.data['times']) and not self.started and self.data['ENDED']):
226 continue
227
228 self.times = self.data['times']
229 self.times.sort()
230 self.throttle_value = self.data['throttle']
231 self.min_time = self.times[0]
232 self.max_time = self.times[-1]
318 self.min_time = self.data.times[0]
319 self.max_time = self.data.times[-1]
233 320
234 321 if self.isConfig is False:
235 print 'setting up'
236 self.setup()
322 self.__setup()
237 323 self.isConfig = True
238 self.__plot()
239 324
240 if self.data['ENDED'] is True:
241 print '********GRAPHIC ENDED********'
242 self.ended = True
243 self.isConfig = False
244 325 self.__plot()
245 self.deleteanotherfiles() #CLPDG
246 elif seconds_passed >= self.data['throttle']:
247 print 'passed', seconds_passed
248 self.__plot()
249 seconds_passed = 0
250 326
251 327 except zmq.Again as e:
252 print 'Waiting for data...'
253 plt.pause(2)
254 seconds_passed += 2
328 log.log('Waiting for data...')
329 if self.data:
330 plt.pause(self.data.throttle)
331 else:
332 time.sleep(2)
255 333
256 334 def close(self):
257 if self.dataOut:
335 if self.data:
258 336 self.__plot()
259 337
260 338
261 339 class PlotSpectraData(PlotData):
340 '''
341 Plot for Spectra data
342 '''
262 343
263 344 CODE = 'spc'
264 345 colormap = 'jro'
265 CONFLATE = False
266 346
267 347 def setup(self):
268
269 ncolspan = 1
270 colspan = 1
271 self.ncols = int(numpy.sqrt(self.dataOut.nChannels)+0.9)
272 self.nrows = int(self.dataOut.nChannels*1./self.ncols + 0.9)
273 self.width = 3.6*self.ncols
274 self.height = 3.2*self.nrows
348 self.nplots = len(self.data.channels)
349 self.ncols = int(numpy.sqrt(self.nplots)+ 0.9)
350 self.nrows = int((1.0*self.nplots/self.ncols) + 0.9)
351 self.width = 3.4*self.ncols
352 self.height = 3*self.nrows
353 self.cb_label = 'dB'
275 354 if self.showprofile:
276 ncolspan = 3
277 colspan = 2
278 self.width += 1.2*self.ncols
355 self.width += 0.8*self.ncols
279 356
280 357 self.ylabel = 'Range [Km]'
281 self.titles = ['Channel {}'.format(x) for x in self.dataOut.channelList]
282
283 if self.figure is None:
284 self.figure = plt.figure(figsize=(self.width, self.height),
285 edgecolor='k',
286 facecolor='w')
287 else:
288 self.figure.clf()
289
290 n = 0
291 for y in range(self.nrows):
292 for x in range(self.ncols):
293 if n >= self.dataOut.nChannels:
294 break
295 ax = plt.subplot2grid((self.nrows, self.ncols*ncolspan), (y, x*ncolspan), 1, colspan)
296 if self.showprofile:
297 ax.ax_profile = plt.subplot2grid((self.nrows, self.ncols*ncolspan), (y, x*ncolspan+colspan), 1, 1)
298
299 ax.firsttime = True
300 self.axes.append(ax)
301 n += 1
302 358
303 359 def plot(self):
304
305 360 if self.xaxis == "frequency":
306 x = self.dataOut.getFreqRange(1)/1000.
307 xlabel = "Frequency (kHz)"
361 x = self.data.xrange[0]
362 self.xlabel = "Frequency (kHz)"
308 363 elif self.xaxis == "time":
309 x = self.dataOut.getAcfRange(1)
310 xlabel = "Time (ms)"
364 x = self.data.xrange[1]
365 self.xlabel = "Time (ms)"
311 366 else:
312 x = self.dataOut.getVelRange(1)
313 xlabel = "Velocity (m/s)"
367 x = self.data.xrange[2]
368 self.xlabel = "Velocity (m/s)"
369
370 if self.CODE == 'spc_mean':
371 x = self.data.xrange[2]
372 self.xlabel = "Velocity (m/s)"
314 373
315 y = self.dataOut.getHeiRange()
316 z = self.data[self.CODE]
374 self.titles = []
375
376 y = self.data.heights
377 self.y = y
378 z = self.data['spc']
317 379
318 380 for n, ax in enumerate(self.axes):
381 noise = self.data['noise'][n][-1]
382 if self.CODE == 'spc_mean':
383 mean = self.data['mean'][n][-1]
319 384 if ax.firsttime:
320 self.xmax = self.xmax if self.xmax else np.nanmax(x)
385 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
321 386 self.xmin = self.xmin if self.xmin else -self.xmax
322 self.ymin = self.ymin if self.ymin else np.nanmin(y)
323 self.ymax = self.ymax if self.ymax else np.nanmax(y)
324 self.zmin = self.zmin if self.zmin else np.nanmin(z)
325 self.zmax = self.zmax if self.zmax else np.nanmax(z)
326 ax.plot = ax.pcolormesh(x, y, z[n].T,
387 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
388 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
389 ax.plt = ax.pcolormesh(x, y, z[n].T,
327 390 vmin=self.zmin,
328 391 vmax=self.zmax,
329 392 cmap=plt.get_cmap(self.colormap)
330 393 )
331 divider = make_axes_locatable(ax)
332 cax = divider.new_horizontal(size='3%', pad=0.05)
333 self.figure.add_axes(cax)
334 plt.colorbar(ax.plot, cax)
335
336 ax.set_xlim(self.xmin, self.xmax)
337 ax.set_ylim(self.ymin, self.ymax)
338
339 ax.set_ylabel(self.ylabel)
340 ax.set_xlabel(xlabel)
341
342 ax.firsttime = False
343 394
344 395 if self.showprofile:
345 ax.plot_profile= ax.ax_profile.plot(self.data['rti'][self.max_time][n], y)[0]
346 ax.ax_profile.set_xlim(self.zmin, self.zmax)
347 ax.ax_profile.set_ylim(self.ymin, self.ymax)
348 ax.ax_profile.set_xlabel('dB')
349 ax.ax_profile.grid(b=True, axis='x')
350 ax.plot_noise = ax.ax_profile.plot(numpy.repeat(self.data['noise'][self.max_time][n], len(y)), y,
351 color="k", linestyle="dashed", lw=2)[0]
352 [tick.set_visible(False) for tick in ax.ax_profile.get_yticklabels()]
396 ax.plt_profile= self.pf_axes[n].plot(self.data['rti'][n][-1], y)[0]
397 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
398 color="k", linestyle="dashed", lw=1)[0]
399 if self.CODE == 'spc_mean':
400 ax.plt_mean = ax.plot(mean, y, color='k')[0]
353 401 else:
354 ax.plot.set_array(z[n].T.ravel())
402 ax.plt.set_array(z[n].T.ravel())
355 403 if self.showprofile:
356 ax.plot_profile.set_data(self.data['rti'][self.max_time][n], y)
357 ax.plot_noise.set_data(numpy.repeat(self.data['noise'][self.max_time][n], len(y)), y)
404 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
405 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
406 if self.CODE == 'spc_mean':
407 ax.plt_mean.set_data(mean, y)
358 408
359 ax.set_title('{} - Noise: {:.2f} dB'.format(self.titles[n], self.data['noise'][self.max_time][n]),
360 size=8)
409 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
361 410 self.saveTime = self.max_time
362 411
363 412
364 413 class PlotCrossSpectraData(PlotData):
365 414
366 415 CODE = 'cspc'
367 416 zmin_coh = None
368 417 zmax_coh = None
369 418 zmin_phase = None
370 419 zmax_phase = None
371 CONFLATE = False
372 420
373 421 def setup(self):
374 422
375 ncolspan = 1
376 colspan = 1
377 self.ncols = 2
378 self.nrows = self.dataOut.nPairs
379 self.width = 3.6*self.ncols
380 self.height = 3.2*self.nrows
381
423 self.ncols = 4
424 self.nrows = len(self.data.pairs)
425 self.nplots = self.nrows*4
426 self.width = 3.4*self.ncols
427 self.height = 3*self.nrows
382 428 self.ylabel = 'Range [Km]'
383 self.titles = ['Channel {}'.format(x) for x in self.dataOut.channelList]
384
385 if self.figure is None:
386 self.figure = plt.figure(figsize=(self.width, self.height),
387 edgecolor='k',
388 facecolor='w')
389 else:
390 self.figure.clf()
391
392 for y in range(self.nrows):
393 for x in range(self.ncols):
394 ax = plt.subplot2grid((self.nrows, self.ncols), (y, x), 1, 1)
395 ax.firsttime = True
396 self.axes.append(ax)
429 self.showprofile = False
397 430
398 431 def plot(self):
399 432
400 433 if self.xaxis == "frequency":
401 x = self.dataOut.getFreqRange(1)/1000.
402 xlabel = "Frequency (kHz)"
434 x = self.data.xrange[0]
435 self.xlabel = "Frequency (kHz)"
403 436 elif self.xaxis == "time":
404 x = self.dataOut.getAcfRange(1)
405 xlabel = "Time (ms)"
437 x = self.data.xrange[1]
438 self.xlabel = "Time (ms)"
406 439 else:
407 x = self.dataOut.getVelRange(1)
408 xlabel = "Velocity (m/s)"
440 x = self.data.xrange[2]
441 self.xlabel = "Velocity (m/s)"
442
443 self.titles = []
409 444
410 y = self.dataOut.getHeiRange()
411 z_coh = self.data['cspc_coh']
412 z_phase = self.data['cspc_phase']
445 y = self.data.heights
446 self.y = y
447 spc = self.data['spc']
448 cspc = self.data['cspc']
413 449
414 450 for n in range(self.nrows):
415 ax = self.axes[2*n]
416 ax1 = self.axes[2*n+1]
451 noise = self.data['noise'][n][-1]
452 pair = self.data.pairs[n]
453 ax = self.axes[4*n]
454 ax3 = self.axes[4*n+3]
417 455 if ax.firsttime:
418 self.xmax = self.xmax if self.xmax else np.nanmax(x)
456 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
419 457 self.xmin = self.xmin if self.xmin else -self.xmax
420 self.ymin = self.ymin if self.ymin else np.nanmin(y)
421 self.ymax = self.ymax if self.ymax else np.nanmax(y)
422 self.zmin_coh = self.zmin_coh if self.zmin_coh else 0.0
423 self.zmax_coh = self.zmax_coh if self.zmax_coh else 1.0
424 self.zmin_phase = self.zmin_phase if self.zmin_phase else -180
425 self.zmax_phase = self.zmax_phase if self.zmax_phase else 180
426
427 ax.plot = ax.pcolormesh(x, y, z_coh[n].T,
428 vmin=self.zmin_coh,
429 vmax=self.zmax_coh,
430 cmap=plt.get_cmap(self.colormap_coh)
431 )
432 divider = make_axes_locatable(ax)
433 cax = divider.new_horizontal(size='3%', pad=0.05)
434 self.figure.add_axes(cax)
435 plt.colorbar(ax.plot, cax)
436
437 ax.set_xlim(self.xmin, self.xmax)
438 ax.set_ylim(self.ymin, self.ymax)
439
440 ax.set_ylabel(self.ylabel)
441 ax.set_xlabel(xlabel)
442 ax.firsttime = False
443
444 ax1.plot = ax1.pcolormesh(x, y, z_phase[n].T,
445 vmin=self.zmin_phase,
446 vmax=self.zmax_phase,
447 cmap=plt.get_cmap(self.colormap_phase)
458 self.zmin = self.zmin if self.zmin else numpy.nanmin(spc)
459 self.zmax = self.zmax if self.zmax else numpy.nanmax(spc)
460 ax.plt = ax.pcolormesh(x, y, spc[pair[0]].T,
461 vmin=self.zmin,
462 vmax=self.zmax,
463 cmap=plt.get_cmap(self.colormap)
448 464 )
449 divider = make_axes_locatable(ax1)
450 cax = divider.new_horizontal(size='3%', pad=0.05)
451 self.figure.add_axes(cax)
452 plt.colorbar(ax1.plot, cax)
453
454 ax1.set_xlim(self.xmin, self.xmax)
455 ax1.set_ylim(self.ymin, self.ymax)
456
457 ax1.set_ylabel(self.ylabel)
458 ax1.set_xlabel(xlabel)
459 ax1.firsttime = False
460 465 else:
461 ax.plot.set_array(z_coh[n].T.ravel())
462 ax1.plot.set_array(z_phase[n].T.ravel())
463
464 ax.set_title('Coherence Ch{} * Ch{}'.format(self.dataOut.pairsList[n][0], self.dataOut.pairsList[n][1]), size=8)
465 ax1.set_title('Phase Ch{} * Ch{}'.format(self.dataOut.pairsList[n][0], self.dataOut.pairsList[n][1]), size=8)
466 self.saveTime = self.max_time
467
468
469 class PlotSpectraMeanData(PlotSpectraData):
470
471 CODE = 'spc_mean'
472 colormap = 'jet'
473
474 def plot(self):
475
476 if self.xaxis == "frequency":
477 x = self.dataOut.getFreqRange(1)/1000.
478 xlabel = "Frequency (kHz)"
479 elif self.xaxis == "time":
480 x = self.dataOut.getAcfRange(1)
481 xlabel = "Time (ms)"
482 else:
483 x = self.dataOut.getVelRange(1)
484 xlabel = "Velocity (m/s)"
485
486 y = self.dataOut.getHeiRange()
487 z = self.data['spc']
488 mean = self.data['mean'][self.max_time]
489
490 for n, ax in enumerate(self.axes):
466 ax.plt.set_array(spc[pair[0]].T.ravel())
467 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
491 468
469 ax = self.axes[4*n+1]
492 470 if ax.firsttime:
493 self.xmax = self.xmax if self.xmax else np.nanmax(x)
494 self.xmin = self.xmin if self.xmin else -self.xmax
495 self.ymin = self.ymin if self.ymin else np.nanmin(y)
496 self.ymax = self.ymax if self.ymax else np.nanmax(y)
497 self.zmin = self.zmin if self.zmin else np.nanmin(z)
498 self.zmax = self.zmax if self.zmax else np.nanmax(z)
499 ax.plt = ax.pcolormesh(x, y, z[n].T,
471 ax.plt = ax.pcolormesh(x, y, spc[pair[1]].T,
500 472 vmin=self.zmin,
501 473 vmax=self.zmax,
502 474 cmap=plt.get_cmap(self.colormap)
503 475 )
504 ax.plt_dop = ax.plot(mean[n], y,
505 color='k')[0]
506
507 divider = make_axes_locatable(ax)
508 cax = divider.new_horizontal(size='3%', pad=0.05)
509 self.figure.add_axes(cax)
510 plt.colorbar(ax.plt, cax)
511
512 ax.set_xlim(self.xmin, self.xmax)
513 ax.set_ylim(self.ymin, self.ymax)
476 else:
477 ax.plt.set_array(spc[pair[1]].T.ravel())
478 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
514 479
515 ax.set_ylabel(self.ylabel)
516 ax.set_xlabel(xlabel)
480 out = cspc[n]/numpy.sqrt(spc[pair[0]]*spc[pair[1]])
481 coh = numpy.abs(out)
482 phase = numpy.arctan2(out.imag, out.real)*180/numpy.pi
517 483
518 ax.firsttime = False
484 ax = self.axes[4*n+2]
485 if ax.firsttime:
486 ax.plt = ax.pcolormesh(x, y, coh.T,
487 vmin=0,
488 vmax=1,
489 cmap=plt.get_cmap(self.colormap_coh)
490 )
491 else:
492 ax.plt.set_array(coh.T.ravel())
493 self.titles.append('Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
519 494
520 if self.showprofile:
521 ax.plt_profile= ax.ax_profile.plot(self.data['rti'][self.max_time][n], y)[0]
522 ax.ax_profile.set_xlim(self.zmin, self.zmax)
523 ax.ax_profile.set_ylim(self.ymin, self.ymax)
524 ax.ax_profile.set_xlabel('dB')
525 ax.ax_profile.grid(b=True, axis='x')
526 ax.plt_noise = ax.ax_profile.plot(numpy.repeat(self.data['noise'][self.max_time][n], len(y)), y,
527 color="k", linestyle="dashed", lw=2)[0]
528 [tick.set_visible(False) for tick in ax.ax_profile.get_yticklabels()]
495 ax = self.axes[4*n+3]
496 if ax.firsttime:
497 ax.plt = ax.pcolormesh(x, y, phase.T,
498 vmin=-180,
499 vmax=180,
500 cmap=plt.get_cmap(self.colormap_phase)
501 )
529 502 else:
530 ax.plt.set_array(z[n].T.ravel())
531 ax.plt_dop.set_data(mean[n], y)
532 if self.showprofile:
533 ax.plt_profile.set_data(self.data['rti'][self.max_time][n], y)
534 ax.plt_noise.set_data(numpy.repeat(self.data['noise'][self.max_time][n], len(y)), y)
503 ax.plt.set_array(phase.T.ravel())
504 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
535 505
536 ax.set_title('{} - Noise: {:.2f} dB'.format(self.titles[n], self.data['noise'][self.max_time][n]),
537 size=8)
538 506 self.saveTime = self.max_time
539 507
540 508
509 class PlotSpectraMeanData(PlotSpectraData):
510 '''
511 Plot for Spectra and Mean
512 '''
513 CODE = 'spc_mean'
514 colormap = 'jro'
515
516
541 517 class PlotRTIData(PlotData):
518 '''
519 Plot for RTI data
520 '''
542 521
543 522 CODE = 'rti'
544 523 colormap = 'jro'
545 524
546 525 def setup(self):
526 self.xaxis = 'time'
547 527 self.ncols = 1
548 self.nrows = self.dataOut.nChannels
549 self.width = 10
550 #TODO : arreglar la altura de la figura, esta hardcodeada.
551 #Se arreglo, testear!
552 if self.ind_plt_ch:
553 self.height = 3.2#*self.nrows if self.nrows<6 else 12
554 else:
555 self.height = 2.2*self.nrows if self.nrows<6 else 12
556
557 '''
558 if self.nrows==1:
559 self.height += 1
560 '''
528 self.nrows = len(self.data.channels)
529 self.nplots = len(self.data.channels)
561 530 self.ylabel = 'Range [Km]'
562 self.titles = ['Channel {}'.format(x) for x in self.dataOut.channelList]
563
564 '''
565 Logica:
566 1) Si la variable ind_plt_ch es True, va a crear mas de 1 figura
567 2) guardamos "Figures" en una lista y "axes" en otra, quizas se deberia guardar el
568 axis dentro de "Figures" como un diccionario.
569 '''
570 if self.ind_plt_ch is False: #standard mode
571
572 if self.figure is None: #solo para la priemra vez
573 self.figure = plt.figure(figsize=(self.width, self.height),
574 edgecolor='k',
575 facecolor='w')
576 else:
577 self.figure.clf()
578 self.axes = []
579
580
581 for n in range(self.nrows):
582 ax = self.figure.add_subplot(self.nrows, self.ncols, n+1)
583 #ax = self.figure(n+1)
584 ax.firsttime = True
585 self.axes.append(ax)
586
587 else : #append one figure foreach channel in channelList
588 if self.figurelist == None:
589 self.figurelist = []
590 for n in range(self.nrows):
591 self.figure = plt.figure(figsize=(self.width, self.height),
592 edgecolor='k',
593 facecolor='w')
594 #add always one subplot
595 self.figurelist.append(self.figure)
596
597 else : # cada dia nuevo limpia el axes, pero mantiene el figure
598 for eachfigure in self.figurelist:
599 eachfigure.clf() # eliminaria todas las figuras de la lista?
600 self.axes = []
601
602 for eachfigure in self.figurelist:
603 ax = eachfigure.add_subplot(1,1,1) #solo 1 axis por figura
604 #ax = self.figure(n+1)
605 ax.firsttime = True
606 #Cada figura tiene un distinto puntero
607 self.axes.append(ax)
608 #plt.close(eachfigure)
609
531 self.cb_label = 'dB'
532 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
610 533
611 534 def plot(self):
535 self.x = self.data.times
536 self.y = self.data.heights
537 self.z = self.data[self.CODE]
538 self.z = numpy.ma.masked_invalid(self.z)
612 539
613 if self.ind_plt_ch is False: #standard mode
614 self.x = np.array(self.times)
615 self.y = self.dataOut.getHeiRange()
616 self.z = []
617
618 for ch in range(self.nrows):
619 self.z.append([self.data[self.CODE][t][ch] for t in self.times])
620
621 self.z = np.array(self.z)
622 540 for n, ax in enumerate(self.axes):
623 541 x, y, z = self.fill_gaps(*self.decimate())
624 if self.xmin is None:
625 xmin = self.min_time
626 else:
627 xmin = fromtimestamp(int(self.xmin), self.min_time)
628 if self.xmax is None:
629 xmax = xmin + self.xrange*60*60
630 else:
631 xmax = xmin + (self.xmax - self.xmin) * 60 * 60
632 self.zmin = self.zmin if self.zmin else np.min(self.z)
633 self.zmax = self.zmax if self.zmax else np.max(self.z)
542 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
543 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
634 544 if ax.firsttime:
635 self.ymin = self.ymin if self.ymin else np.nanmin(self.y)
636 self.ymax = self.ymax if self.ymax else np.nanmax(self.y)
637 plot = ax.pcolormesh(x, y, z[n].T,
545 ax.plt = ax.pcolormesh(x, y, z[n].T,
638 546 vmin=self.zmin,
639 547 vmax=self.zmax,
640 548 cmap=plt.get_cmap(self.colormap)
641 549 )
642 divider = make_axes_locatable(ax)
643 cax = divider.new_horizontal(size='2%', pad=0.05)
644 self.figure.add_axes(cax)
645 plt.colorbar(plot, cax)
646 ax.set_ylim(self.ymin, self.ymax)
647 ax.xaxis.set_major_formatter(FuncFormatter(func))
648 ax.xaxis.set_major_locator(LinearLocator(6))
649 ax.set_ylabel(self.ylabel)
650 # if self.xmin is None:
651 # xmin = self.min_time
652 # else:
653 # xmin = (datetime.datetime.combine(self.dataOut.datatime.date(),
654 # datetime.time(self.xmin, 0, 0))-d1970).total_seconds()
655
656 ax.set_xlim(xmin, xmax)
657 ax.firsttime = False
550 if self.showprofile:
551 ax.plot_profile= self.pf_axes[n].plot(self.data['rti'][n][-1], self.y)[0]
552 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
553 color="k", linestyle="dashed", lw=1)[0]
658 554 else:
659 555 ax.collections.remove(ax.collections[0])
660 ax.set_xlim(xmin, xmax)
661 plot = ax.pcolormesh(x, y, z[n].T,
662 vmin=self.zmin,
663 vmax=self.zmax,
664 cmap=plt.get_cmap(self.colormap)
665 )
666 ax.set_title('{} {}'.format(self.titles[n],
667 datetime.datetime.fromtimestamp(self.max_time).strftime('%y/%m/%d %H:%M:%S')),
668 size=8)
669
670 self.saveTime = self.min_time
671 else :
672 self.x = np.array(self.times)
673 self.y = self.dataOut.getHeiRange()
674 self.z = []
675
676 for ch in range(self.nrows):
677 self.z.append([self.data[self.CODE][t][ch] for t in self.times])
678
679 self.z = np.array(self.z)
680 for n, eachfigure in enumerate(self.figurelist): #estaba ax in axes
681
682 x, y, z = self.fill_gaps(*self.decimate())
683 xmin = self.min_time
684 xmax = xmin+self.xrange*60*60
685 self.zmin = self.zmin if self.zmin else np.min(self.z)
686 self.zmax = self.zmax if self.zmax else np.max(self.z)
687 if self.axes[n].firsttime:
688 self.ymin = self.ymin if self.ymin else np.nanmin(self.y)
689 self.ymax = self.ymax if self.ymax else np.nanmax(self.y)
690 plot = self.axes[n].pcolormesh(x, y, z[n].T,
691 vmin=self.zmin,
692 vmax=self.zmax,
693 cmap=plt.get_cmap(self.colormap)
694 )
695 divider = make_axes_locatable(self.axes[n])
696 cax = divider.new_horizontal(size='2%', pad=0.05)
697 eachfigure.add_axes(cax)
698 #self.figure2.add_axes(cax)
699 plt.colorbar(plot, cax)
700 self.axes[n].set_ylim(self.ymin, self.ymax)
701
702 self.axes[n].xaxis.set_major_formatter(FuncFormatter(func))
703 self.axes[n].xaxis.set_major_locator(LinearLocator(6))
704
705 self.axes[n].set_ylabel(self.ylabel)
706
707 if self.xmin is None:
708 xmin = self.min_time
709 else:
710 xmin = (datetime.datetime.combine(self.dataOut.datatime.date(),
711 datetime.time(self.xmin, 0, 0))-d1970).total_seconds()
712
713 self.axes[n].set_xlim(xmin, xmax)
714 self.axes[n].firsttime = False
715 else:
716 self.axes[n].collections.remove(self.axes[n].collections[0])
717 self.axes[n].set_xlim(xmin, xmax)
718 plot = self.axes[n].pcolormesh(x, y, z[n].T,
556 ax.plt = ax.pcolormesh(x, y, z[n].T,
719 557 vmin=self.zmin,
720 558 vmax=self.zmax,
721 559 cmap=plt.get_cmap(self.colormap)
722 560 )
723 self.axes[n].set_title('{} {}'.format(self.titles[n],
724 datetime.datetime.fromtimestamp(self.max_time).strftime('%y/%m/%d %H:%M:%S')),
725 size=8)
561 if self.showprofile:
562 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
563 ax.plot_noise.set_data(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y)
726 564
727 565 self.saveTime = self.min_time
728 566
729 567
730 568 class PlotCOHData(PlotRTIData):
569 '''
570 Plot for Coherence data
571 '''
731 572
732 573 CODE = 'coh'
733 574
734 575 def setup(self):
735
576 self.xaxis = 'time'
736 577 self.ncols = 1
737 self.nrows = self.dataOut.nPairs
738 self.width = 10
739 self.height = 2.2*self.nrows if self.nrows<6 else 12
740 self.ind_plt_ch = False #just for coherence and phase
741 if self.nrows==1:
742 self.height += 1
578 self.nrows = len(self.data.pairs)
579 self.nplots = len(self.data.pairs)
743 580 self.ylabel = 'Range [Km]'
744 self.titles = ['{} Ch{} * Ch{}'.format(self.CODE.upper(), x[0], x[1]) for x in self.dataOut.pairsList]
745
746 if self.figure is None:
747 self.figure = plt.figure(figsize=(self.width, self.height),
748 edgecolor='k',
749 facecolor='w')
581 if self.CODE == 'coh':
582 self.cb_label = ''
583 self.titles = ['Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
750 584 else:
751 self.figure.clf()
752 self.axes = []
585 self.cb_label = 'Degrees'
586 self.titles = ['Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
753 587
754 for n in range(self.nrows):
755 ax = self.figure.add_subplot(self.nrows, self.ncols, n+1)
756 ax.firsttime = True
757 self.axes.append(ax)
588
589 class PlotPHASEData(PlotCOHData):
590 '''
591 Plot for Phase map data
592 '''
593
594 CODE = 'phase'
595 colormap = 'seismic'
758 596
759 597
760 598 class PlotNoiseData(PlotData):
599 '''
600 Plot for noise
601 '''
602
761 603 CODE = 'noise'
762 604
763 605 def setup(self):
764
606 self.xaxis = 'time'
765 607 self.ncols = 1
766 608 self.nrows = 1
767 self.width = 10
768 self.height = 3.2
609 self.nplots = 1
769 610 self.ylabel = 'Intensity [dB]'
770 611 self.titles = ['Noise']
771
772 if self.figure is None:
773 self.figure = plt.figure(figsize=(self.width, self.height),
774 edgecolor='k',
775 facecolor='w')
776 else:
777 self.figure.clf()
778 self.axes = []
779
780 self.ax = self.figure.add_subplot(self.nrows, self.ncols, 1)
781 self.ax.firsttime = True
782
783 def plot(self):
784
785 x = self.times
786 xmin = self.min_time
787 xmax = xmin+self.xrange*60*60
788 if self.ax.firsttime:
789 for ch in self.dataOut.channelList:
790 y = [self.data[self.CODE][t][ch] for t in self.times]
791 self.ax.plot(x, y, lw=1, label='Ch{}'.format(ch))
792 self.ax.firsttime = False
793 self.ax.xaxis.set_major_formatter(FuncFormatter(func))
794 self.ax.xaxis.set_major_locator(LinearLocator(6))
795 self.ax.set_ylabel(self.ylabel)
796 plt.legend()
797 else:
798 for ch in self.dataOut.channelList:
799 y = [self.data[self.CODE][t][ch] for t in self.times]
800 self.ax.lines[ch].set_data(x, y)
801
802 self.ax.set_xlim(xmin, xmax)
803 self.ax.set_ylim(min(y)-5, max(y)+5)
804 self.saveTime = self.min_time
805
806
807 class PlotWindProfilerData(PlotRTIData):
808
809 CODE = 'wind'
810 colormap = 'seismic'
811
812 def setup(self):
813 self.ncols = 1
814 self.nrows = self.dataOut.data_output.shape[0]
815 self.width = 10
816 self.height = 2.2*self.nrows
817 self.ylabel = 'Height [Km]'
818 self.titles = ['Zonal Wind' ,'Meridional Wind', 'Vertical Wind']
819 self.clabels = ['Velocity (m/s)','Velocity (m/s)','Velocity (cm/s)']
820 self.windFactor = [1, 1, 100]
821
822 if self.figure is None:
823 self.figure = plt.figure(figsize=(self.width, self.height),
824 edgecolor='k',
825 facecolor='w')
826 else:
827 self.figure.clf()
828 self.axes = []
829
830 for n in range(self.nrows):
831 ax = self.figure.add_subplot(self.nrows, self.ncols, n+1)
832 ax.firsttime = True
833 self.axes.append(ax)
612 self.colorbar = False
834 613
835 614 def plot(self):
836 615
837 self.x = np.array(self.times)
838 self.y = self.dataOut.heightList
839 self.z = []
840
841 for ch in range(self.nrows):
842 self.z.append([self.data['output'][t][ch] for t in self.times])
843
844 self.z = np.array(self.z)
845 self.z = numpy.ma.masked_invalid(self.z)
846
847 cmap=plt.get_cmap(self.colormap)
848 cmap.set_bad('black', 1.)
849
850 for n, ax in enumerate(self.axes):
851 x, y, z = self.fill_gaps(*self.decimate())
616 x = self.data.times
852 617 xmin = self.min_time
853 618 xmax = xmin+self.xrange*60*60
854 if ax.firsttime:
855 self.ymin = self.ymin if self.ymin else np.nanmin(self.y)
856 self.ymax = self.ymax if self.ymax else np.nanmax(self.y)
857 self.zmax = self.zmax if self.zmax else numpy.nanmax(abs(self.z[:-1, :]))
858 self.zmin = self.zmin if self.zmin else -self.zmax
859
860 plot = ax.pcolormesh(x, y, z[n].T*self.windFactor[n],
861 vmin=self.zmin,
862 vmax=self.zmax,
863 cmap=cmap
864 )
865 divider = make_axes_locatable(ax)
866 cax = divider.new_horizontal(size='2%', pad=0.05)
867 self.figure.add_axes(cax)
868 cb = plt.colorbar(plot, cax)
869 cb.set_label(self.clabels[n])
870 ax.set_ylim(self.ymin, self.ymax)
619 Y = self.data[self.CODE]
871 620
872 ax.xaxis.set_major_formatter(FuncFormatter(func))
873 ax.xaxis.set_major_locator(LinearLocator(6))
874
875 ax.set_ylabel(self.ylabel)
876
877 ax.set_xlim(xmin, xmax)
878 ax.firsttime = False
621 if self.axes[0].firsttime:
622 for ch in self.data.channels:
623 y = Y[ch]
624 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
625 plt.legend()
879 626 else:
880 ax.collections.remove(ax.collections[0])
881 ax.set_xlim(xmin, xmax)
882 plot = ax.pcolormesh(x, y, z[n].T*self.windFactor[n],
883 vmin=self.zmin,
884 vmax=self.zmax,
885 cmap=plt.get_cmap(self.colormap)
886 )
887 ax.set_title('{} {}'.format(self.titles[n],
888 datetime.datetime.fromtimestamp(self.max_time).strftime('%y/%m/%d %H:%M:%S')),
889 size=8)
627 for ch in self.data.channels:
628 y = Y[ch]
629 self.axes[0].lines[ch].set_data(x, y)
890 630
631 self.ymin = numpy.nanmin(Y) - 5
632 self.ymax = numpy.nanmax(Y) + 5
891 633 self.saveTime = self.min_time
892 634
893 635
894 636 class PlotSNRData(PlotRTIData):
637 '''
638 Plot for SNR Data
639 '''
640
895 641 CODE = 'snr'
896 642 colormap = 'jet'
897 643
644
898 645 class PlotDOPData(PlotRTIData):
646 '''
647 Plot for DOPPLER Data
648 '''
649
899 650 CODE = 'dop'
900 651 colormap = 'jet'
901 652
902 653
903 class PlotPHASEData(PlotCOHData):
904 CODE = 'phase'
905 colormap = 'seismic'
906
907
908 654 class PlotSkyMapData(PlotData):
655 '''
656 Plot for meteors detection data
657 '''
909 658
910 659 CODE = 'met'
911 660
912 661 def setup(self):
913 662
914 663 self.ncols = 1
915 664 self.nrows = 1
916 665 self.width = 7.2
917 666 self.height = 7.2
918 667
919 668 self.xlabel = 'Zonal Zenith Angle (deg)'
920 669 self.ylabel = 'Meridional Zenith Angle (deg)'
921 670
922 671 if self.figure is None:
923 672 self.figure = plt.figure(figsize=(self.width, self.height),
924 673 edgecolor='k',
925 674 facecolor='w')
926 675 else:
927 676 self.figure.clf()
928 677
929 678 self.ax = plt.subplot2grid((self.nrows, self.ncols), (0, 0), 1, 1, polar=True)
930 679 self.ax.firsttime = True
931 680
932 681
933 682 def plot(self):
934 683
935 arrayParameters = np.concatenate([self.data['param'][t] for t in self.times])
684 arrayParameters = numpy.concatenate([self.data['param'][t] for t in self.data.times])
936 685 error = arrayParameters[:,-1]
937 686 indValid = numpy.where(error == 0)[0]
938 687 finalMeteor = arrayParameters[indValid,:]
939 688 finalAzimuth = finalMeteor[:,3]
940 689 finalZenith = finalMeteor[:,4]
941 690
942 691 x = finalAzimuth*numpy.pi/180
943 692 y = finalZenith
944 693
945 694 if self.ax.firsttime:
946 695 self.ax.plot = self.ax.plot(x, y, 'bo', markersize=5)[0]
947 696 self.ax.set_ylim(0,90)
948 697 self.ax.set_yticks(numpy.arange(0,90,20))
949 698 self.ax.set_xlabel(self.xlabel)
950 699 self.ax.set_ylabel(self.ylabel)
951 700 self.ax.yaxis.labelpad = 40
952 701 self.ax.firsttime = False
953 702 else:
954 703 self.ax.plot.set_data(x, y)
955 704
956 705
957 706 dt1 = datetime.datetime.fromtimestamp(self.min_time).strftime('%y/%m/%d %H:%M:%S')
958 707 dt2 = datetime.datetime.fromtimestamp(self.max_time).strftime('%y/%m/%d %H:%M:%S')
959 708 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
960 709 dt2,
961 710 len(x))
962 711 self.ax.set_title(title, size=8)
963 712
964 713 self.saveTime = self.max_time
714
715 class PlotParamData(PlotRTIData):
716 '''
717 Plot for data_param object
718 '''
719
720 CODE = 'param'
721 colormap = 'seismic'
722
723 def setup(self):
724 self.xaxis = 'time'
725 self.ncols = 1
726 self.nrows = self.data.shape(self.CODE)[0]
727 self.nplots = self.nrows
728 if self.showSNR:
729 self.nrows += 1
730
731 self.ylabel = 'Height [Km]'
732 self.titles = self.data.parameters \
733 if self.data.parameters else ['Param {}'.format(x) for x in xrange(self.nrows)]
734 if self.showSNR:
735 self.titles.append('SNR')
736
737 def plot(self):
738 self.data.normalize_heights()
739 self.x = self.data.times
740 self.y = self.data.heights
741 if self.showSNR:
742 self.z = numpy.concatenate(
743 (self.data[self.CODE], self.data['snr'])
744 )
745 else:
746 self.z = self.data[self.CODE]
747
748 self.z = numpy.ma.masked_invalid(self.z)
749
750 for n, ax in enumerate(self.axes):
751
752 x, y, z = self.fill_gaps(*self.decimate())
753
754 if ax.firsttime:
755 if self.zlimits is not None:
756 self.zmin, self.zmax = self.zlimits[n]
757 self.zmax = self.zmax if self.zmax is not None else numpy.nanmax(abs(self.z[:-1, :]))
758 self.zmin = self.zmin if self.zmin is not None else -self.zmax
759 ax.plt = ax.pcolormesh(x, y, z[n, :, :].T*self.factors[n],
760 vmin=self.zmin,
761 vmax=self.zmax,
762 cmap=self.cmaps[n]
763 )
764 else:
765 if self.zlimits is not None:
766 self.zmin, self.zmax = self.zlimits[n]
767 ax.collections.remove(ax.collections[0])
768 ax.plt = ax.pcolormesh(x, y, z[n, :, :].T*self.factors[n],
769 vmin=self.zmin,
770 vmax=self.zmax,
771 cmap=self.cmaps[n]
772 )
773
774 self.saveTime = self.min_time
775
776 class PlotOuputData(PlotParamData):
777 '''
778 Plot data_output object
779 '''
780
781 CODE = 'output'
782 colormap = 'seismic' No newline at end of file
Show file before | Show file after | Hide comments
1 NO CONTENT: modified file
The requested commit or file is too big and content was truncated. Show full diff
Show file before | Show file after | Hide comments
@@ -1,904 +1,903
1 import itertools
2
1 3 import numpy
2 4
3 5 from jroproc_base import ProcessingUnit, Operation
4 6 from schainpy.model.data.jrodata import Spectra
5 7 from schainpy.model.data.jrodata import hildebrand_sekhon
6 8
7 9 class SpectraProc(ProcessingUnit):
8 10
9 11 def __init__(self, **kwargs):
10 12
11 13 ProcessingUnit.__init__(self, **kwargs)
12 14
13 15 self.buffer = None
14 16 self.firstdatatime = None
15 17 self.profIndex = 0
16 18 self.dataOut = Spectra()
17 19 self.id_min = None
18 20 self.id_max = None
19 21
20 22 def __updateSpecFromVoltage(self):
21 23
22 24 self.dataOut.timeZone = self.dataIn.timeZone
23 25 self.dataOut.dstFlag = self.dataIn.dstFlag
24 26 self.dataOut.errorCount = self.dataIn.errorCount
25 27 self.dataOut.useLocalTime = self.dataIn.useLocalTime
26 28
27 29 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
28 30 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
29 31 self.dataOut.channelList = self.dataIn.channelList
30 32 self.dataOut.heightList = self.dataIn.heightList
31 33 self.dataOut.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
32 34
33 35 self.dataOut.nBaud = self.dataIn.nBaud
34 36 self.dataOut.nCode = self.dataIn.nCode
35 37 self.dataOut.code = self.dataIn.code
36 38 self.dataOut.nProfiles = self.dataOut.nFFTPoints
37 39
38 40 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
39 41 self.dataOut.utctime = self.firstdatatime
40 42 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
41 43 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
42 44 self.dataOut.flagShiftFFT = False
43 45
44 46 self.dataOut.nCohInt = self.dataIn.nCohInt
45 47 self.dataOut.nIncohInt = 1
46 48
47 49 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
48 50
49 51 self.dataOut.frequency = self.dataIn.frequency
50 52 self.dataOut.realtime = self.dataIn.realtime
51 53
52 54 self.dataOut.azimuth = self.dataIn.azimuth
53 55 self.dataOut.zenith = self.dataIn.zenith
54 56
55 57 self.dataOut.beam.codeList = self.dataIn.beam.codeList
56 58 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
57 59 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
58 60
59 61 def __getFft(self):
60 62 """
61 63 Convierte valores de Voltaje a Spectra
62 64
63 65 Affected:
64 66 self.dataOut.data_spc
65 67 self.dataOut.data_cspc
66 68 self.dataOut.data_dc
67 69 self.dataOut.heightList
68 70 self.profIndex
69 71 self.buffer
70 72 self.dataOut.flagNoData
71 73 """
72 74 fft_volt = numpy.fft.fft(self.buffer,n=self.dataOut.nFFTPoints,axis=1)
73 75 fft_volt = fft_volt.astype(numpy.dtype('complex'))
74 76 dc = fft_volt[:,0,:]
75 77
76 78 #calculo de self-spectra
77 79 fft_volt = numpy.fft.fftshift(fft_volt,axes=(1,))
78 80 spc = fft_volt * numpy.conjugate(fft_volt)
79 81 spc = spc.real
80 82
81 83 blocksize = 0
82 84 blocksize += dc.size
83 85 blocksize += spc.size
84 86
85 87 cspc = None
86 88 pairIndex = 0
87 89 if self.dataOut.pairsList != None:
88 90 #calculo de cross-spectra
89 91 cspc = numpy.zeros((self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
90 92 for pair in self.dataOut.pairsList:
91 93 if pair[0] not in self.dataOut.channelList:
92 94 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
93 95 if pair[1] not in self.dataOut.channelList:
94 96 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
95 97
96 98 cspc[pairIndex,:,:] = fft_volt[pair[0],:,:] * numpy.conjugate(fft_volt[pair[1],:,:])
97 99 pairIndex += 1
98 100 blocksize += cspc.size
99 101
100 102 self.dataOut.data_spc = spc
101 103 self.dataOut.data_cspc = cspc
102 104 self.dataOut.data_dc = dc
103 105 self.dataOut.blockSize = blocksize
104 106 self.dataOut.flagShiftFFT = True
105 107
106 108 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], ippFactor=None):
107 109
108 110 self.dataOut.flagNoData = True
109 111
110 112 if self.dataIn.type == "Spectra":
111 113 self.dataOut.copy(self.dataIn)
112 # self.__selectPairs(pairsList)
114 if not pairsList:
115 pairsList = itertools.combinations(self.dataOut.channelList, 2)
116 if self.dataOut.data_cspc is not None:
117 self.__selectPairs(pairsList)
113 118 return True
114 119
115 120 if self.dataIn.type == "Voltage":
116 121
117 122 if nFFTPoints == None:
118 123 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
119 124
120 125 if nProfiles == None:
121 126 nProfiles = nFFTPoints
122 127
123 128 if ippFactor == None:
124 129 ippFactor = 1
125 130
126 131 self.dataOut.ippFactor = ippFactor
127 132
128 133 self.dataOut.nFFTPoints = nFFTPoints
129 134 self.dataOut.pairsList = pairsList
130 135
131 136 if self.buffer is None:
132 137 self.buffer = numpy.zeros( (self.dataIn.nChannels,
133 138 nProfiles,
134 139 self.dataIn.nHeights),
135 140 dtype='complex')
136 141
137 142 if self.dataIn.flagDataAsBlock:
138 143 #data dimension: [nChannels, nProfiles, nSamples]
139 144 nVoltProfiles = self.dataIn.data.shape[1]
140 145 # nVoltProfiles = self.dataIn.nProfiles
141 146
142 147 if nVoltProfiles == nProfiles:
143 148 self.buffer = self.dataIn.data.copy()
144 149 self.profIndex = nVoltProfiles
145 150
146 151 elif nVoltProfiles < nProfiles:
147 152
148 153 if self.profIndex == 0:
149 154 self.id_min = 0
150 155 self.id_max = nVoltProfiles
151 156
152 157 self.buffer[:,self.id_min:self.id_max,:] = self.dataIn.data
153 158 self.profIndex += nVoltProfiles
154 159 self.id_min += nVoltProfiles
155 160 self.id_max += nVoltProfiles
156 161 else:
157 162 raise ValueError, "The type object %s has %d profiles, it should just has %d profiles"%(self.dataIn.type,self.dataIn.data.shape[1],nProfiles)
158 163 self.dataOut.flagNoData = True
159 164 return 0
160 165 else:
161 166 self.buffer[:,self.profIndex,:] = self.dataIn.data.copy()
162 167 self.profIndex += 1
163 168
164 169 if self.firstdatatime == None:
165 170 self.firstdatatime = self.dataIn.utctime
166 171
167 172 if self.profIndex == nProfiles:
168 173 self.__updateSpecFromVoltage()
169 174 self.__getFft()
170 175
171 176 self.dataOut.flagNoData = False
172 177 self.firstdatatime = None
173 178 self.profIndex = 0
174 179
175 180 return True
176 181
177 182 raise ValueError, "The type of input object '%s' is not valid"%(self.dataIn.type)
178 183
179 184 def __selectPairs(self, pairsList):
180 185
181 if channelList == None:
186 if not pairsList:
182 187 return
183 188
184 pairsIndexListSelected = []
185
186 for thisPair in pairsList:
189 pairs = []
190 pairsIndex = []
187 191
188 if thisPair not in self.dataOut.pairsList:
192 for pair in pairsList:
193 if pair[0] not in self.dataOut.channelList or pair[1] not in self.dataOut.channelList:
189 194 continue
195 pairs.append(pair)
196 pairsIndex.append(pairs.index(pair))
190 197
191 pairIndex = self.dataOut.pairsList.index(thisPair)
192
193 pairsIndexListSelected.append(pairIndex)
194
195 if not pairsIndexListSelected:
196 self.dataOut.data_cspc = None
197 self.dataOut.pairsList = []
198 return
199
200 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
201 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
198 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndex]
199 self.dataOut.pairsList = pairs
200 self.dataOut.pairsIndexList = pairsIndex
202 201
203 202 return
204 203
205 204 def __selectPairsByChannel(self, channelList=None):
206 205
207 206 if channelList == None:
208 207 return
209 208
210 209 pairsIndexListSelected = []
211 210 for pairIndex in self.dataOut.pairsIndexList:
212 211 #First pair
213 212 if self.dataOut.pairsList[pairIndex][0] not in channelList:
214 213 continue
215 214 #Second pair
216 215 if self.dataOut.pairsList[pairIndex][1] not in channelList:
217 216 continue
218 217
219 218 pairsIndexListSelected.append(pairIndex)
220 219
221 220 if not pairsIndexListSelected:
222 221 self.dataOut.data_cspc = None
223 222 self.dataOut.pairsList = []
224 223 return
225 224
226 225 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
227 226 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
228 227
229 228 return
230 229
231 230 def selectChannels(self, channelList):
232 231
233 232 channelIndexList = []
234 233
235 234 for channel in channelList:
236 235 if channel not in self.dataOut.channelList:
237 236 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" %(channel, str(self.dataOut.channelList))
238 237
239 238 index = self.dataOut.channelList.index(channel)
240 239 channelIndexList.append(index)
241 240
242 241 self.selectChannelsByIndex(channelIndexList)
243 242
244 243 def selectChannelsByIndex(self, channelIndexList):
245 244 """
246 245 Selecciona un bloque de datos en base a canales segun el channelIndexList
247 246
248 247 Input:
249 248 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
250 249
251 250 Affected:
252 251 self.dataOut.data_spc
253 252 self.dataOut.channelIndexList
254 253 self.dataOut.nChannels
255 254
256 255 Return:
257 256 None
258 257 """
259 258
260 259 for channelIndex in channelIndexList:
261 260 if channelIndex not in self.dataOut.channelIndexList:
262 261 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " %(channelIndex, self.dataOut.channelIndexList)
263 262
264 263 # nChannels = len(channelIndexList)
265 264
266 265 data_spc = self.dataOut.data_spc[channelIndexList,:]
267 266 data_dc = self.dataOut.data_dc[channelIndexList,:]
268 267
269 268 self.dataOut.data_spc = data_spc
270 269 self.dataOut.data_dc = data_dc
271 270
272 271 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
273 272 # self.dataOut.nChannels = nChannels
274 273
275 274 self.__selectPairsByChannel(self.dataOut.channelList)
276 275
277 276 return 1
278 277
279 278 def selectHeights(self, minHei, maxHei):
280 279 """
281 280 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
282 281 minHei <= height <= maxHei
283 282
284 283 Input:
285 284 minHei : valor minimo de altura a considerar
286 285 maxHei : valor maximo de altura a considerar
287 286
288 287 Affected:
289 288 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
290 289
291 290 Return:
292 291 1 si el metodo se ejecuto con exito caso contrario devuelve 0
293 292 """
294 293
295 294 if (minHei > maxHei):
296 295 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (minHei, maxHei)
297 296
298 297 if (minHei < self.dataOut.heightList[0]):
299 298 minHei = self.dataOut.heightList[0]
300 299
301 300 if (maxHei > self.dataOut.heightList[-1]):
302 301 maxHei = self.dataOut.heightList[-1]
303 302
304 303 minIndex = 0
305 304 maxIndex = 0
306 305 heights = self.dataOut.heightList
307 306
308 307 inda = numpy.where(heights >= minHei)
309 308 indb = numpy.where(heights <= maxHei)
310 309
311 310 try:
312 311 minIndex = inda[0][0]
313 312 except:
314 313 minIndex = 0
315 314
316 315 try:
317 316 maxIndex = indb[0][-1]
318 317 except:
319 318 maxIndex = len(heights)
320 319
321 320 self.selectHeightsByIndex(minIndex, maxIndex)
322 321
323 322 return 1
324 323
325 324 def getBeaconSignal(self, tauindex = 0, channelindex = 0, hei_ref=None):
326 325 newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
327 326
328 327 if hei_ref != None:
329 328 newheis = numpy.where(self.dataOut.heightList>hei_ref)
330 329
331 330 minIndex = min(newheis[0])
332 331 maxIndex = max(newheis[0])
333 332 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
334 333 heightList = self.dataOut.heightList[minIndex:maxIndex+1]
335 334
336 335 # determina indices
337 336 nheis = int(self.dataOut.radarControllerHeaderObj.txB/(self.dataOut.heightList[1]-self.dataOut.heightList[0]))
338 337 avg_dB = 10*numpy.log10(numpy.sum(data_spc[channelindex,:,:],axis=0))
339 338 beacon_dB = numpy.sort(avg_dB)[-nheis:]
340 339 beacon_heiIndexList = []
341 340 for val in avg_dB.tolist():
342 341 if val >= beacon_dB[0]:
343 342 beacon_heiIndexList.append(avg_dB.tolist().index(val))
344 343
345 344 #data_spc = data_spc[:,:,beacon_heiIndexList]
346 345 data_cspc = None
347 346 if self.dataOut.data_cspc is not None:
348 347 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
349 348 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
350 349
351 350 data_dc = None
352 351 if self.dataOut.data_dc is not None:
353 352 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
354 353 #data_dc = data_dc[:,beacon_heiIndexList]
355 354
356 355 self.dataOut.data_spc = data_spc
357 356 self.dataOut.data_cspc = data_cspc
358 357 self.dataOut.data_dc = data_dc
359 358 self.dataOut.heightList = heightList
360 359 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
361 360
362 361 return 1
363 362
364 363
365 364 def selectHeightsByIndex(self, minIndex, maxIndex):
366 365 """
367 366 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
368 367 minIndex <= index <= maxIndex
369 368
370 369 Input:
371 370 minIndex : valor de indice minimo de altura a considerar
372 371 maxIndex : valor de indice maximo de altura a considerar
373 372
374 373 Affected:
375 374 self.dataOut.data_spc
376 375 self.dataOut.data_cspc
377 376 self.dataOut.data_dc
378 377 self.dataOut.heightList
379 378
380 379 Return:
381 380 1 si el metodo se ejecuto con exito caso contrario devuelve 0
382 381 """
383 382
384 383 if (minIndex < 0) or (minIndex > maxIndex):
385 384 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex)
386 385
387 386 if (maxIndex >= self.dataOut.nHeights):
388 387 maxIndex = self.dataOut.nHeights-1
389 388
390 389 #Spectra
391 390 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
392 391
393 392 data_cspc = None
394 393 if self.dataOut.data_cspc is not None:
395 394 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
396 395
397 396 data_dc = None
398 397 if self.dataOut.data_dc is not None:
399 398 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
400 399
401 400 self.dataOut.data_spc = data_spc
402 401 self.dataOut.data_cspc = data_cspc
403 402 self.dataOut.data_dc = data_dc
404 403
405 404 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex+1]
406 405
407 406 return 1
408 407
409 408 def removeDC(self, mode = 2):
410 409 jspectra = self.dataOut.data_spc
411 410 jcspectra = self.dataOut.data_cspc
412 411
413 412
414 413 num_chan = jspectra.shape[0]
415 414 num_hei = jspectra.shape[2]
416 415
417 416 if jcspectra is not None:
418 417 jcspectraExist = True
419 418 num_pairs = jcspectra.shape[0]
420 419 else: jcspectraExist = False
421 420
422 421 freq_dc = jspectra.shape[1]/2
423 422 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
424 423
425 424 if ind_vel[0]<0:
426 425 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
427 426
428 427 if mode == 1:
429 428 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
430 429
431 430 if jcspectraExist:
432 431 jcspectra[:,freq_dc,:] = (jcspectra[:,ind_vel[1],:] + jcspectra[:,ind_vel[2],:])/2
433 432
434 433 if mode == 2:
435 434
436 435 vel = numpy.array([-2,-1,1,2])
437 436 xx = numpy.zeros([4,4])
438 437
439 438 for fil in range(4):
440 439 xx[fil,:] = vel[fil]**numpy.asarray(range(4))
441 440
442 441 xx_inv = numpy.linalg.inv(xx)
443 442 xx_aux = xx_inv[0,:]
444 443
445 444 for ich in range(num_chan):
446 445 yy = jspectra[ich,ind_vel,:]
447 446 jspectra[ich,freq_dc,:] = numpy.dot(xx_aux,yy)
448 447
449 448 junkid = jspectra[ich,freq_dc,:]<=0
450 449 cjunkid = sum(junkid)
451 450
452 451 if cjunkid.any():
453 452 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
454 453
455 454 if jcspectraExist:
456 455 for ip in range(num_pairs):
457 456 yy = jcspectra[ip,ind_vel,:]
458 457 jcspectra[ip,freq_dc,:] = numpy.dot(xx_aux,yy)
459 458
460 459
461 460 self.dataOut.data_spc = jspectra
462 461 self.dataOut.data_cspc = jcspectra
463 462
464 463 return 1
465 464
466 465 def removeInterference(self, interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None):
467 466
468 467 jspectra = self.dataOut.data_spc
469 468 jcspectra = self.dataOut.data_cspc
470 469 jnoise = self.dataOut.getNoise()
471 470 num_incoh = self.dataOut.nIncohInt
472 471
473 472 num_channel = jspectra.shape[0]
474 473 num_prof = jspectra.shape[1]
475 474 num_hei = jspectra.shape[2]
476 475
477 476 #hei_interf
478 477 if hei_interf is None:
479 478 count_hei = num_hei/2 #Como es entero no importa
480 479 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
481 480 hei_interf = numpy.asarray(hei_interf)[0]
482 481 #nhei_interf
483 482 if (nhei_interf == None):
484 483 nhei_interf = 5
485 484 if (nhei_interf < 1):
486 485 nhei_interf = 1
487 486 if (nhei_interf > count_hei):
488 487 nhei_interf = count_hei
489 488 if (offhei_interf == None):
490 489 offhei_interf = 0
491 490
492 491 ind_hei = range(num_hei)
493 492 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
494 493 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
495 494 mask_prof = numpy.asarray(range(num_prof))
496 495 num_mask_prof = mask_prof.size
497 496 comp_mask_prof = [0, num_prof/2]
498 497
499 498
500 499 #noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
501 500 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
502 501 jnoise = numpy.nan
503 502 noise_exist = jnoise[0] < numpy.Inf
504 503
505 504 #Subrutina de Remocion de la Interferencia
506 505 for ich in range(num_channel):
507 506 #Se ordena los espectros segun su potencia (menor a mayor)
508 507 power = jspectra[ich,mask_prof,:]
509 508 power = power[:,hei_interf]
510 509 power = power.sum(axis = 0)
511 510 psort = power.ravel().argsort()
512 511
513 512 #Se estima la interferencia promedio en los Espectros de Potencia empleando
514 513 junkspc_interf = jspectra[ich,:,hei_interf[psort[range(offhei_interf, nhei_interf + offhei_interf)]]]
515 514
516 515 if noise_exist:
517 516 # tmp_noise = jnoise[ich] / num_prof
518 517 tmp_noise = jnoise[ich]
519 518 junkspc_interf = junkspc_interf - tmp_noise
520 519 #junkspc_interf[:,comp_mask_prof] = 0
521 520
522 521 jspc_interf = junkspc_interf.sum(axis = 0) / nhei_interf
523 522 jspc_interf = jspc_interf.transpose()
524 523 #Calculando el espectro de interferencia promedio
525 524 noiseid = numpy.where(jspc_interf <= tmp_noise/ numpy.sqrt(num_incoh))
526 525 noiseid = noiseid[0]
527 526 cnoiseid = noiseid.size
528 527 interfid = numpy.where(jspc_interf > tmp_noise/ numpy.sqrt(num_incoh))
529 528 interfid = interfid[0]
530 529 cinterfid = interfid.size
531 530
532 531 if (cnoiseid > 0): jspc_interf[noiseid] = 0
533 532
534 533 #Expandiendo los perfiles a limpiar
535 534 if (cinterfid > 0):
536 535 new_interfid = (numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof)%num_prof
537 536 new_interfid = numpy.asarray(new_interfid)
538 537 new_interfid = {x for x in new_interfid}
539 538 new_interfid = numpy.array(list(new_interfid))
540 539 new_cinterfid = new_interfid.size
541 540 else: new_cinterfid = 0
542 541
543 542 for ip in range(new_cinterfid):
544 543 ind = junkspc_interf[:,new_interfid[ip]].ravel().argsort()
545 544 jspc_interf[new_interfid[ip]] = junkspc_interf[ind[nhei_interf/2],new_interfid[ip]]
546 545
547 546
548 547 jspectra[ich,:,ind_hei] = jspectra[ich,:,ind_hei] - jspc_interf #Corregir indices
549 548
550 549 #Removiendo la interferencia del punto de mayor interferencia
551 550 ListAux = jspc_interf[mask_prof].tolist()
552 551 maxid = ListAux.index(max(ListAux))
553 552
554 553
555 554 if cinterfid > 0:
556 555 for ip in range(cinterfid*(interf == 2) - 1):
557 556 ind = (jspectra[ich,interfid[ip],:] < tmp_noise*(1 + 1/numpy.sqrt(num_incoh))).nonzero()
558 557 cind = len(ind)
559 558
560 559 if (cind > 0):
561 560 jspectra[ich,interfid[ip],ind] = tmp_noise*(1 + (numpy.random.uniform(cind) - 0.5)/numpy.sqrt(num_incoh))
562 561
563 562 ind = numpy.array([-2,-1,1,2])
564 563 xx = numpy.zeros([4,4])
565 564
566 565 for id1 in range(4):
567 566 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
568 567
569 568 xx_inv = numpy.linalg.inv(xx)
570 569 xx = xx_inv[:,0]
571 570 ind = (ind + maxid + num_mask_prof)%num_mask_prof
572 571 yy = jspectra[ich,mask_prof[ind],:]
573 572 jspectra[ich,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
574 573
575 574
576 575 indAux = (jspectra[ich,:,:] < tmp_noise*(1-1/numpy.sqrt(num_incoh))).nonzero()
577 576 jspectra[ich,indAux[0],indAux[1]] = tmp_noise * (1 - 1/numpy.sqrt(num_incoh))
578 577
579 578 #Remocion de Interferencia en el Cross Spectra
580 579 if jcspectra is None: return jspectra, jcspectra
581 580 num_pairs = jcspectra.size/(num_prof*num_hei)
582 581 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
583 582
584 583 for ip in range(num_pairs):
585 584
586 585 #-------------------------------------------
587 586
588 587 cspower = numpy.abs(jcspectra[ip,mask_prof,:])
589 588 cspower = cspower[:,hei_interf]
590 589 cspower = cspower.sum(axis = 0)
591 590
592 591 cspsort = cspower.ravel().argsort()
593 592 junkcspc_interf = jcspectra[ip,:,hei_interf[cspsort[range(offhei_interf, nhei_interf + offhei_interf)]]]
594 593 junkcspc_interf = junkcspc_interf.transpose()
595 594 jcspc_interf = junkcspc_interf.sum(axis = 1)/nhei_interf
596 595
597 596 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
598 597
599 598 median_real = numpy.median(numpy.real(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
600 599 median_imag = numpy.median(numpy.imag(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
601 600 junkcspc_interf[comp_mask_prof,:] = numpy.complex(median_real, median_imag)
602 601
603 602 for iprof in range(num_prof):
604 603 ind = numpy.abs(junkcspc_interf[iprof,:]).ravel().argsort()
605 604 jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf/2]]
606 605
607 606 #Removiendo la Interferencia
608 607 jcspectra[ip,:,ind_hei] = jcspectra[ip,:,ind_hei] - jcspc_interf
609 608
610 609 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
611 610 maxid = ListAux.index(max(ListAux))
612 611
613 612 ind = numpy.array([-2,-1,1,2])
614 613 xx = numpy.zeros([4,4])
615 614
616 615 for id1 in range(4):
617 616 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
618 617
619 618 xx_inv = numpy.linalg.inv(xx)
620 619 xx = xx_inv[:,0]
621 620
622 621 ind = (ind + maxid + num_mask_prof)%num_mask_prof
623 622 yy = jcspectra[ip,mask_prof[ind],:]
624 623 jcspectra[ip,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
625 624
626 625 #Guardar Resultados
627 626 self.dataOut.data_spc = jspectra
628 627 self.dataOut.data_cspc = jcspectra
629 628
630 629 return 1
631 630
632 631 def setRadarFrequency(self, frequency=None):
633 632
634 633 if frequency != None:
635 634 self.dataOut.frequency = frequency
636 635
637 636 return 1
638 637
639 638 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
640 639 #validacion de rango
641 640 if minHei == None:
642 641 minHei = self.dataOut.heightList[0]
643 642
644 643 if maxHei == None:
645 644 maxHei = self.dataOut.heightList[-1]
646 645
647 646 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
648 647 print 'minHei: %.2f is out of the heights range'%(minHei)
649 648 print 'minHei is setting to %.2f'%(self.dataOut.heightList[0])
650 649 minHei = self.dataOut.heightList[0]
651 650
652 651 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
653 652 print 'maxHei: %.2f is out of the heights range'%(maxHei)
654 653 print 'maxHei is setting to %.2f'%(self.dataOut.heightList[-1])
655 654 maxHei = self.dataOut.heightList[-1]
656 655
657 656 # validacion de velocidades
658 657 velrange = self.dataOut.getVelRange(1)
659 658
660 659 if minVel == None:
661 660 minVel = velrange[0]
662 661
663 662 if maxVel == None:
664 663 maxVel = velrange[-1]
665 664
666 665 if (minVel < velrange[0]) or (minVel > maxVel):
667 666 print 'minVel: %.2f is out of the velocity range'%(minVel)
668 667 print 'minVel is setting to %.2f'%(velrange[0])
669 668 minVel = velrange[0]
670 669
671 670 if (maxVel > velrange[-1]) or (maxVel < minVel):
672 671 print 'maxVel: %.2f is out of the velocity range'%(maxVel)
673 672 print 'maxVel is setting to %.2f'%(velrange[-1])
674 673 maxVel = velrange[-1]
675 674
676 675 # seleccion de indices para rango
677 676 minIndex = 0
678 677 maxIndex = 0
679 678 heights = self.dataOut.heightList
680 679
681 680 inda = numpy.where(heights >= minHei)
682 681 indb = numpy.where(heights <= maxHei)
683 682
684 683 try:
685 684 minIndex = inda[0][0]
686 685 except:
687 686 minIndex = 0
688 687
689 688 try:
690 689 maxIndex = indb[0][-1]
691 690 except:
692 691 maxIndex = len(heights)
693 692
694 693 if (minIndex < 0) or (minIndex > maxIndex):
695 694 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
696 695
697 696 if (maxIndex >= self.dataOut.nHeights):
698 697 maxIndex = self.dataOut.nHeights-1
699 698
700 699 # seleccion de indices para velocidades
701 700 indminvel = numpy.where(velrange >= minVel)
702 701 indmaxvel = numpy.where(velrange <= maxVel)
703 702 try:
704 703 minIndexVel = indminvel[0][0]
705 704 except:
706 705 minIndexVel = 0
707 706
708 707 try:
709 708 maxIndexVel = indmaxvel[0][-1]
710 709 except:
711 710 maxIndexVel = len(velrange)
712 711
713 712 #seleccion del espectro
714 713 data_spc = self.dataOut.data_spc[:,minIndexVel:maxIndexVel+1,minIndex:maxIndex+1]
715 714 #estimacion de ruido
716 715 noise = numpy.zeros(self.dataOut.nChannels)
717 716
718 717 for channel in range(self.dataOut.nChannels):
719 718 daux = data_spc[channel,:,:]
720 719 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
721 720
722 721 self.dataOut.noise_estimation = noise.copy()
723 722
724 723 return 1
725 724
726 725 class IncohInt(Operation):
727 726
728 727
729 728 __profIndex = 0
730 729 __withOverapping = False
731 730
732 731 __byTime = False
733 732 __initime = None
734 733 __lastdatatime = None
735 734 __integrationtime = None
736 735
737 736 __buffer_spc = None
738 737 __buffer_cspc = None
739 738 __buffer_dc = None
740 739
741 740 __dataReady = False
742 741
743 742 __timeInterval = None
744 743
745 744 n = None
746 745
747 746
748 747
749 748 def __init__(self, **kwargs):
750 749
751 750 Operation.__init__(self, **kwargs)
752 751 # self.isConfig = False
753 752
754 753 def setup(self, n=None, timeInterval=None, overlapping=False):
755 754 """
756 755 Set the parameters of the integration class.
757 756
758 757 Inputs:
759 758
760 759 n : Number of coherent integrations
761 760 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
762 761 overlapping :
763 762
764 763 """
765 764
766 765 self.__initime = None
767 766 self.__lastdatatime = 0
768 767
769 768 self.__buffer_spc = 0
770 769 self.__buffer_cspc = 0
771 770 self.__buffer_dc = 0
772 771
773 772 self.__profIndex = 0
774 773 self.__dataReady = False
775 774 self.__byTime = False
776 775
777 776 if n is None and timeInterval is None:
778 777 raise ValueError, "n or timeInterval should be specified ..."
779 778
780 779 if n is not None:
781 780 self.n = int(n)
782 781 else:
783 782 self.__integrationtime = int(timeInterval) #if (type(timeInterval)!=integer) -> change this line
784 783 self.n = None
785 784 self.__byTime = True
786 785
787 786 def putData(self, data_spc, data_cspc, data_dc):
788 787
789 788 """
790 789 Add a profile to the __buffer_spc and increase in one the __profileIndex
791 790
792 791 """
793 792
794 793 self.__buffer_spc += data_spc
795 794
796 795 if data_cspc is None:
797 796 self.__buffer_cspc = None
798 797 else:
799 798 self.__buffer_cspc += data_cspc
800 799
801 800 if data_dc is None:
802 801 self.__buffer_dc = None
803 802 else:
804 803 self.__buffer_dc += data_dc
805 804
806 805 self.__profIndex += 1
807 806
808 807 return
809 808
810 809 def pushData(self):
811 810 """
812 811 Return the sum of the last profiles and the profiles used in the sum.
813 812
814 813 Affected:
815 814
816 815 self.__profileIndex
817 816
818 817 """
819 818
820 819 data_spc = self.__buffer_spc
821 820 data_cspc = self.__buffer_cspc
822 821 data_dc = self.__buffer_dc
823 822 n = self.__profIndex
824 823
825 824 self.__buffer_spc = 0
826 825 self.__buffer_cspc = 0
827 826 self.__buffer_dc = 0
828 827 self.__profIndex = 0
829 828
830 829 return data_spc, data_cspc, data_dc, n
831 830
832 831 def byProfiles(self, *args):
833 832
834 833 self.__dataReady = False
835 834 avgdata_spc = None
836 835 avgdata_cspc = None
837 836 avgdata_dc = None
838 837
839 838 self.putData(*args)
840 839
841 840 if self.__profIndex == self.n:
842 841
843 842 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
844 843 self.n = n
845 844 self.__dataReady = True
846 845
847 846 return avgdata_spc, avgdata_cspc, avgdata_dc
848 847
849 848 def byTime(self, datatime, *args):
850 849
851 850 self.__dataReady = False
852 851 avgdata_spc = None
853 852 avgdata_cspc = None
854 853 avgdata_dc = None
855 854
856 855 self.putData(*args)
857 856
858 857 if (datatime - self.__initime) >= self.__integrationtime:
859 858 avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
860 859 self.n = n
861 860 self.__dataReady = True
862 861
863 862 return avgdata_spc, avgdata_cspc, avgdata_dc
864 863
865 864 def integrate(self, datatime, *args):
866 865
867 866 if self.__profIndex == 0:
868 867 self.__initime = datatime
869 868
870 869 if self.__byTime:
871 870 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(datatime, *args)
872 871 else:
873 872 avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
874 873
875 874 if not self.__dataReady:
876 875 return None, None, None, None
877 876
878 877 return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
879 878
880 879 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
881 880
882 881 if n==1:
883 882 return
884 883
885 884 dataOut.flagNoData = True
886 885
887 886 if not self.isConfig:
888 887 self.setup(n, timeInterval, overlapping)
889 888 self.isConfig = True
890 889
891 890 avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
892 891 dataOut.data_spc,
893 892 dataOut.data_cspc,
894 893 dataOut.data_dc)
895 894
896 895 if self.__dataReady:
897 896
898 897 dataOut.data_spc = avgdata_spc
899 898 dataOut.data_cspc = avgdata_cspc
900 899 dataOut.data_dc = avgdata_dc
901 900
902 901 dataOut.nIncohInt *= self.n
903 902 dataOut.utctime = avgdatatime
904 903 dataOut.flagNoData = False
Show file before | Show file after | Hide comments
@@ -1,501 +1,604
1 1 '''
2 2 @author: Juan C. Espinoza
3 3 '''
4 4
5 5 import time
6 6 import json
7 7 import numpy
8 8 import paho.mqtt.client as mqtt
9 9 import zmq
10 10 import datetime
11 11 from zmq.utils.monitor import recv_monitor_message
12 12 from functools import wraps
13 13 from threading import Thread
14 14 from multiprocessing import Process
15 15
16 16 from schainpy.model.proc.jroproc_base import Operation, ProcessingUnit
17 17 from schainpy.model.data.jrodata import JROData
18 from schainpy.utils import log
18 19
19 20 MAXNUMX = 100
20 21 MAXNUMY = 100
21 22
22 23 class PrettyFloat(float):
23 24 def __repr__(self):
24 25 return '%.2f' % self
25 26
26 27 def roundFloats(obj):
27 28 if isinstance(obj, list):
28 29 return map(roundFloats, obj)
29 30 elif isinstance(obj, float):
30 31 return round(obj, 2)
31 32
32 33 def decimate(z, MAXNUMY):
33 # dx = int(len(self.x)/self.__MAXNUMX) + 1
34
35 34 dy = int(len(z[0])/MAXNUMY) + 1
36 35
37 36 return z[::, ::dy]
38 37
39 38 class throttle(object):
40 """Decorator that prevents a function from being called more than once every
39 '''
40 Decorator that prevents a function from being called more than once every
41 41 time period.
42 42 To create a function that cannot be called more than once a minute, but
43 43 will sleep until it can be called:
44 44 @throttle(minutes=1)
45 45 def foo():
46 46 pass
47 47
48 48 for i in range(10):
49 49 foo()
50 50 print "This function has run %s times." % i
51 """
51 '''
52 52
53 53 def __init__(self, seconds=0, minutes=0, hours=0):
54 54 self.throttle_period = datetime.timedelta(
55 55 seconds=seconds, minutes=minutes, hours=hours
56 56 )
57 57
58 58 self.time_of_last_call = datetime.datetime.min
59 59
60 60 def __call__(self, fn):
61 61 @wraps(fn)
62 62 def wrapper(*args, **kwargs):
63 63 now = datetime.datetime.now()
64 64 time_since_last_call = now - self.time_of_last_call
65 65 time_left = self.throttle_period - time_since_last_call
66 66
67 67 if time_left > datetime.timedelta(seconds=0):
68 68 return
69 69
70 70 self.time_of_last_call = datetime.datetime.now()
71 71 return fn(*args, **kwargs)
72 72
73 73 return wrapper
74 74
75 class Data(object):
76 '''
77 Object to hold data to be plotted
78 '''
79
80 def __init__(self, plottypes, throttle_value):
81 self.plottypes = plottypes
82 self.throttle = throttle_value
83 self.ended = False
84 self.__times = []
85
86 def __str__(self):
87 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
88 return 'Data[{}][{}]'.format(';'.join(dum), len(self.__times))
89
90 def __len__(self):
91 return len(self.__times)
92
93 def __getitem__(self, key):
94 if key not in self.data:
95 raise KeyError(log.error('Missing key: {}'.format(key)))
96
97 if 'spc' in key:
98 ret = self.data[key]
99 else:
100 ret = numpy.array([self.data[key][x] for x in self.times])
101 if ret.ndim > 1:
102 ret = numpy.swapaxes(ret, 0, 1)
103 return ret
104
105 def setup(self):
106 '''
107 Configure object
108 '''
109
110 self.ended = False
111 self.data = {}
112 self.__times = []
113 self.__heights = []
114 self.__all_heights = set()
115 for plot in self.plottypes:
116 self.data[plot] = {}
117
118 def shape(self, key):
119 '''
120 Get the shape of the one-element data for the given key
121 '''
122
123 if len(self.data[key]):
124 if 'spc' in key:
125 return self.data[key].shape
126 return self.data[key][self.__times[0]].shape
127 return (0,)
128
129 def update(self, dataOut):
130 '''
131 Update data object with new dataOut
132 '''
133
134 tm = dataOut.utctime
135 if tm in self.__times:
136 return
137
138 self.parameters = getattr(dataOut, 'parameters', [])
139 self.pairs = dataOut.pairsList
140 self.channels = dataOut.channelList
141 self.xrange = (dataOut.getFreqRange(1)/1000. , dataOut.getAcfRange(1) , dataOut.getVelRange(1))
142 self.interval = dataOut.getTimeInterval()
143 self.__heights.append(dataOut.heightList)
144 self.__all_heights.update(dataOut.heightList)
145 self.__times.append(tm)
146
147 for plot in self.plottypes:
148 if plot == 'spc':
149 z = dataOut.data_spc/dataOut.normFactor
150 self.data[plot] = 10*numpy.log10(z)
151 if plot == 'cspc':
152 self.data[plot] = dataOut.data_cspc
153 if plot == 'noise':
154 self.data[plot][tm] = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
155 if plot == 'rti':
156 self.data[plot][tm] = dataOut.getPower()
157 if plot == 'snr_db':
158 self.data['snr'][tm] = dataOut.data_SNR
159 if plot == 'snr':
160 self.data[plot][tm] = 10*numpy.log10(dataOut.data_SNR)
161 if plot == 'dop':
162 self.data[plot][tm] = 10*numpy.log10(dataOut.data_DOP)
163 if plot == 'mean':
164 self.data[plot][tm] = dataOut.data_MEAN
165 if plot == 'std':
166 self.data[plot][tm] = dataOut.data_STD
167 if plot == 'coh':
168 self.data[plot][tm] = dataOut.getCoherence()
169 if plot == 'phase':
170 self.data[plot][tm] = dataOut.getCoherence(phase=True)
171 if plot == 'output':
172 self.data[plot][tm] = dataOut.data_output
173 if plot == 'param':
174 self.data[plot][tm] = dataOut.data_param
175
176 def normalize_heights(self):
177 '''
178 Ensure same-dimension of the data for different heighList
179 '''
180
181 H = numpy.array(list(self.__all_heights))
182 H.sort()
183 for key in self.data:
184 shape = self.shape(key)[:-1] + H.shape
185 for tm, obj in self.data[key].items():
186 h = self.__heights[self.__times.index(tm)]
187 if H.size == h.size:
188 continue
189 index = numpy.where(numpy.in1d(H, h))[0]
190 dummy = numpy.zeros(shape) + numpy.nan
191 if len(shape) == 2:
192 dummy[:, index] = obj
193 else:
194 dummy[index] = obj
195 self.data[key][tm] = dummy
196
197 self.__heights = [H for tm in self.__times]
198
199 def jsonify(self, decimate=False):
200 '''
201 Convert data to json
202 '''
203
204 ret = {}
205 tm = self.times[-1]
206
207 for key, value in self.data:
208 if key in ('spc', 'cspc'):
209 ret[key] = roundFloats(self.data[key].to_list())
210 else:
211 ret[key] = roundFloats(self.data[key][tm].to_list())
212
213 ret['timestamp'] = tm
214 ret['interval'] = self.interval
215
216 @property
217 def times(self):
218 '''
219 Return the list of times of the current data
220 '''
221
222 ret = numpy.array(self.__times)
223 ret.sort()
224 return ret
225
226 @property
227 def heights(self):
228 '''
229 Return the list of heights of the current data
230 '''
231
232 return numpy.array(self.__heights[-1])
75 233
76 234 class PublishData(Operation):
77 """Clase publish."""
235 '''
236 Operation to send data over zmq.
237 '''
78 238
79 239 def __init__(self, **kwargs):
80 240 """Inicio."""
81 241 Operation.__init__(self, **kwargs)
82 242 self.isConfig = False
83 243 self.client = None
84 244 self.zeromq = None
85 245 self.mqtt = None
86 246
87 247 def on_disconnect(self, client, userdata, rc):
88 248 if rc != 0:
89 print("Unexpected disconnection.")
249 log.warning('Unexpected disconnection.')
90 250 self.connect()
91 251
92 252 def connect(self):
93 print 'trying to connect'
253 log.warning('trying to connect')
94 254 try:
95 255 self.client.connect(
96 256 host=self.host,
97 257 port=self.port,
98 258 keepalive=60*10,
99 259 bind_address='')
100 260 self.client.loop_start()
101 261 # self.client.publish(
102 262 # self.topic + 'SETUP',
103 263 # json.dumps(setup),
104 264 # retain=True
105 265 # )
106 266 except:
107 print "MQTT Conection error."
267 log.error('MQTT Conection error.')
108 268 self.client = False
109 269
110 270 def setup(self, port=1883, username=None, password=None, clientId="user", zeromq=1, verbose=True, **kwargs):
111 271 self.counter = 0
112 272 self.topic = kwargs.get('topic', 'schain')
113 273 self.delay = kwargs.get('delay', 0)
114 274 self.plottype = kwargs.get('plottype', 'spectra')
115 275 self.host = kwargs.get('host', "10.10.10.82")
116 276 self.port = kwargs.get('port', 3000)
117 277 self.clientId = clientId
118 278 self.cnt = 0
119 279 self.zeromq = zeromq
120 280 self.mqtt = kwargs.get('plottype', 0)
121 281 self.client = None
122 282 self.verbose = verbose
123 self.dataOut.firstdata = True
124 283 setup = []
125 284 if mqtt is 1:
126 285 self.client = mqtt.Client(
127 286 client_id=self.clientId + self.topic + 'SCHAIN',
128 287 clean_session=True)
129 288 self.client.on_disconnect = self.on_disconnect
130 289 self.connect()
131 290 for plot in self.plottype:
132 291 setup.append({
133 292 'plot': plot,
134 293 'topic': self.topic + plot,
135 294 'title': getattr(self, plot + '_' + 'title', False),
136 295 'xlabel': getattr(self, plot + '_' + 'xlabel', False),
137 296 'ylabel': getattr(self, plot + '_' + 'ylabel', False),
138 297 'xrange': getattr(self, plot + '_' + 'xrange', False),
139 298 'yrange': getattr(self, plot + '_' + 'yrange', False),
140 299 'zrange': getattr(self, plot + '_' + 'zrange', False),
141 300 })
142 301 if zeromq is 1:
143 302 context = zmq.Context()
144 303 self.zmq_socket = context.socket(zmq.PUSH)
145 304 server = kwargs.get('server', 'zmq.pipe')
146 305
147 306 if 'tcp://' in server:
148 307 address = server
149 308 else:
150 309 address = 'ipc:///tmp/%s' % server
151 310
152 311 self.zmq_socket.connect(address)
153 312 time.sleep(1)
154 313
155 314
156 315 def publish_data(self):
157 316 self.dataOut.finished = False
158 317 if self.mqtt is 1:
159 318 yData = self.dataOut.heightList[:2].tolist()
160 319 if self.plottype == 'spectra':
161 320 data = getattr(self.dataOut, 'data_spc')
162 321 z = data/self.dataOut.normFactor
163 322 zdB = 10*numpy.log10(z)
164 323 xlen, ylen = zdB[0].shape
165 324 dx = int(xlen/MAXNUMX) + 1
166 325 dy = int(ylen/MAXNUMY) + 1
167 326 Z = [0 for i in self.dataOut.channelList]
168 327 for i in self.dataOut.channelList:
169 328 Z[i] = zdB[i][::dx, ::dy].tolist()
170 329 payload = {
171 330 'timestamp': self.dataOut.utctime,
172 331 'data': roundFloats(Z),
173 332 'channels': ['Ch %s' % ch for ch in self.dataOut.channelList],
174 333 'interval': self.dataOut.getTimeInterval(),
175 334 'type': self.plottype,
176 335 'yData': yData
177 336 }
178 # print payload
179 337
180 338 elif self.plottype in ('rti', 'power'):
181 339 data = getattr(self.dataOut, 'data_spc')
182 340 z = data/self.dataOut.normFactor
183 341 avg = numpy.average(z, axis=1)
184 342 avgdB = 10*numpy.log10(avg)
185 343 xlen, ylen = z[0].shape
186 344 dy = numpy.floor(ylen/self.__MAXNUMY) + 1
187 345 AVG = [0 for i in self.dataOut.channelList]
188 346 for i in self.dataOut.channelList:
189 347 AVG[i] = avgdB[i][::dy].tolist()
190 348 payload = {
191 349 'timestamp': self.dataOut.utctime,
192 350 'data': roundFloats(AVG),
193 351 'channels': ['Ch %s' % ch for ch in self.dataOut.channelList],
194 352 'interval': self.dataOut.getTimeInterval(),
195 353 'type': self.plottype,
196 354 'yData': yData
197 355 }
198 356 elif self.plottype == 'noise':
199 357 noise = self.dataOut.getNoise()/self.dataOut.normFactor
200 358 noisedB = 10*numpy.log10(noise)
201 359 payload = {
202 360 'timestamp': self.dataOut.utctime,
203 361 'data': roundFloats(noisedB.reshape(-1, 1).tolist()),
204 362 'channels': ['Ch %s' % ch for ch in self.dataOut.channelList],
205 363 'interval': self.dataOut.getTimeInterval(),
206 364 'type': self.plottype,
207 365 'yData': yData
208 366 }
209 367 elif self.plottype == 'snr':
210 368 data = getattr(self.dataOut, 'data_SNR')
211 369 avgdB = 10*numpy.log10(data)
212 370
213 371 ylen = data[0].size
214 372 dy = numpy.floor(ylen/self.__MAXNUMY) + 1
215 373 AVG = [0 for i in self.dataOut.channelList]
216 374 for i in self.dataOut.channelList:
217 375 AVG[i] = avgdB[i][::dy].tolist()
218 376 payload = {
219 377 'timestamp': self.dataOut.utctime,
220 378 'data': roundFloats(AVG),
221 379 'channels': ['Ch %s' % ch for ch in self.dataOut.channelList],
222 380 'type': self.plottype,
223 381 'yData': yData
224 382 }
225 383 else:
226 384 print "Tipo de grafico invalido"
227 385 payload = {
228 386 'data': 'None',
229 387 'timestamp': 'None',
230 388 'type': None
231 389 }
232 # print 'Publishing data to {}'.format(self.host)
390
233 391 self.client.publish(self.topic + self.plottype, json.dumps(payload), qos=0)
234 392
235 393 if self.zeromq is 1:
236 394 if self.verbose:
237 print '[Sending] {} - {}'.format(self.dataOut.type, self.dataOut.datatime)
395 log.log(
396 '{} - {}'.format(self.dataOut.type, self.dataOut.datatime),
397 'Sending'
398 )
238 399 self.zmq_socket.send_pyobj(self.dataOut)
239 self.dataOut.firstdata = False
240
241 400
242 401 def run(self, dataOut, **kwargs):
243 402 self.dataOut = dataOut
244 403 if not self.isConfig:
245 404 self.setup(**kwargs)
246 405 self.isConfig = True
247 406
248 407 self.publish_data()
249 408 time.sleep(self.delay)
250 409
251 410 def close(self):
252 411 if self.zeromq is 1:
253 412 self.dataOut.finished = True
254 413 self.zmq_socket.send_pyobj(self.dataOut)
414 time.sleep(0.1)
255 415 self.zmq_socket.close()
256 416 if self.client:
257 417 self.client.loop_stop()
258 418 self.client.disconnect()
259 419
260 420
261 421 class ReceiverData(ProcessingUnit):
262 422
263 423 def __init__(self, **kwargs):
264 424
265 425 ProcessingUnit.__init__(self, **kwargs)
266 426
267 427 self.isConfig = False
268 428 server = kwargs.get('server', 'zmq.pipe')
269 429 if 'tcp://' in server:
270 430 address = server
271 431 else:
272 432 address = 'ipc:///tmp/%s' % server
273 433
274 434 self.address = address
275 435 self.dataOut = JROData()
276 436
277 437 def setup(self):
278 438
279 439 self.context = zmq.Context()
280 440 self.receiver = self.context.socket(zmq.PULL)
281 441 self.receiver.bind(self.address)
282 442 time.sleep(0.5)
283 print '[Starting] ReceiverData from {}'.format(self.address)
443 log.success('ReceiverData from {}'.format(self.address))
284 444
285 445
286 446 def run(self):
287 447
288 448 if not self.isConfig:
289 449 self.setup()
290 450 self.isConfig = True
291 451
292 452 self.dataOut = self.receiver.recv_pyobj()
293 print '[Receiving] {} - {}'.format(self.dataOut.type,
294 self.dataOut.datatime.ctime())
453 log.log('{} - {}'.format(self.dataOut.type,
454 self.dataOut.datatime.ctime(),),
455 'Receiving')
295 456
296 457
297 458 class PlotterReceiver(ProcessingUnit, Process):
298 459
299 460 throttle_value = 5
300 461
301 462 def __init__(self, **kwargs):
302 463
303 464 ProcessingUnit.__init__(self, **kwargs)
304 465 Process.__init__(self)
305 466 self.mp = False
306 467 self.isConfig = False
307 468 self.isWebConfig = False
308 self.plottypes = []
309 469 self.connections = 0
310 470 server = kwargs.get('server', 'zmq.pipe')
311 471 plot_server = kwargs.get('plot_server', 'zmq.web')
312 472 if 'tcp://' in server:
313 473 address = server
314 474 else:
315 475 address = 'ipc:///tmp/%s' % server
316 476
317 477 if 'tcp://' in plot_server:
318 478 plot_address = plot_server
319 479 else:
320 480 plot_address = 'ipc:///tmp/%s' % plot_server
321 481
322 482 self.address = address
323 483 self.plot_address = plot_address
324 484 self.plottypes = [s.strip() for s in kwargs.get('plottypes', 'rti').split(',')]
325 485 self.realtime = kwargs.get('realtime', False)
326 486 self.throttle_value = kwargs.get('throttle', 5)
327 487 self.sendData = self.initThrottle(self.throttle_value)
488 self.dates = []
328 489 self.setup()
329 490
330 491 def setup(self):
331 492
332 self.data = {}
333 self.data['times'] = []
334 for plottype in self.plottypes:
335 self.data[plottype] = {}
336 self.data['noise'] = {}
337 self.data['throttle'] = self.throttle_value
338 self.data['ENDED'] = False
493 self.data = Data(self.plottypes, self.throttle_value)
339 494 self.isConfig = True
340 self.data_web = {}
341 495
342 496 def event_monitor(self, monitor):
343 497
344 498 events = {}
345 499
346 500 for name in dir(zmq):
347 501 if name.startswith('EVENT_'):
348 502 value = getattr(zmq, name)
349 503 events[value] = name
350 504
351 505 while monitor.poll():
352 506 evt = recv_monitor_message(monitor)
353 507 if evt['event'] == 32:
354 508 self.connections += 1
355 509 if evt['event'] == 512:
356 510 pass
357 if self.connections == 0 and self.started is True:
358 self.ended = True
359 511
360 512 evt.update({'description': events[evt['event']]})
361 513
362 514 if evt['event'] == zmq.EVENT_MONITOR_STOPPED:
363 515 break
364 516 monitor.close()
365 print("event monitor thread done!")
517 print('event monitor thread done!')
366 518
367 519 def initThrottle(self, throttle_value):
368 520
369 521 @throttle(seconds=throttle_value)
370 522 def sendDataThrottled(fn_sender, data):
371 523 fn_sender(data)
372 524
373 525 return sendDataThrottled
374 526
375
376 527 def send(self, data):
377 # print '[sending] data=%s size=%s' % (data.keys(), len(data['times']))
528 log.success('Sending {}'.format(data), self.name)
378 529 self.sender.send_pyobj(data)
379 530
380
381 def update(self):
382 t = self.dataOut.utctime
383
384 if t in self.data['times']:
385 return
386
387 self.data['times'].append(t)
388 self.data['dataOut'] = self.dataOut
389
390 for plottype in self.plottypes:
391 if plottype == 'spc':
392 z = self.dataOut.data_spc/self.dataOut.normFactor
393 self.data[plottype] = 10*numpy.log10(z)
394 self.data['noise'][t] = 10*numpy.log10(self.dataOut.getNoise()/self.dataOut.normFactor)
395 if plottype == 'cspc':
396 jcoherence = self.dataOut.data_cspc/numpy.sqrt(self.dataOut.data_spc*self.dataOut.data_spc)
397 self.data['cspc_coh'] = numpy.abs(jcoherence)
398 self.data['cspc_phase'] = numpy.arctan2(jcoherence.imag, jcoherence.real)*180/numpy.pi
399 if plottype == 'rti':
400 self.data[plottype][t] = self.dataOut.getPower()
401 if plottype == 'snr':
402 self.data[plottype][t] = 10*numpy.log10(self.dataOut.data_SNR)
403 if plottype == 'dop':
404 self.data[plottype][t] = 10*numpy.log10(self.dataOut.data_DOP)
405 if plottype == 'mean':
406 self.data[plottype][t] = self.dataOut.data_MEAN
407 if plottype == 'std':
408 self.data[plottype][t] = self.dataOut.data_STD
409 if plottype == 'coh':
410 self.data[plottype][t] = self.dataOut.getCoherence()
411 if plottype == 'phase':
412 self.data[plottype][t] = self.dataOut.getCoherence(phase=True)
413 if plottype == 'output':
414 self.data[plottype][t] = self.dataOut.data_output
415 if plottype == 'param':
416 self.data[plottype][t] = self.dataOut.data_param
417 if self.realtime:
418 self.data_web['timestamp'] = t
419 if plottype == 'spc':
420 self.data_web[plottype] = roundFloats(decimate(self.data[plottype]).tolist())
421 elif plottype == 'cspc':
422 self.data_web['cspc_coh'] = roundFloats(decimate(self.data['cspc_coh']).tolist())
423 self.data_web['cspc_phase'] = roundFloats(decimate(self.data['cspc_phase']).tolist())
424 elif plottype == 'noise':
425 self.data_web['noise'] = roundFloats(self.data['noise'][t].tolist())
426 else:
427 self.data_web[plottype] = roundFloats(decimate(self.data[plottype][t]).tolist())
428 self.data_web['interval'] = self.dataOut.getTimeInterval()
429 self.data_web['type'] = plottype
430
431 531 def run(self):
432 532
433 print '[Starting] {} from {}'.format(self.name, self.address)
533 log.success(
534 'Starting from {}'.format(self.address),
535 self.name
536 )
434 537
435 538 self.context = zmq.Context()
436 539 self.receiver = self.context.socket(zmq.PULL)
437 540 self.receiver.bind(self.address)
438 541 monitor = self.receiver.get_monitor_socket()
439 542 self.sender = self.context.socket(zmq.PUB)
440 543 if self.realtime:
441 544 self.sender_web = self.context.socket(zmq.PUB)
442 545 self.sender_web.connect(self.plot_address)
443 546 time.sleep(1)
444 547
445 548 if 'server' in self.kwargs:
446 549 self.sender.bind("ipc:///tmp/{}.plots".format(self.kwargs['server']))
447 550 else:
448 551 self.sender.bind("ipc:///tmp/zmq.plots")
449 552
450 time.sleep(3)
553 time.sleep(2)
451 554
452 555 t = Thread(target=self.event_monitor, args=(monitor,))
453 556 t.start()
454 557
455 558 while True:
456 self.dataOut = self.receiver.recv_pyobj()
457 # print '[Receiving] {} - {}'.format(self.dataOut.type,
458 # self.dataOut.datatime.ctime())
459
460 self.update()
559 dataOut = self.receiver.recv_pyobj()
560 dt = datetime.datetime.fromtimestamp(dataOut.utctime).date()
561 sended = False
562 if dt not in self.dates:
563 if self.data:
564 self.data.ended = True
565 self.send(self.data)
566 sended = True
567 self.data.setup()
568 self.dates.append(dt)
461 569
462 if self.dataOut.firstdata is True:
463 self.data['STARTED'] = True
570 self.data.update(dataOut)
464 571
465 if self.dataOut.finished is True:
466 self.send(self.data)
572 if dataOut.finished is True:
467 573 self.connections -= 1
468 if self.connections == 0 and self.started:
469 self.ended = True
470 self.data['ENDED'] = True
574 if self.connections == 0 and dt in self.dates:
575 self.data.ended = True
471 576 self.send(self.data)
472 self.setup()
473 self.started = False
577 self.data.setup()
474 578 else:
475 579 if self.realtime:
476 580 self.send(self.data)
477 self.sender_web.send_string(json.dumps(self.data_web))
581 # self.sender_web.send_string(self.data.jsonify())
478 582 else:
583 if not sended:
479 584 self.sendData(self.send, self.data)
480 self.started = True
481 585
482 self.data['STARTED'] = False
483 586 return
484 587
485 588 def sendToWeb(self):
486 589
487 590 if not self.isWebConfig:
488 591 context = zmq.Context()
489 592 sender_web_config = context.socket(zmq.PUB)
490 593 if 'tcp://' in self.plot_address:
491 594 dum, address, port = self.plot_address.split(':')
492 595 conf_address = '{}:{}:{}'.format(dum, address, int(port)+1)
493 596 else:
494 597 conf_address = self.plot_address + '.config'
495 598 sender_web_config.bind(conf_address)
496 599 time.sleep(1)
497 600 for kwargs in self.operationKwargs.values():
498 601 if 'plot' in kwargs:
499 print '[Sending] Config data to web for {}'.format(kwargs['code'].upper())
602 log.success('[Sending] Config data to web for {}'.format(kwargs['code'].upper()))
500 603 sender_web_config.send_string(json.dumps(kwargs))
501 604 self.isWebConfig = True No newline at end of file
General Comments 0
You need to be logged in to leave comments. Login now