#!/usr/bin/env python import numpy as np import matplotlib.pyplot as plt from filval.result_set import ResultSet from filval.histogram import hist, hist2d, hist_add, hist_norm, hist_scale, hist2d_norm from filval.plotting import (decl_plot, render_plots, hist_plot, hist_plot_stack, hist2d_plot, Plot, generate_dashboard, hists_to_table) an_tttt = ([0.47, 0.33, 0.18, 0.78, 0.49, 0.52, 0.33, 0.49], [0, 0, 0, 0, 0, 0, 0, 0], [0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5]) an_ttw = ([2.29663, 0.508494, 0.161166, 1.03811, 0.256401, 0.127582, 0.181522, 0.141659], [0, 0, 0, 0, 0, 0, 0, 0], [0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5]) an_ttz = ([0.974751, 0.269195, 1e-06, 0.395831, 0.0264703, 0.06816, 0.8804, 0.274265], [0, 0, 0, 0, 0, 0, 0, 0], [0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5]) an_tth = ([1.13826, 0.361824, 0.162123, 0.683917, 0.137608, 0.0632719, 0.554491, 0.197864], [0, 0, 0, 0, 0, 0, 0, 0], [0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5]) @decl_plot def plot_yield_grid(rss, tau_category=-1): r"""## Event Yield The event yield for the eight signal regions defined in AN-17-115. Data is normalized to the integrated luminosity of $35.9\textrm{fb}^{-1}$. Ignoring taus means both that there is no requirement on number of good taus *and* the taus, if present, are not considered for the SS pair. If taus are not ignored, any good tau with $p_T$>20Gev is considered in constructing the SS lepton pair. The yields are then further broken down by the number of good taus in the event. A "good" tau in the above means any tau candidate passing the `byTightIsolationMVArun2v1DBoldDMwLT` ID w/ pt>20GeV. It is also required to pass tau-lepton cross cleaning where it must not match any electron or muon within $\delta R < 0.4$. Truth-matched taus are those that match within $\delta R < 0.3$ with gen-level taus that pass the flag `fromHardProcessDecayed`. """ def get_sr(rs, tm=False): if tm: if tau_category == 0: return hist(rs.SRs_0tmtau) if tau_category == 1: return hist(rs.SRs_1tmtau) elif tau_category == 2: return hist(rs.SRs_2tmtau) else: if tau_category == -1: return hist(rs.ignore_tau_SRs) elif tau_category == 0: return hist(rs.SRs_0tau) elif tau_category == 1: return hist(rs.SRs_1tau) elif tau_category == 2: return hist(rs.SRs_2tau) _, ((ax_tttt, ax_ttw), (ax_ttz, ax_tth)) = plt.subplots(2, 2) tttt, ttw, ttz, tth = [get_sr(rs) for rs in rss] tm_tttt, tm_ttw, tm_ttz, tm_tth = [get_sr(rs, True) for rs in rss] plt.sca(ax_tttt) hist_plot(tttt, title='TTTT', stats=False, label='Mock', include_errors=True) if tau_category == -1 and len(tttt[0]) == len(an_tttt[0]): hist_plot(an_tttt, title='TTTT', stats=False, label='AN') elif tau_category >= 0: hist_plot(tm_tttt, title='TTTT', stats=False, label='Truth-Matched Taus', include_errors=True) plt.ylim((0, None)) plt.sca(ax_ttw) hist_plot(ttw, title='TTW', stats=False, label='Mock', include_errors=True) if tau_category == -1 and len(tttt[0]) == len(an_tttt[0]): hist_plot(an_ttw, title='TTW', stats=False, label='AN') elif tau_category >= 0: hist_plot(tm_ttw, title='TTW', stats=False, label='Truth-Matched Taus', include_errors=True) plt.ylim((0, None)) plt.legend() plt.sca(ax_ttz) hist_plot(ttz, title='TTZ', stats=False, label='Mock', include_errors=True) if tau_category == -1 and len(tttt[0]) == len(an_tttt[0]): hist_plot(an_ttz, title='TTZ', stats=False, label='AN') elif tau_category >= 0: hist_plot(tm_ttz, title='TTZ', stats=False, label='Truth-Matched Taus', include_errors=True) plt.ylim((0, None)) plt.xlabel('Signal Region') plt.sca(ax_tth) hist_plot(tth, title='TTH', stats=False, label='Mock', include_errors=True) if tau_category == -1 and len(tttt[0]) == len(an_tttt[0]): hist_plot(an_tth, title='TTH', stats=False, label='AN') elif tau_category >= 0: hist_plot(tm_tth, title='TTH', stats=False, label='Truth-Matched Taus', include_errors=True) plt.ylim((0, None)) plt.xlabel('Signal Region') def to_table(hists): return hists_to_table(hists, row_labels=['TTTT', 'TTW', 'TTZ', 'TTH'], column_labels=[f'SR{n}' for n in range(1, len(tttt[0])+1)]) tables = '

