1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19 """
20 This module provides some wrapper functions for plotting TimeSeries and FrequencySeries using classes defined in pylal.plotutils.
21 """
22
23
24
25
26
27 from __future__ import division
28
29 import datetime
30 import re
31 import itertools
32 import numpy
33 import numbers
34
35 from scipy import interpolate
36
37 import lal
38
39 from pylal import plotutils
40 from pylal import seriesutils
41 from pylal import git_version
42 from lal import rate
43
44 __author__ = "Duncan M. Macleod <duncan.macleod@ligo.org>"
45 __version__ = git_version.id
46 __date__ = git_version.date
47
48
49
50
51
53 """
54 Plot a REALXTimeSeries.
55 """
56
57 if hasattr(series, "__contains__"):
58 serieslist = series
59 else:
60 serieslist = [series]
61
62
63 xlim = kwargs.pop("xlim", None)
64 ylim = kwargs.pop("ylim", None)
65 if xlim:
66 start,end = xlim
67 else:
68 start = min(float(s.epoch) for s in serieslist)
69 end = max(float(s.epoch + s.data.length*s.deltaT) for s in serieslist)
70
71
72 logx = kwargs.pop("logx", False)
73 logy = kwargs.pop("logy", False)
74
75
76 loc = kwargs.pop("loc", 0)
77 alpha = kwargs.pop("alpha", 0.8)
78
79
80 hidden_colorbar = kwargs.pop("hidden_colorbar", False)
81
82
83 bbox_inches = kwargs.pop("bbox_inches", None)
84
85
86
87
88
89 xlabel = kwargs.pop("xlabel", None)
90 if xlabel:
91 unit = 1
92 if not xlabel:
93 unit, timestr = plotutils.time_axis_unit(end-start)
94 if not t0:
95 t0 = start
96 t0 = lal.LIGOTimeGPS(t0)
97 if int(t0.gpsNanoSeconds) == 0:
98 xlabel = datetime.datetime(*lal.GPSToUTC(int(t0))[:6])\
99 .strftime("%B %d %Y, %H:%M:%S %ZUTC")
100 xlabel = "Time (%s) since %s (%s)" % (timestr, xlabel, int(t0))
101 else:
102 xlabel = datetime.datetime(*lal.GPSToUTC(t0.gpsSeconds)[:6])\
103 .strftime("%B %d %Y, %H:%M:%S %ZUTC")
104 xlabel = "Time (%s) since %s (%s)"\
105 % (timestr,
106 xlabel.replace(" UTC",".%.3s UTC" % t0.gpsNanoSeconds),\
107 t0)
108 t0 = float(t0)
109 ylabel = kwargs.pop("ylabel", "Amplitude (counts)")
110 title = kwargs.pop("title", "")
111 subtitle = kwargs.pop("subtitle", "")
112
113
114
115
116
117 color = kwargs.pop("color", None)
118 if isinstance(color, str):
119 color = color.split(',')
120 if len(color) == 1:
121 color = [color[0]]*len(serieslist)
122 kwargs2 = dict()
123 kwargs2.update(kwargs)
124 if kwargs.has_key("marker"):
125 kwargs.setdefault("markersize", 5)
126 kwargs2.setdefault("markersize", 2)
127 else:
128 kwargs.setdefault("linestyle", "-")
129 kwargs.setdefault("linewidth", "1")
130 kwargs2.setdefault("linewidth", "0.1")
131
132
133
134
135
136 allnames = [s.name for s in serieslist]
137 namedseries = [s for s in serieslist if not re.search("(min|max)\Z",s.name)]
138
139 plot = plotutils.SimplePlot(xlabel, ylabel, title, subtitle)
140 for i,(series,c) in enumerate(itertools.izip(namedseries,\
141 plotutils.default_colors())):
142 if color:
143 c = color[serieslist.index(series)]
144 x = numpy.arange(series.data.length) * series.deltaT +\
145 float(series.epoch) - float(t0)
146 x = x.astype(float)
147 x /= unit
148 d = series.data.data
149 if logy and ylim:
150 numpy.putmask(d, d==0, ylim[0]-abs(ylim[0])*0.01)
151 plot.add_content(x, d, color=c,\
152 label=plotutils.display_name(series.name), **kwargs)
153
154 for i,name in enumerate(allnames):
155 for ext in ["min", "max"]:
156 if re.match("%s[- _]%s" % (re.escape(series.name), ext), name):
157 if color:
158 c = color[i]
159 series2 = serieslist[i]
160 x2 = numpy.arange(series2.data.length) * series2.deltaT\
161 + float(series2.epoch) - float(t0)
162 x2 /= unit
163 d2 = series2.data.data
164 if logy and ylim:
165 numpy.putmask(d2, d2==0, ylim[0]-abs(ylim[0])*0.01)
166 plot.ax.plot(x2, d2, color=c, **kwargs2)
167 plot.ax.fill_between(x2, d, d2, alpha=0.1,\
168 color=c)
169
170
171 plot.finalize(loc=loc, alpha=alpha)
172 if hidden_colorbar:
173 plotutils.add_colorbar(plot.ax, visible=False)
174
175
176 axis_lims = plot.ax.get_ylim()
177 if zeroline:
178 plot.ax.plot([0, 0], [axis_lims[0], axis_lims[1]], 'r--', linewidth=2)
179 plot.ax.set_ylim([ axis_lims[0], axis_lims[1] ])
180
181
182 if logx:
183 plot.ax.xaxis.set_scale("log")
184 if logy:
185 plot.ax.yaxis.set_scale("log")
186 plot.ax._update_transScale()
187
188
189 if xlim:
190 xlim = (numpy.asarray(xlim).astype(float)-t0)/unit
191 plot.ax.set_xlim(xlim)
192 if ylim:
193 plot.ax.set_ylim(ylim)
194
195 plotutils.set_time_ticks(plot.ax)
196 plotutils.set_minor_ticks(plot.ax, x=False)
197
198
199 plot.savefig(outfile, bbox_inches=bbox_inches,\
200 bbox_extra_artists=plot.ax.texts)
201 plot.close()
202
203
204
205
206
208 """
209 Plot a (swig)LAL FrequencySeries.
210 """
211
212 if hasattr(series, "__contains__"):
213 serieslist = series
214 else:
215 serieslist = [series]
216
217
218 xlim = kwargs.pop("xlim", None)
219 ylim = kwargs.pop("ylim", None)
220 if not xlim:
221 fmin = min(float(s.f0) for s in serieslist)
222 fmax = max(float(s.f0 + s.data.length*s.deltaF) for s in serieslist)
223 xlim = [fmin, fmax]
224
225
226 logx = kwargs.pop("logx", False)
227 logy = kwargs.pop("logy", False)
228
229
230 loc = kwargs.pop("loc", 0)
231 alpha = kwargs.pop("alpha", 0.8)
232
233
234 hidden_colorbar = kwargs.pop("hidden_colorbar", False)
235
236
237 bbox_inches = kwargs.pop("bbox_inches", None)
238
239
240
241
242
243 xlabel = kwargs.pop("xlabel", "Frequency (Hz)")
244 ylabel = kwargs.pop("ylabel", "Amplitude")
245 title = kwargs.pop("title", "")
246 subtitle = kwargs.pop("subtitle", "")
247
248
249
250
251
252 color = kwargs.pop("color", None)
253 if isinstance(color, str):
254 color = [color]*len(serieslist)
255 kwargs2 = dict()
256 kwargs2.update(kwargs)
257 if kwargs.has_key("marker"):
258 kwargs.setdefault("markersize", 5)
259 kwargs2.setdefault("markersize", 2)
260 else:
261 kwargs.setdefault("linestyle", "-")
262 kwargs.setdefault("linewidth", "1")
263 kwargs2.setdefault("linewidth", "0.1")
264
265
266
267
268
269 allnames = [s.name for s in serieslist]
270 namedseries = [s for s in serieslist if not re.search("(min|max)\Z",s.name)]
271
272 plot = plotutils.SimplePlot(xlabel, ylabel, title, subtitle)
273 for i,(series,c) in enumerate(itertools.izip(namedseries,\
274 plotutils.default_colors())):
275 if color:
276 c = color[serieslist.index(series)]
277 if series.f_array is not None:
278 x = series.f_array
279 else:
280 x = numpy.arange(series.data.length) * series.deltaF + series.f0
281 x = x.astype(float)
282 plot.add_content(x, series.data.data, color=c,
283 label=plotutils.display_name(series.name), **kwargs)
284
285 for i,name in enumerate(allnames):
286 for ext in ["min", "max"]:
287 if re.match("%s[- _]%s" % (re.escape(series.name), ext), name):
288 if color:
289 c = color[i]
290 series2 = serieslist[i]
291 if series2.f_array is not None:
292 x2 = series2.f_array
293 else:
294 x2 = numpy.arange(series2.data.length) * series2.deltaF\
295 + series2.f0
296 plot.ax.plot(x2, series2.data.data, color=c, **kwargs2)
297
298 if series.data.data.shape == series2.data.data.shape:
299 plot.ax.fill_between(x2, series.data.data,\
300 series2.data.data, alpha=0.1,\
301 color=c)
302
303
304 plot.finalize(loc=loc, alpha=alpha)
305 if hidden_colorbar:
306 plotutils.add_colorbar(plot.ax, visible=False)
307
308
309 if logx:
310 plot.ax.xaxis.set_scale("log")
311 if logy:
312 plot.ax.yaxis.set_scale("log")
313 plot.ax._update_transScale()
314
315
316 if xlim:
317 plot.ax.set_xlim(xlim)
318 if ylim:
319 plot.ax.set_ylim(ylim)
320
321 plot.ax.grid(True, which="both")
322 plotutils.set_minor_ticks(plot.ax)
323
324
325 plot.savefig(outfile, bbox_inches=bbox_inches,\
326 bbox_extra_artists=plot.ax.texts)
327 plot.close()
328
329
330
331
332
333 -def plotspectrogram(sequencelist, outfile, epoch=0, deltaT=1, f0=0, deltaF=1,\
334 t0=0, ydata=None, **kwargs):
335 """
336 Plots a list of REAL?VectorSequences on a time-frequency-amplitude colour
337 map. The epochand deltaT arguments define the spacing in the x direction
338 for the given VectorSequence, and similarly f0 and deltaF in the
339 y-direction. If a list of VectorSequences is given, any of these arguments
340 can be in list form, one for each sequence.
341
342 ydata can be given as to explicitly define the frequencies at which the
343 sequences are sampled.
344 """
345
346 if not hasattr(sequencelist, "__contains__"):
347 sequencelist = [sequencelist]
348 numseq = len(sequencelist)
349
350
351 if isinstance(epoch, numbers.Number) or isinstance(epoch, lal.LIGOTimeGPS):
352 epoch = [epoch]*numseq
353 epoch = map(float, epoch)
354 if not len(epoch) == numseq:
355 raise ValueError("Wrong number of epoch arguments given.")
356 if isinstance(deltaT, numbers.Number):
357 deltaT = [deltaT]*numseq
358 deltaT = map(float, deltaT)
359 if not len(deltaT) == numseq:
360 raise ValueError("Wrong number of deltaT arguments given.")
361 if not ydata is None:
362 if isinstance(f0, numbers.Number):
363 f0 = [f0]*numseq
364 f0 = map(float, f0)
365 if not len(f0) == numseq:
366 raise ValueError("Wrong number of f0 arguments given.")
367 if isinstance(deltaF, numbers.Number):
368 deltaF = [deltaF]*numseq
369 deltaF = map(float, deltaF)
370 if not len(deltaF) == numseq:
371 raise ValueError("Wrong number of deltaF arguments given.")
372
373
374 xlim = kwargs.pop("xlim", None)
375 ylim = kwargs.pop("ylim", None)
376 colorlim = kwargs.pop("colorlim", None)
377 if xlim:
378 start,end = xlim
379 else:
380 start = min(epoch)
381 end = max(e + l.length * dt\
382 for e,l,dt in zip(epoch, sequencelist, deltaT))
383 if not ydata is None and not ylim:
384 ylim = [ydata.min(), ydata.max()]
385
386
387 logx = kwargs.pop("logx", False)
388 logy = kwargs.pop("logy", False)
389 logcolor = kwargs.pop("logcolor", False)
390
391
392 loc = kwargs.pop("loc", 0)
393 alpha = kwargs.pop("alpha", 0.8)
394
395
396 hidden_colorbar = kwargs.pop("hidden_colorbar", False)
397
398
399 bbox_inches = kwargs.pop("bbox_inches", None)
400
401
402
403
404
405 xlabel = kwargs.pop("xlabel", None)
406 if xlabel:
407 unit = 1
408 if not xlabel:
409 unit, timestr = plotutils.time_axis_unit(end-start)
410 if not t0:
411 t0 = start
412 t0 = lal.LIGOTimeGPS(t0)
413 if int(t0.gpsNanoSeconds) == 0:
414 xlabel = datetime.datetime(*lal.GPSToUTC(int(t0))[:6])\
415 .strftime("%B %d %Y, %H:%M:%S %ZUTC")
416 xlabel = "Time (%s) since %s (%s)" % (timestr, xlabel, int(t0))
417 else:
418 xlabel = datetime.datetime(*lal.GPSToUTC(t0.gpsSeconds)[:6])\
419 .strftime("%B %d %Y, %H:%M:%S %ZUTC")
420 xlabel = "Time (%s) since %s (%s)"\
421 % (timestr,
422 xlabel.replace(" UTC",".%.3s UTC" % t0.gpsNanoSeconds),\
423 t0)
424 t0 = float(t0)
425 ylabel = kwargs.pop("ylabel", "Frequency (Hz)")
426 colorlabel = kwargs.pop("colorlabel", "Amplitude")
427 title = kwargs.pop("title", "")
428 subtitle = kwargs.pop("subtitle", "")
429
430
431
432
433
434 interpolate = logy and ydata is None
435
436 for i,sequence in enumerate(sequencelist):
437 if interpolate:
438
439 sequence,ydata = loginterpolate(sequence, f0[i], deltaF[i])
440 if logy and ylim:
441 plotted = (ydata > ylim[0]) & (ydata <= ylim[1])
442 newVectorLength = int(plotted.sum())
443 newsequence = lal.CreateREAL8VectorSequence(sequence.length,\
444 newVectorLength)
445 for j in range(sequence.length):
446 newsequence.data[j,:] = sequence.data[j,:][plotted]
447 del sequence
448 sequencelist[i] = newsequence
449 if len(sequencelist) and logy and ylim:
450 ydata = ydata[plotted]
451
452
453
454
455
456 xbins = []
457 for i in range(numseq):
458 xmin = epoch[i]
459 xmax = epoch[i] + sequencelist[i].length * deltaT[i]
460 xbins.append(rate.LinearBins(float(xmin-t0)/unit, float(xmax-t0)/unit,\
461 2))
462
463 ybins = []
464 for i in range(numseq):
465 if ydata is not None:
466 ydata = numpy.asarray(ydata)
467 ymin = ydata.min()
468 ymax = ydata.max()
469 else:
470 ymin = f0[i]
471 ymax = f0[i] + sequencelist[i].vectorLength * deltaF[i]
472 if logy:
473 if ymin == 0:
474 ymin = deltaF[i]
475 ybins.append(rate.LogarithmicBins(ymin, ymax, 2))
476 else:
477 ybins.append(rate.LinearBins(ymin, ymax, 2))
478
479
480
481
482
483 kwargs.setdefault("interpolation", "kaiser")
484
485 plot = plotutils.ImagePlot(xlabel=xlabel, ylabel=ylabel, title=title,\
486 subtitle=subtitle, colorlabel=colorlabel)
487
488 for sequence,x,y in zip(sequencelist, xbins, ybins):
489 data = numpy.ma.masked_where(numpy.isnan(sequence.data), sequence.data,\
490 copy=False)
491 plot.add_content(data.T, x, y, **kwargs)
492
493
494 plot.finalize(colorbar=True, logcolor=logcolor, minorticks=True,\
495 clim=colorlim)
496 if hidden_colorbar:
497 plotutils.add_colorbar(plot.ax, visible=False)
498
499
500 if logx:
501 plot.ax.xaxis.set_scale("log")
502 if logy:
503 plot.ax.yaxis.set_scale("log")
504 plot.ax._update_transScale()
505
506
507 if xlim:
508 xlim = (numpy.asarray(xlim).astype(float)-t0)/unit
509 plot.ax.set_xlim(xlim)
510 if ylim:
511 plot.ax.set_ylim(ylim)
512
513
514 plot.ax.grid(True, which="both")
515 plotutils.set_time_ticks(plot.ax)
516 plotutils.set_minor_ticks(plot.ax)
517
518
519 plot.savefig(outfile, bbox_inches=bbox_inches,\
520 bbox_extra_artists=plot.ax.texts)
521 plot.close()
522
524 """
525 Interpolate a REAL?VectorSequence into logarithmic spacing in the second
526 dimension (y-axis, or frequency).
527 """
528
529 if not N:
530 N = sequence.vectorLength
531 ylin = numpy.arange(sequence.vectorLength)*deltaY + y0
532 ymin = y0 and y0 or deltaY
533 ymax = y0 + sequence.vectorLength * deltaY
534 ylog = numpy.logspace(numpy.log10(ymin*1.01), numpy.log10(ymax*0.99), num=N,\
535 endpoint=False)
536
537
538 datatype = seriesutils.typecode(type(sequence))
539 TYPESTR = seriesutils._typestr[datatype]
540 func = getattr(lal, "Create%sVectorSequence" % TYPESTR)
541 out = func(sequence.length, N)
542
543
544 for i in range(sequence.length):
545 intplt = interpolate.interp1d(ylin, sequence.data[i,:])
546 out.data[i,:] = intplt(ylog)
547
548 return out,ylog
549
550
551
552
553
554 -def plothistogram(serieslist, outfile, nonzero=False, num_bins=100,\
555 cumulative=False, bar=False, **kwargs):
556 """
557 Plots a line histogram of the data contained in the series, or list of
558 series'.
559 """
560
561 if not hasattr(serieslist, "__contains__"):
562 serieslist = [serieslist]
563
564
565 xlim = kwargs.pop("xlim", None)
566 ylim = kwargs.pop("ylim", None)
567
568
569 xlabel = kwargs.pop("xlabel", "Amplitude")
570 ylabel = kwargs.pop("ylabel", "Fraction of data")
571 title = kwargs.pop("title", "")
572 subtitle = kwargs.pop("subtitle","")
573
574
575 logx = kwargs.pop("logx", False)
576 logy = kwargs.pop("logy", False)
577
578
579 hidden_colorbar = kwargs.pop("hidden_colorbar", False)
580
581
582 bbox_inches = kwargs.pop("bbox_inches", "tight")
583
584
585 fill = kwargs.pop("fill", False)
586
587
588 loc = kwargs.pop("loc", 0)
589 alpha = kwargs.pop("alpha", 0.8)
590
591
592 color = kwargs.pop("color", None)
593 if isinstance(color, str):
594 color = color.split(',')
595 if len(color) == 1:
596 color = [color[0]]*len(serieslist)
597
598
599
600
601
602
603 plot = plotutils.LineHistogram(xlabel, ylabel, title, subtitle)
604
605 for i, series in enumerate(serieslist):
606 data = series.data.data
607 if nonzero:
608 data = data[data!=0]
609 if color:
610 kwargs["color"] = color[i]
611 plot.add_content(data, normalize=data.size, label=series.name, **kwargs)
612 plot.finalize(loc=loc, alpha=alpha, logx=logx, logy=logy, bar=bar,\
613 xlim=xlim, num_bins=num_bins, cumulative=cumulative)
614 if hidden_colorbar:
615 plotutils.add_colorbar(plot.ax, visible=False)
616
617 plot.ax.autoscale_view(tight=True, scalex=True, scaley=True)
618
619 if xlim:
620 plot.ax.set_xlim(xlim)
621 if ylim:
622 plot.ax.set_ylim(ylim)
623
624 plotutils.set_minor_ticks(plot.ax)
625
626
627 plot.savefig(outfile, bbox_inches=bbox_inches,\
628 bbox_extra_artists=plot.ax.texts)
629 plot.close()
630