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