Mock

' tables += to_table([tttt, ttw, ttz, tth]) if tau_category == -1: tables += '

AN

' tables += to_table([an_tttt, an_ttw, an_ttz, an_tth]) elif tau_category >= 0: tables += '

TM Taus

' tables += to_table([tm_tttt, tm_ttw, tm_ttz, tm_tth]) return tables @decl_plot def plot_yield_v_gen(rss, tau_category=-1): r"""## Event Yield Vs. # of Generated Taus """ def get_sr(rs): h = None if tau_category == 0: h = rs.SRs_0tmtau_diff_nGenTau if tau_category == 1: h = rs.SRs_1tmtau_diff_nGenTau elif tau_category == 2: h = rs.SRs_2tmtau_diff_nGenTau return hist2d_norm(hist2d(h), norm=1, axis=0) _, ((ax_tttt, ax_ttw), (ax_ttz, ax_tth)) = plt.subplots(2, 2) tttt, ttw, ttz, tth = [get_sr(rs) for rs in rss] def do_plot(h, title): hist2d_plot(h, title=title, txt_format='{:.2f}') plt.sca(ax_tttt) do_plot(tttt, 'TTTT') plt.ylabel("\# Gen Taus") plt.sca(ax_ttw) do_plot(ttw, 'TTW') plt.legend() plt.sca(ax_ttz) do_plot(ttz, 'TTZ') plt.xlabel('Signal Region') plt.ylabel("\# Gen Taus") plt.sca(ax_tth) do_plot(tth, 'TTH') plt.xlabel('Signal Region') @decl_plot def plot_nGen_v_nSel(rss): _, ((ax_tttt, ax_ttw), (ax_ttz, ax_tth)) = plt.subplots(2, 2) tttt, ttw, ttz, tth = [hist2d(rs.nGen_v_RecoTaus_in_SR) for rs in rss] def do_plot(h, title): hist2d_plot(h, title=title, txt_format='{:.2f}') plt.sca(ax_tttt) do_plot(tttt, 'TTTT') plt.ylabel("\# Selected Taus") plt.sca(ax_ttw) do_plot(ttw, 'TTW') plt.legend() plt.sca(ax_ttz) do_plot(ttz, 'TTZ') plt.xlabel('\# Gen Taus') plt.ylabel("\# Selected Taus") plt.sca(ax_tth) do_plot(tth, 'TTH') plt.xlabel('\# Gen Taus') @decl_plot def plot_yield_stack(rss): r"""## Event Yield - Stacked The event yield for the eight signal regions defined in AN-17-115. Data is normalized to the Moriond 2018 integrated luminosity ($35.9\textrm{fb}^{-1}$). Code for the histogram generation is here: """ tttt, ttw, ttz, tth = map(lambda rs: hist(rs.SRs), rss) hist_plot_stack([ttw, ttz, tth], labels=['TTW', 'TTZ', 'TTH']) tttt = tttt[0]*10, tttt[1], tttt[2] hist_plot(tttt, label='TTTT (x10)', stats=False, color='k') plt.ylim((0, 60)) plt.xlabel('Signal Region') plt.legend() @decl_plot def plot_lep_multi(rss, dataset): _, (ax_els, ax_mus, ax_taus) = plt.subplots(3, 1) els = list(map(lambda rs: hist_norm(hist(rs.nEls)), rss)) mus = list(map(lambda rs: hist_norm(hist(rs.nMus)), rss)) taus = list(map(lambda rs: hist_norm(hist(rs.nTaus)), rss)) def _plot(ax, procs): plt.sca(ax) tttt, ttw, ttz, tth = procs h = {'TTTT': tttt, 'TTW': ttw, 'TTZ': ttz, 'TTH': tth}[dataset] hist_plot(h, stats=False, label=dataset) _plot(ax_els, els) plt.xlabel('\\# Good Electrons') plt.legend() _plot(ax_mus, mus) plt.xlabel('\\# Good Muons') _plot(ax_taus, taus) plt.xlabel('\\# Good Taus') @decl_plot def plot_sig_strength(rss): r""" The signal strength of the TTTT signal defined as $\frac{S}{\sqrt{S+B}}$ """ tttt, ttw, ttz, tth = map(lambda rs: hist(rs.SRs), rss) bg = hist_add(ttw, ttz, tth) strength = tttt[0] / np.sqrt(tttt[0] + bg[0]) hist_plot((strength, tttt[1], tttt[2]), stats=False) @decl_plot def plot_event_obs(rss, dataset, in_signal_region=True): r"""The distribution of $N_{jet}$, $N_{Bjet}$, MET, and $H_T$ in either all events or only signal region (igoring taus) events. """ _, ((ax_njet, ax_nbjet), (ax_ht, ax_met)) = plt.subplots(2, 2) # tttt, ttw, ttz, tth = map(lambda rs: hist(rs.SRs), rss) tttt, ttw, ttz, tth = rss rs = {'TTTT': tttt, 'TTW': ttw, 'TTZ': ttz, 'TTH': tth}[dataset] def _plot(ax, obs): plt.sca(ax) if in_signal_region: h = {'MET': rs.met_in_SR, 'HT': rs.ht_in_SR, 'NJET': rs.njet_in_SR, 'NBJET': rs.nbjet_in_SR}[obs] else: h = {'MET': rs.met, 'HT': rs.ht, 'NJET': rs.njet, 'NBJET': rs.nbjet}[obs] hist_plot(hist(h), stats=False, label=dataset, xlabel=obs) _plot(ax_njet, 'NJET') _plot(ax_nbjet, 'NBJET') _plot(ax_ht, 'HT') _plot(ax_met, 'MET') @decl_plot def plot_event_obs_stack(rss, in_signal_region=True): r""" """ _, ((ax_njet, ax_nbjet), (ax_ht, ax_met)) = plt.subplots(2, 2) def _plot(ax, obs): plt.sca(ax) if in_signal_region: attr = {'MET': 'met_in_SR', 'HT': 'ht_in_SR', 'NJET': 'njet_in_SR', 'NBJET': 'nbjet_in_SR'}[obs] else: attr = {'MET': 'met', 'HT': 'ht', 'NJET': 'njet', 'NBJET': 'nbjet'}[obs] tttt, ttw, ttz, tth = map(lambda rs: hist(getattr(rs, attr)), rss) hist_plot_stack([ttw, ttz, tth], labels=["TTW", "TTZ", "TTH"]) hist_plot(hist_scale(tttt, 5), label="TTTT (x5)", color='k') plt.xlabel(obs) _plot(ax_njet, 'NJET') _plot(ax_nbjet, 'NBJET') plt.legend() _plot(ax_ht, 'HT') _plot(ax_met, 'MET') # @decl_plot # def plot_s_over_b(rss): # def get_sr(rs): # h = None # if tau_category == 0: # h = rs.SRs_0tmtau_diff_nGenTau # if tau_category == 1: # h = rs.SRs_1tmtau_diff_nGenTau # elif tau_category == 2: # h = rs.SRs_2tmtau_diff_nGenTau # return hist2d_norm(hist2d(h), norm=1, axis=0) # # _, ((ax_tttt, ax_ttw), (ax_ttz, ax_tth)) = plt.subplots(2, 2) # tttt = # ttw, ttz, tth = [get_sr(rs) for rs in rss] # pass @decl_plot def plot_tau_purity(rss): _, ((ax_tttt, ax_ttw), (ax_ttz, ax_tth)) = plt.subplots(2, 2) tttt, ttw, ttz, tth = list(map(lambda rs: hist(rs.tau_purity_v_pt), rss)) def _plot(ax, dataset): plt.sca(ax) h = {'TTTT': tttt, 'TTW': ttw, 'TTZ': ttz, 'TTH': tth}[dataset] hist_plot(h, stats=False, label=dataset) plt.text(200, 0.05, dataset) plt.xlabel(r"$P_T$(GeV)") _plot(ax_tttt, 'TTTT') _plot(ax_ttw, 'TTW') _plot(ax_ttz, 'TTZ') _plot(ax_tth, 'TTH') if __name__ == '__main__': data_path = 'data/output_new_sr_new_id_binning/' save_plots = True # First create a ResultSet object which loads all of the objects from root file # into memory and makes them available as attributes rss = (ResultSet("tttt", data_path+'yield_tttt.root'), ResultSet("ttw", data_path+'yield_ttw.root'), ResultSet("ttz", data_path+'yield_ttz.root'), ResultSet("tth", data_path+'yield_tth.root')) # Next, declare all of the (sub)plots that will be assembled into full # figures later yield_tau_ignore_tau = plot_yield_grid, (rss, -1) yield_tau_0tau = plot_yield_grid, (rss, 0) yield_tau_1tau = plot_yield_grid, (rss, 1) yield_tau_2tau = plot_yield_grid, (rss, 2) tttt_event_obs_in_sr = plot_event_obs, (rss, 'TTTT'), {'in_signal_region': True} ttw_event_obs_in_sr = plot_event_obs, (rss, 'TTW'), {'in_signal_region': True} ttz_event_obs_in_sr = plot_event_obs, (rss, 'TTZ'), {'in_signal_region': True} tth_event_obs_in_sr = plot_event_obs, (rss, 'TTH'), {'in_signal_region': True} tttt_event_obs = plot_event_obs, (rss, 'TTTT'), {'in_signal_region': False} ttw_event_obs = plot_event_obs, (rss, 'TTW'), {'in_signal_region': False} ttz_event_obs = plot_event_obs, (rss, 'TTZ'), {'in_signal_region': False} tth_event_obs = plot_event_obs, (rss, 'TTH'), {'in_signal_region': False} tttt_lep_multi = plot_lep_multi, (rss, 'TTTT') ttw_lep_multi = plot_lep_multi, (rss, 'TTW') ttz_lep_multi = plot_lep_multi, (rss, 'TTZ') tth_lep_multi = plot_lep_multi, (rss, 'TTH') yield_v_gen_0tm = plot_yield_v_gen, (rss, 0) yield_v_gen_1tm = plot_yield_v_gen, (rss, 1) yield_v_gen_2tm = plot_yield_v_gen, (rss, 2) # tau_purity = plot_tau_purity, (rss) nGen_v_nSel = plot_nGen_v_nSel, (rss,) # Now assemble the plots into figures. plots = [ Plot([[yield_tau_ignore_tau]], 'Yield Ignoring Taus'), Plot([[yield_tau_0tau]], 'Yield For events with 0 Tau'), Plot([[yield_tau_1tau]], 'Yield For events with 1 Tau'), Plot([[yield_tau_2tau]], 'Yield For events with 2 or more Tau'), # Plot([[yield_v_gen_0tm]], # 'Yield For events with 0 TM Tau Vs Gen Tau'), # Plot([[yield_v_gen_1tm]], # 'Yield For events with 1 TM Tau Vs Gen Tau'), # Plot([[yield_v_gen_2tm]], # 'Yield For events with 2 TM Tau Vs Gen Tau'), Plot([[nGen_v_nSel]], r'#Generated tau vs. #Selected tau'), Plot([[tttt_lep_multi]], 'Lepton Multiplicity - TTTT'), Plot([[ttw_lep_multi]], 'Lepton Multiplicity - TTW'), Plot([[ttz_lep_multi]], 'Lepton Multiplicity - TTZ'), Plot([[tth_lep_multi]], 'Lepton Multiplicity - TTH'), Plot([[tttt_event_obs_in_sr]], 'TTTT - Event Observables (In SR)'), Plot([[ttw_event_obs_in_sr]], 'TTW - Event Observables (In SR)'), Plot([[ttz_event_obs_in_sr]], 'TTZ - Event Observables (In SR)'), Plot([[tth_event_obs_in_sr]], 'TTH - Event Observables (In SR)'), # Plot([[tttt_event_obs]], # 'TTTT - Event Observables (All Events)'), # Plot([[ttw_event_obs]], # 'TTW - Event Observables (All Events)'), # Plot([[ttz_event_obs]], # 'TTZ - Event Observables (All Events)'), # Plot([[tth_event_obs]], # 'TTH - Event Observables (All Events)'), # Plot([[yield_notau]], # 'Yield Without Tau'), # Plot([[yield_tau_stack]], # 'Yield With Tau Stacked'), # Plot([[yield_notau_stack]], # 'Yield Without Tau Stacked'), # Plot([[yield_tau_stack], # [yield_notau_stack]], # 'Event Yield, top: with tau, bottom: no tau'), # Plot([[sig_strength_tau], # [sig_strength_notau]], # 'Signal Strength'), # Plot([[event_obs_stack]], # 'Event Observables'), # Plot([[tau_purity]], # 'Tau Purity'), ] # Finally, render and save the plots and generate the html+bootstrap # dashboard to view them render_plots(plots, to_disk=save_plots) if not save_plots: generate_dashboard(plots, 'TTTT Yields', output='yields.html', source=__file__, ana_source=("https://github.com/cfangmeier/FTAnalysis/commit/" "0cbdac4509391fffb9fff87d0521b7dd0a30a55c"), config=data_path+'config.yaml' )