From 6ce14f2252e99d8abe5fa7e3fb46d84acc10420f Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Thu, 10 Apr 2025 21:22:20 -0500 Subject: [PATCH 01/10] add signal --- .../darkhiggs_postfit_all_year.ipynb | 98 ++++++++++--------- 1 file changed, 54 insertions(+), 44 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index 648285643..af20b0fe3 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -18,7 +18,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -59,10 +59,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+'_postfit']\n", + " dir120to300 = input_file[identifier120to300+'_postfit']\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", @@ -81,8 +112,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", @@ -224,34 +255,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()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "indiv_mc.keys()" + "indiv_mc, indiv_mc_variance = merge_individual_bkg(f, processes40to120, processes120to300)\n", + "sig_dict = merge_sig(f)" ] }, { @@ -314,6 +319,15 @@ " ### 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, $\" +\n", + " r\"$M_{h_{D}} = 90,$\" + \"\\n\" +\n", + " r\"$M_{\\chi} = 500\\ \\mathrm{GeV}$\"\n", + " )\n", + " # Mz2000_mhs90_Mdm500\n", + " hep.histplot(sig_dict[f'{k}'], edges, ax=ax, label=signal_label, color='#8e4b00', linewidth=6, linestyle='dashdot')\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", @@ -354,9 +368,14 @@ " ax.set_xticks([50,100,150,200,250])\n", " rax.tick_params(axis='x', labelsize=30)\n", "\n", + " if k == 1:\n", + " order = [11] # Signal\n", + " ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n", + " loc='upper center', fontsize=30, bbox_to_anchor=(0.5, 0.94))\n", + "\n", "\n", " if k == 2:\n", - " order = [11, 9, 10, 8] # Data, Postfit, Prefit, Stat. Unc\n", + " order = [12, 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", @@ -377,16 +396,7 @@ " 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", - "\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", - "\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.legend([rax_handles[0]], [rax_labels[0]], loc='lower right', fontsize=30, bbox_to_anchor=(1.03, -0.1))" ] }, { From 696c1ef20ec761d883b1e07a068a06a389b2d60b Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Sat, 12 Apr 2025 18:12:33 -0500 Subject: [PATCH 02/10] update style --- analysis/notebooks/darkhiggs_postfit_all_year.ipynb | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index af20b0fe3..a5b85d167 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -322,7 +322,7 @@ " ### Call signal ###\n", " signal_label = (\n", " r\"$M_{Z'} = 2000, $\" +\n", - " r\"$M_{h_{D}} = 90,$\" + \"\\n\" +\n", + " r\"$M_{H_{D}} = 90,$\" + \"\\n\" +\n", " r\"$M_{\\chi} = 500\\ \\mathrm{GeV}$\"\n", " )\n", " # Mz2000_mhs90_Mdm500\n", @@ -355,7 +355,7 @@ "\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", + " ax.set_ylabel('Events', fontsize=40)\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", From 5faf51818fda2db95b5eb2dbdb339cadb2b9be82 Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Tue, 15 Apr 2025 07:21:13 -0500 Subject: [PATCH 03/10] switch to use prefit for signal histogram --- analysis/notebooks/darkhiggs_postfit_all_year.ipynb | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index a5b85d167..0aca222d9 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -77,8 +77,8 @@ " identifier40to120 = f'{region}{year}passmass40to120recoil{itag}'\n", " identifier120to300 = f'{region}{year}passmass120to300recoil{itag}'\n", "\n", - " dir40to120 = input_file[identifier40to120+'_postfit']\n", - " dir120to300 = input_file[identifier120to300+'_postfit']\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", From bb8059f1045379573c96cce3fbb91a0ef008305d Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Fri, 18 Apr 2025 10:35:26 -0500 Subject: [PATCH 04/10] update --- .../darkhiggs_postfit_all_year.ipynb | 115 ++++++++++-------- 1 file changed, 67 insertions(+), 48 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index 0aca222d9..52d5dcebe 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -23,20 +23,24 @@ "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", @@ -259,6 +263,39 @@ "sig_dict = merge_sig(f)" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# 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)" + ] + }, { "cell_type": "code", "execution_count": null, @@ -277,10 +314,10 @@ " 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'$U$ $\\in$ {recoilDict[str(k)]}',x=0.225,y=0.613, fontsize=30,transform=ax.transAxes)\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", @@ -313,20 +350,18 @@ " **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, $\" +\n", - " r\"$M_{H_{D}} = 90,$\" + \"\\n\" +\n", - " r\"$M_{\\chi} = 500\\ \\mathrm{GeV}$\"\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='#8e4b00', linewidth=6, linestyle='dashdot')\n", + " hep.histplot(sig_dict[f'{k}'], edges, ax=ax, label=signal_label, color='#11f0e6', linewidth=6, linestyle='dashed')\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", @@ -358,45 +393,29 @@ " ax.set_ylabel('Events', fontsize=40)\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", " else:\n", " ax.set_xticks([50,100,150,200,250])\n", " rax.tick_params(axis='x', labelsize=30)\n", "\n", - " if k == 1:\n", - " order = [11] # Signal\n", - " ax.legend([handles[idx] for idx in order], [labels[idx] for idx in order],\n", - " loc='upper center', fontsize=30, bbox_to_anchor=(0.5, 0.94))\n", - "\n", - "\n", - " if k == 2:\n", - " order = [12, 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", "\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))" + "orders = [10, 8, 9, 0, 1, 2, 3, 4, 5, 6, 7]\n", + "legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=5, loc='center', fontsize=27,\n", + " bbox_to_anchor=(0.63, 0.79), frameon=True, facecolor='white')\n", + "legend.get_frame().set_alpha(1)\n", + "\n", + "rax_orders = [1, 2, 0]\n", + "rax_handles, rax_labels = rax.get_legend_handles_labels()\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=27)\n", + "rax_leg.get_frame().set_alpha(1)" ] }, { From a55681086571446dc42afc72ea97138c7ca80fd8 Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Sat, 19 Apr 2025 21:40:47 -0500 Subject: [PATCH 05/10] style update --- analysis/notebooks/darkhiggs_postfit_all_year.ipynb | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index 52d5dcebe..20c4f0db4 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -314,7 +314,7 @@ " ax=axs[0][k]\n", " rax=axs[1][k]\n", "\n", - " ax.text(s=f'$U$ $\\in$ {recoilDict[str(k)]}',x=0.225,y=0.613, fontsize=30,transform=ax.transAxes)\n", + " ax.text(s=f'$U$ $\\in$ {recoilDict[str(k)]}',x=0.225,y=0.65, fontsize=30,transform=ax.transAxes)\n", " ax._get_lines.prop_cycler = ax._get_patches_for_fill.prop_cycler\n", " ax.set_yscale('log')\n", " ax.set_ylim(1e-2, 9e+7)\n", @@ -404,11 +404,11 @@ "\n", " if k == 4:\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", - "legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=5, loc='center', fontsize=27,\n", - " bbox_to_anchor=(0.63, 0.79), frameon=True, facecolor='white')\n", + "legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=6, loc='center', fontsize=27,\n", + " bbox_to_anchor=(0.5, 0.81), frameon=True, facecolor='white')\n", "legend.get_frame().set_alpha(1)\n", "\n", "rax_orders = [1, 2, 0]\n", From 4516acd231ae37fb5a46f60266f7956e494d6955 Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Tue, 5 Aug 2025 20:16:51 -0500 Subject: [PATCH 06/10] nb for 2D limit plotting --- analysis/notebooks/limit_plot.ipynb | 252 ++++++++++++++++++++++++++++ 1 file changed, 252 insertions(+) create mode 100644 analysis/notebooks/limit_plot.ipynb diff --git a/analysis/notebooks/limit_plot.ipynb b/analysis/notebooks/limit_plot.ipynb new file mode 100644 index 000000000..cd5a0c5ee --- /dev/null +++ b/analysis/notebooks/limit_plot.ipynb @@ -0,0 +1,252 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "id": "54b5dc95", + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import re\n", + "from pathlib import Path\n", + "\n", + "value = 'mhs50'\n", + "\n", + "def parse_limit_logs(value: str):\n", + " # Use pathlib for modern, object-oriented path handling\n", + " log_folder = Path('./Input_for_Limits_Plots')\n", + "\n", + " # Use glob to find all matching files directly, which is more efficient\n", + " # than listing the whole directory and filtering manually.\n", + " log_pattern = f'darkhiggs_*{value}*_obsAsymptoticLimits*'\n", + " files = list(log_folder.glob(log_pattern))\n", + "\n", + " # Compile regular expressions once for efficiency\n", + " filename_regex = re.compile(r'Mz(\\d+\\.?\\d*)_mhs([\\w\\d.-]+)_Mdm(\\d+\\.?\\d*)')\n", + "\n", + " # A mapping from text in the file to the desired column name\n", + " # This avoids a long series of 'if' statements.\n", + " limit_parsers = {\n", + " 'Observed Limit:': 'obs',\n", + " 'Expected 50.0%:': 'exp',\n", + " 'Expected 84.0%:': 'p1s',\n", + " 'Expected 16.0%:': 'm1s',\n", + " 'Expected 97.5%:': 'p2s',\n", + " 'Expected 2.5%:': 'm2s',\n", + " }\n", + "\n", + " # Collect data in a list of dictionaries, a standard way to build a DataFrame\n", + " all_results = []\n", + "\n", + " for ifile in files:\n", + " filename = ifile.name\n", + "\n", + " # Use regex to robustly extract parameters from the filename\n", + " match = filename_regex.search(filename)\n", + "\n", + " # Start a dictionary to hold data for this specific file\n", + " # This is much cleaner than managing many separate lists.\n", + " point_data = {\n", + " 'Mz': float(match.group(1)),\n", + " 'mhs': match.group(2), # mhs seems to be treated as a string/category\n", + " 'Mdm': float(match.group(3)),\n", + " }\n", + "\n", + " found_limits = {}\n", + "\n", + " with ifile.open('r') as f:\n", + " for line in f.readlines():\n", + " for key, col_name in limit_parsers.items():\n", + " if key in line:\n", + " found_limits[col_name] = float(line.split('<')[1].strip())\n", + " break\n", + "\n", + " if len(found_limits) == len(limit_parsers):\n", + " point_data.update(found_limits)\n", + " all_results.append(point_data)\n", + "\n", + " # Create the DataFrame once from the list of dictionaries\n", + " df = pd.DataFrame(all_results)\n", + " # Sort by multiple columns in a single, efficient operation\n", + " # reset_index is good practice after sorting to get a clean 0-based index.\n", + " df = df.sort_values(by=['Mz', 'Mdm']).reset_index(drop=True)\n", + "\n", + " return df" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "af9fc88b", + "metadata": {}, + "outputs": [], + "source": [ + "mhs50_df = parse_limit_logs('mhs50')\n", + "mhs70_df = parse_limit_logs('mhs70')\n", + "mhs90_df = parse_limit_logs('mhs90')\n", + "mhs110_df = parse_limit_logs('mhs110')\n", + "mhs130_df = parse_limit_logs('mhs130')\n", + "mhs150_df = parse_limit_logs('mhs150')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a3a5cd8d", + "metadata": {}, + "outputs": [], + "source": [ + "import sys\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "import matplotlib.colors as mcolors\n", + "import matplotlib.ticker as mticker\n", + "from matplotlib.lines import Line2D\n", + "import mplhep as hep\n", + "hep.style.use('CMS')\n", + "\n", + "sys.path.append('../')\n", + "\n", + "# --- Using your specific library import ---\n", + "from libs.limitlib import interpolate_rbf\n", + "\n", + "cmap = mcolors.LinearSegmentedColormap.from_list(\"n\", list(reversed([\n", + " # '#fff5f0',\n", + " '#fee0d2',\n", + " '#ffffff',\n", + " '#fcbba1',\n", + " '#fc9272',\n", + " '#fb6a4a',\n", + " '#ef3b2c',\n", + " '#cb181d',\n", + " '#a50f15',\n", + " '#67000d',\n", + " # '#000000',\n", + " ])))\n", + "\n", + "def plot2d(df):\n", + " fig, ax = plt.subplots(figsize=(14,10))\n", + " hep.cms.label(data=True, year='2016-2018', lumi=138)#, loc=1)#, label=\"Preliminary\")\n", + "\n", + " x = df['Mz']\n", + " y = df['Mdm']\n", + " obs = df['obs']\n", + "\n", + " contours_filled = np.log10(np.logspace(-1,1,7))\n", + " contours_line = [0]\n", + "\n", + " def get_x_y_z(x,y,z):\n", + " ix, iy, iz = interpolate_rbf(x,y,z,maxval=5000)\n", + " iz [iy>ix] = 1e9 #* np.exp(-(iy/ix))\n", + " if True:\n", + " iz = np.log10(iz)\n", + " iz[iz Date: Tue, 5 Aug 2025 20:18:31 -0500 Subject: [PATCH 07/10] update --- .../darkhiggs_postfit_all_year.ipynb | 40 ++++++++++++++----- 1 file changed, 29 insertions(+), 11 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index 20c4f0db4..a44e45308 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -314,7 +314,7 @@ " ax=axs[0][k]\n", " rax=axs[1][k]\n", "\n", - " ax.text(s=f'$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=30,transform=ax.transAxes)\n", " ax._get_lines.prop_cycler = ax._get_patches_for_fill.prop_cycler\n", " ax.set_yscale('log')\n", " ax.set_ylim(1e-2, 9e+7)\n", @@ -335,7 +335,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", @@ -355,13 +355,13 @@ "\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", + " 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')\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", @@ -385,7 +385,7 @@ " 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", @@ -407,17 +407,35 @@ " 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", - "legend = fig.legend([handles[idx] for idx in orders], [labels[idx] for idx in orders], ncol=6, loc='center', fontsize=27,\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", "legend.get_frame().set_alpha(1)\n", + "legend.get_frame().set_linewidth(0)\n", "\n", "rax_orders = [1, 2, 0]\n", "rax_handles, rax_labels = rax.get_legend_handles_labels()\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=27)\n", - "rax_leg.get_frame().set_alpha(1)" + " bbox_to_anchor=(0.5, 0.1), frameon=True, facecolor='white', fontsize=30)\n", + "rax_leg.get_frame().set_alpha(1)\n", + "rax_leg.get_frame().set_linewidth(0)\n", + "\n", + "fig.savefig('./sr.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, @@ -428,7 +446,7 @@ ], "metadata": { "kernelspec": { - "display_name": "packages", + "display_name": "packages (3.10.12)", "language": "python", "name": "python3" }, From a2acbafe4b20e41503d385af2aab649ef4dcec0e Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Wed, 6 Aug 2025 08:44:18 -0500 Subject: [PATCH 08/10] final version for PAS --- analysis/notebooks/limit_plot.ipynb | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/analysis/notebooks/limit_plot.ipynb b/analysis/notebooks/limit_plot.ipynb index cd5a0c5ee..5d8ee2013 100644 --- a/analysis/notebooks/limit_plot.ipynb +++ b/analysis/notebooks/limit_plot.ipynb @@ -127,7 +127,7 @@ "\n", "def plot2d(df):\n", " fig, ax = plt.subplots(figsize=(14,10))\n", - " hep.cms.label(data=True, year='2016-2018', lumi=138)#, loc=1)#, label=\"Preliminary\")\n", + " hep.cms.label(data=True, lumi=138, label=\"Preliminary\")#, loc=1)#)\n", "\n", " x = df['Mz']\n", " y = df['Mdm']\n", From dab19f62f082f264d7a3e2c965a7e3e4c6cec15b Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Wed, 13 Aug 2025 08:54:37 -0500 Subject: [PATCH 09/10] updated for bigger plot --- .../darkhiggs_postfit_all_year.ipynb | 29 ++++++++++++++----- 1 file changed, 21 insertions(+), 8 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index a44e45308..bef17555c 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -303,18 +303,20 @@ "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'$\\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=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, 9e+7)\n", @@ -390,7 +392,7 @@ "\n", " if k == 0:\n", " hep.cms.text(ax=ax, loc=0, text='Preliminary',fontsize=45)\n", - " ax.set_ylabel('Events', fontsize=40)\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('Data/Pred.', fontsize=36)\n", @@ -407,19 +409,30 @@ " 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", - "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", + "### 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", + "### 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.5, 0.85), frameon=True, facecolor='white')\n", "legend.get_frame().set_alpha(1)\n", "legend.get_frame().set_linewidth(0)\n", "\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.1), frameon=True, facecolor='white', fontsize=30)\n", + " bbox_to_anchor=(0.5, 0.093), frameon=True, facecolor='white', fontsize=33)\n", + "\n", "rax_leg.get_frame().set_alpha(1)\n", "rax_leg.get_frame().set_linewidth(0)\n", "\n", - "fig.savefig('./sr.pdf')" + "# fig.savefig('/home/jongho/Figure_002_v2.pdf')" ] }, { From ed9427429469935bd010dd889b5ec3db3c073936 Mon Sep 17 00:00:00 2001 From: Quantumapple Date: Sun, 7 Sep 2025 13:18:17 -0500 Subject: [PATCH 10/10] final version for paper --- analysis/notebooks/darkhiggs_postfit_all_year.ipynb | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb index bef17555c..952249759 100644 --- a/analysis/notebooks/darkhiggs_postfit_all_year.ipynb +++ b/analysis/notebooks/darkhiggs_postfit_all_year.ipynb @@ -391,7 +391,7 @@ " rax.set_xlim(40, 300)\n", "\n", " if k == 0:\n", - " hep.cms.text(ax=ax, loc=0, text='Preliminary',fontsize=45)\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", @@ -400,6 +400,8 @@ " rax.tick_params(axis='x', labelsize=30)\n", " handles, labels = ax.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", @@ -415,7 +417,7 @@ "\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.5, 0.85), frameon=True, facecolor='white')\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", @@ -432,7 +434,7 @@ "rax_leg.get_frame().set_alpha(1)\n", "rax_leg.get_frame().set_linewidth(0)\n", "\n", - "# fig.savefig('/home/jongho/Figure_002_v2.pdf')" + "fig.savefig('/home/jongho/Figure_002_v3.pdf')" ] }, {