plotting.py 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423
  1. """
  2. plotting.py
  3. The functions in this module are meant for plotting the histogram objects created via
  4. filval.histogram
  5. """
  6. from collections import defaultdict
  7. from itertools import zip_longest
  8. from io import BytesIO
  9. from base64 import b64encode
  10. import numpy as np
  11. import matplotlib.pyplot as plt
  12. from markdown import Markdown
  13. import latexipy as lp
  14. from filval.histogram import (hist, hist2d, hist_bin_centers, hist_fit,
  15. hist_norm, hist_stats)
  16. __all__ = ['Plot',
  17. 'decl_plot',
  18. 'grid_plot',
  19. 'render_plots',
  20. 'generate_dashboard',
  21. 'hist_plot',
  22. 'hist_plot_stack',
  23. 'hist2d_plot',
  24. 'hists_to_table',
  25. 'simple_plot']
  26. class Plot:
  27. def __init__(self, subplots, name, title=None, docs="N/A", arg_dicts=None):
  28. if type(subplots) is not list:
  29. subplots = [[subplots]]
  30. elif len(subplots) > 0 and type(subplots[0]) is not list:
  31. subplots = [subplots]
  32. self.subplots = subplots
  33. self.name = name
  34. self.title = title
  35. self.docs = docs
  36. self.arg_dicts = arg_dicts if arg_dicts is not None else {}
  37. MD = Markdown(extensions=['mdx_math', 'tables'],
  38. extension_configs={'mdx_math': {'enable_dollar_delimiter': True}})
  39. lp.latexify(params={'pgf.texsystem': 'pdflatex',
  40. 'text.usetex': True,
  41. 'font.family': 'serif',
  42. 'pgf.preamble': [],
  43. 'font.size': 15,
  44. 'axes.labelsize': 15,
  45. 'axes.titlesize': 13,
  46. 'legend.fontsize': 13,
  47. 'xtick.labelsize': 11,
  48. 'ytick.labelsize': 11,
  49. 'figure.dpi': 150,
  50. 'savefig.transparent': False,
  51. },
  52. new_backend='TkAgg')
  53. def _fn_call_to_dict(fn, *args, **kwargs):
  54. from inspect import signature
  55. from html import escape
  56. pnames = list(signature(fn).parameters)
  57. pvals = list(args) + list(kwargs.values())
  58. return {escape(str(k)): escape(str(v)) for k, v in zip(pnames, pvals)}
  59. def _process_docs(fn):
  60. from inspect import getdoc
  61. raw = getdoc(fn)
  62. if raw:
  63. return MD.convert(raw)
  64. else:
  65. return None
  66. def decl_plot(fn):
  67. from functools import wraps
  68. @wraps(fn)
  69. def f(*args, **kwargs):
  70. txt = fn(*args, **kwargs)
  71. argdict = _fn_call_to_dict(fn, *args, **kwargs)
  72. docs = _process_docs(fn)
  73. if not txt:
  74. txt = ''
  75. txt = MD.convert(txt)
  76. return argdict, docs, txt
  77. return f
  78. def simple_plot(thx, *args, log=None, **kwargs):
  79. import ROOT
  80. if isinstance(thx, ROOT.TH2):
  81. def f(h):
  82. hist2d_plot(hist2d(h), *args, **kwargs)
  83. plt.xlabel(h.GetXaxis().GetTitle())
  84. plt.ylabel(h.GetYaxis().GetTitle())
  85. if log == 'x':
  86. plt.semilogx()
  87. elif log == 'y':
  88. plt.semilogy()
  89. elif log == 'xy':
  90. plt.loglog()
  91. return dict(), "", ""
  92. return Plot([[(f, (thx,), {})]], thx.GetName())
  93. elif isinstance(thx, ROOT.TH1):
  94. def f(h):
  95. hist_plot(hist(h), *args, **kwargs)
  96. plt.xlabel(h.GetXaxis().GetTitle())
  97. plt.ylabel(h.GetYaxis().GetTitle())
  98. if log == 'x':
  99. plt.semilogx()
  100. elif log == 'y':
  101. plt.semilogy()
  102. elif log == 'xy':
  103. plt.loglog()
  104. return dict(), "", ""
  105. return Plot([[(f, (thx,), {})]], thx.GetName())
  106. else:
  107. raise ValueError("must call simple_plot with a ROOT TH1 or TH2 object")
  108. def generate_dashboard(plots, title, output='dashboard.html', template='dashboard.j2',
  109. source=None, ana_source=None, config=None, header=None):
  110. from jinja2 import Environment, PackageLoader, select_autoescape
  111. from os.path import join, isdir
  112. from os import mkdir
  113. from urllib.parse import quote
  114. env = Environment(
  115. loader=PackageLoader('filval', 'templates'),
  116. autoescape=select_autoescape(['htm', 'html', 'xml']),
  117. )
  118. env.globals.update({'quote': quote,
  119. 'enumerate': enumerate,
  120. 'zip': zip,
  121. })
  122. def get_by_n(objects, n=2):
  123. objects = list(objects)
  124. while objects:
  125. yield objects[:n]
  126. objects = objects[n:]
  127. if source is not None:
  128. with open(source, 'r') as f:
  129. source = f.read()
  130. if not isdir('output'):
  131. mkdir('output')
  132. if header is not None:
  133. header = MD.convert(header)
  134. dashboard_path = join('output', output)
  135. with open(dashboard_path, 'w') as tempout:
  136. templ = env.get_template(template)
  137. tempout.write(templ.render(
  138. plots=get_by_n(plots, 3),
  139. title=title,
  140. source=source,
  141. ana_source=ana_source,
  142. config=config,
  143. header=header,
  144. ))
  145. return dashboard_path
  146. def _add_stats(hist, title=''):
  147. fmt = r'''\begin{{eqnarray*}}
  148. \sum{{x_i}} &=& {sum:5.3f} \\
  149. \sum{{\Delta x_i \cdot x_i}} &=& {int:5.3G} \\
  150. \mu &=& {mean:5.3G} \\
  151. \sigma^2 &=& {var:5.3G} \\
  152. \sigma &=& {std:5.3G}
  153. \end{{eqnarray*}}'''
  154. txt = fmt.format(**hist_stats(hist), title=title)
  155. txt = txt.replace('\n', ' ')
  156. plt.text(0.7, 0.9, txt,
  157. bbox={'facecolor': 'white',
  158. 'alpha': 0.7,
  159. 'boxstyle': 'square,pad=0.8'},
  160. transform=plt.gca().transAxes,
  161. verticalalignment='top',
  162. horizontalalignment='left',
  163. size='small')
  164. if title:
  165. plt.text(0.72, 0.97, title,
  166. bbox={'facecolor': 'white',
  167. 'alpha': 0.8},
  168. transform=plt.gca().transAxes,
  169. verticalalignment='top',
  170. horizontalalignment='left')
  171. def grid_plot(subplots):
  172. if any(len(row) != len(subplots[0]) for row in subplots):
  173. raise ValueError('make_plot requires a rectangular list-of-lists as '
  174. 'input. Fill empty slots with None')
  175. def calc_row_span(fig, row, col):
  176. span = 1
  177. for r in range(row + 1, len(fig)):
  178. if fig[r][col] == 'FU':
  179. span += 1
  180. else:
  181. break
  182. return span
  183. def calc_column_span(fig, row, col):
  184. span = 1
  185. for c in range(col + 1, len(fig[row])):
  186. if fig[row][c] == 'FL':
  187. span += 1
  188. else:
  189. break
  190. return span
  191. rows = len(subplots)
  192. cols = len(subplots[0])
  193. argdicts = defaultdict(list)
  194. docs = defaultdict(list)
  195. txts = defaultdict(list)
  196. for i in range(rows):
  197. for j in range(cols):
  198. cell = subplots[i][j]
  199. if cell in ('FL', 'FU', None):
  200. continue
  201. if not isinstance(cell, list):
  202. cell = [cell]
  203. column_span = calc_column_span(subplots, i, j)
  204. row_span = calc_row_span(subplots, i, j)
  205. plt.subplot2grid((rows, cols), (i, j),
  206. colspan=column_span, rowspan=row_span)
  207. for plot in cell:
  208. if not isinstance(plot, tuple):
  209. plot_fn, args, kwargs = plot, (), {}
  210. elif len(plot) == 1:
  211. plot_fn, args, kwargs = plot[0], (), {}
  212. elif len(plot) == 2:
  213. plot_fn, args, kwargs = plot[0], plot[1], {}
  214. elif len(plot) == 3:
  215. plot_fn, args, kwargs = plot[0], plot[1], plot[2]
  216. else:
  217. raise ValueError('Plot tuple must be of format (func), '
  218. f'or (func, tuple), or (func, tuple, dict). Got {plot}')
  219. this_args, this_docs, txt = plot_fn(*args, **kwargs)
  220. argdicts[(i, j)].append(this_args)
  221. docs[(i, j)].append(this_docs)
  222. txts[(i, j)].append(txt)
  223. return argdicts, docs, txts
  224. def render_plots(plots, exts=('png',), directory='output/figures/', scale=1.0, to_disk=True):
  225. for plot in plots:
  226. print(f'Building plot {plot.name}')
  227. plot.data = None
  228. if to_disk:
  229. with lp.figure(plot.name.replace(' ', '_'), directory=directory,
  230. exts=exts,
  231. size=(scale * 10, scale * 10)):
  232. argdicts, docs, txts = grid_plot(plot.subplots)
  233. else:
  234. out = BytesIO()
  235. with lp.mem_figure(out,
  236. ext=exts[0],
  237. size=(scale * 10, scale * 10)):
  238. argdicts, docs, txts = grid_plot(plot.subplots)
  239. out.seek(0)
  240. plot.data = b64encode(out.read()).decode()
  241. plot.argdicts = argdicts
  242. plot.docs = docs
  243. plot.txts = txts
  244. def add_decorations(axes, luminosity, energy):
  245. cms_prelim = r'{\raggedright{}\textsf{\textbf{CMS}}\\ \emph{Preliminary}}'
  246. axes.text(0.01, 0.98, cms_prelim,
  247. horizontalalignment='left',
  248. verticalalignment='top',
  249. transform=axes.transAxes)
  250. lumi = ""
  251. energy_str = ""
  252. if luminosity is not None:
  253. lumi = r'${} \mathrm{{fb}}^{{-1}}$'.format(luminosity)
  254. if energy is not None:
  255. energy_str = r'({} TeV)'.format(energy)
  256. axes.text(1, 1, ' '.join([lumi, energy_str]),
  257. horizontalalignment='right',
  258. verticalalignment='bottom',
  259. transform=axes.transAxes)
  260. def hist_plot(h, *args, include_errors=False, fit=None, stats=False, **kwargs):
  261. """ Plots a 1D ROOT histogram object using matplotlib """
  262. from inspect import signature
  263. values, errors, edges = h
  264. left, right = np.array(edges[:-1]), np.array(edges[1:])
  265. x = np.array([left, right]).T.flatten()
  266. y = np.array([values, values]).T.flatten()
  267. title = kwargs.pop('title', '')
  268. plt.plot(x, y, *args, linewidth=1, **kwargs)
  269. if include_errors:
  270. plt.errorbar(hist_bin_centers(h), values, yerr=errors,
  271. color='k', marker=None, linestyle='None',
  272. barsabove=True, elinewidth=.7, capsize=1)
  273. if fit:
  274. f, p0 = fit
  275. popt, pcov = hist_fit(h, f, p0)
  276. fit_xs = np.linspace(x[0], x[-1], 100)
  277. fit_ys = f(fit_xs, *popt)
  278. plt.plot(fit_xs, fit_ys, '--g')
  279. arglabels = list(signature(f).parameters)[1:]
  280. label_txt = "\n".join('{:7s}={: 0.2G}'.format(label, value)
  281. for label, value in zip(arglabels, popt))
  282. plt.text(0.60, 0.95, label_txt, va='top', transform=plt.gca().transAxes,
  283. fontsize='medium', family='monospace', usetex=False)
  284. if stats:
  285. _add_stats(h, title)
  286. def hist2d_plot(h, txt_format=None, colorbar=False, **kwargs):
  287. """ Plots a 2D ROOT histogram object using matplotlib """
  288. try:
  289. values, errors, xs, ys = h
  290. except (TypeError, ValueError):
  291. values, errors, xs, ys = hist2d(h)
  292. plt.xlabel(kwargs.pop('xlabel', ''))
  293. plt.ylabel(kwargs.pop('ylabel', ''))
  294. plt.title(kwargs.pop('title', ''))
  295. plt.pcolormesh(xs, ys, values, **kwargs)
  296. if txt_format is not None:
  297. cmap = plt.get_cmap()
  298. min_, max_ = float(np.min(values)), float(np.max(values))
  299. def get_intensity(val):
  300. cmap_idx = int((cmap.N-1) * (val - min_) / (max_-min_))
  301. color = cmap.colors[cmap_idx]
  302. return color[0]*0.25 + color[1]*0.5 + color[2]*0.25
  303. for idx_row in range(values.shape[0]):
  304. for idx_col in range(values.shape[1]):
  305. x_mid = (xs[idx_row, idx_col] + xs[idx_row, idx_col+1]) / 2
  306. y_mid = (ys[idx_row, idx_col] + ys[idx_row+1, idx_col]) / 2
  307. val = txt_format.format(values[idx_row, idx_col])
  308. txt_color = 'w' if get_intensity(values[idx_row, idx_col]) < 0.5 else 'k'
  309. plt.text(x_mid, y_mid, val, verticalalignment='center', horizontalalignment='center',
  310. color=txt_color, fontsize=12)
  311. if colorbar:
  312. plt.colorbar()
  313. def hist_plot_stack(hists: list, labels: list = None):
  314. """
  315. Creates a stacked histogram in the current axes.
  316. :param hists: list of histogram
  317. :param labels:
  318. :return:
  319. """
  320. if len(hists) == 0:
  321. return
  322. if len(set([len(hist[0]) for hist in hists])) != 1:
  323. raise ValueError("all histograms must have the same number of bins")
  324. if labels is None:
  325. labels = [None for _ in hists]
  326. if len(labels) != len(hists):
  327. raise ValueError("Label mismatch")
  328. bottoms = [0 for _ in hists[0][0]]
  329. for hist, label in zip(hists, labels):
  330. centers = []
  331. widths = []
  332. heights = []
  333. for left, right, content in zip(hist[2][:-1], hist[2][1:], hist[0]):
  334. centers.append((right + left) / 2)
  335. widths.append(right - left)
  336. heights.append(content)
  337. plt.bar(centers, heights, widths, bottoms, label=label)
  338. for i, content in enumerate(hist[0]):
  339. bottoms[i] += content
  340. def hists_to_table(hists, row_labels=(), column_labels=(), format="{:.2f}"):
  341. table = ['<table class="table table-condensed">']
  342. if column_labels:
  343. table.append('<thead><tr>')
  344. if row_labels:
  345. table.append('<th></th>')
  346. table.extend(f'<th>{label}</th>' for label in column_labels)
  347. table.append('</tr></thead>')
  348. table.append('<tbody>\n')
  349. for row_label, (vals, *_) in zip_longest(row_labels, hists):
  350. table.append('<tr>')
  351. if row_label:
  352. table.append(f'<td><strong>{row_label}</strong></td>')
  353. table.extend(('<td>'+format.format(val)+'</td>') for val in vals)
  354. table.append('</tr>\n')
  355. table.append('</tbody></table>')
  356. return ''.join(table)