##// END OF EJS Templates
schain
RSS
Fix xmin & xmax in plots, fix SendToFTP
jespinoza -
r1334:96199ee477a0
parent child
Show More
Side by Side Unified
Context file:
r1334:96199ee477a0 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,102 +1,102
1 1 # Signal Chain
2 2
3 3 Signal Chain is a radar data processing library wich includes modules to read,
4 4 and write different files formats, besides modules to process and visualize the
5 5 data.
6 6
7 7 ## Dependencies
8 8
9 9 - GCC (gcc or gfortran)
10 10 - Python.h (python-dev or python-devel)
11 11 - Python-TK (python-tk)
12 12 - HDF5 libraries (libhdf5-dev)
13 13
14 14 ## Installation
15 15
16 16 To get started the easiest way to install it is through [PyPI](https://pypi.org/project/schainpy/) with pip:
17 17
18 18 ```bash
19 19 pip install schainpy
20 20 ```
21 21
22 22 ### From source
23 23
24 24 First, ensure that you have the above-listed dependencies installed, then clone
25 25 the repository and install as normal python package:
26 26
27 27 ```bash
28 28 git clone https://github.com/JRO-Peru/schain.git
29 29 cd schain
30 30 git checkout `branch-name` (optional)
31 31 sudo pip install ./
32 32 ```
33 33
34 34 ### Using Docker
35 35
36 36 Download Dockerfile from the repository, and create a docker image:
37 37
38 38 ```bash
39 39 docker build -t schain .
40 40 ```
41 41
42 42 You can run a container using an xml file or a schain script also you need to
43 43 mount a volume for the data input and for the output files/plots:
44 44
45 45 ```bash
46 46 docker run -it --rm --volume /path/to/host/data:/data schain xml /data/test.xml
47 47 docker run -it --rm --volume /path/to/host/data:/data --entrypoint /urs/local/bin/python schain /data/test.py
48 48 ```
49 49
50 50 ## CLI (command line interface)
51 51
52 52 Signal Chain provides the following commands:
53 53
54 54 - schainGUI: Open the GUI
55 55 - schain: Signal chain command line
56 56
57 57 ## Example
58 58
59 59 Here you can find an script to read Spectra data (.pdata), remove dc and plot
60 60 self-spectra & RTI:
61 61
62 62 ```python
63 63 #!/usr/bin/python
64 64
65 65 from schainpy.controller import Project
66 66
67 67 prj = Project()
68 68
69 69 read_unit = prj.addReadUnit(
70 70 datatype='Spectra',
71 71 path='/path/to/pdata/',
72 72 startDate='2014/01/31',
73 73 endDate='2014/03/31',
74 74 startTime='00:00:00',
75 75 endTime='23:59:59',
76 76 online=0,
77 77 walk=0
78 78 )
79 79
80 80 proc_unit = prj.addProcUnit(
81 81 datatype='Spectra',
82 82 inputId=read_unit.getId()
83 83 )
84 84
85 85 op = proc_unit.addOperation(name='selectChannels')
86 86 op.addParameter(name='channelList', value='0,1')
87 87
88 88 op = proc_unit.addOperation(name='selectHeights')
89 89 op.addParameter(name='minHei', value='80')
90 90 op.addParameter(name='maxHei', value='200')
91 91
92 92 op = proc_unit.addOperation(name='removeDC')
93 93
94 94 op = proc_unit.addOperation(name='SpectraPlot')
95 95 op.addParameter(name='wintitle', value='Spectra', format='str')
96 96
97 op = procUnitConfObj1.addOperation(name='RTIPlot')
97 op = proc_unit.addOperation(name='RTIPlot')
98 98 op.addParameter(name='wintitle', value='RTI', format='str')
99 99
100 100 prj.start()
101 101
102 102 ```
Show file before | Show file after | Hide comments
@@ -1,716 +1,665
1 # Copyright (c) 2012-2020 Jicamarca Radio Observatory
2 # All rights reserved.
3 #
4 # Distributed under the terms of the BSD 3-clause license.
5 """Base class to create plot operations
6
7 """
1 8
2 9 import os
3 10 import sys
4 11 import zmq
5 12 import time
6 13 import numpy
7 14 import datetime
8 try:
9 from queue import Queue
10 except:
11 from Queue import Queue
15 from multiprocessing import Queue
12 16 from functools import wraps
13 17 from threading import Thread
14 18 import matplotlib
15 19
16 20 if 'BACKEND' in os.environ:
17 21 matplotlib.use(os.environ['BACKEND'])
18 22 elif 'linux' in sys.platform:
19 23 matplotlib.use("TkAgg")
20 24 elif 'darwin' in sys.platform:
21 25 matplotlib.use('WxAgg')
22 26 else:
23 27 from schainpy.utils import log
24 28 log.warning('Using default Backend="Agg"', 'INFO')
25 29 matplotlib.use('Agg')
26 30
27 31 import matplotlib.pyplot as plt
28 32 from matplotlib.patches import Polygon
29 33 from mpl_toolkits.axes_grid1 import make_axes_locatable
30 34 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
31 35
32 36 from schainpy.model.data.jrodata import PlotterData
33 37 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
34 38 from schainpy.utils import log
35 39
36 40 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
37 41 blu_values = matplotlib.pyplot.get_cmap(
38 42 'seismic_r', 20)(numpy.arange(20))[10:15]
39 43 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
40 44 'jro', numpy.vstack((blu_values, jet_values)))
41 45 matplotlib.pyplot.register_cmap(cmap=ncmap)
42 46
43 47 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis',
44 48 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm')]
45 49
46 50 EARTH_RADIUS = 6.3710e3
47 51
48 52 def ll2xy(lat1, lon1, lat2, lon2):
49 53
50 54 p = 0.017453292519943295
51 55 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
52 56 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
53 57 r = 12742 * numpy.arcsin(numpy.sqrt(a))
54 58 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
55 59 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
56 60 theta = -theta + numpy.pi/2
57 61 return r*numpy.cos(theta), r*numpy.sin(theta)
58 62
59 63
60 64 def km2deg(km):
61 65 '''
62 66 Convert distance in km to degrees
63 67 '''
64 68
65 69 return numpy.rad2deg(km/EARTH_RADIUS)
66 70
67 71
68 72 def figpause(interval):
69 73 backend = plt.rcParams['backend']
70 74 if backend in matplotlib.rcsetup.interactive_bk:
71 75 figManager = matplotlib._pylab_helpers.Gcf.get_active()
72 76 if figManager is not None:
73 77 canvas = figManager.canvas
74 78 if canvas.figure.stale:
75 79 canvas.draw()
76 80 try:
77 81 canvas.start_event_loop(interval)
78 82 except:
79 83 pass
80 84 return
81 85
82 86
83 87 def popup(message):
84 88 '''
85 89 '''
86 90
87 91 fig = plt.figure(figsize=(12, 8), facecolor='r')
88 92 text = '\n'.join([s.strip() for s in message.split(':')])
89 93 fig.text(0.01, 0.5, text, ha='left', va='center',
90 94 size='20', weight='heavy', color='w')
91 95 fig.show()
92 96 figpause(1000)
93 97
94 98
95 99 class Throttle(object):
96 100 '''
97 101 Decorator that prevents a function from being called more than once every
98 102 time period.
99 103 To create a function that cannot be called more than once a minute, but
100 104 will sleep until it can be called:
101 105 @Throttle(minutes=1)
102 106 def foo():
103 107 pass
104 108
105 109 for i in range(10):
106 110 foo()
107 111 print "This function has run %s times." % i
108 112 '''
109 113
110 114 def __init__(self, seconds=0, minutes=0, hours=0):
111 115 self.throttle_period = datetime.timedelta(
112 116 seconds=seconds, minutes=minutes, hours=hours
113 117 )
114 118
115 119 self.time_of_last_call = datetime.datetime.min
116 120
117 121 def __call__(self, fn):
118 122 @wraps(fn)
119 123 def wrapper(*args, **kwargs):
120 124 coerce = kwargs.pop('coerce', None)
121 125 if coerce:
122 126 self.time_of_last_call = datetime.datetime.now()
123 127 return fn(*args, **kwargs)
124 128 else:
125 129 now = datetime.datetime.now()
126 130 time_since_last_call = now - self.time_of_last_call
127 131 time_left = self.throttle_period - time_since_last_call
128 132
129 133 if time_left > datetime.timedelta(seconds=0):
130 134 return
131 135
132 136 self.time_of_last_call = datetime.datetime.now()
133 137 return fn(*args, **kwargs)
134 138
135 139 return wrapper
136 140
137 141 def apply_throttle(value):
138 142
139 143 @Throttle(seconds=value)
140 144 def fnThrottled(fn):
141 145 fn()
142 146
143 147 return fnThrottled
144 148
145 149
146 150 @MPDecorator
147 151 class Plot(Operation):
148 '''
149 Base class for Schain plotting operations
150 '''
152 """Base class for Schain plotting operations
153
154 This class should never be use directtly you must subclass a new operation,
155 children classes must be defined as follow:
156
157 ExamplePlot(Plot):
158
159 CODE = 'code'
160 colormap = 'jet'
161 plot_type = 'pcolor' # options are ('pcolor', 'pcolorbuffer', 'scatter', 'scatterbuffer')
162
163 def setup(self):
164 pass
165
166 def plot(self):
167 pass
168
169 """
151 170
152 171 CODE = 'Figure'
153 172 colormap = 'jet'
154 173 bgcolor = 'white'
155 174 buffering = True
156 175 __missing = 1E30
157 176
158 177 __attrs__ = ['show', 'save', 'ymin', 'ymax', 'zmin', 'zmax', 'title',
159 178 'showprofile']
160 179
161 180 def __init__(self):
162 181
163 182 Operation.__init__(self)
164 183 self.isConfig = False
165 184 self.isPlotConfig = False
166 self.save_counter = 1
185 self.save_time = 0
167 186 self.sender_time = 0
168 187 self.data = None
169 188 self.firsttime = True
170 189 self.sender_queue = Queue(maxsize=60)
171 190 self.plots_adjust = {'left': 0.125, 'right': 0.9, 'bottom': 0.15, 'top': 0.9, 'wspace': 0.2, 'hspace': 0.2}
172 191
173 192 def __fmtTime(self, x, pos):
174 193 '''
175 194 '''
176 195
177 196 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
178 197
179 198 def __setup(self, **kwargs):
180 199 '''
181 200 Initialize variables
182 201 '''
183 202
184 203 self.figures = []
185 204 self.axes = []
186 205 self.cb_axes = []
187 206 self.localtime = kwargs.pop('localtime', True)
188 207 self.show = kwargs.get('show', True)
189 208 self.save = kwargs.get('save', False)
190 self.save_period = kwargs.get('save_period', 1)
209 self.save_period = kwargs.get('save_period', 60)
191 210 self.colormap = kwargs.get('colormap', self.colormap)
192 211 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
193 212 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
194 213 self.colormaps = kwargs.get('colormaps', None)
195 214 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
196 215 self.showprofile = kwargs.get('showprofile', False)
197 216 self.title = kwargs.get('wintitle', self.CODE.upper())
198 217 self.cb_label = kwargs.get('cb_label', None)
199 218 self.cb_labels = kwargs.get('cb_labels', None)
200 219 self.labels = kwargs.get('labels', None)
201 220 self.xaxis = kwargs.get('xaxis', 'frequency')
202 221 self.zmin = kwargs.get('zmin', None)
203 222 self.zmax = kwargs.get('zmax', None)
204 223 self.zlimits = kwargs.get('zlimits', None)
205 224 self.xmin = kwargs.get('xmin', None)
206 225 self.xmax = kwargs.get('xmax', None)
207 226 self.xrange = kwargs.get('xrange', 12)
208 227 self.xscale = kwargs.get('xscale', None)
209 228 self.ymin = kwargs.get('ymin', None)
210 229 self.ymax = kwargs.get('ymax', None)
211 230 self.yscale = kwargs.get('yscale', None)
212 231 self.xlabel = kwargs.get('xlabel', None)
213 232 self.attr_time = kwargs.get('attr_time', 'utctime')
214 233 self.decimation = kwargs.get('decimation', None)
215 234 self.showSNR = kwargs.get('showSNR', False)
216 235 self.oneFigure = kwargs.get('oneFigure', True)
217 236 self.width = kwargs.get('width', None)
218 237 self.height = kwargs.get('height', None)
219 238 self.colorbar = kwargs.get('colorbar', True)
220 239 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
221 240 self.channels = kwargs.get('channels', None)
222 241 self.titles = kwargs.get('titles', [])
223 242 self.polar = False
224 243 self.type = kwargs.get('type', 'iq')
225 244 self.grid = kwargs.get('grid', False)
226 245 self.pause = kwargs.get('pause', False)
227 self.save_code = kwargs.get('save_code', None)
246 self.save_code = kwargs.get('save_code', self.CODE)
228 247 self.throttle = kwargs.get('throttle', 0)
229 248 self.exp_code = kwargs.get('exp_code', None)
230 249 self.plot_server = kwargs.get('plot_server', False)
231 250 self.sender_period = kwargs.get('sender_period', 60)
232 251 self.tag = kwargs.get('tag', '')
233 252 self.height_index = kwargs.get('height_index', None)
234 253 self.__throttle_plot = apply_throttle(self.throttle)
235 254 self.data = PlotterData(
236 255 self.CODE, self.throttle, self.exp_code, self.localtime, self.buffering, snr=self.showSNR)
237 256
238 257 if self.plot_server:
239 258 if not self.plot_server.startswith('tcp://'):
240 259 self.plot_server = 'tcp://{}'.format(self.plot_server)
241 260 log.success(
242 261 'Sending to server: {}'.format(self.plot_server),
243 262 self.name
244 263 )
245 if 'plot_name' in kwargs:
246 self.plot_name = kwargs['plot_name']
247 264
248 265 def __setup_plot(self):
249 266 '''
250 267 Common setup for all figures, here figures and axes are created
251 268 '''
252 269
253 270 self.setup()
254 271
255 272 self.time_label = 'LT' if self.localtime else 'UTC'
256 273
257 274 if self.width is None:
258 275 self.width = 8
259 276
260 277 self.figures = []
261 278 self.axes = []
262 279 self.cb_axes = []
263 280 self.pf_axes = []
264 281 self.cmaps = []
265 282
266 283 size = '15%' if self.ncols == 1 else '30%'
267 284 pad = '4%' if self.ncols == 1 else '8%'
268 285
269 286 if self.oneFigure:
270 287 if self.height is None:
271 288 self.height = 1.4 * self.nrows + 1
272 289 fig = plt.figure(figsize=(self.width, self.height),
273 290 edgecolor='k',
274 291 facecolor='w')
275 292 self.figures.append(fig)
276 293 for n in range(self.nplots):
277 294 ax = fig.add_subplot(self.nrows, self.ncols,
278 295 n + 1, polar=self.polar)
279 296 ax.tick_params(labelsize=8)
280 297 ax.firsttime = True
281 298 ax.index = 0
282 299 ax.press = None
283 300 self.axes.append(ax)
284 301 if self.showprofile:
285 302 cax = self.__add_axes(ax, size=size, pad=pad)
286 303 cax.tick_params(labelsize=8)
287 304 self.pf_axes.append(cax)
288 305 else:
289 306 if self.height is None:
290 307 self.height = 3
291 308 for n in range(self.nplots):
292 309 fig = plt.figure(figsize=(self.width, self.height),
293 310 edgecolor='k',
294 311 facecolor='w')
295 312 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
296 313 ax.tick_params(labelsize=8)
297 314 ax.firsttime = True
298 315 ax.index = 0
299 316 ax.press = None
300 317 self.figures.append(fig)
301 318 self.axes.append(ax)
302 319 if self.showprofile:
303 320 cax = self.__add_axes(ax, size=size, pad=pad)
304 321 cax.tick_params(labelsize=8)
305 322 self.pf_axes.append(cax)
306 323
307 324 for n in range(self.nrows):
308 325 if self.colormaps is not None:
309 326 cmap = plt.get_cmap(self.colormaps[n])
310 327 else:
311 328 cmap = plt.get_cmap(self.colormap)
312 329 cmap.set_bad(self.bgcolor, 1.)
313 330 self.cmaps.append(cmap)
314 331
315 332 def __add_axes(self, ax, size='30%', pad='8%'):
316 333 '''
317 334 Add new axes to the given figure
318 335 '''
319 336 divider = make_axes_locatable(ax)
320 337 nax = divider.new_horizontal(size=size, pad=pad)
321 338 ax.figure.add_axes(nax)
322 339 return nax
323 340
324 341 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
325 342 '''
326 343 Create a masked array for missing data
327 344 '''
328 345 if x_buffer.shape[0] < 2:
329 346 return x_buffer, y_buffer, z_buffer
330 347
331 348 deltas = x_buffer[1:] - x_buffer[0:-1]
332 349 x_median = numpy.median(deltas)
333 350
334 351 index = numpy.where(deltas > 5 * x_median)
335 352
336 353 if len(index[0]) != 0:
337 354 z_buffer[::, index[0], ::] = self.__missing
338 355 z_buffer = numpy.ma.masked_inside(z_buffer,
339 356 0.99 * self.__missing,
340 357 1.01 * self.__missing)
341 358
342 359 return x_buffer, y_buffer, z_buffer
343 360
344 361 def decimate(self):
345 362
346 363 # dx = int(len(self.x)/self.__MAXNUMX) + 1
347 364 dy = int(len(self.y) / self.decimation) + 1
348 365
349 366 # x = self.x[::dx]
350 367 x = self.x
351 368 y = self.y[::dy]
352 369 z = self.z[::, ::, ::dy]
353 370
354 371 return x, y, z
355 372
356 373 def format(self):
357 374 '''
358 375 Set min and max values, labels, ticks and titles
359 376 '''
360
361 if self.xmin is None:
362 xmin = self.data.min_time
363 else:
364 if self.xaxis is 'time':
365 dt = self.getDateTime(self.data.min_time)
366 xmin = (dt.replace(hour=int(self.xmin), minute=0, second=0) -
367 datetime.datetime(1970, 1, 1)).total_seconds()
368 if self.data.localtime:
369 xmin += time.timezone
370 else:
371 xmin = self.xmin
372
373 if self.xmax is None:
374 xmax = xmin + self.xrange * 60 * 60
375 else:
376 if self.xaxis is 'time':
377 dt = self.getDateTime(self.data.max_time)
378 xmax = self.xmax - 1
379 xmax = (dt.replace(hour=int(xmax), minute=59, second=59) -
380 datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds()
381 if self.data.localtime:
382 xmax += time.timezone
383 else:
384 xmax = self.xmax
385
386 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
387 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
388 377
389 378 for n, ax in enumerate(self.axes):
390 379 if ax.firsttime:
391
392 dig = int(numpy.log10(ymax))
393 if dig == 0:
394 digD = len(str(ymax)) - 2
395 ydec = ymax*(10**digD)
396
397 dig = int(numpy.log10(ydec))
398 ystep = ((ydec + (10**(dig)))//10**(dig))*(10**(dig))
399 ystep = ystep/5
400 ystep = ystep/(10**digD)
401
402 else:
403 ystep = ((ymax + (10**(dig)))//10**(dig))*(10**(dig))
404 ystep = ystep/5
405
406 if self.xaxis is not 'time':
407
408 dig = int(numpy.log10(xmax))
409
410 if dig <= 0:
411 digD = len(str(xmax)) - 2
412 xdec = xmax*(10**digD)
413
414 dig = int(numpy.log10(xdec))
415 xstep = ((xdec + (10**(dig)))//10**(dig))*(10**(dig))
416 xstep = xstep*0.5
417 xstep = xstep/(10**digD)
418
419 else:
420 xstep = ((xmax + (10**(dig)))//10**(dig))*(10**(dig))
421 xstep = xstep/5
422
380 if self.xaxis != 'time':
381 xmin = self.xmin
382 xmax = self.xmax
383 else:
384 xmin = self.tmin
385 xmax = self.tmin + self.xrange*60*60
386 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
387 ax.xaxis.set_major_locator(LinearLocator(9))
388 ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
389 ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
423 390 ax.set_facecolor(self.bgcolor)
424 ax.yaxis.set_major_locator(MultipleLocator(ystep))
425 391 if self.xscale:
426 392 ax.xaxis.set_major_formatter(FuncFormatter(
427 393 lambda x, pos: '{0:g}'.format(x*self.xscale)))
428 if self.xscale:
394 if self.yscale:
429 395 ax.yaxis.set_major_formatter(FuncFormatter(
430 396 lambda x, pos: '{0:g}'.format(x*self.yscale)))
431 if self.xaxis is 'time':
432 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
433 ax.xaxis.set_major_locator(LinearLocator(9))
434 else:
435 ax.xaxis.set_major_locator(MultipleLocator(xstep))
436 397 if self.xlabel is not None:
437 398 ax.set_xlabel(self.xlabel)
438 ax.set_ylabel(self.ylabel)
439 ax.firsttime = False
399 if self.ylabel is not None:
400 ax.set_ylabel(self.ylabel)
440 401 if self.showprofile:
441 402 self.pf_axes[n].set_ylim(ymin, ymax)
442 403 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
443 404 self.pf_axes[n].set_xlabel('dB')
444 405 self.pf_axes[n].grid(b=True, axis='x')
445 406 [tick.set_visible(False)
446 407 for tick in self.pf_axes[n].get_yticklabels()]
447 408 if self.colorbar:
448 409 ax.cbar = plt.colorbar(
449 410 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
450 411 ax.cbar.ax.tick_params(labelsize=8)
451 412 ax.cbar.ax.press = None
452 413 if self.cb_label:
453 414 ax.cbar.set_label(self.cb_label, size=8)
454 415 elif self.cb_labels:
455 416 ax.cbar.set_label(self.cb_labels[n], size=8)
456 417 else:
457 418 ax.cbar = None
458 if self.grid:
459 ax.grid(True)
460
461 if not self.polar:
462 419 ax.set_xlim(xmin, xmax)
463 420 ax.set_ylim(ymin, ymax)
421 ax.firsttime = False
422 if self.grid:
423 ax.grid(True)
424 if not self.polar:
464 425 ax.set_title('{} {} {}'.format(
465 426 self.titles[n],
466 427 self.getDateTime(self.data.max_time).strftime(
467 428 '%Y-%m-%d %H:%M:%S'),
468 429 self.time_label),
469 430 size=8)
470 431 else:
471 432 ax.set_title('{}'.format(self.titles[n]), size=8)
472 433 ax.set_ylim(0, 90)
473 434 ax.set_yticks(numpy.arange(0, 90, 20))
474 435 ax.yaxis.labelpad = 40
475 436
476 437 if self.firsttime:
477 438 for n, fig in enumerate(self.figures):
478 439 fig.subplots_adjust(**self.plots_adjust)
479 440 self.firsttime = False
480 441
481 442 def clear_figures(self):
482 443 '''
483 444 Reset axes for redraw plots
484 445 '''
485 446
486 447 for ax in self.axes+self.pf_axes+self.cb_axes:
487 448 ax.clear()
488 449 ax.firsttime = True
489 450 if hasattr(ax, 'cbar') and ax.cbar:
490 451 ax.cbar.remove()
491 452
492 453 def __plot(self):
493 454 '''
494 455 Main function to plot, format and save figures
495 456 '''
496 457
497 458 self.plot()
498 459 self.format()
499 460
500 461 for n, fig in enumerate(self.figures):
501 462 if self.nrows == 0 or self.nplots == 0:
502 463 log.warning('No data', self.name)
503 464 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
504 465 fig.canvas.manager.set_window_title(self.CODE)
505 466 continue
506 467
507 468 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
508 469 self.getDateTime(self.data.max_time).strftime('%Y/%m/%d')))
509 470 fig.canvas.draw()
510 471 if self.show:
511 472 fig.show()
512 473 figpause(0.01)
513 474
514 475 if self.save:
515 476 self.save_figure(n)
516 477
517 478 if self.plot_server:
518 479 self.send_to_server()
519 480
520 481 def save_figure(self, n):
521 482 '''
522 483 '''
523 484
524 if self.save_counter < self.save_period:
525 self.save_counter += 1
485 if (self.data.tm - self.save_time) < self.sender_period:
526 486 return
527 487
528 self.save_counter = 1
488 self.save_time = self.data.tm
529 489
530 490 fig = self.figures[n]
531 491
532 if self.save_code:
533 if isinstance(self.save_code, str):
534 labels = [self.save_code for x in self.figures]
535 else:
536 labels = self.save_code
537 else:
538 labels = [self.CODE for x in self.figures]
539
540 492 figname = os.path.join(
541 493 self.save,
542 labels[n],
494 self.save_code,
543 495 '{}_{}.png'.format(
544 labels[n],
496 self.save_code,
545 497 self.getDateTime(self.data.max_time).strftime(
546 498 '%Y%m%d_%H%M%S'
547 499 ),
548 500 )
549 501 )
550 502 log.log('Saving figure: {}'.format(figname), self.name)
551 503 if not os.path.isdir(os.path.dirname(figname)):
552 504 os.makedirs(os.path.dirname(figname))
553 505 fig.savefig(figname)
554 506
555 507 if self.throttle == 0:
556 508 figname = os.path.join(
557 509 self.save,
558 510 '{}_{}.png'.format(
559 labels[n],
511 self.save_code,
560 512 self.getDateTime(self.data.min_time).strftime(
561 513 '%Y%m%d'
562 514 ),
563 515 )
564 516 )
565 517 fig.savefig(figname)
566 518
567 519 def send_to_server(self):
568 520 '''
569 521 '''
570 522
571 523 interval = self.data.tm - self.sender_time
572 524 if interval < self.sender_period:
573 525 return
574 526
575 527 self.sender_time = self.data.tm
576 528
577 529 attrs = ['titles', 'zmin', 'zmax', 'tag', 'ymin', 'ymax']
578 530 for attr in attrs:
579 531 value = getattr(self, attr)
580 532 if value:
581 533 if isinstance(value, (numpy.float32, numpy.float64)):
582 534 value = round(float(value), 2)
583 535 self.data.meta[attr] = value
584 536 if self.colormap == 'jet':
585 537 self.data.meta['colormap'] = 'Jet'
586 538 elif 'RdBu' in self.colormap:
587 539 self.data.meta['colormap'] = 'RdBu'
588 540 else:
589 541 self.data.meta['colormap'] = 'Viridis'
590 542 self.data.meta['interval'] = int(interval)
591 # msg = self.data.jsonify(self.data.tm, self.plot_name, self.plot_type)
543
592 544 try:
593 545 self.sender_queue.put(self.data.tm, block=False)
594 546 except:
595 547 tm = self.sender_queue.get()
596 548 self.sender_queue.put(self.data.tm)
597 549
598 550 while True:
599 551 if self.sender_queue.empty():
600 552 break
601 553 tm = self.sender_queue.get()
602 554 try:
603 msg = self.data.jsonify(tm, self.plot_name, self.plot_type)
555 msg = self.data.jsonify(tm, self.save_code, self.plot_type)
604 556 except:
605 557 continue
606 558 self.socket.send_string(msg)
607 559 socks = dict(self.poll.poll(5000))
608 560 if socks.get(self.socket) == zmq.POLLIN:
609 561 reply = self.socket.recv_string()
610 562 if reply == 'ok':
611 563 log.log("Response from server ok", self.name)
612 564 time.sleep(0.2)
613 565 continue
614 566 else:
615 567 log.warning(
616 568 "Malformed reply from server: {}".format(reply), self.name)
617 569 else:
618 570 log.warning(
619 571 "No response from server, retrying...", self.name)
620 572 self.sender_queue.put(self.data.tm)
621 573 self.socket.setsockopt(zmq.LINGER, 0)
622 574 self.socket.close()
623 575 self.poll.unregister(self.socket)
624 576 time.sleep(0.1)
625 577 self.socket = self.context.socket(zmq.REQ)
626 578 self.socket.connect(self.plot_server)
627 579 self.poll.register(self.socket, zmq.POLLIN)
628 580 break
629 581
630 582 def setup(self):
631 583 '''
632 584 This method should be implemented in the child class, the following
633 585 attributes should be set:
634 586
635 587 self.nrows: number of rows
636 588 self.ncols: number of cols
637 589 self.nplots: number of plots (channels or pairs)
638 590 self.ylabel: label for Y axes
639 591 self.titles: list of axes title
640 592
641 593 '''
642 594 raise NotImplementedError
643 595
644 596 def plot(self):
645 597 '''
646 598 Must be defined in the child class
647 599 '''
648 600 raise NotImplementedError
649 601
650 602 def run(self, dataOut, **kwargs):
651 603 '''
652 604 Main plotting routine
653 605 '''
654 606
655 607 if self.isConfig is False:
656 608 self.__setup(**kwargs)
657
658 t = getattr(dataOut, self.attr_time)
659 609
660 610 if self.localtime:
661 611 self.getDateTime = datetime.datetime.fromtimestamp
662 612 else:
663 613 self.getDateTime = datetime.datetime.utcfromtimestamp
664
665 if self.xmin is None:
666 self.tmin = t
667 if 'buffer' in self.plot_type:
668 self.xmin = self.getDateTime(t).hour
669 else:
670 self.tmin = (
671 self.getDateTime(t).replace(
672 hour=int(self.xmin),
673 minute=0,
674 second=0) - self.getDateTime(0)).total_seconds()
675 614
676 615 self.data.setup()
677 616 self.isConfig = True
678 617 if self.plot_server:
679 618 self.context = zmq.Context()
680 619 self.socket = self.context.socket(zmq.REQ)
681 620 self.socket.connect(self.plot_server)
682 621 self.poll = zmq.Poller()
683 622 self.poll.register(self.socket, zmq.POLLIN)
684 623
685 624 tm = getattr(dataOut, self.attr_time)
686 625
687 if self.data and (tm - self.tmin) >= self.xrange*60*60:
626 if self.data and 'time' in self.xaxis and (tm - self.tmin) >= self.xrange*60*60:
688 627 self.save_counter = self.save_period
689 628 self.__plot()
690 if 'time' in self.xaxis:
691 self.xmin += self.xrange
692 if self.xmin >= 24:
693 self.xmin -= 24
694 629 self.tmin += self.xrange*60*60
695 630 self.data.setup()
696 631 self.clear_figures()
697 632
698 633 self.data.update(dataOut, tm)
699 634
700 635 if self.isPlotConfig is False:
701 636 self.__setup_plot()
702 637 self.isPlotConfig = True
638 if self.xaxis == 'time':
639 dt = self.getDateTime(tm)
640 if self.xmin is None:
641 self.tmin = tm
642 self.xmin = dt.hour
643 minutes = (self.xmin-int(self.xmin)) * 60
644 seconds = (minutes - int(minutes)) * 60
645 self.tmin = (dt.replace(hour=int(self.xmin), minute=int(minutes), second=int(seconds)) -
646 datetime.datetime(1970, 1, 1)).total_seconds()
647 if self.localtime:
648 self.tmin += time.timezone
649
650 if self.xmin is not None and self.xmax is not None:
651 self.xrange = self.xmax - self.xmin
703 652
704 653 if self.throttle == 0:
705 654 self.__plot()
706 655 else:
707 656 self.__throttle_plot(self.__plot)#, coerce=coerce)
708 657
709 658 def close(self):
710 659
711 660 if self.data and not self.data.flagNoData:
712 661 self.save_counter = self.save_period
713 662 self.__plot()
714 663 if self.data and not self.data.flagNoData and self.pause:
715 664 figpause(10)
716 665
Show file before | Show file after | Hide comments
@@ -1,346 +1,339
1 1 import os
2 2 import datetime
3 3 import numpy
4 4
5 5 from schainpy.model.graphics.jroplot_base import Plot, plt
6 6 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot
7 7 from schainpy.utils import log
8 8
9 9 EARTH_RADIUS = 6.3710e3
10 10
11 11
12 12 def ll2xy(lat1, lon1, lat2, lon2):
13 13
14 14 p = 0.017453292519943295
15 15 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
16 16 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
17 17 r = 12742 * numpy.arcsin(numpy.sqrt(a))
18 18 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
19 19 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
20 20 theta = -theta + numpy.pi/2
21 21 return r*numpy.cos(theta), r*numpy.sin(theta)
22 22
23 23
24 24 def km2deg(km):
25 25 '''
26 26 Convert distance in km to degrees
27 27 '''
28 28
29 29 return numpy.rad2deg(km/EARTH_RADIUS)
30 30
31 31
32 32
33 33 class SpectralMomentsPlot(SpectraPlot):
34 34 '''
35 35 Plot for Spectral Moments
36 36 '''
37 37 CODE = 'spc_moments'
38 38 colormap = 'jet'
39 plot_name = 'SpectralMoments'
40 39 plot_type = 'pcolor'
41 40
42 41
43 42 class SnrPlot(RTIPlot):
44 43 '''
45 44 Plot for SNR Data
46 45 '''
47 46
48 47 CODE = 'snr'
49 48 colormap = 'jet'
50 plot_name = 'SNR'
51 49
52 50
53 51 class DopplerPlot(RTIPlot):
54 52 '''
55 53 Plot for DOPPLER Data (1st moment)
56 54 '''
57 55
58 56 CODE = 'dop'
59 57 colormap = 'jet'
60 plot_name = 'DopplerShift'
61 58
62 59
63 60 class PowerPlot(RTIPlot):
64 61 '''
65 62 Plot for Power Data (0 moment)
66 63 '''
67 64
68 65 CODE = 'pow'
69 66 colormap = 'jet'
70 plot_name = 'TotalPower'
71 67
72 68
73 69 class SpectralWidthPlot(RTIPlot):
74 70 '''
75 71 Plot for Spectral Width Data (2nd moment)
76 72 '''
77 73
78 74 CODE = 'width'
79 75 colormap = 'jet'
80 plot_name = 'SpectralWidth'
81 76
82 77
83 78 class SkyMapPlot(Plot):
84 79 '''
85 80 Plot for meteors detection data
86 81 '''
87 82
88 83 CODE = 'param'
89 84
90 85 def setup(self):
91 86
92 87 self.ncols = 1
93 88 self.nrows = 1
94 89 self.width = 7.2
95 90 self.height = 7.2
96 91 self.nplots = 1
97 92 self.xlabel = 'Zonal Zenith Angle (deg)'
98 93 self.ylabel = 'Meridional Zenith Angle (deg)'
99 94 self.polar = True
100 95 self.ymin = -180
101 96 self.ymax = 180
102 97 self.colorbar = False
103 98
104 99 def plot(self):
105 100
106 101 arrayParameters = numpy.concatenate(self.data['param'])
107 102 error = arrayParameters[:, -1]
108 103 indValid = numpy.where(error == 0)[0]
109 104 finalMeteor = arrayParameters[indValid, :]
110 105 finalAzimuth = finalMeteor[:, 3]
111 106 finalZenith = finalMeteor[:, 4]
112 107
113 108 x = finalAzimuth * numpy.pi / 180
114 109 y = finalZenith
115 110
116 111 ax = self.axes[0]
117 112
118 113 if ax.firsttime:
119 114 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
120 115 else:
121 116 ax.plot.set_data(x, y)
122 117
123 118 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
124 119 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
125 120 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
126 121 dt2,
127 122 len(x))
128 123 self.titles[0] = title
129 124
130 125
131 126 class ParametersPlot(RTIPlot):
132 127 '''
133 128 Plot for data_param object
134 129 '''
135 130
136 131 CODE = 'param'
137 132 colormap = 'seismic'
138 plot_name = 'Parameters'
139 133
140 134 def setup(self):
141 135 self.xaxis = 'time'
142 136 self.ncols = 1
143 137 self.nrows = self.data.shape(self.CODE)[0]
144 138 self.nplots = self.nrows
145 139 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
146 140
147 141 if not self.xlabel:
148 142 self.xlabel = 'Time'
149 143
150 144 if self.showSNR:
151 145 self.nrows += 1
152 146 self.nplots += 1
153 147
154 148 self.ylabel = 'Height [km]'
155 149 if not self.titles:
156 150 self.titles = self.data.parameters \
157 151 if self.data.parameters else ['Param {}'.format(x) for x in range(self.nrows)]
158 152 if self.showSNR:
159 153 self.titles.append('SNR')
160 154
161 155 def plot(self):
162 156 self.data.normalize_heights()
163 157 self.x = self.data.times
164 158 self.y = self.data.heights
165 159 if self.showSNR:
166 160 self.z = numpy.concatenate(
167 161 (self.data[self.CODE], self.data['snr'])
168 162 )
169 163 else:
170 164 self.z = self.data[self.CODE]
171 165
172 166 self.z = numpy.ma.masked_invalid(self.z)
173 167
174 168 if self.decimation is None:
175 169 x, y, z = self.fill_gaps(self.x, self.y, self.z)
176 170 else:
177 171 x, y, z = self.fill_gaps(*self.decimate())
178 172
179 173 for n, ax in enumerate(self.axes):
180 174
181 175 self.zmax = self.zmax if self.zmax is not None else numpy.max(
182 176 self.z[n])
183 177 self.zmin = self.zmin if self.zmin is not None else numpy.min(
184 178 self.z[n])
185 179
186 180 if ax.firsttime:
187 181 if self.zlimits is not None:
188 182 self.zmin, self.zmax = self.zlimits[n]
189 183
190 184 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
191 185 vmin=self.zmin,
192 186 vmax=self.zmax,
193 187 cmap=self.cmaps[n]
194 188 )
195 189 else:
196 190 if self.zlimits is not None:
197 191 self.zmin, self.zmax = self.zlimits[n]
198 192 ax.collections.remove(ax.collections[0])
199 193 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
200 194 vmin=self.zmin,
201 195 vmax=self.zmax,
202 196 cmap=self.cmaps[n]
203 197 )
204 198
205 199
206 200 class OutputPlot(ParametersPlot):
207 201 '''
208 202 Plot data_output object
209 203 '''
210 204
211 205 CODE = 'output'
212 206 colormap = 'seismic'
213 plot_name = 'Output'
214 207
215 208
216 209 class PolarMapPlot(Plot):
217 210 '''
218 211 Plot for weather radar
219 212 '''
220 213
221 214 CODE = 'param'
222 215 colormap = 'seismic'
223 216
224 217 def setup(self):
225 218 self.ncols = 1
226 219 self.nrows = 1
227 220 self.width = 9
228 221 self.height = 8
229 222 self.mode = self.data.meta['mode']
230 223 if self.channels is not None:
231 224 self.nplots = len(self.channels)
232 225 self.nrows = len(self.channels)
233 226 else:
234 227 self.nplots = self.data.shape(self.CODE)[0]
235 228 self.nrows = self.nplots
236 229 self.channels = list(range(self.nplots))
237 230 if self.mode == 'E':
238 231 self.xlabel = 'Longitude'
239 232 self.ylabel = 'Latitude'
240 233 else:
241 234 self.xlabel = 'Range (km)'
242 235 self.ylabel = 'Height (km)'
243 236 self.bgcolor = 'white'
244 237 self.cb_labels = self.data.meta['units']
245 238 self.lat = self.data.meta['latitude']
246 239 self.lon = self.data.meta['longitude']
247 240 self.xmin, self.xmax = float(
248 241 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
249 242 self.ymin, self.ymax = float(
250 243 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
251 244 # self.polar = True
252 245
253 246 def plot(self):
254 247
255 248 for n, ax in enumerate(self.axes):
256 249 data = self.data['param'][self.channels[n]]
257 250
258 251 zeniths = numpy.linspace(
259 252 0, self.data.meta['max_range'], data.shape[1])
260 253 if self.mode == 'E':
261 254 azimuths = -numpy.radians(self.data.heights)+numpy.pi/2
262 255 r, theta = numpy.meshgrid(zeniths, azimuths)
263 256 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
264 257 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
265 258 x = km2deg(x) + self.lon
266 259 y = km2deg(y) + self.lat
267 260 else:
268 261 azimuths = numpy.radians(self.data.heights)
269 262 r, theta = numpy.meshgrid(zeniths, azimuths)
270 263 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
271 264 self.y = zeniths
272 265
273 266 if ax.firsttime:
274 267 if self.zlimits is not None:
275 268 self.zmin, self.zmax = self.zlimits[n]
276 269 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
277 270 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
278 271 vmin=self.zmin,
279 272 vmax=self.zmax,
280 273 cmap=self.cmaps[n])
281 274 else:
282 275 if self.zlimits is not None:
283 276 self.zmin, self.zmax = self.zlimits[n]
284 277 ax.collections.remove(ax.collections[0])
285 278 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
286 279 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
287 280 vmin=self.zmin,
288 281 vmax=self.zmax,
289 282 cmap=self.cmaps[n])
290 283
291 284 if self.mode == 'A':
292 285 continue
293 286
294 287 # plot district names
295 288 f = open('/data/workspace/schain_scripts/distrito.csv')
296 289 for line in f:
297 290 label, lon, lat = [s.strip() for s in line.split(',') if s]
298 291 lat = float(lat)
299 292 lon = float(lon)
300 293 # ax.plot(lon, lat, '.b', ms=2)
301 294 ax.text(lon, lat, label.decode('utf8'), ha='center',
302 295 va='bottom', size='8', color='black')
303 296
304 297 # plot limites
305 298 limites = []
306 299 tmp = []
307 300 for line in open('/data/workspace/schain_scripts/lima.csv'):
308 301 if '#' in line:
309 302 if tmp:
310 303 limites.append(tmp)
311 304 tmp = []
312 305 continue
313 306 values = line.strip().split(',')
314 307 tmp.append((float(values[0]), float(values[1])))
315 308 for points in limites:
316 309 ax.add_patch(
317 310 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
318 311
319 312 # plot Cuencas
320 313 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
321 314 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
322 315 values = [line.strip().split(',') for line in f]
323 316 points = [(float(s[0]), float(s[1])) for s in values]
324 317 ax.add_patch(Polygon(points, ec='b', fc='none'))
325 318
326 319 # plot grid
327 320 for r in (15, 30, 45, 60):
328 321 ax.add_artist(plt.Circle((self.lon, self.lat),
329 322 km2deg(r), color='0.6', fill=False, lw=0.2))
330 323 ax.text(
331 324 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
332 325 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
333 326 '{}km'.format(r),
334 327 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
335 328
336 329 if self.mode == 'E':
337 330 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
338 331 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
339 332 else:
340 333 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
341 334 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
342 335
343 336 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
344 337 self.titles = ['{} {}'.format(
345 338 self.data.parameters[x], title) for x in self.channels]
346 339
Show file before | Show file after | Hide comments
@@ -1,651 +1,643
1 1 '''
2 2 Created on Jul 9, 2014
3 3 Modified on May 10, 2020
4 4
5 5 @author: Juan C. Espinoza
6 6 '''
7 7
8 8 import os
9 9 import datetime
10 10 import numpy
11 11
12 12 from schainpy.model.graphics.jroplot_base import Plot, plt
13 13
14 14
15 15 class SpectraPlot(Plot):
16 16 '''
17 17 Plot for Spectra data
18 18 '''
19 19
20 20 CODE = 'spc'
21 21 colormap = 'jet'
22 plot_name = 'Spectra'
23 22 plot_type = 'pcolor'
24 23
25 24 def setup(self):
26 25 self.nplots = len(self.data.channels)
27 26 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
28 27 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
29 self.height = 3 * self.nrows
28 self.height = 2.6 * self.nrows
30 29 self.cb_label = 'dB'
31 30 if self.showprofile:
32 31 self.width = 4 * self.ncols
33 32 else:
34 33 self.width = 3.5 * self.ncols
35 34 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
36 35 self.ylabel = 'Range [km]'
37 36
38 37 def plot(self):
39 38 if self.xaxis == "frequency":
40 39 x = self.data.xrange[0]
41 40 self.xlabel = "Frequency (kHz)"
42 41 elif self.xaxis == "time":
43 42 x = self.data.xrange[1]
44 43 self.xlabel = "Time (ms)"
45 44 else:
46 45 x = self.data.xrange[2]
47 46 self.xlabel = "Velocity (m/s)"
48 47
49 48 if self.CODE == 'spc_moments':
50 49 x = self.data.xrange[2]
51 50 self.xlabel = "Velocity (m/s)"
52 51
53 52 self.titles = []
54 53
55 54 y = self.data.heights
56 55 self.y = y
57 56 z = self.data['spc']
58 57
59 58 for n, ax in enumerate(self.axes):
60 59 noise = self.data['noise'][n][-1]
61 60 if self.CODE == 'spc_moments':
62 61 mean = self.data['moments'][n, :, 1, :][-1]
63 62 if ax.firsttime:
64 63 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
65 64 self.xmin = self.xmin if self.xmin else -self.xmax
66 65 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
67 66 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
68 67 ax.plt = ax.pcolormesh(x, y, z[n].T,
69 68 vmin=self.zmin,
70 69 vmax=self.zmax,
71 70 cmap=plt.get_cmap(self.colormap)
72 71 )
73 72
74 73 if self.showprofile:
75 74 ax.plt_profile = self.pf_axes[n].plot(
76 75 self.data['rti'][n][-1], y)[0]
77 76 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
78 77 color="k", linestyle="dashed", lw=1)[0]
79 78 if self.CODE == 'spc_moments':
80 79 ax.plt_mean = ax.plot(mean, y, color='k')[0]
81 80 else:
82 81 ax.plt.set_array(z[n].T.ravel())
83 82 if self.showprofile:
84 83 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
85 84 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
86 85 if self.CODE == 'spc_moments':
87 86 ax.plt_mean.set_data(mean, y)
88 87 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
89 88
90 89
91 90 class CrossSpectraPlot(Plot):
92 91
93 92 CODE = 'cspc'
94 93 colormap = 'jet'
95 plot_name = 'CrossSpectra'
96 94 plot_type = 'pcolor'
97 95 zmin_coh = None
98 96 zmax_coh = None
99 97 zmin_phase = None
100 98 zmax_phase = None
101 99
102 100 def setup(self):
103 101
104 102 self.ncols = 4
105 103 self.nrows = len(self.data.pairs)
106 104 self.nplots = self.nrows * 4
107 self.width = 3.4 * self.ncols
108 self.height = 3 * self.nrows
105 self.width = 3.1 * self.ncols
106 self.height = 2.6 * self.nrows
109 107 self.ylabel = 'Range [km]'
110 108 self.showprofile = False
111 self.plots_adjust.update({'bottom': 0.08})
109 self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
112 110
113 111 def plot(self):
114 112
115 113 if self.xaxis == "frequency":
116 114 x = self.data.xrange[0]
117 115 self.xlabel = "Frequency (kHz)"
118 116 elif self.xaxis == "time":
119 117 x = self.data.xrange[1]
120 118 self.xlabel = "Time (ms)"
121 119 else:
122 120 x = self.data.xrange[2]
123 121 self.xlabel = "Velocity (m/s)"
124 122
125 123 self.titles = []
126 124
127 125 y = self.data.heights
128 126 self.y = y
129 127 nspc = self.data['spc']
130 128 spc = self.data['cspc'][0]
131 129 cspc = self.data['cspc'][1]
132 130
133 131 for n in range(self.nrows):
134 132 noise = self.data['noise'][:,-1]
135 133 pair = self.data.pairs[n]
136 134 ax = self.axes[4 * n]
137 135 if ax.firsttime:
138 136 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
139 137 self.xmin = self.xmin if self.xmin else -self.xmax
140 138 self.zmin = self.zmin if self.zmin else numpy.nanmin(nspc)
141 139 self.zmax = self.zmax if self.zmax else numpy.nanmax(nspc)
142 140 ax.plt = ax.pcolormesh(x , y , nspc[pair[0]].T,
143 141 vmin=self.zmin,
144 142 vmax=self.zmax,
145 143 cmap=plt.get_cmap(self.colormap)
146 144 )
147 145 else:
148 146 ax.plt.set_array(nspc[pair[0]].T.ravel())
149 147 self.titles.append('CH {}: {:3.2f}dB'.format(pair[0], noise[pair[0]]))
150 148
151 149 ax = self.axes[4 * n + 1]
152 150 if ax.firsttime:
153 151 ax.plt = ax.pcolormesh(x , y, nspc[pair[1]].T,
154 152 vmin=self.zmin,
155 153 vmax=self.zmax,
156 154 cmap=plt.get_cmap(self.colormap)
157 155 )
158 156 else:
159 157 ax.plt.set_array(nspc[pair[1]].T.ravel())
160 158 self.titles.append('CH {}: {:3.2f}dB'.format(pair[1], noise[pair[1]]))
161 159
162 160 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
163 161 coh = numpy.abs(out)
164 162 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
165 163
166 164 ax = self.axes[4 * n + 2]
167 165 if ax.firsttime:
168 166 ax.plt = ax.pcolormesh(x, y, coh.T,
169 167 vmin=0,
170 168 vmax=1,
171 169 cmap=plt.get_cmap(self.colormap_coh)
172 170 )
173 171 else:
174 172 ax.plt.set_array(coh.T.ravel())
175 173 self.titles.append(
176 174 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
177 175
178 176 ax = self.axes[4 * n + 3]
179 177 if ax.firsttime:
180 178 ax.plt = ax.pcolormesh(x, y, phase.T,
181 179 vmin=-180,
182 180 vmax=180,
183 181 cmap=plt.get_cmap(self.colormap_phase)
184 182 )
185 183 else:
186 184 ax.plt.set_array(phase.T.ravel())
187 185 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
188 186
189 187
190 188 class RTIPlot(Plot):
191 189 '''
192 190 Plot for RTI data
193 191 '''
194 192
195 193 CODE = 'rti'
196 194 colormap = 'jet'
197 plot_name = 'RTI'
198 195 plot_type = 'pcolorbuffer'
199 196
200 197 def setup(self):
201 198 self.xaxis = 'time'
202 199 self.ncols = 1
203 200 self.nrows = len(self.data.channels)
204 201 self.nplots = len(self.data.channels)
205 202 self.ylabel = 'Range [km]'
206 203 self.xlabel = 'Time'
207 204 self.cb_label = 'dB'
208 205 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95})
209 206 self.titles = ['{} Channel {}'.format(
210 207 self.CODE.upper(), x) for x in range(self.nrows)]
211 208
212 209 def plot(self):
213 210 self.x = self.data.times
214 211 self.y = self.data.heights
215 212 self.z = self.data[self.CODE]
216 213 self.z = numpy.ma.masked_invalid(self.z)
217 214
218 215 if self.decimation is None:
219 216 x, y, z = self.fill_gaps(self.x, self.y, self.z)
220 217 else:
221 218 x, y, z = self.fill_gaps(*self.decimate())
222 219
223 220 for n, ax in enumerate(self.axes):
224 221 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
225 222 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
226 223 if ax.firsttime:
227 224 ax.plt = ax.pcolormesh(x, y, z[n].T,
228 225 vmin=self.zmin,
229 226 vmax=self.zmax,
230 227 cmap=plt.get_cmap(self.colormap)
231 228 )
232 229 if self.showprofile:
233 230 ax.plot_profile = self.pf_axes[n].plot(
234 231 self.data['rti'][n][-1], self.y)[0]
235 232 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
236 233 color="k", linestyle="dashed", lw=1)[0]
237 234 else:
238 235 ax.collections.remove(ax.collections[0])
239 236 ax.plt = ax.pcolormesh(x, y, z[n].T,
240 237 vmin=self.zmin,
241 238 vmax=self.zmax,
242 239 cmap=plt.get_cmap(self.colormap)
243 240 )
244 241 if self.showprofile:
245 242 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
246 243 ax.plot_noise.set_data(numpy.repeat(
247 244 self.data['noise'][n][-1], len(self.y)), self.y)
248 245
249 246
250 247 class CoherencePlot(RTIPlot):
251 248 '''
252 249 Plot for Coherence data
253 250 '''
254 251
255 252 CODE = 'coh'
256 plot_name = 'Coherence'
257 253
258 254 def setup(self):
259 255 self.xaxis = 'time'
260 256 self.ncols = 1
261 257 self.nrows = len(self.data.pairs)
262 258 self.nplots = len(self.data.pairs)
263 259 self.ylabel = 'Range [km]'
264 260 self.xlabel = 'Time'
265 261 self.plots_adjust.update({'hspace':0.6, 'left': 0.1, 'bottom': 0.1,'right':0.95})
266 262 if self.CODE == 'coh':
267 263 self.cb_label = ''
268 264 self.titles = [
269 265 'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
270 266 else:
271 267 self.cb_label = 'Degrees'
272 268 self.titles = [
273 269 'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
274 270
275 271
276 272 class PhasePlot(CoherencePlot):
277 273 '''
278 274 Plot for Phase map data
279 275 '''
280 276
281 277 CODE = 'phase'
282 278 colormap = 'seismic'
283 plot_name = 'Phase'
284 279
285 280
286 281 class NoisePlot(Plot):
287 282 '''
288 283 Plot for noise
289 284 '''
290 285
291 286 CODE = 'noise'
292 plot_name = 'Noise'
293 287 plot_type = 'scatterbuffer'
294 288
295 289
296 290 def setup(self):
297 291 self.xaxis = 'time'
298 292 self.ncols = 1
299 293 self.nrows = 1
300 294 self.nplots = 1
301 295 self.ylabel = 'Intensity [dB]'
302 296 self.xlabel = 'Time'
303 297 self.titles = ['Noise']
304 298 self.colorbar = False
305 299
306 300 def plot(self):
307 301
308 302 x = self.data.times
309 303 xmin = self.data.min_time
310 304 xmax = xmin + self.xrange * 60 * 60
311 305 Y = self.data[self.CODE]
312 306
313 307 if self.axes[0].firsttime:
314 308 for ch in self.data.channels:
315 309 y = Y[ch]
316 310 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
317 311 plt.legend()
318 312 else:
319 313 for ch in self.data.channels:
320 314 y = Y[ch]
321 315 self.axes[0].lines[ch].set_data(x, y)
322 316
323 317 self.ymin = numpy.nanmin(Y) - 5
324 318 self.ymax = numpy.nanmax(Y) + 5
325 319
326 320
327 321 class PowerProfilePlot(Plot):
328 322
329 323 CODE = 'spcprofile'
330 plot_name = 'Power Profile'
331 324 plot_type = 'scatter'
332 325 buffering = False
333 326
334 327 def setup(self):
335 328
336 329 self.ncols = 1
337 330 self.nrows = 1
338 331 self.nplots = 1
339 332 self.height = 4
340 333 self.width = 3
341 334 self.ylabel = 'Range [km]'
342 335 self.xlabel = 'Intensity [dB]'
343 336 self.titles = ['Power Profile']
344 337 self.colorbar = False
345 338
346 339 def plot(self):
347 340
348 341 y = self.data.heights
349 342 self.y = y
350 343
351 344 x = self.data['spcprofile']
352 345
353 346 if self.xmin is None: self.xmin = numpy.nanmin(x)*0.9
354 347 if self.xmax is None: self.xmax = numpy.nanmax(x)*1.1
355 348
356 349 if self.axes[0].firsttime:
357 350 for ch in self.data.channels:
358 351 self.axes[0].plot(x[ch], y, lw=1, label='Ch{}'.format(ch))
359 352 plt.legend()
360 353 else:
361 354 for ch in self.data.channels:
362 355 self.axes[0].lines[ch].set_data(x[ch], y)
363 356
364 357
365 358 class SpectraCutPlot(Plot):
366 359
367 360 CODE = 'spc_cut'
368 plot_name = 'Spectra Cut'
369 361 plot_type = 'scatter'
370 362 buffering = False
371 363
372 364 def setup(self):
373 365
374 366 self.nplots = len(self.data.channels)
375 367 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
376 368 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
377 369 self.width = 3.4 * self.ncols + 1.5
378 370 self.height = 3 * self.nrows
379 371 self.ylabel = 'Power [dB]'
380 372 self.colorbar = False
381 373 self.plots_adjust.update({'left':0.1, 'hspace':0.3, 'right': 0.75, 'bottom':0.08})
382 374
383 375 def plot(self):
384 376 if self.xaxis == "frequency":
385 377 x = self.data.xrange[0][1:]
386 378 self.xlabel = "Frequency (kHz)"
387 379 elif self.xaxis == "time":
388 380 x = self.data.xrange[1]
389 381 self.xlabel = "Time (ms)"
390 382 else:
391 383 x = self.data.xrange[2]
392 384 self.xlabel = "Velocity (m/s)"
393 385
394 386 self.titles = []
395 387
396 388 y = self.data.heights
397 389 #self.y = y
398 390 z = self.data['spc_cut']
399 391
400 392 if self.height_index:
401 393 index = numpy.array(self.height_index)
402 394 else:
403 395 index = numpy.arange(0, len(y), int((len(y))/9))
404 396
405 397 for n, ax in enumerate(self.axes):
406 398 if ax.firsttime:
407 399 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
408 400 self.xmin = self.xmin if self.xmin else -self.xmax
409 401 self.ymin = self.ymin if self.ymin else numpy.nanmin(z)
410 402 self.ymax = self.ymax if self.ymax else numpy.nanmax(z)
411 403 ax.plt = ax.plot(x, z[n, :, index].T)
412 404 labels = ['Range = {:2.1f}km'.format(y[i]) for i in index]
413 405 self.figures[0].legend(ax.plt, labels, loc='center right')
414 406 else:
415 407 for i, line in enumerate(ax.plt):
416 408 line.set_data(x, z[n, :, i])
417 409 self.titles.append('CH {}'.format(n))
418 410
419 411
420 412 class BeaconPhase(Plot):
421 413
422 414 __isConfig = None
423 415 __nsubplots = None
424 416
425 417 PREFIX = 'beacon_phase'
426 418
427 419 def __init__(self):
428 420 Plot.__init__(self)
429 421 self.timerange = 24*60*60
430 422 self.isConfig = False
431 423 self.__nsubplots = 1
432 424 self.counter_imagwr = 0
433 425 self.WIDTH = 800
434 426 self.HEIGHT = 400
435 427 self.WIDTHPROF = 120
436 428 self.HEIGHTPROF = 0
437 429 self.xdata = None
438 430 self.ydata = None
439 431
440 432 self.PLOT_CODE = BEACON_CODE
441 433
442 434 self.FTP_WEI = None
443 435 self.EXP_CODE = None
444 436 self.SUB_EXP_CODE = None
445 437 self.PLOT_POS = None
446 438
447 439 self.filename_phase = None
448 440
449 441 self.figfile = None
450 442
451 443 self.xmin = None
452 444 self.xmax = None
453 445
454 446 def getSubplots(self):
455 447
456 448 ncol = 1
457 449 nrow = 1
458 450
459 451 return nrow, ncol
460 452
461 453 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
462 454
463 455 self.__showprofile = showprofile
464 456 self.nplots = nplots
465 457
466 458 ncolspan = 7
467 459 colspan = 6
468 460 self.__nsubplots = 2
469 461
470 462 self.createFigure(id = id,
471 463 wintitle = wintitle,
472 464 widthplot = self.WIDTH+self.WIDTHPROF,
473 465 heightplot = self.HEIGHT+self.HEIGHTPROF,
474 466 show=show)
475 467
476 468 nrow, ncol = self.getSubplots()
477 469
478 470 self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
479 471
480 472 def save_phase(self, filename_phase):
481 473 f = open(filename_phase,'w+')
482 474 f.write('\n\n')
483 475 f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
484 476 f.write('DD MM YYYY HH MM SS pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n' )
485 477 f.close()
486 478
487 479 def save_data(self, filename_phase, data, data_datetime):
488 480 f=open(filename_phase,'a')
489 481 timetuple_data = data_datetime.timetuple()
490 482 day = str(timetuple_data.tm_mday)
491 483 month = str(timetuple_data.tm_mon)
492 484 year = str(timetuple_data.tm_year)
493 485 hour = str(timetuple_data.tm_hour)
494 486 minute = str(timetuple_data.tm_min)
495 487 second = str(timetuple_data.tm_sec)
496 488 f.write(day+' '+month+' '+year+' '+hour+' '+minute+' '+second+' '+str(data[0])+' '+str(data[1])+' '+str(data[2])+' '+str(data[3])+'\n')
497 489 f.close()
498 490
499 491 def plot(self):
500 492 log.warning('TODO: Not yet implemented...')
501 493
502 494 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
503 495 xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
504 496 timerange=None,
505 497 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
506 498 server=None, folder=None, username=None, password=None,
507 499 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
508 500
509 501 if dataOut.flagNoData:
510 502 return dataOut
511 503
512 504 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
513 505 return
514 506
515 507 if pairsList == None:
516 508 pairsIndexList = dataOut.pairsIndexList[:10]
517 509 else:
518 510 pairsIndexList = []
519 511 for pair in pairsList:
520 512 if pair not in dataOut.pairsList:
521 513 raise ValueError("Pair %s is not in dataOut.pairsList" %(pair))
522 514 pairsIndexList.append(dataOut.pairsList.index(pair))
523 515
524 516 if pairsIndexList == []:
525 517 return
526 518
527 519 # if len(pairsIndexList) > 4:
528 520 # pairsIndexList = pairsIndexList[0:4]
529 521
530 522 hmin_index = None
531 523 hmax_index = None
532 524
533 525 if hmin != None and hmax != None:
534 526 indexes = numpy.arange(dataOut.nHeights)
535 527 hmin_list = indexes[dataOut.heightList >= hmin]
536 528 hmax_list = indexes[dataOut.heightList <= hmax]
537 529
538 530 if hmin_list.any():
539 531 hmin_index = hmin_list[0]
540 532
541 533 if hmax_list.any():
542 534 hmax_index = hmax_list[-1]+1
543 535
544 536 x = dataOut.getTimeRange()
545 537 #y = dataOut.getHeiRange()
546 538
547 539
548 540 thisDatetime = dataOut.datatime
549 541
550 542 title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
551 543 xlabel = "Local Time"
552 544 ylabel = "Phase (degrees)"
553 545
554 546 update_figfile = False
555 547
556 548 nplots = len(pairsIndexList)
557 549 #phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
558 550 phase_beacon = numpy.zeros(len(pairsIndexList))
559 551 for i in range(nplots):
560 552 pair = dataOut.pairsList[pairsIndexList[i]]
561 553 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
562 554 powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
563 555 powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
564 556 avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
565 557 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
566 558
567 559 if dataOut.beacon_heiIndexList:
568 560 phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
569 561 else:
570 562 phase_beacon[i] = numpy.average(phase)
571 563
572 564 if not self.isConfig:
573 565
574 566 nplots = len(pairsIndexList)
575 567
576 568 self.setup(id=id,
577 569 nplots=nplots,
578 570 wintitle=wintitle,
579 571 showprofile=showprofile,
580 572 show=show)
581 573
582 574 if timerange != None:
583 575 self.timerange = timerange
584 576
585 577 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
586 578
587 579 if ymin == None: ymin = 0
588 580 if ymax == None: ymax = 360
589 581
590 582 self.FTP_WEI = ftp_wei
591 583 self.EXP_CODE = exp_code
592 584 self.SUB_EXP_CODE = sub_exp_code
593 585 self.PLOT_POS = plot_pos
594 586
595 587 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
596 588 self.isConfig = True
597 589 self.figfile = figfile
598 590 self.xdata = numpy.array([])
599 591 self.ydata = numpy.array([])
600 592
601 593 update_figfile = True
602 594
603 595 #open file beacon phase
604 596 path = '%s%03d' %(self.PREFIX, self.id)
605 597 beacon_file = os.path.join(path,'%s.txt'%self.name)
606 598 self.filename_phase = os.path.join(figpath,beacon_file)
607 599 #self.save_phase(self.filename_phase)
608 600
609 601
610 602 #store data beacon phase
611 603 #self.save_data(self.filename_phase, phase_beacon, thisDatetime)
612 604
613 605 self.setWinTitle(title)
614 606
615 607
616 608 title = "Phase Plot %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
617 609
618 610 legendlabels = ["Pair (%d,%d)"%(pair[0], pair[1]) for pair in dataOut.pairsList]
619 611
620 612 axes = self.axesList[0]
621 613
622 614 self.xdata = numpy.hstack((self.xdata, x[0:1]))
623 615
624 616 if len(self.ydata)==0:
625 617 self.ydata = phase_beacon.reshape(-1,1)
626 618 else:
627 619 self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1,1)))
628 620
629 621
630 622 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
631 623 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
632 624 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
633 625 XAxisAsTime=True, grid='both'
634 626 )
635 627
636 628 self.draw()
637 629
638 630 if dataOut.ltctime >= self.xmax:
639 631 self.counter_imagwr = wr_period
640 632 self.isConfig = False
641 633 update_figfile = True
642 634
643 635 self.save(figpath=figpath,
644 636 figfile=figfile,
645 637 save=save,
646 638 ftp=ftp,
647 639 wr_period=wr_period,
648 640 thisDatetime=thisDatetime,
649 641 update_figfile=update_figfile)
650 642
651 643 return dataOut No newline at end of file
Show file before | Show file after | Hide comments
@@ -1,302 +1,297
1 1 '''
2 2 Created on Jul 9, 2014
3 3
4 4 @author: roj-idl71
5 5 '''
6 6 import os
7 7 import datetime
8 8 import numpy
9 9
10 10 from schainpy.model.graphics.jroplot_base import Plot, plt
11 11
12 12
13 13 class ScopePlot(Plot):
14 14
15 15 '''
16 16 Plot for Scope
17 17 '''
18 18
19 19 CODE = 'scope'
20 plot_name = 'Scope'
21 20 plot_type = 'scatter'
22 21
23 22 def setup(self):
24 23
25 24 self.xaxis = 'Range (Km)'
26 25 self.ncols = 1
27 26 self.nrows = 1
28 27 self.nplots = 1
29 28 self.ylabel = 'Intensity [dB]'
30 29 self.titles = ['Scope']
31 30 self.colorbar = False
32 31 self.width = 6
33 32 self.height = 4
34 33
35 34 def plot_iq(self, x, y, channelIndexList, thisDatetime, wintitle):
36 35
37 36 yreal = y[channelIndexList,:].real
38 37 yimag = y[channelIndexList,:].imag
39 38 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
40 39 self.xlabel = "Range (Km)"
41 40 self.ylabel = "Intensity - IQ"
42 41
43 42 self.y = yreal
44 43 self.x = x
45 44 self.xmin = min(x)
46 45 self.xmax = max(x)
47 46
48 47
49 48 self.titles[0] = title
50 49
51 50 for i,ax in enumerate(self.axes):
52 51 title = "Channel %d" %(i)
53 52 if ax.firsttime:
54 53 ax.plt_r = ax.plot(x, yreal[i,:], color='b')[0]
55 54 ax.plt_i = ax.plot(x, yimag[i,:], color='r')[0]
56 55 else:
57 56 ax.plt_r.set_data(x, yreal[i,:])
58 57 ax.plt_i.set_data(x, yimag[i,:])
59 58
60 59 def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle):
61 60 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
62 61 yreal = y.real
63 62 yreal = 10*numpy.log10(yreal)
64 63 self.y = yreal
65 64 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
66 65 self.xlabel = "Range (Km)"
67 66 self.ylabel = "Intensity"
68 67 self.xmin = min(x)
69 68 self.xmax = max(x)
70 69
71 70
72 71 self.titles[0] = title
73 72
74 73 for i,ax in enumerate(self.axes):
75 74 title = "Channel %d" %(i)
76 75
77 76 ychannel = yreal[i,:]
78 77
79 78 if ax.firsttime:
80 79 ax.plt_r = ax.plot(x, ychannel)[0]
81 80 else:
82 81 #pass
83 82 ax.plt_r.set_data(x, ychannel)
84 83
85 84 def plot_weatherpower(self, x, y, channelIndexList, thisDatetime, wintitle):
86 85
87 86
88 87 y = y[channelIndexList,:]
89 88 yreal = y.real
90 89 yreal = 10*numpy.log10(yreal)
91 90 self.y = yreal
92 91 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
93 92 self.xlabel = "Range (Km)"
94 93 self.ylabel = "Intensity"
95 94 self.xmin = min(x)
96 95 self.xmax = max(x)
97 96
98 97 self.titles[0] =title
99 98 for i,ax in enumerate(self.axes):
100 99 title = "Channel %d" %(i)
101 100
102 101 ychannel = yreal[i,:]
103 102
104 103 if ax.firsttime:
105 104 ax.plt_r = ax.plot(x, ychannel)[0]
106 105 else:
107 106 #pass
108 107 ax.plt_r.set_data(x, ychannel)
109 108
110 109 def plot_weathervelocity(self, x, y, channelIndexList, thisDatetime, wintitle):
111 110
112 111 x = x[channelIndexList,:]
113 112 yreal = y
114 113 self.y = yreal
115 114 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
116 115 self.xlabel = "Velocity (m/s)"
117 116 self.ylabel = "Range (Km)"
118 117 self.xmin = numpy.min(x)
119 118 self.xmax = numpy.max(x)
120 119 self.titles[0] =title
121 120 for i,ax in enumerate(self.axes):
122 121 title = "Channel %d" %(i)
123 122 xchannel = x[i,:]
124 123 if ax.firsttime:
125 124 ax.plt_r = ax.plot(xchannel, yreal)[0]
126 125 else:
127 126 #pass
128 127 ax.plt_r.set_data(xchannel, yreal)
129 128
130 129 def plot_weatherspecwidth(self, x, y, channelIndexList, thisDatetime, wintitle):
131 130
132 131 x = x[channelIndexList,:]
133 132 yreal = y
134 133 self.y = yreal
135 134 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
136 135 self.xlabel = "width "
137 136 self.ylabel = "Range (Km)"
138 137 self.xmin = numpy.min(x)
139 138 self.xmax = numpy.max(x)
140 139 self.titles[0] =title
141 140 for i,ax in enumerate(self.axes):
142 141 title = "Channel %d" %(i)
143 142 xchannel = x[i,:]
144 143 if ax.firsttime:
145 144 ax.plt_r = ax.plot(xchannel, yreal)[0]
146 145 else:
147 146 #pass
148 147 ax.plt_r.set_data(xchannel, yreal)
149 148
150 149 def plot(self):
151 150 if self.channels:
152 151 channels = self.channels
153 152 else:
154 153 channels = self.data.channels
155 154
156 155 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1])
157 156 if self.CODE == "pp_power":
158 157 scope = self.data['pp_power']
159 158 elif self.CODE == "pp_signal":
160 159 scope = self.data["pp_signal"]
161 160 elif self.CODE == "pp_velocity":
162 161 scope = self.data["pp_velocity"]
163 162 elif self.CODE == "pp_specwidth":
164 163 scope = self.data["pp_specwidth"]
165 164 else:
166 165 scope =self.data["scope"]
167 166
168 167 if self.data.flagDataAsBlock:
169 168
170 169 for i in range(self.data.nProfiles):
171 170
172 171 wintitle1 = " [Profile = %d] " %i
173 172 if self.CODE =="scope":
174 173 if self.type == "power":
175 174 self.plot_power(self.data.heights,
176 175 scope[:,i,:],
177 176 channels,
178 177 thisDatetime,
179 178 wintitle1
180 179 )
181 180
182 181 if self.type == "iq":
183 182 self.plot_iq(self.data.heights,
184 183 scope[:,i,:],
185 184 channels,
186 185 thisDatetime,
187 186 wintitle1
188 187 )
189 188 if self.CODE=="pp_power":
190 189 self.plot_weatherpower(self.data.heights,
191 190 scope[:,i,:],
192 191 channels,
193 192 thisDatetime,
194 193 wintitle
195 194 )
196 195 if self.CODE=="pp_signal":
197 196 self.plot_weatherpower(self.data.heights,
198 197 scope[:,i,:],
199 198 channels,
200 199 thisDatetime,
201 200 wintitle
202 201 )
203 202 if self.CODE=="pp_velocity":
204 203 self.plot_weathervelocity(scope[:,i,:],
205 204 self.data.heights,
206 205 channels,
207 206 thisDatetime,
208 207 wintitle
209 208 )
210 209 if self.CODE=="pp_spcwidth":
211 210 self.plot_weatherspecwidth(scope[:,i,:],
212 211 self.data.heights,
213 212 channels,
214 213 thisDatetime,
215 214 wintitle
216 215 )
217 216 else:
218 217 wintitle = " [Profile = %d] " %self.data.profileIndex
219 218 if self.CODE== "scope":
220 219 if self.type == "power":
221 220 self.plot_power(self.data.heights,
222 221 scope,
223 222 channels,
224 223 thisDatetime,
225 224 wintitle
226 225 )
227 226
228 227 if self.type == "iq":
229 228 self.plot_iq(self.data.heights,
230 229 scope,
231 230 channels,
232 231 thisDatetime,
233 232 wintitle
234 233 )
235 234 if self.CODE=="pp_power":
236 235 self.plot_weatherpower(self.data.heights,
237 236 scope,
238 237 channels,
239 238 thisDatetime,
240 239 wintitle
241 240 )
242 241 if self.CODE=="pp_signal":
243 242 self.plot_weatherpower(self.data.heights,
244 243 scope,
245 244 channels,
246 245 thisDatetime,
247 246 wintitle
248 247 )
249 248 if self.CODE=="pp_velocity":
250 249 self.plot_weathervelocity(scope,
251 250 self.data.heights,
252 251 channels,
253 252 thisDatetime,
254 253 wintitle
255 254 )
256 255 if self.CODE=="pp_specwidth":
257 256 self.plot_weatherspecwidth(scope,
258 257 self.data.heights,
259 258 channels,
260 259 thisDatetime,
261 260 wintitle
262 261 )
263 262
264 263
265 264
266 265 class PulsepairPowerPlot(ScopePlot):
267 266 '''
268 267 Plot for P= S+N
269 268 '''
270 269
271 270 CODE = 'pp_power'
272 plot_name = 'PulsepairPower'
273 271 plot_type = 'scatter'
274 272 buffering = False
275 273
276 274 class PulsepairVelocityPlot(ScopePlot):
277 275 '''
278 276 Plot for VELOCITY
279 277 '''
280 278 CODE = 'pp_velocity'
281 plot_name = 'PulsepairVelocity'
282 279 plot_type = 'scatter'
283 280 buffering = False
284 281
285 282 class PulsepairSpecwidthPlot(ScopePlot):
286 283 '''
287 284 Plot for WIDTH
288 285 '''
289 286 CODE = 'pp_specwidth'
290 plot_name = 'PulsepairSpecwidth'
291 287 plot_type = 'scatter'
292 288 buffering = False
293 289
294 290 class PulsepairSignalPlot(ScopePlot):
295 291 '''
296 292 Plot for S
297 293 '''
298 294
299 295 CODE = 'pp_signal'
300 plot_name = 'PulsepairSignal'
301 296 plot_type = 'scatter'
302 297 buffering = False
Show file before | Show file after | Hide comments
@@ -1,326 +1,351
1 '''
2 @author: Juan C. Espinoza
3 '''
1 # Copyright (c) 2012-2020 Jicamarca Radio Observatory
2 # All rights reserved.
3 #
4 # Distributed under the terms of the BSD 3-clause license.
5 """Utilities for publish/send data, files & plots over different protocols
6 """
4 7
5 8 import os
6 9 import glob
7 10 import time
8 11 import json
9 12 import numpy
10 13 import zmq
11 14 import datetime
12 15 import ftplib
13 16 from functools import wraps
14 17 from threading import Thread
15 18 from multiprocessing import Process
16 19
17 20 from schainpy.model.proc.jroproc_base import Operation, ProcessingUnit, MPDecorator
18 21 from schainpy.model.data.jrodata import JROData
19 22 from schainpy.utils import log
20 23
21 MAXNUMX = 500
22 MAXNUMY = 500
23 24
24 25 PLOT_CODES = {
25 26 'rti': 0, # Range time intensity (RTI).
26 27 'spc': 1, # Spectra (and Cross-spectra) information.
27 28 'cspc': 2, # Cross-Correlation information.
28 29 'coh': 3, # Coherence map.
29 30 'base': 4, # Base lines graphic.
30 31 'row': 5, # Row Spectra.
31 32 'total': 6, # Total Power.
32 33 'drift': 7, # Drifts graphics.
33 34 'height': 8, # Height profile.
34 35 'phase': 9, # Signal Phase.
35 36 'power': 16,
36 37 'noise': 17,
37 38 'beacon': 18,
38 39 'wind': 22,
39 40 'skymap': 23,
40 41 'Unknown': 24,
41 42 'V-E': 25, # PIP Velocity.
42 43 'Z-E': 26, # PIP Reflectivity.
43 44 'V-A': 27, # RHI Velocity.
44 45 'Z-A': 28, # RHI Reflectivity.
45 46 }
46 47
47 48 def get_plot_code(s):
48 49 label = s.split('_')[0]
49 50 codes = [key for key in PLOT_CODES if key in label]
50 51 if codes:
51 52 return PLOT_CODES[codes[0]]
52 53 else:
53 54 return 24
54 55
55 def decimate(z, MAXNUMY):
56 dy = int(len(z[0])/MAXNUMY) + 1
57
58 return z[::, ::dy]
59
60 56
61 57 class PublishData(Operation):
62 58 '''
63 59 Operation to send data over zmq.
64 60 '''
65 61
66 62 __attrs__ = ['host', 'port', 'delay', 'verbose']
67 63
68 64 def setup(self, server='zmq.pipe', delay=0, verbose=True, **kwargs):
69 65 self.counter = 0
70 66 self.delay = kwargs.get('delay', 0)
71 67 self.cnt = 0
72 68 self.verbose = verbose
73 69 context = zmq.Context()
74 70 self.zmq_socket = context.socket(zmq.PUSH)
75 71 server = kwargs.get('server', 'zmq.pipe')
76 72
77 73 if 'tcp://' in server:
78 74 address = server
79 75 else:
80 76 address = 'ipc:///tmp/%s' % server
81 77
82 78 self.zmq_socket.connect(address)
83 79 time.sleep(1)
84 80
85 81
86 82 def publish_data(self):
87 83 self.dataOut.finished = False
88 84
89 85 if self.verbose:
90 86 log.log(
91 87 'Sending {} - {}'.format(self.dataOut.type, self.dataOut.datatime),
92 88 self.name
93 89 )
94 90 self.zmq_socket.send_pyobj(self.dataOut)
95 91
96 92 def run(self, dataOut, **kwargs):
97 93 self.dataOut = dataOut
98 94 if not self.isConfig:
99 95 self.setup(**kwargs)
100 96 self.isConfig = True
101 97
102 98 self.publish_data()
103 99 time.sleep(self.delay)
104 100
105 101 def close(self):
106 102
107 103 self.dataOut.finished = True
108 104 self.zmq_socket.send_pyobj(self.dataOut)
109 105 time.sleep(0.1)
110 106 self.zmq_socket.close()
111 107
112 108
113 109 class ReceiverData(ProcessingUnit):
114 110
115 111 __attrs__ = ['server']
116 112
117 113 def __init__(self, **kwargs):
118 114
119 115 ProcessingUnit.__init__(self, **kwargs)
120 116
121 117 self.isConfig = False
122 118 server = kwargs.get('server', 'zmq.pipe')
123 119 if 'tcp://' in server:
124 120 address = server
125 121 else:
126 122 address = 'ipc:///tmp/%s' % server
127 123
128 124 self.address = address
129 125 self.dataOut = JROData()
130 126
131 127 def setup(self):
132 128
133 129 self.context = zmq.Context()
134 130 self.receiver = self.context.socket(zmq.PULL)
135 131 self.receiver.bind(self.address)
136 132 time.sleep(0.5)
137 133 log.success('ReceiverData from {}'.format(self.address))
138 134
139 135
140 136 def run(self):
141 137
142 138 if not self.isConfig:
143 139 self.setup()
144 140 self.isConfig = True
145 141
146 142 self.dataOut = self.receiver.recv_pyobj()
147 143 log.log('{} - {}'.format(self.dataOut.type,
148 144 self.dataOut.datatime.ctime(),),
149 145 'Receiving')
150 146
151 147 @MPDecorator
152 148 class SendToFTP(Operation):
153
154 '''
155 Operation to send data over FTP.
156 patternX = 'local, remote, ext, period, exp_code, sub_exp_code'
157 '''
158
149 """Operation for send files over FTP
150
151 This operation is used to send files over FTP, you can send different files
152 from different folders by adding as many `pattern` as you wish.
153
154 Parameters:
155 -----------
156 server : str
157 FTP server address.
158 username : str
159 FTP username
160 password : str
161 FTP password
162 timeout : int
163 timeout to restart the connection
164 patternX : list
165 detail of files to be send must have the following order: local, remote
166 ext, period, exp_code, sub_exp_code
167
168 Example:
169 --------
170
171 ftp = proc_unit.addOperation(name='SendToFTP', optype='external')
172 ftp.addParameter(name='server', value='jro-app.igp.gob.pe')
173 ftp.addParameter(name='username', value='wmaster')
174 ftp.addParameter(name='password', value='mst2010vhf')
175 ftp.addParameter(
176 name='pattern1',
177 value='/local/path/rti,/remote/path,png,300,11,0'
178 )
179 ftp.addParameter(
180 name='pattern2',
181 value='/local/path/spc,/remote/path,png,300,11,0'
182 )
183 ftp.addParameter(
184 name='pattern3',
185 value='/local/path/param,/remote/path,hdf5,300,,'
186 )
187
188 """
189
159 190 __attrs__ = ['server', 'username', 'password', 'timeout', 'patternX']
160 191
161 192 def __init__(self):
162 193 '''
163 194 '''
164 195 Operation.__init__(self)
165 196 self.ftp = None
166 197 self.ready = False
167 198
168 199 def setup(self, server, username, password, timeout, **kwargs):
169 200 '''
170 201 '''
171 202
172 203 self.server = server
173 204 self.username = username
174 205 self.password = password
175 206 self.timeout = timeout
176 207 self.patterns = []
177 208 self.times = []
178 209 self.latest = []
179 210 for arg, value in kwargs.items():
180 211 if 'pattern' in arg:
181 212 self.patterns.append(value)
182 self.times.append(time.time())
213 self.times.append(0)
183 214 self.latest.append('')
184 215
185 216 def connect(self):
186 217 '''
187 218 '''
188 219
189 220 log.log('Connecting to ftp://{}'.format(self.server), self.name)
190 221 try:
191 222 self.ftp = ftplib.FTP(self.server, timeout=self.timeout)
192 223 except ftplib.all_errors:
193 224 log.error('Server connection fail: {}'.format(self.server), self.name)
194 225 if self.ftp is not None:
195 226 self.ftp.close()
196 227 self.ftp = None
197 228 self.ready = False
198 229 return
199 230
200 231 try:
201 232 self.ftp.login(self.username, self.password)
202 233 except ftplib.all_errors:
203 234 log.error('The given username y/o password are incorrect', self.name)
204 235 if self.ftp is not None:
205 236 self.ftp.close()
206 237 self.ftp = None
207 238 self.ready = False
208 239 return
209 240
210 241 log.success('Connection success', self.name)
211 242 self.ready = True
212 243 return
213 244
214 245 def check(self):
215 246
216 247 try:
217 248 self.ftp.voidcmd("NOOP")
218 249 except:
219 250 log.warning('Connection lost... trying to reconnect', self.name)
220 251 if self.ftp is not None:
221 252 self.ftp.close()
222 253 self.ftp = None
223 254 self.connect()
224 255
225 256 def find_files(self, path, ext):
226 257
227 files = glob.glob1(path, '*{}'.format(ext))
258 files = glob.glob1(path.strip(), '*{}'.format(ext.strip()))
228 259 files.sort()
229 260 if files:
230 261 return files[-1]
231 262 return None
232 263
233 264 def getftpname(self, filename, exp_code, sub_exp_code):
234 265
235 266 thisDatetime = datetime.datetime.strptime(filename.split('_')[1], '%Y%m%d')
236 267 YEAR_STR = '%4.4d' % thisDatetime.timetuple().tm_year
237 268 DOY_STR = '%3.3d' % thisDatetime.timetuple().tm_yday
238 269 exp_code = '%3.3d' % exp_code
239 270 sub_exp_code = '%2.2d' % sub_exp_code
240 271 plot_code = '%2.2d' % get_plot_code(filename)
241 272 name = YEAR_STR + DOY_STR + '00' + exp_code + sub_exp_code + plot_code + '00.png'
242 273 return name
243 274
244 275 def upload(self, src, dst):
245 276
246 277 log.log('Uploading {} -> {} '.format(
247 278 src.split('/')[-1], dst.split('/')[-1]),
248 279 self.name,
249 280 nl=False
250 281 )
251 282
252 283 fp = open(src, 'rb')
253 284 command = 'STOR {}'.format(dst)
254 285
255 286 try:
256 287 self.ftp.storbinary(command, fp, blocksize=1024)
257 288 except Exception as e:
258 289 log.error('{}'.format(e), self.name)
259 if self.ftp is not None:
260 self.ftp.close()
261 self.ftp = None
262 290 return 0
263 291
264 292 try:
265 293 self.ftp.sendcmd('SITE CHMOD 755 {}'.format(dst))
266 294 except Exception as e:
267 295 log.error('{}'.format(e), self.name)
268 if self.ftp is not None:
269 self.ftp.close()
270 self.ftp = None
271 296 return 0
272 297
273 298 fp.close()
274 299 log.success('OK', tag='')
275 300 return 1
276 301
277 302 def send_files(self):
278 303
279 304 for x, pattern in enumerate(self.patterns):
280 305 local, remote, ext, period, exp_code, sub_exp_code = pattern
281 if time.time()-self.times[x] >= int(period):
282 srcname = self.find_files(local, ext)
283 src = os.path.join(local, srcname)
284 if os.path.getmtime(src) < time.time() - 30*60:
285 log.warning('Skipping old file {}'.format(srcname))
286 continue
287
288 if srcname is None or srcname == self.latest[x]:
289 log.warning('File alreday uploaded {}'.format(srcname))
290 continue
291
292 if 'png' in ext:
293 dstname = self.getftpname(srcname, int(exp_code), int(sub_exp_code))
294 else:
295 dstname = srcname
296
297 dst = os.path.join(remote, dstname)
298
299 if self.upload(src, dst):
300 self.times[x] = time.time()
301 self.latest[x] = srcname
302 else:
303 self.ready = False
304 break
306
307 if (self.dataOut.utctime - self.times[x]) < int(period):
308 continue
309
310 srcname = self.find_files(local, ext)
311
312 if srcname is None:
313 continue
314
315 if srcname == self.latest[x]:
316 log.warning('File alreday uploaded {}'.format(srcname))
317 continue
318
319 if exp_code.strip():
320 dstname = self.getftpname(srcname, int(exp_code), int(sub_exp_code))
321 else:
322 dstname = srcname
323
324 src = os.path.join(local, srcname)
325 dst = os.path.join(remote.strip(), dstname)
326
327 if self.upload(src, dst):
328 self.times[x] = self.dataOut.utctime
329 self.latest[x] = srcname
305 330
306 331 def run(self, dataOut, server, username, password, timeout=10, **kwargs):
307 332
308 333 if not self.isConfig:
309 334 self.setup(
310 335 server=server,
311 336 username=username,
312 337 password=password,
313 338 timeout=timeout,
314 339 **kwargs
315 340 )
316 341 self.isConfig = True
317 if not self.ready:
318 342 self.connect()
319 if self.ftp is not None:
320 self.check()
321 self.send_files()
343
344 self.dataOut = dataOut
345 self.check()
346 self.send_files()
322 347
323 348 def close(self):
324 349
325 350 if self.ftp is not None:
326 351 self.ftp.close()
General Comments 0
You need to be logged in to leave comments. Login now