Skip to content

Add signal in postfit stack #128

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
Quantumapple wants to merge 12 commits into
base: master
Choose a base branch
from
Open

Add signal in postfit stack #128

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
6ce14f2
add signal
Quantumapple Apr 11, 2025
696c1ef
update style
Quantumapple Apr 12, 2025
5faf518
switch to use prefit for signal histogram
Quantumapple Apr 15, 2025
bb8059f
update
Quantumapple Apr 18, 2025
a556810
style update
Quantumapple Apr 20, 2025
4516acd
nb for 2D limit plotting
Quantumapple Aug 6, 2025
67fdce6
Merge branch 'mcremone:master' into master
Quantumapple Aug 6, 2025
86a3b88
update
Quantumapple Aug 6, 2025
7c7078d
Merge branch 'master' of github.com:Quantumapple/decaf
Quantumapple Aug 6, 2025
a2acbaf
final version for PAS
Quantumapple Aug 6, 2025
dab19f6
updated for bigger plot
Quantumapple Aug 13, 2025
ed94274
final version for paper
Quantumapple Sep 7, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
218 changes: 140 additions & 78 deletions analysis/notebooks/darkhiggs_postfit_all_year.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
"cells": [
{
"cell_type": "code",
"execution_count": 2,
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -18,25 +18,29 @@
},
{
"cell_type": "code",
"execution_count": 18,
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"config_dict = {\n",
" 'dyjetsMC': {'color': '#e5ccff', 'name': 'DY+jets'},\n",
" 'hbbMC': {'color': '#ccccff', 'name': r'$H\\rightarrow b\\bar{b}$'},\n",
" 'qcdMC': {'color': '#ffccff', 'name': 'QCD'},\n",
" 'vvMC': {'color': '#ffff99', 'name': 'VV'},\n",
" 'st': {'color': '#ff9999', 'name': 'Single t'},\n",
" 'stMC': {'color': '#ff9999', 'name': 'Single t'},\n",
" 'tt': {'color': '#ffcc66', 'name': r'$t\\bar{t}$'},\n",
" 'ttMC': {'color': '#ffcc66', 'name': r'$t\\bar{t}$'},\n",
" 'wjets': {'color': '#ccffcc', 'name': 'W+jets'},\n",
" 'wjetsMC': {'color': '#ccffcc', 'name': 'W+jets'},\n",
" 'zjets': {'color': '#99ffff', 'name': 'Z+jets'},\n",
" 'zjetsMC': {'color': '#99ffff', 'name': 'Z+jets'},\n",
" 'dyjetsMC': {'color': '#b9ac70', 'name': 'DY+jets'},\n",
" 'hbbMC': {'color': '#e76300', 'name': r'$H\\rightarrow b\\bar{b}$'},\n",
" 'qcdMC': {'color': '#a96b59', 'name': 'QCD'},\n",
" 'vvMC': {'color': '#832db6', 'name': 'VV'},\n",
" 'st': {'color': '#94a4a2', 'name': 'Single top'},\n",
" 'stMC': {'color': '#94a4a2', 'name': 'Single top'},\n",
" 'tt': {'color': '#bd1f01', 'name': r'$t\\bar{t}$'},\n",
" 'ttMC': {'color': '#bd1f01', 'name': r'$t\\bar{t}$'},\n",
" 'wjets': {'color': '#ffa90e', 'name': 'W+jets'},\n",
" 'wjetsMC': {'color': '#ffa90e', 'name': 'W+jets'},\n",
" 'zjets': {'color': '#3f90da', 'name': 'Z+jets'},\n",
" 'zjetsMC': {'color': '#3f90da', 'name': 'Z+jets'},\n",
"}\n",
"\n",
"\n",
"# color_cycle = ['#3f90da','#ffa90e','#bd1f01','#94a4a2','#832db6','#a96b59','#e76300','#b9ac70','#717581','#92dadd']\n",
"\n",
"\n",
"recoilDict = {\n",
" '0':'[250, 310) GeV',\n",
" '1':'[310, 370) GeV',\n",
Expand All @@ -59,10 +63,41 @@
},
{
"cell_type": "code",
"execution_count": 12,
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def merge_sig(input_file, sig_mass='Mz2000_mhs90_Mdm500'):\n",
" years = ['2016', '2017', '2018']\n",
" region = 'sr'\n",
" recoil_tags = [str(i) for i in range(5)]\n",
"\n",
" outputs = {}\n",
"\n",
" for itag in recoil_tags:\n",
" merged_sig = None\n",
"\n",
" for year in years:\n",
" identifier40to120 = f'{region}{year}passmass40to120recoil{itag}'\n",
" identifier120to300 = f'{region}{year}passmass120to300recoil{itag}'\n",
"\n",
" dir40to120 = input_file[identifier40to120+'_prefit']\n",
" dir120to300 = input_file[identifier120to300+'_prefit']\n",
"\n",
" sig40to120 = dir40to120[sig_mass].values\n",
" sig120to300 = dir120to300[sig_mass].values\n",
"\n",
" sig = np.concatenate((sig40to120, sig120to300), axis=None)\n",
"\n",
" if merged_sig is None:\n",
" merged_sig = np.array(sig) # Initialize as a NumPy array\n",
" else:\n",
" merged_sig += sig # Element-wise addition\n",
"\n",
" outputs[f'{itag}'] = merged_sig\n",
"\n",
" return outputs\n",
"\n",
"def merge_data(input_file):\n",
"\n",
" years = ['2016', '2017', '2018']\n",
Expand All @@ -81,8 +116,8 @@
" identifier40to120 = f'{region}{year}{category}mass40to120recoil{itag}'\n",
" identifier120to300 = f'{region}{year}{category}mass120to300recoil{itag}'\n",
"\n",
" prefit_dir40to120 = input_file[identifier40to120+'_prefit']\n",
" prefit_dir120to300 = input_file[identifier120to300+'_prefit']\n",
" prefit_dir40to120 = input_file[identifier40to120+'_postfit']\n",
" prefit_dir120to300 = input_file[identifier120to300+'_postfit']\n",
"\n",
" data40to120 = prefit_dir40to120[\"data_obs\"].values\n",
" data120to300 = prefit_dir120to300[\"data_obs\"].values\n",
Expand Down Expand Up @@ -224,25 +259,8 @@
"processes120to300 = ['hbbMC', 'dyjetsMC', 'qcdMC', 'vvMC', 'stMC', 'tt', 'ttMC', 'wjets', 'wjetsMC', 'zjets', 'zjetsMC']\n",
"\n",
"total_bkg_dict = merge_bkg(f, processes40to120, processes120to300)\n",
"indiv_mc, indiv_mc_variance = merge_individual_bkg(f, processes40to120, processes120to300)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"data_values_dict.keys()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"total_bkg_dict.keys()"
"indiv_mc, indiv_mc_variance = merge_individual_bkg(f, processes40to120, processes120to300)\n",
"sig_dict = merge_sig(f)"
]
},
{
Expand All @@ -251,7 +269,31 @@
"metadata": {},
"outputs": [],
"source": [
"indiv_mc.keys()"
"# import pandas as pd\n",
"# import sqlite3\n",
"\n",
"# ### Dump into csv format\n",
"# edges = np.array([40., 50., 60., 70., 80., 90., 100., 120., 150., 180., 240., 300.])\n",
"# centers = edges[1:]+edges[:-1]\n",
"\n",
"# # df = pd.DataFrame(total_bkg_dict)\n",
"# # df.insert(0, 'bin_center', centers/2)\n",
"\n",
"# # df = pd.DataFrame(indiv_mc)\n",
"# # df.insert(0, 'bin_center', centers/2)\n",
"\n",
"# # df = pd.DataFrame(indiv_mc_variance)\n",
"# # df.insert(0, 'bin_center', centers/2)\n",
"\n",
"# # df = pd.DataFrame(data_values_dict)\n",
"# # df.insert(0, 'bin_center', centers/2)\n",
"\n",
"# df = pd.DataFrame(sig_dict)\n",
"# df.insert(0, 'bin_center', centers/2)\n",
"\n",
"# outfile = 'darkhiggs_histogram_values.sqlite3'\n",
"# with sqlite3.connect(outfile) as output:\n",
"# df.to_sql('mz2000_mhd90_mchi500', output, if_exists='append', index=False)"
]
},
{
Expand All @@ -261,21 +303,23 @@
"outputs": [],
"source": [
"fig, axs = plt.subplots(2, 5, figsize=(40,12), sharex='col', sharey='row',\n",
" gridspec_kw=dict(height_ratios=[3, 1], hspace=0.15, wspace=0.))\n",
" gridspec_kw=dict(height_ratios=[3, 1], hspace=0.05, wspace=0.), layout=\"constrained\")\n",
"errps = {'hatch':'////', 'facecolor':'none', 'lw': 0, 'color': 'k', 'alpha': 0.3}\n",
"edges = [40., 50., 60., 70., 80., 90., 100., 120., 150., 180., 240., 300.]\n",
"\n",
"fig.supxlabel('AK15 Soft-Drop Mass [GeV]', fontsize=40)\n",
"fig.supxlabel('AK15 Soft-Drop Mass [GeV]', fontsize=45)\n",
"\n",
"for k in range(5):\n",
"\n",
" ax=axs[0][k]\n",
" rax=axs[1][k]\n",
"\n",
" ax.text(s=f'$U$ $\\in$ {recoilDict[str(k)]}',x=0.02,y=0.92, fontsize=35,transform=ax.transAxes)\n",
" # ax.text(s=f'$\\mathit{{U}}$ $\\in$ {recoilDict[str(k)]}',x=0.225,y=0.65, fontsize=30,transform=ax.transAxes)\n",
" ax.text(s=f'$\\mathit{{U}}$ $\\in$ {recoilDict[str(k)]}',x=0.225,y=0.65, fontsize=32,transform=ax.transAxes)\n",
"\n",
" ax._get_lines.prop_cycler = ax._get_patches_for_fill.prop_cycler\n",
" ax.set_yscale('log')\n",
" ax.set_ylim(1e-2, 8e+5)\n",
" ax.set_ylim(1e-2, 9e+7)\n",
"\n",
" colors=[]\n",
" mc_process_data = []\n",
Expand All @@ -293,7 +337,7 @@
"\n",
" ### Try to draw stack plots\n",
" hep.histplot(mc_process_data, edges, ax=ax, stack=True,\n",
" histtype='fill', edgecolor = 'k', linewidth=3, label=mc_names)\n",
" histtype='fill', edgecolor = 'k', linewidth=1, label=mc_names)\n",
"\n",
" ### Draw stat uncs.\n",
" y1 = total_bkg_dict[f'pass_{k}_postfit'] - np.sqrt(total_bkg_dict[f'pass_{k}_postfit_variance'])\n",
Expand All @@ -308,12 +352,19 @@
" **errps, label='Stat. Unc.'\n",
" )\n",
"\n",
" hep.histplot(total_bkg_dict[f'pass_{k}_postfit'], edges, ax=ax, label=[\"Postift\"], color='b', linewidth=6)\n",
" hep.histplot(total_bkg_dict[f'pass_{k}_prefit'], edges, ax=ax, label=[\"Prefit\"], color='r', linestyle='dashed', linewidth=6)\n",
"\n",
" ### Call data ###\n",
" hep.histplot(data_values_dict[f'pass_{k}'], edges, ax=ax, histtype='errorbar', label=\"Data\", color='k', markersize=20)\n",
"\n",
" ### Call signal ###\n",
" signal_label = (\n",
" r\"$m_{Z'} = 2000\\, GeV, $\" +\n",
" r\"$m_{H_{D}} = 90\\, GeV,$\" +\n",
" r\"$m_{\\chi} = 500\\, GeV$\"\n",
" )\n",
"\n",
" # Mz2000_mhs90_Mdm500\n",
" hep.histplot(sig_dict[f'{k}'], edges, ax=ax, label=signal_label, color='#11f0e6', linewidth=6, linestyle='dashed', edges=False)\n",
"\n",
" #### Draw ratio ####\n",
" hep.histplot(data_values_dict[f'pass_{k}']/total_bkg_dict[f'pass_{k}_prefit'], edges, yerr=np.sqrt(data_values_dict[f'pass_{k}'])/total_bkg_dict[f'pass_{k}_prefit'],\n",
" ax=rax, histtype='errorbar', color='r', capsize=10, label=\"Prefit\", markersize=20)\n",
Expand All @@ -336,59 +387,70 @@
" rax.axhline(1, ls='--', color='k')\n",
" ymax=abs(rax.get_ylim()[1])*1.1\n",
" ymin=1.-(ymax-1.)\n",
" rax.set_ylim(max(ymin,0),min(ymax,2))\n",
" rax.set_ylim(0.3, 1.55)\n",
" rax.set_xlim(40, 300)\n",
"\n",
" if k == 0:\n",
" hep.cms.text(ax=ax, loc=0, text='Preliminary',fontsize=45)\n",
" ax.set_ylabel('Events/GeV', fontsize=40)\n",
" hep.cms.text(ax=ax, loc=0, fontsize=45)\n",
" ax.set_ylabel('Events', fontsize=45)\n",
" ax.tick_params(axis='y', labelsize=30)\n",
" ax.set_xticks([50,100,150,200,250])\n",
" rax.set_ylabel('Obs/Exp', fontsize=40)\n",
" rax.set_ylabel('Data/Pred.', fontsize=36)\n",
" rax.tick_params(axis='y', labelsize=30)\n",
" rax.tick_params(axis='x', labelsize=30)\n",
" handles, labels = ax.get_legend_handles_labels()\n",
" rax_handles, rax_labels = rax.get_legend_handles_labels()\n",
"\n",
" ax.text(s=f'$\\mathit{{U}}$: recoil', x=0.2, y=0.85, fontsize=38, transform=ax.transAxes)\n",
"\n",
" else:\n",
" ax.set_xticks([50,100,150,200,250])\n",
" rax.tick_params(axis='x', labelsize=30)\n",
"\n",
"\n",
" if k == 2:\n",
" order = [11, 9, 10, 8] # Data, Postfit, Prefit, Stat. Unc\n",
" ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n",
" loc='center left', fontsize=30, ncol=2, bbox_to_anchor=(-0.03, 0.78))\n",
" rax.legend([rax_handles[1]], [rax_labels[1]], loc='lower right', fontsize=30, bbox_to_anchor=(1.03, -0.1))\n",
"\n",
" if k == 3:\n",
" order = [7, 6, 5, 4] # Z+jets, W+jets, tt, st\n",
" ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n",
" loc='center left', fontsize=30, ncol=2, bbox_to_anchor=(-0.01, 0.78))\n",
" rax.legend([rax_handles[2]], [rax_labels[2]], loc='lower right', fontsize=30, bbox_to_anchor=(1.03, -0.1))\n",
"\n",
"\n",
" if k == 4:\n",
"\n",
" ### Lumi text\n",
" hep.cms.lumitext(ax=ax, text=r\"138 fb$^{-1}$, 2016-2018 (13 TeV)\", fontsize=45)\n",
" hep.cms.lumitext(ax=ax, text=r\"138 fb$^{-1}$ (13 TeV)\", fontsize=45)\n",
"\n",
"orders = [10, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7]\n",
"### Default version\n",
"# legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=6, loc='center', fontsize=30,\n",
"# bbox_to_anchor=(0.5, 0.81), frameon=True, facecolor='white')\n",
"\n",
" #### Legend\n",
" order = [3, 2, 1, 0] ## VV, QCD, DY+jets, Htobb\n",
" ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n",
" loc='center left', fontsize=30, ncol=2, bbox_to_anchor=(-0.01, 0.78))\n",
" rax.legend([rax_handles[0]], [rax_labels[0]], loc='lower right', fontsize=30, bbox_to_anchor=(1.03, -0.1))\n",
"### If layout=\"constrained\" is used\n",
"legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=6, loc='center', fontsize=33,\n",
" bbox_to_anchor=(0.55, 0.85), frameon=True, facecolor='white')\n",
"legend.get_frame().set_alpha(1)\n",
"legend.get_frame().set_linewidth(0)\n",
"\n",
" # order = [11, 9, 10, 7, 6, 5, 4, 3, 2, 1, 0, 8]\n",
" # ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n",
" # loc='center left', fontsize=20, ncol=3, bbox_to_anchor=(-0.02, 0.72))\n",
"rax_orders = [1, 2, 0]\n",
"rax_handles, rax_labels = rax.get_legend_handles_labels()\n",
"### Default version\n",
"# rax_leg = fig.legend([rax_handles[idx] for idx in rax_orders], [rax_labels[idx] for idx in rax_orders], ncol=3, loc='lower center',\n",
"# bbox_to_anchor=(0.5, 0.1), frameon=True, facecolor='white', fontsize=30)\n",
"\n",
"### If layout=\"constrained\" is used\n",
"rax_leg = fig.legend([rax_handles[idx] for idx in rax_orders], [rax_labels[idx] for idx in rax_orders], ncol=3, loc='lower center',\n",
" bbox_to_anchor=(0.5, 0.093), frameon=True, facecolor='white', fontsize=33)\n",
"\n",
" # order = [1, 2, 0]\n",
" # rax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n",
" # loc='center left', fontsize=20, ncol=3, bbox_to_anchor=(0.0, 0.12))"
"rax_leg.get_frame().set_alpha(1)\n",
"rax_leg.get_frame().set_linewidth(0)\n",
"\n",
"fig.savefig('/home/jongho/Figure_002_v3.pdf')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
Expand All @@ -399,7 +461,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "packages",
"display_name": "packages (3.10.12)",
"language": "python",
"name": "python3"
},
Expand Down
Loading