##// END OF EJS Templates
schain
RSS
plot extra point depending on matplotlib version
rflores -
r1732:8e0a97a787ce
parent child
Show More
Side by Side Unified
Context file:
r1732:8e0a97a787ce Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,1295 +1,1302
1 1 # Copyright (c) 2012-2021 Jicamarca Radio Observatory
2 2 # All rights reserved.
3 3 #
4 4 # Distributed under the terms of the BSD 3-clause license.
5 5 """Classes to plot Spectra data
6 6
7 7 """
8 8
9 9 import os
10 10 import numpy
11 11
12 12 from schainpy.model.graphics.jroplot_base import Plot, plt, log
13 13
14 from matplotlib import __version__ as plt_version
15
16 if plt_version >='3.3.4':
17 EXTRA_POINTS = 0
18 else:
19 EXTRA_POINTS = 1
20
14 21
15 22 class SpectraPlot(Plot):
16 23 '''
17 24 Plot for Spectra data
18 25 '''
19 26
20 27 CODE = 'spc'
21 28 colormap = 'jet'
22 29 plot_type = 'pcolor'
23 30 buffering = False
24 31
25 32 def setup(self):
26 33
27 34 self.nplots = len(self.data.channels)
28 35 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
29 36 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
30 37 self.height = 2.6 * self.nrows
31 38 self.cb_label = 'dB'
32 39 if self.showprofile:
33 40 self.width = 4 * self.ncols
34 41 else:
35 42 self.width = 3.5 * self.ncols
36 43 self.plots_adjust.update({'wspace': 0.8, 'hspace':0.2, 'left': 0.2, 'right': 0.9, 'bottom': 0.18})
37 44 self.ylabel = 'Range [km]'
38 45
39 46 def update(self, dataOut):
40 47
41 48 data = {}
42 49 meta = {}
43 50 spc = 10 * numpy.log10(dataOut.data_spc / dataOut.normFactor)
44 51 data['spc'] = spc
45 52 data['rti'] = dataOut.getPower()
46 53 if hasattr(dataOut, 'LagPlot'): #Double Pulse
47 54 max_hei_id = dataOut.nHeights - 2*dataOut.LagPlot
48 55 data['noise'] = 10*numpy.log10(dataOut.getNoise(ymin_index=53,ymax_index=max_hei_id)/dataOut.normFactor)
49 56 data['noise'][0] = 10*numpy.log10(dataOut.getNoise(ymin_index=53)[0]/dataOut.normFactor)
50 57 else:
51 58 data['noise'] = 10 * numpy.log10(dataOut.getNoise() / dataOut.normFactor)
52 extrapoints = spc.shape[1] % dataOut.nFFTPoints
53 extrapoints=1
54 meta['xrange'] = (dataOut.getFreqRange(extrapoints) / 1000., dataOut.getAcfRange(extrapoints), dataOut.getVelRange(extrapoints))
59 extrapoints = spc.shape[1] % dataOut.nFFTPoints
60 extrapoints=1
61 meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS) / 1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
55 62 if self.CODE == 'spc_moments':
56 63 data['moments'] = dataOut.moments
57 64 if self.CODE == 'gaussian_fit':
58 65 data['gaussfit'] = dataOut.DGauFitParams
59 66
60 return data, meta
61
67 return data, meta
68
62 69 def plot(self):
63 70
64 71 if self.xaxis == "frequency":
65 72 x = self.data.xrange[0]
66 73 self.xlabel = "Frequency (kHz)"
67 74 elif self.xaxis == "time":
68 75 x = self.data.xrange[1]
69 76 self.xlabel = "Time (ms)"
70 77 else:
71 78 x = self.data.xrange[2]
72 79 self.xlabel = "Velocity (m/s)"
73 80
74 81 if (self.CODE == 'spc_moments') | (self.CODE == 'gaussian_fit'):
75 82 x = self.data.xrange[2]
76 83 self.xlabel = "Velocity (m/s)"
77 84
78 85 self.titles = []
79
86
80 87 y = self.data.yrange
81 88 self.y = y
82 89
83 90 data = self.data[-1]
84 91 z = data['spc']
85 92
86 93 for n, ax in enumerate(self.axes):
87 94 noise = data['noise'][n]
88 95
89 96 if self.CODE == 'spc_moments':
90 97 mean = data['moments'][n, 1]
91 if self.CODE == 'gaussian_fit':
98 if self.CODE == 'gaussian_fit':
92 99 gau0 = data['gaussfit'][n][2,:,0]
93 100 gau1 = data['gaussfit'][n][2,:,1]
94 101 if ax.firsttime:
95 102 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
96 103 self.xmin = self.xmin if self.xmin else numpy.nanmin(x)#-self.xmax
97 104 #self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
98 105 #self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
99 106 if self.zlimits is not None:
100 107 self.zmin, self.zmax = self.zlimits[n]
101 108 ax.plt = ax.pcolormesh(x, y, z[n].T,
102 109 vmin=self.zmin,
103 110 vmax=self.zmax,
104 111 cmap=plt.get_cmap(self.colormap)
105 112 )
106 113
107 114 if self.showprofile:
108 115 ax.plt_profile = self.pf_axes[n].plot(
109 116 data['rti'][n], y)[0]
110 117 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
111 118 color="k", linestyle="dashed", lw=1)[0]
112 119 if self.CODE == 'spc_moments':
113 120 ax.plt_mean = ax.plot(mean, y, color='k', lw=1)[0]
114 121 if self.CODE == 'gaussian_fit':
115 122 ax.plt_gau0 = ax.plot(gau0, y, color='r', lw=1)[0]
116 123 ax.plt_gau1 = ax.plot(gau1, y, color='y', lw=1)[0]
117 124 else:
118 125 if self.zlimits is not None:
119 126 self.zmin, self.zmax = self.zlimits[n]
120 127 ax.plt.set_array(z[n].T.ravel())
121 128 if self.showprofile:
122 129 ax.plt_profile.set_data(data['rti'][n], y)
123 130 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
124 131 if self.CODE == 'spc_moments':
125 132 ax.plt_mean.set_data(mean, y)
126 133 if self.CODE == 'gaussian_fit':
127 134 ax.plt_gau0.set_data(gau0, y)
128 135 ax.plt_gau1.set_data(gau1, y)
129 136 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
130 137
131 138 class SpectraObliquePlot(Plot):
132 139 '''
133 140 Plot for Spectra data
134 141
135 142 Written by R. Flores
136 143 '''
137 144
138 145 CODE = 'spc_oblique'
139 146 colormap = 'jet'
140 147 plot_type = 'pcolor'
141 148
142 149 def setup(self):
143 150 self.xaxis = "oblique"
144 151 self.nplots = len(self.data.channels)
145 152 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
146 153 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
147 154 self.height = 2.6 * self.nrows
148 155 self.cb_label = 'dB'
149 156 if self.showprofile:
150 157 self.width = 4 * self.ncols
151 158 else:
152 159 self.width = 3.5 * self.ncols
153 160 self.plots_adjust.update({'wspace': 0.8, 'hspace':0.2, 'left': 0.2, 'right': 0.9, 'bottom': 0.18})
154 161 self.ylabel = 'Range [km]'
155 162
156 163 def update(self, dataOut):
157
164
158 165 data = {}
159 166 meta = {}
160 167 spc = 10*numpy.log10(dataOut.data_spc/dataOut.normFactor)
161 168 data['spc'] = spc
162 169 data['rti'] = dataOut.getPower()
163 170 data['noise'] = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
164 meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
165
171 meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
172
166 173 data['shift1'] = dataOut.Dop_EEJ_T1[0]
167 174 data['shift2'] = dataOut.Dop_EEJ_T2[0]
168 175 data['max_val_2'] = dataOut.Oblique_params[0,-1,:]
169 176 data['shift1_error'] = dataOut.Err_Dop_EEJ_T1[0]
170 177 data['shift2_error'] = dataOut.Err_Dop_EEJ_T2[0]
171
172 return data, meta
178
179 return data, meta
173 180
174 181 def plot(self):
175 182
176 183 if self.xaxis == "frequency":
177 184 x = self.data.xrange[0]
178 185 self.xlabel = "Frequency (kHz)"
179 186 elif self.xaxis == "time":
180 187 x = self.data.xrange[1]
181 188 self.xlabel = "Time (ms)"
182 189 else:
183 190 x = self.data.xrange[2]
184 191 self.xlabel = "Velocity (m/s)"
185 192
186 193 self.titles = []
187 194
188 195 y = self.data.yrange
189 196 self.y = y
190 197 data = self.data[-1]
191 198 z = data['spc']
192 199
193 200 for n, ax in enumerate(self.axes):
194 201 noise = self.data['noise'][n][-1]
195 202 shift1 = data['shift1']
196 203 shift2 = data['shift2']
197 204 max_val_2 = data['max_val_2']
198 205 err1 = data['shift1_error']
199 206 err2 = data['shift2_error']
200 207 if ax.firsttime:
201 208 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
202 209 self.xmin = self.xmin if self.xmin else -self.xmax
203 210 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
204 211 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
205 212 ax.plt = ax.pcolormesh(x, y, z[n].T,
206 213 vmin=self.zmin,
207 214 vmax=self.zmax,
208 215 cmap=plt.get_cmap(self.colormap)
209 216 )
210 217
211 218 if self.showprofile:
212 219 ax.plt_profile = self.pf_axes[n].plot(
213 220 self.data['rti'][n][-1], y)[0]
214 221 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
215 222 color="k", linestyle="dashed", lw=1)[0]
216
223
217 224 self.ploterr1 = ax.errorbar(shift1, y, xerr=err1, fmt='k^', elinewidth=2.2, marker='o', linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
218 225 self.ploterr2 = ax.errorbar(shift2, y, xerr=err2, fmt='m^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
219 226 self.ploterr3 = ax.errorbar(max_val_2, y, xerr=0, fmt='g^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
220 227 else:
221 228 self.ploterr1.remove()
222 229 self.ploterr2.remove()
223 230 self.ploterr3.remove()
224 231 ax.plt.set_array(z[n].T.ravel())
225 232 if self.showprofile:
226 233 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
227 234 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
228 235 self.ploterr1 = ax.errorbar(shift1, y, xerr=err1, fmt='k^', elinewidth=2.2, marker='o', linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
229 236 self.ploterr2 = ax.errorbar(shift2, y, xerr=err2, fmt='m^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
230 237 self.ploterr3 = ax.errorbar(max_val_2, y, xerr=0, fmt='g^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
231
238
232 239 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
233 240
234 241
235 242 class CrossSpectraPlot(Plot):
236 243
237 244 CODE = 'cspc'
238 245 colormap = 'jet'
239 246 plot_type = 'pcolor'
240 247 zmin_coh = None
241 248 zmax_coh = None
242 249 zmin_phase = None
243 250 zmax_phase = None
244 251
245 252 def setup(self):
246 253
247 254 self.ncols = 4
248 255 self.nplots = len(self.data.pairs) * 2
249 256 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
250 257 self.width = 3.1 * self.ncols
251 258 self.height = 5 * self.nrows
252 259 self.ylabel = 'Range [km]'
253 260 self.showprofile = False
254 261 self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
255 262
256 263 def update(self, dataOut):
257 264
258 265 data = {}
259 266 meta = {}
260 267
261 268 spc = dataOut.data_spc
262 269 cspc = dataOut.data_cspc
263 270 extrapoints = spc.shape[1] % dataOut.nFFTPoints
264 meta['xrange'] = (dataOut.getFreqRange(extrapoints) / 1000., dataOut.getAcfRange(extrapoints), dataOut.getVelRange(extrapoints))
271 meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS) / 1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
265 272 meta['pairs'] = dataOut.pairsList
266 273
267 274 tmp = []
268 275
269 276 for n, pair in enumerate(meta['pairs']):
270 277 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
271 278 coh = numpy.abs(out)
272 279 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
273 280 tmp.append(coh)
274 281 tmp.append(phase)
275 282
276 283 data['cspc'] = numpy.array(tmp)
277 284
278 return data, meta
279
285 return data, meta
286
280 287 def plot(self):
281 288
282 289 if self.xaxis == "frequency":
283 290 x = self.data.xrange[0]
284 291 self.xlabel = "Frequency (kHz)"
285 292 elif self.xaxis == "time":
286 293 x = self.data.xrange[1]
287 294 self.xlabel = "Time (ms)"
288 295 else:
289 296 x = self.data.xrange[2]
290 297 self.xlabel = "Velocity (m/s)"
291 298
292 299 self.titles = []
293 300
294 301 y = self.data.yrange
295 302 self.y = y
296 303
297 304 data = self.data[-1]
298 305 cspc = data['cspc']
299 306
300 307 for n in range(len(self.data.pairs)):
301 308 pair = self.data.pairs[n]
302 309 coh = cspc[n * 2]
303 310 phase = cspc[n * 2 + 1]
304 311 ax = self.axes[2 * n]
305 312 if ax.firsttime:
306 313 ax.plt = ax.pcolormesh(x, y, coh.T,
307 314 vmin=0,
308 315 vmax=1,
309 316 cmap=plt.get_cmap(self.colormap_coh)
310 317 )
311 318 else:
312 319 ax.plt.set_array(coh.T.ravel())
313 320 self.titles.append(
314 321 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
315
322
316 323 ax = self.axes[2 * n + 1]
317 324 if ax.firsttime:
318 325 ax.plt = ax.pcolormesh(x, y, phase.T,
319 326 vmin=-180,
320 327 vmax=180,
321 328 cmap=plt.get_cmap(self.colormap_phase)
322 329 )
323 330 else:
324 331 ax.plt.set_array(phase.T.ravel())
325 332 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
326 333
327 334
328 335 class CrossSpectra4Plot(Plot):
329 336
330 337 CODE = 'cspc'
331 338 colormap = 'jet'
332 339 plot_type = 'pcolor'
333 340 zmin_coh = None
334 341 zmax_coh = None
335 342 zmin_phase = None
336 343 zmax_phase = None
337 344
338 345 def setup(self):
339 346
340 347 self.ncols = 4
341 348 self.nrows = len(self.data.pairs)
342 349 self.nplots = self.nrows * 4
343 350 self.width = 3.1 * self.ncols
344 351 self.height = 5 * self.nrows
345 352 self.ylabel = 'Range [km]'
346 353 self.showprofile = False
347 354 self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
348 355
349 356 def plot(self):
350 357
351 358 if self.xaxis == "frequency":
352 359 x = self.data.xrange[0]
353 360 self.xlabel = "Frequency (kHz)"
354 361 elif self.xaxis == "time":
355 362 x = self.data.xrange[1]
356 363 self.xlabel = "Time (ms)"
357 364 else:
358 365 x = self.data.xrange[2]
359 366 self.xlabel = "Velocity (m/s)"
360 367
361 368 self.titles = []
362 369
363 370
364 371 y = self.data.heights
365 372 self.y = y
366 373 nspc = self.data['spc']
367 374 #print(numpy.shape(self.data['spc']))
368 375 spc = self.data['cspc'][0]
369 376 #print(numpy.shape(nspc))
370 377 #exit()
371 378 #nspc[1,:,:] = numpy.flip(nspc[1,:,:],axis=0)
372 379 #print(numpy.shape(spc))
373 380 #exit()
374 381 cspc = self.data['cspc'][1]
375 382
376 383 #xflip=numpy.flip(x)
377 384 #print(numpy.shape(cspc))
378 385 #exit()
379 386
380 387 for n in range(self.nrows):
381 388 noise = self.data['noise'][:,-1]
382 389 pair = self.data.pairs[n]
383 390 #print(pair)
384 391 #exit()
385 392 ax = self.axes[4 * n]
386 393 if ax.firsttime:
387 394 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
388 395 self.xmin = self.xmin if self.xmin else -self.xmax
389 396 self.zmin = self.zmin if self.zmin else numpy.nanmin(nspc)
390 397 self.zmax = self.zmax if self.zmax else numpy.nanmax(nspc)
391 398 ax.plt = ax.pcolormesh(x , y , nspc[pair[0]].T,
392 399 vmin=self.zmin,
393 400 vmax=self.zmax,
394 401 cmap=plt.get_cmap(self.colormap)
395 402 )
396 403 else:
397 404 #print(numpy.shape(nspc[pair[0]].T))
398 405 #exit()
399 406 ax.plt.set_array(nspc[pair[0]].T.ravel())
400 407 self.titles.append('CH {}: {:3.2f}dB'.format(pair[0], noise[pair[0]]))
401 408
402 409 ax = self.axes[4 * n + 1]
403 410
404 411 if ax.firsttime:
405 412 ax.plt = ax.pcolormesh(x , y, numpy.flip(nspc[pair[1]],axis=0).T,
406 413 vmin=self.zmin,
407 414 vmax=self.zmax,
408 415 cmap=plt.get_cmap(self.colormap)
409 416 )
410 417 else:
411 418
412 419 ax.plt.set_array(numpy.flip(nspc[pair[1]],axis=0).T.ravel())
413 420 self.titles.append('CH {}: {:3.2f}dB'.format(pair[1], noise[pair[1]]))
414 421
415 422 out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
416 423 coh = numpy.abs(out)
417 424 phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
418 425
419 426 ax = self.axes[4 * n + 2]
420 427 if ax.firsttime:
421 428 ax.plt = ax.pcolormesh(x, y, numpy.flip(coh,axis=0).T,
422 429 vmin=0,
423 430 vmax=1,
424 431 cmap=plt.get_cmap(self.colormap_coh)
425 432 )
426 433 else:
427 434 ax.plt.set_array(numpy.flip(coh,axis=0).T.ravel())
428 435 self.titles.append(
429 436 'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
430 437
431 438 ax = self.axes[4 * n + 3]
432 439 if ax.firsttime:
433 440 ax.plt = ax.pcolormesh(x, y, numpy.flip(phase,axis=0).T,
434 441 vmin=-180,
435 442 vmax=180,
436 443 cmap=plt.get_cmap(self.colormap_phase)
437 444 )
438 445 else:
439 446 ax.plt.set_array(numpy.flip(phase,axis=0).T.ravel())
440 447 self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
441 448
442 449
443 450 class CrossSpectra2Plot(Plot):
444 451
445 452 CODE = 'cspc'
446 453 colormap = 'jet'
447 454 plot_type = 'pcolor'
448 455 zmin_coh = None
449 456 zmax_coh = None
450 457 zmin_phase = None
451 458 zmax_phase = None
452 459
453 460 def setup(self):
454 461
455 462 self.ncols = 1
456 463 self.nrows = len(self.data.pairs)
457 464 self.nplots = self.nrows * 1
458 465 self.width = 3.1 * self.ncols
459 466 self.height = 5 * self.nrows
460 467 self.ylabel = 'Range [km]'
461 468 self.showprofile = False
462 469 self.plots_adjust.update({'left': 0.22, 'right': .90, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
463 470
464 471 def plot(self):
465 472
466 473 if self.xaxis == "frequency":
467 474 x = self.data.xrange[0]
468 475 self.xlabel = "Frequency (kHz)"
469 476 elif self.xaxis == "time":
470 477 x = self.data.xrange[1]
471 478 self.xlabel = "Time (ms)"
472 479 else:
473 480 x = self.data.xrange[2]
474 481 self.xlabel = "Velocity (m/s)"
475 482
476 483 self.titles = []
477 484
478 485
479 486 y = self.data.heights
480 487 self.y = y
481 488 #nspc = self.data['spc']
482 489 #print(numpy.shape(self.data['spc']))
483 490 #spc = self.data['cspc'][0]
484 491 #print(numpy.shape(spc))
485 492 #exit()
486 493 cspc = self.data['cspc'][1]
487 494 #print(numpy.shape(cspc))
488 495 #exit()
489 496
490 497 for n in range(self.nrows):
491 498 noise = self.data['noise'][:,-1]
492 499 pair = self.data.pairs[n]
493 500 #print(pair) #exit()
494 501
495 502
496 503
497 504 out = cspc[n]# / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
498 505
499 506 #print(out[:,53])
500 507 #exit()
501 508 cross = numpy.abs(out)
502 509 z = cross/self.data.nFactor
503 510 #print("here")
504 511 #print(dataOut.data_spc[0,0,0])
505 512 #exit()
506 513
507 514 cross = 10*numpy.log10(z)
508 515 #print(numpy.shape(cross))
509 516 #print(cross[0,:])
510 517 #print(self.data.nFactor)
511 518 #exit()
512 519 #phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
513 520
514 521 ax = self.axes[1 * n]
515 522 if ax.firsttime:
516 523 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
517 524 self.xmin = self.xmin if self.xmin else -self.xmax
518 525 self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
519 526 self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
520 527 ax.plt = ax.pcolormesh(x, y, cross.T,
521 528 vmin=self.zmin,
522 529 vmax=self.zmax,
523 530 cmap=plt.get_cmap(self.colormap)
524 531 )
525 532 else:
526 533 ax.plt.set_array(cross.T.ravel())
527 534 self.titles.append(
528 535 'Cross Spectra Power Ch{} * Ch{}'.format(pair[0], pair[1]))
529 536
530 537
531 538 class CrossSpectra3Plot(Plot):
532 539
533 540 CODE = 'cspc'
534 541 colormap = 'jet'
535 542 plot_type = 'pcolor'
536 543 zmin_coh = None
537 544 zmax_coh = None
538 545 zmin_phase = None
539 546 zmax_phase = None
540 547
541 548 def setup(self):
542 549
543 550 self.ncols = 3
544 551 self.nrows = len(self.data.pairs)
545 552 self.nplots = self.nrows * 3
546 553 self.width = 3.1 * self.ncols
547 554 self.height = 5 * self.nrows
548 555 self.ylabel = 'Range [km]'
549 556 self.showprofile = False
550 557 self.plots_adjust.update({'left': 0.22, 'right': .90, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
551 558
552 559 def plot(self):
553 560
554 561 if self.xaxis == "frequency":
555 562 x = self.data.xrange[0]
556 563 self.xlabel = "Frequency (kHz)"
557 564 elif self.xaxis == "time":
558 565 x = self.data.xrange[1]
559 566 self.xlabel = "Time (ms)"
560 567 else:
561 568 x = self.data.xrange[2]
562 569 self.xlabel = "Velocity (m/s)"
563 570
564 571 self.titles = []
565 572
566 573
567 574 y = self.data.heights
568 575 self.y = y
569 576 #nspc = self.data['spc']
570 577 #print(numpy.shape(self.data['spc']))
571 578 #spc = self.data['cspc'][0]
572 579 #print(numpy.shape(spc))
573 580 #exit()
574 581 cspc = self.data['cspc'][1]
575 582 #print(numpy.shape(cspc))
576 583 #exit()
577 584
578 585 for n in range(self.nrows):
579 586 noise = self.data['noise'][:,-1]
580 587 pair = self.data.pairs[n]
581 588 #print(pair) #exit()
582 589
583 590
584 591
585 592 out = cspc[n]# / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
586 593
587 594 #print(out[:,53])
588 595 #exit()
589 596 cross = numpy.abs(out)
590 597 z = cross/self.data.nFactor
591 598 cross = 10*numpy.log10(z)
592 599
593 600 out_r= out.real/self.data.nFactor
594 601 #out_r = 10*numpy.log10(out_r)
595 602
596 603 out_i= out.imag/self.data.nFactor
597 604 #out_i = 10*numpy.log10(out_i)
598 605 #print(numpy.shape(cross))
599 606 #print(cross[0,:])
600 607 #print(self.data.nFactor)
601 608 #exit()
602 609 #phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
603 610
604 611 ax = self.axes[3 * n]
605 612 if ax.firsttime:
606 613 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
607 614 self.xmin = self.xmin if self.xmin else -self.xmax
608 615 self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
609 616 self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
610 617 ax.plt = ax.pcolormesh(x, y, cross.T,
611 618 vmin=self.zmin,
612 619 vmax=self.zmax,
613 620 cmap=plt.get_cmap(self.colormap)
614 621 )
615 622 else:
616 623 ax.plt.set_array(cross.T.ravel())
617 624 self.titles.append(
618 625 'Cross Spectra Power Ch{} * Ch{}'.format(pair[0], pair[1]))
619 626
620 627 ax = self.axes[3 * n + 1]
621 628 if ax.firsttime:
622 629 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
623 630 self.xmin = self.xmin if self.xmin else -self.xmax
624 631 self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
625 632 self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
626 633 ax.plt = ax.pcolormesh(x, y, out_r.T,
627 634 vmin=-1.e6,
628 635 vmax=0,
629 636 cmap=plt.get_cmap(self.colormap)
630 637 )
631 638 else:
632 639 ax.plt.set_array(out_r.T.ravel())
633 640 self.titles.append(
634 641 'Cross Spectra Real Ch{} * Ch{}'.format(pair[0], pair[1]))
635 642
636 643 ax = self.axes[3 * n + 2]
637 644
638 645
639 646 if ax.firsttime:
640 647 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
641 648 self.xmin = self.xmin if self.xmin else -self.xmax
642 649 self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
643 650 self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
644 651 ax.plt = ax.pcolormesh(x, y, out_i.T,
645 652 vmin=-1.e6,
646 653 vmax=1.e6,
647 654 cmap=plt.get_cmap(self.colormap)
648 655 )
649 656 else:
650 657 ax.plt.set_array(out_i.T.ravel())
651 658 self.titles.append(
652 659 'Cross Spectra Imag Ch{} * Ch{}'.format(pair[0], pair[1]))
653 660
654 661 class RTIPlot(Plot):
655 662 '''
656 663 Plot for RTI data
657 664 '''
658 665
659 666 CODE = 'rti'
660 667 colormap = 'jet'
661 668 plot_type = 'pcolorbuffer'
662 669
663 670 def setup(self):
664 671 self.xaxis = 'time'
665 672 self.ncols = 1
666 673 self.nrows = len(self.data.channels)
667 674 self.nplots = len(self.data.channels)
668 675 self.ylabel = 'Range [km]'
669 676 self.xlabel = 'Time'
670 677 self.cb_label = 'dB'
671 678 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
672 679 self.titles = ['{} Channel {}'.format(
673 680 self.CODE.upper(), x) for x in range(self.nrows)]
674 681
675 682 def update(self, dataOut):
676 683
677 684 data = {}
678 685 meta = {}
679 686 data['rti'] = dataOut.getPower()
680 687 data['noise'] = 10 * numpy.log10(dataOut.getNoise() / dataOut.normFactor)
681 688
682 689 return data, meta
683 690
684 691 def plot(self):
685 692 self.x = self.data.times
686 693 self.y = self.data.yrange
687 694 self.z = self.data[self.CODE]
688 695
689 696 self.z = numpy.ma.masked_invalid(self.z)
690 697
691 698 if self.decimation is None:
692 699 x, y, z = self.fill_gaps(self.x, self.y, self.z)
693 700 else:
694 701 x, y, z = self.fill_gaps(*self.decimate())
695 702
696 703 for n, ax in enumerate(self.axes):
697 704 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
698 705 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
699 706 if ax.firsttime:
700 707 if self.zlimits is not None:
701 708 self.zmin, self.zmax = self.zlimits[n]
702 709 ax.plt = ax.pcolormesh(x, y, z[n].T,
703 710 vmin=self.zmin,
704 711 vmax=self.zmax,
705 712 cmap=plt.get_cmap(self.colormap)
706 713 )
707 714 if self.showprofile:
708 715 ax.plot_profile = self.pf_axes[n].plot(
709 716 self.data['rti'][n][-1], self.y)[0]
710 717 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(self.data['noise'][n][-1], len(self.y)), self.y,
711 718 color="k", linestyle="dashed", lw=1)[0]
712 719 else:
713 720 if self.zlimits is not None:
714 721 self.zmin, self.zmax = self.zlimits[n]
715 722 ax.plt.remove()
716 723 ax.plt = ax.pcolormesh(x, y, z[n].T,
717 724 vmin=self.zmin,
718 725 vmax=self.zmax,
719 726 cmap=plt.get_cmap(self.colormap)
720 727 )
721 728 if self.showprofile:
722 729 ax.plot_profile.set_data(self.data['rti'][n][-1], self.y)
723 730 ax.plot_noise.set_data(numpy.repeat(
724 731 self.data['noise'][n][-1], len(self.y)), self.y)
725 732
726 733
727 734 class SpectrogramPlot(Plot):
728 735 '''
729 736 Plot for Spectrogram data
730 737 '''
731 738
732 739 CODE = 'Spectrogram_Profile'
733 740 colormap = 'binary'
734 741 plot_type = 'pcolorbuffer'
735 742
736 743 def setup(self):
737 744 self.xaxis = 'time'
738 745 self.ncols = 1
739 746 self.nrows = len(self.data.channels)
740 747 self.nplots = len(self.data.channels)
741 748 self.xlabel = 'Time'
742 749 #self.cb_label = 'dB'
743 750 self.plots_adjust.update({'hspace':1.2, 'left': 0.1, 'bottom': 0.12, 'right':0.95})
744 751 self.titles = []
745 752
746 753 #self.titles = ['{} Channel {} \n H = {} km ({} - {})'.format(
747 754 #self.CODE.upper(), x, self.data.heightList[self.data.hei], self.data.heightList[self.data.hei],self.data.heightList[self.data.hei]+(self.data.DH*self.data.nProfiles)) for x in range(self.nrows)]
748 755
749 756 self.titles = ['{} Channel {}'.format(
750 757 self.CODE.upper(), x) for x in range(self.nrows)]
751 758
752 759
753 760 def update(self, dataOut):
754 761 data = {}
755 762 meta = {}
756 763
757 764 maxHei = 1620#+12000
758 765 maxHei = 1180
759 766 indb = numpy.where(dataOut.heightList <= maxHei)
760 767 hei = indb[0][-1]
761 768 #print(dataOut.heightList)
762 769
763 770 factor = dataOut.nIncohInt
764 771 z = dataOut.data_spc[:,:,hei] / factor
765 772 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
766 773 #buffer = 10 * numpy.log10(z)
767 774
768 meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
775 meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
769 776
770 777
771 778 #self.hei = hei
772 779 #self.heightList = dataOut.heightList
773 780 #self.DH = (dataOut.heightList[1] - dataOut.heightList[0])/dataOut.step
774 781 #self.nProfiles = dataOut.nProfiles
775 782
776 783 data['Spectrogram_Profile'] = 10 * numpy.log10(z)
777 784
778 785 data['hei'] = hei
779 786 data['DH'] = (dataOut.heightList[1] - dataOut.heightList[0])/dataOut.step
780 787 data['nProfiles'] = dataOut.nProfiles
781 788 #meta['yrange'] = dataOut.heightList[0:dataOut.NSHTS]
782 789 '''
783 790 import matplotlib.pyplot as plt
784 791 plt.plot(10 * numpy.log10(z[0,:]))
785 792 plt.show()
786 793
787 794 from time import sleep
788 795 sleep(10)
789 796 '''
790 797 return data, meta
791 798
792 799 def plot(self):
793 800
794 801 self.x = self.data.times
795 802 self.z = self.data[self.CODE]
796 803 self.y = self.data.xrange[0]
797 804
798 805 hei = self.data['hei'][-1]
799 806 DH = self.data['DH'][-1]
800 807 nProfiles = self.data['nProfiles'][-1]
801 808
802 809 self.ylabel = "Frequency (kHz)"
803 810
804 811 self.z = numpy.ma.masked_invalid(self.z)
805 812
806 813 if self.decimation is None:
807 814 x, y, z = self.fill_gaps(self.x, self.y, self.z)
808 815 else:
809 816 x, y, z = self.fill_gaps(*self.decimate())
810 817
811 818 for n, ax in enumerate(self.axes):
812 819 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
813 820 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
814 821 data = self.data[-1]
815 822 if ax.firsttime:
816 823 ax.plt = ax.pcolormesh(x, y, z[n].T,
817 824 vmin=self.zmin,
818 825 vmax=self.zmax,
819 826 cmap=plt.get_cmap(self.colormap)
820 827 )
821 828 else:
822 829 ax.collections.remove(ax.collections[0])
823 830 ax.plt = ax.pcolormesh(x, y, z[n].T,
824 831 vmin=self.zmin,
825 832 vmax=self.zmax,
826 833 cmap=plt.get_cmap(self.colormap)
827 834 )
828 835
829 836 #self.titles.append('Spectrogram')
830 837
831 838 #self.titles.append('{} Channel {} \n H = {} km ({} - {})'.format(
832 839 #self.CODE.upper(), x, y[hei], y[hei],y[hei]+(DH*nProfiles)))
833 840
834 841 class CoherencePlot(RTIPlot):
835 842 '''
836 843 Plot for Coherence data
837 844 '''
838 845
839 846 CODE = 'coh'
840 847
841 848 def setup(self):
842 849 self.xaxis = 'time'
843 850 self.ncols = 1
844 851 self.nrows = len(self.data.pairs)
845 852 self.nplots = len(self.data.pairs)
846 853 self.ylabel = 'Range [km]'
847 854 self.xlabel = 'Time'
848 855 self.plots_adjust.update({'hspace':0.6, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
849 856 if self.CODE == 'coh':
850 857 self.cb_label = ''
851 858 self.titles = [
852 859 'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
853 860 else:
854 861 self.cb_label = 'Degrees'
855 862 self.titles = [
856 863 'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
857 864
858 865 def update(self, dataOut):
859 866
860 867 data = {}
861 868 meta = {}
862 869 data['coh'] = dataOut.getCoherence()
863 870 meta['pairs'] = dataOut.pairsList
864 871
865 872 return data, meta
866 873
867 874 class PhasePlot(CoherencePlot):
868 875 '''
869 876 Plot for Phase map data
870 877 '''
871 878
872 879 CODE = 'phase'
873 880 colormap = 'seismic'
874 881
875 882 def update(self, dataOut):
876 883
877 884 data = {}
878 885 meta = {}
879 886 data['phase'] = dataOut.getCoherence(phase=True)
880 887 meta['pairs'] = dataOut.pairsList
881 888
882 889 return data, meta
883 890
884 891 class NoisePlot(Plot):
885 892 '''
886 893 Plot for noise
887 894 '''
888 895
889 896 CODE = 'noise'
890 897 plot_type = 'scatterbuffer'
891 898
892 899 def setup(self):
893 900 self.xaxis = 'time'
894 901 self.ncols = 1
895 902 self.nrows = 1
896 903 self.nplots = 1
897 904 self.ylabel = 'Intensity [dB]'
898 905 self.xlabel = 'Time'
899 906 self.titles = ['Noise']
900 907 self.colorbar = False
901 908 self.plots_adjust.update({'right': 0.85 })
902 909
903 910 def update(self, dataOut):
904 911
905 912 data = {}
906 913 meta = {}
907 914 data['noise'] = 10 * numpy.log10(dataOut.getNoise() / dataOut.normFactor).reshape(dataOut.nChannels, 1)
908 915 meta['yrange'] = numpy.array([])
909 916
910 917 return data, meta
911 918
912 919 def plot(self):
913 920
914 921 x = self.data.times
915 922 xmin = self.data.min_time
916 923 xmax = xmin + self.xrange * 60 * 60
917 924 Y = self.data['noise']
918 925
919 926 if self.axes[0].firsttime:
920 927 self.ymin = numpy.nanmin(Y) - 5
921 928 self.ymax = numpy.nanmax(Y) + 5
922 929 for ch in self.data.channels:
923 930 y = Y[ch]
924 931 self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
925 932 plt.legend(bbox_to_anchor=(1.18, 1.0))
926 933 else:
927 934 for ch in self.data.channels:
928 935 y = Y[ch]
929 936 self.axes[0].lines[ch].set_data(x, y)
930 937
931 938 self.ymin = numpy.nanmin(Y) - 5
932 939 self.ymax = numpy.nanmax(Y) + 10
933 940
934 941
935 942 class PowerProfilePlot(Plot):
936 943
937 944 CODE = 'pow_profile'
938 945 plot_type = 'scatter'
939 946
940 947 def setup(self):
941 948
942 949 self.ncols = 1
943 950 self.nrows = 1
944 951 self.nplots = 1
945 952 self.height = 4
946 953 self.width = 3
947 954 self.ylabel = 'Range [km]'
948 955 self.xlabel = 'Intensity [dB]'
949 956 self.titles = ['Power Profile']
950 957 self.colorbar = False
951 958
952 959 def update(self, dataOut):
953 960
954 961 data = {}
955 962 meta = {}
956 963 data[self.CODE] = dataOut.getPower()
957 964
958 965 return data, meta
959 966
960 967 def plot(self):
961 968
962 969 y = self.data.yrange
963 970 self.y = y
964 971
965 972 x = self.data[-1][self.CODE]
966
973
967 974 if self.xmin is None: self.xmin = numpy.nanmin(x) * 0.9
968 975 if self.xmax is None: self.xmax = numpy.nanmax(x) * 1.1
969
976
970 977 if self.axes[0].firsttime:
971 978 for ch in self.data.channels:
972 979 self.axes[0].plot(x[ch], y, lw=1, label='Ch{}'.format(ch))
973 980 plt.legend()
974 981 else:
975 982 for ch in self.data.channels:
976 983 self.axes[0].lines[ch].set_data(x[ch], y)
977 984
978 985
979 986 class SpectraCutPlot(Plot):
980 987
981 988 CODE = 'spc_cut'
982 989 plot_type = 'scatter'
983 990 buffering = False
984 991
985 992 def setup(self):
986 993
987 994 self.nplots = len(self.data.channels)
988 995 self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
989 996 self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
990 997 self.width = 3.4 * self.ncols + 1.5
991 998 self.height = 3 * self.nrows
992 999 self.ylabel = 'Power [dB]'
993 1000 self.colorbar = False
994 1001 self.plots_adjust.update({'left':0.1, 'hspace':0.3, 'right': 0.75, 'bottom':0.08})
995 1002
996 1003 def update(self, dataOut):
997 1004
998 1005 data = {}
999 1006 meta = {}
1000 1007 spc = 10 * numpy.log10(dataOut.data_pre[0] / dataOut.normFactor)
1001 1008 data['spc'] = spc
1002 meta['xrange'] = (dataOut.getFreqRange(1) / 1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
1009 meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS) / 1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
1003 1010 if self.CODE == 'cut_gaussian_fit':
1004 1011 data['gauss_fit0'] = 10 * numpy.log10(dataOut.GaussFit0 / dataOut.normFactor)
1005 1012 data['gauss_fit1'] = 10 * numpy.log10(dataOut.GaussFit1 / dataOut.normFactor)
1006 1013 return data, meta
1007 1014
1008 1015 def plot(self):
1009 1016 if self.xaxis == "frequency":
1010 1017 x = self.data.xrange[0][1:]
1011 1018 self.xlabel = "Frequency (kHz)"
1012 1019 elif self.xaxis == "time":
1013 1020 x = self.data.xrange[1]
1014 1021 self.xlabel = "Time (ms)"
1015 1022 else:
1016 1023 x = self.data.xrange[2][:-1]
1017 1024 self.xlabel = "Velocity (m/s)"
1018
1025
1019 1026 if self.CODE == 'cut_gaussian_fit':
1020 1027 x = self.data.xrange[2][:-1]
1021 1028 self.xlabel = "Velocity (m/s)"
1022 1029
1023 1030 self.titles = []
1024 1031
1025 1032 y = self.data.yrange
1026 1033 data = self.data[-1]
1027 1034 z = data['spc']
1028 1035
1029 1036 if self.height_index:
1030 1037 index = numpy.array(self.height_index)
1031 1038 else:
1032 1039 index = numpy.arange(0, len(y), int((len(y)) / 9))
1033 1040
1034 1041 for n, ax in enumerate(self.axes):
1035 if self.CODE == 'cut_gaussian_fit':
1042 if self.CODE == 'cut_gaussian_fit':
1036 1043 gau0 = data['gauss_fit0']
1037 1044 gau1 = data['gauss_fit1']
1038 1045 if ax.firsttime:
1039 1046 self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
1040 1047 self.xmin = self.xmin if self.xmin else -self.xmax
1041 1048 self.ymin = self.ymin if self.ymin else numpy.nanmin(z[:,:,index])
1042 1049 self.ymax = self.ymax if self.ymax else numpy.nanmax(z[:,:,index])
1043
1050
1044 1051 ax.plt = ax.plot(x, z[n, :, index].T, lw=0.25)
1045 1052 if self.CODE == 'cut_gaussian_fit':
1046 1053 ax.plt_gau0 = ax.plot(x, gau0[n, :, index].T, lw=1, linestyle='-.')
1047 1054 for i, line in enumerate(ax.plt_gau0):
1048 line.set_color(ax.plt[i].get_color())
1055 line.set_color(ax.plt[i].get_color())
1049 1056 ax.plt_gau1 = ax.plot(x, gau1[n, :, index].T, lw=1, linestyle='--')
1050 1057 for i, line in enumerate(ax.plt_gau1):
1051 line.set_color(ax.plt[i].get_color())
1058 line.set_color(ax.plt[i].get_color())
1052 1059 labels = ['Range = {:2.1f}km'.format(y[i]) for i in index]
1053 1060 self.figures[0].legend(ax.plt, labels, loc='center right')
1054 1061 else:
1055 1062 for i, line in enumerate(ax.plt):
1056 1063 line.set_data(x, z[n, :, index[i]].T)
1057 1064 for i, line in enumerate(ax.plt_gau0):
1058 1065 line.set_data(x, gau0[n, :, index[i]].T)
1059 1066 line.set_color(ax.plt[i].get_color())
1060 1067 for i, line in enumerate(ax.plt_gau1):
1061 1068 line.set_data(x, gau1[n, :, index[i]].T)
1062 1069 line.set_color(ax.plt[i].get_color())
1063 1070 self.titles.append('CH {}'.format(n))
1064 1071
1065 1072
1066 1073 class BeaconPhase(Plot):
1067 1074
1068 1075 __isConfig = None
1069 1076 __nsubplots = None
1070 1077
1071 1078 PREFIX = 'beacon_phase'
1072 1079
1073 1080 def __init__(self):
1074 1081 Plot.__init__(self)
1075 1082 self.timerange = 24 * 60 * 60
1076 1083 self.isConfig = False
1077 1084 self.__nsubplots = 1
1078 1085 self.counter_imagwr = 0
1079 1086 self.WIDTH = 800
1080 1087 self.HEIGHT = 400
1081 1088 self.WIDTHPROF = 120
1082 1089 self.HEIGHTPROF = 0
1083 1090 self.xdata = None
1084 1091 self.ydata = None
1085 1092
1086 1093 self.PLOT_CODE = BEACON_CODE
1087 1094
1088 1095 self.FTP_WEI = None
1089 1096 self.EXP_CODE = None
1090 1097 self.SUB_EXP_CODE = None
1091 1098 self.PLOT_POS = None
1092 1099
1093 1100 self.filename_phase = None
1094 1101
1095 1102 self.figfile = None
1096 1103
1097 1104 self.xmin = None
1098 1105 self.xmax = None
1099 1106
1100 1107 def getSubplots(self):
1101 1108
1102 1109 ncol = 1
1103 1110 nrow = 1
1104 1111
1105 1112 return nrow, ncol
1106 1113
1107 1114 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
1108 1115
1109 1116 self.__showprofile = showprofile
1110 1117 self.nplots = nplots
1111 1118
1112 1119 ncolspan = 7
1113 1120 colspan = 6
1114 1121 self.__nsubplots = 2
1115 1122
1116 1123 self.createFigure(id=id,
1117 1124 wintitle=wintitle,
1118 1125 widthplot=self.WIDTH + self.WIDTHPROF,
1119 1126 heightplot=self.HEIGHT + self.HEIGHTPROF,
1120 1127 show=show)
1121 1128
1122 1129 nrow, ncol = self.getSubplots()
1123 1130
1124 1131 self.addAxes(nrow, ncol * ncolspan, 0, 0, colspan, 1)
1125 1132
1126 1133 def save_phase(self, filename_phase):
1127 1134 f = open(filename_phase, 'w+')
1128 1135 f.write('\n\n')
1129 1136 f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
1130 1137 f.write('DD MM YYYY HH MM SS pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n')
1131 1138 f.close()
1132 1139
1133 1140 def save_data(self, filename_phase, data, data_datetime):
1134 1141 f = open(filename_phase, 'a')
1135 1142 timetuple_data = data_datetime.timetuple()
1136 1143 day = str(timetuple_data.tm_mday)
1137 1144 month = str(timetuple_data.tm_mon)
1138 1145 year = str(timetuple_data.tm_year)
1139 1146 hour = str(timetuple_data.tm_hour)
1140 1147 minute = str(timetuple_data.tm_min)
1141 1148 second = str(timetuple_data.tm_sec)
1142 1149 f.write(day + ' ' + month + ' ' + year + ' ' + hour + ' ' + minute + ' ' + second + ' ' + str(data[0]) + ' ' + str(data[1]) + ' ' + str(data[2]) + ' ' + str(data[3]) + '\n')
1143 1150 f.close()
1144 1151
1145 1152 def plot(self):
1146 1153 log.warning('TODO: Not yet implemented...')
1147 1154
1148 1155 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
1149 1156 xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
1150 1157 timerange=None,
1151 1158 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
1152 1159 server=None, folder=None, username=None, password=None,
1153 1160 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
1154 1161
1155 1162 if dataOut.flagNoData:
1156 1163 return dataOut
1157 1164
1158 1165 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
1159 1166 return
1160 1167
1161 1168 if pairsList == None:
1162 1169 pairsIndexList = dataOut.pairsIndexList[:10]
1163 1170 else:
1164 1171 pairsIndexList = []
1165 1172 for pair in pairsList:
1166 1173 if pair not in dataOut.pairsList:
1167 1174 raise ValueError("Pair %s is not in dataOut.pairsList" % (pair))
1168 1175 pairsIndexList.append(dataOut.pairsList.index(pair))
1169 1176
1170 1177 if pairsIndexList == []:
1171 1178 return
1172 1179
1173 1180 # if len(pairsIndexList) > 4:
1174 1181 # pairsIndexList = pairsIndexList[0:4]
1175 1182
1176 1183 hmin_index = None
1177 1184 hmax_index = None
1178 1185
1179 1186 if hmin != None and hmax != None:
1180 1187 indexes = numpy.arange(dataOut.nHeights)
1181 1188 hmin_list = indexes[dataOut.heightList >= hmin]
1182 1189 hmax_list = indexes[dataOut.heightList <= hmax]
1183 1190
1184 1191 if hmin_list.any():
1185 1192 hmin_index = hmin_list[0]
1186 1193
1187 1194 if hmax_list.any():
1188 1195 hmax_index = hmax_list[-1] + 1
1189 1196
1190 1197 x = dataOut.getTimeRange()
1191 1198
1192 1199 thisDatetime = dataOut.datatime
1193 1200
1194 1201 title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
1195 1202 xlabel = "Local Time"
1196 1203 ylabel = "Phase (degrees)"
1197 1204
1198 1205 update_figfile = False
1199 1206
1200 1207 nplots = len(pairsIndexList)
1201 1208 # phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
1202 1209 phase_beacon = numpy.zeros(len(pairsIndexList))
1203 1210 for i in range(nplots):
1204 1211 pair = dataOut.pairsList[pairsIndexList[i]]
1205 1212 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
1206 1213 powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
1207 1214 powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
1208 1215 avgcoherenceComplex = ccf / numpy.sqrt(powa * powb)
1209 1216 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real) * 180 / numpy.pi
1210 1217
1211 1218 if dataOut.beacon_heiIndexList:
1212 1219 phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
1213 1220 else:
1214 1221 phase_beacon[i] = numpy.average(phase)
1215 1222
1216 1223 if not self.isConfig:
1217 1224
1218 1225 nplots = len(pairsIndexList)
1219 1226
1220 1227 self.setup(id=id,
1221 1228 nplots=nplots,
1222 1229 wintitle=wintitle,
1223 1230 showprofile=showprofile,
1224 1231 show=show)
1225 1232
1226 1233 if timerange != None:
1227 1234 self.timerange = timerange
1228 1235
1229 1236 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
1230 1237
1231 1238 if ymin == None: ymin = 0
1232 1239 if ymax == None: ymax = 360
1233 1240
1234 1241 self.FTP_WEI = ftp_wei
1235 1242 self.EXP_CODE = exp_code
1236 1243 self.SUB_EXP_CODE = sub_exp_code
1237 1244 self.PLOT_POS = plot_pos
1238 1245
1239 1246 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
1240 1247 self.isConfig = True
1241 1248 self.figfile = figfile
1242 1249 self.xdata = numpy.array([])
1243 1250 self.ydata = numpy.array([])
1244 1251
1245 1252 update_figfile = True
1246 1253
1247 1254 # open file beacon phase
1248 1255 path = '%s%03d' % (self.PREFIX, self.id)
1249 1256 beacon_file = os.path.join(path, '%s.txt' % self.name)
1250 1257 self.filename_phase = os.path.join(figpath, beacon_file)
1251 1258 # self.save_phase(self.filename_phase)
1252 1259
1253 1260
1254 1261 # store data beacon phase
1255 1262 # self.save_data(self.filename_phase, phase_beacon, thisDatetime)
1256 1263
1257 1264 self.setWinTitle(title)
1258 1265
1259 1266
1260 1267 title = "Phase Plot %s" % (thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
1261 1268
1262 1269 legendlabels = ["Pair (%d,%d)" % (pair[0], pair[1]) for pair in dataOut.pairsList]
1263 1270
1264 1271 axes = self.axesList[0]
1265 1272
1266 1273 self.xdata = numpy.hstack((self.xdata, x[0:1]))
1267 1274
1268 1275 if len(self.ydata) == 0:
1269 1276 self.ydata = phase_beacon.reshape(-1, 1)
1270 1277 else:
1271 1278 self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1, 1)))
1272 1279
1273 1280
1274 1281 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
1275 1282 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
1276 1283 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
1277 1284 XAxisAsTime=True, grid='both'
1278 1285 )
1279 1286
1280 1287 self.draw()
1281 1288
1282 1289 if dataOut.ltctime >= self.xmax:
1283 1290 self.counter_imagwr = wr_period
1284 1291 self.isConfig = False
1285 1292 update_figfile = True
1286 1293
1287 1294 self.save(figpath=figpath,
1288 1295 figfile=figfile,
1289 1296 save=save,
1290 1297 ftp=ftp,
1291 1298 wr_period=wr_period,
1292 1299 thisDatetime=thisDatetime,
1293 1300 update_figfile=update_figfile)
1294 1301
1295 1302 return dataOut
General Comments 0
You need to be logged in to leave comments. Login now