Skip to content

Commit

Permalink
divide bin widths
Browse files Browse the repository at this point in the history
  • Loading branch information
@rkansal47
rkansal47 committed Sep 1, 2023
1 parent d059487 commit 653681e
Show file tree
Hide file tree
Showing 9 changed files with 233 additions and 292 deletions.
45 changes: 35 additions & 10 deletions src/HHbbVV/combine/binder/BiasTest.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,7 @@
"source": [
"MAIN_DIR = \"../../../../\"\n",
"\n",
"plot_dir = f\"{MAIN_DIR}/plots/BiasTest/23Jul18Res\"\n",
"plot_dir = f\"{MAIN_DIR}/plots/BiasTest/23Sep1\"\n",
"_ = os.system(f\"mkdir -p {plot_dir}\")"
]
},
Expand Down Expand Up @@ -57,18 +57,34 @@
" )\n",
"\n",
" r = np.array(file.limit)[::4]\n",
" r_negerr = r - np.array(file.limit)[1::4]\n",
" r_poserr = np.array(file.limit)[2::4] - r\n",
" neg_lim = np.array(file.limit)[1::4]\n",
" pos_lim = np.array(file.limit)[2::4]\n",
" r_negerr = r - neg_lim\n",
" r_poserr = pos_lim - r\n",
" reldiff = r - bias\n",
" reldiff[reldiff < 0] = (reldiff / r_poserr)[reldiff < 0]\n",
" reldiff[reldiff > 0] = (reldiff / r_negerr)[reldiff > 0]\n",
"\n",
" r_dict[sample][bias] = {\n",
" \"r\": r,\n",
" \"reldiff\": reldiff,\n",
" \"neg_lim\": neg_lim,\n",
" \"pos_lim\": pos_lim,\n",
" }"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"for sample in samples:\n",
" print(sample)\n",
" for i, bias in enumerate(biases):\n",
" print(f\"For bias = {bias}, {np.sum(r_dict[sample][bias]['pos_lim'] == 20)} failed\")"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand All @@ -77,22 +93,24 @@
"source": [
"from scipy import stats\n",
"\n",
"xrange = 4\n",
"xrange = 3\n",
"x = np.linspace(-xrange, xrange, 101)\n",
"\n",
"for sample in samples:\n",
" fig, ax = plt.subplots(len(biases), 1, figsize=(12, 40))\n",
"\n",
" for i, bias in enumerate(biases):\n",
" r = r_dict[sample][bias][\"r\"]\n",
" reldiff = r_dict[sample][bias][\"reldiff\"]\n",
" reldiff = reldiff[(reldiff > -4) * (reldiff < 4)]\n",
" fit_fail = r_dict[sample][bias][\"pos_lim\"] == 20\n",
" # fit_fail = r_dict[sample][bias]['pos_lim'] == 2000000\n",
" r = r_dict[sample][bias][\"r\"][~fit_fail]\n",
" reldiff = r_dict[sample][bias][\"reldiff\"][~fit_fail]\n",
" reldiff = reldiff[(reldiff > -xrange) * (reldiff < xrange)]\n",
"\n",
" mu, sigma = np.mean(reldiff), np.std(reldiff)\n",
" ax[i].hist(reldiff, np.linspace(-xrange, xrange, 21), histtype=\"step\")\n",
" ax[i].hist(reldiff, np.linspace(-xrange, xrange, 20), histtype=\"step\")\n",
" ax[i].plot(\n",
" x,\n",
" stats.norm.pdf(x, loc=mu, scale=sigma) * len(r) * xrange / 10,\n",
" stats.norm.pdf(x, loc=mu, scale=sigma) * len(r) * xrange / 10, # scale by bin width\n",
" label=rf\"$\\mu = {mu:.2f}, \\sigma = {sigma:.2f}$\",\n",
" )\n",
" ax[i].set_xlabel(rf\"$\\frac{{\\hat{{r}} - {bias}}}{{\\Delta \\hat r}}$\")\n",
Expand All @@ -111,6 +129,13 @@
" plt.savefig(f\"{plot_dir}/{sample}.pdf\", bbox_inches=\"tight\")\n",
" plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
Expand All @@ -129,7 +154,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.15"
"version": "3.10.4"
},
"orig_nbformat": 4
},
Expand Down
9 changes: 5 additions & 4 deletions src/HHbbVV/combine/binder/PlotScan.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,11 +45,12 @@
"scan_thww_wps = [0.4, 0.6, 0.8, 0.9, 0.94, 0.96, 0.98]\n",
"scan_bdt_wps = [0.6, 0.9, 0.96, 0.99, 0.997, 0.998, 0.999]\n",
"\n",
"plot_dir = \"../../../../plots/Scans/23May13/\"\n",
"plot_dir = \"../../../../plots/Scans/23Aug30NonResIntLumi/\"\n",
"_ = os.system(f\"mkdir -p {plot_dir}\")\n",
"\n",
"# cards_dir = \"/eos/uscms/store/user/rkansal/bbVV/cards/23May13NonresScan/\"\n",
"cards_dir = \"/uscms/home/rkansal/hhcombine/cards/23May14NonresScan/\""
"cards_dir = \"/uscms/home/rkansal/hhcombine/cards/23May14NonresScan/\"\n",
"# cards_dir = \"/uscms/home/rkansal/hhcombine/eoscards/23May1Scan/\" # resonant"
]
},
{
Expand Down Expand Up @@ -178,7 +179,7 @@
" mx, my = mxmy(sample)\n",
" plt.title(f\"$(M_X, M_Y) = ({mx}, {my})$ GeV\", y=1.12)\n",
"\n",
" hep.cms.label(\"Work in Progress\", data=True, lumi=\"138.0\")\n",
" hep.cms.label(\"Work in Progress\", data=True, lumi=\"138\")\n",
"\n",
" # plt.title(\"Signifance Scan\")\n",
" plt.savefig(f\"{plot_dir}/{sample}_limscan.pdf\", bbox_inches=\"tight\")"
Expand Down Expand Up @@ -301,7 +302,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.4"
"version": "3.9.15"
},
"orig_nbformat": 4,
"vscode": {
Expand Down
4 changes: 3 additions & 1 deletion src/HHbbVV/combine/run_blinded.sh
View file
Original file line number Diff line number Diff line change
Expand Up @@ -229,7 +229,9 @@ if [ $workspace = 1 ]; then
combineCards.py $ccargs > $ws.txt

echo "Running text2workspace"
text2workspace.py -D $dataset $ws.txt --channel-masks -o $wsm.root 2>&1 | tee $outsdir/text2workspace.txt
# text2workspace.py -D $dataset $ws.txt --channel-masks -o $wsm.root 2>&1 | tee $outsdir/text2workspace.txt
# new version got rid of -D arg??
text2workspace.py $ws.txt --channel-masks -o $wsm.root 2>&1 | tee $outsdir/text2workspace.txt
else
if [ ! -f "$wsm.root" ]; then
echo "Workspace doesn't exist! Use the -w|--workspace option to make workspace first"
Expand Down
41 changes: 36 additions & 5 deletions src/HHbbVV/postprocessing/CreateDatacard.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@

# from utils import add_bool_arg
from hh_vars import LUMI, res_sig_keys, qcd_key, data_key, years, jecs, jmsr, res_mps

import utils

import argparse

Expand Down Expand Up @@ -108,6 +108,14 @@ def mxmy(sample):
add_bool_arg(parser, "sig-separate", "separate templates for signals and bgs", default=False)
add_bool_arg(parser, "do-jshifts", "Do JEC/JMC corrections.", default=True)

add_bool_arg(parser, "only-sm", "Only add SM HH samples for (for debugging nonres)", default=False)
add_bool_arg(
parser,
"combine-lasttwo",
"Check combining last two bins in nonres for HWW review",
default=False,
)

parser.add_argument("--cards-dir", default="cards", type=str, help="output card directory")

parser.add_argument("--mcstats-threshold", default=100, type=float, help="mcstats threshold n_eff")
Expand All @@ -128,8 +136,8 @@ def mxmy(sample):
default=None,
nargs="*",
type=int,
help="order of polynomial for TF in [dim 1, dim 2] = [mH(bb), -] for nonresonant or [mY, mX] for resonant"
"",
help="order of polynomial for TF in [dim 1, dim 2] = [mH(bb), -] for nonresonant or [mY, mX] for resonant."
"Default is 0 for nonresonant and (1, 2) for resonant.",
)

parser.add_argument("--model-name", default=None, type=str, help="output model name")
Expand Down Expand Up @@ -178,6 +186,10 @@ def mxmy(sample):
"qqHH_CV_1_C2V_1_kl_0_HHbbVV",
"qqHH_CV_0p5_C2V_1_kl_1_HHbbVV",
]

if args.only_sm:
nonres_sig_keys_ggf, nonres_sig_keys_vbf = ["HHbbVV"], []

nonres_sig_keys = nonres_sig_keys_ggf + nonres_sig_keys_vbf
sig_keys = []
hist_names = {} # names of hist files for the samples
Expand Down Expand Up @@ -348,7 +360,7 @@ def main(args):

# templates per region per year, templates per region summed across years
templates_dict, templates_summed = get_templates(
args.templates_dir, years, args.sig_separate, args.scale_templates
args.templates_dir, years, args.sig_separate, args.scale_templates, args.combine_lasttwo
)

# TODO: check if / how to include signal trig eff uncs. (rn only using bg uncs.)
Expand Down Expand Up @@ -510,7 +522,23 @@ def combine_templates(
return ctemplates


def get_templates(templates_dir: str, years: List[str], sig_separate: bool, scale: float = None):
def _combine_last_two_bins(templates_dict):
for year in years:
for region in templates_dict[year]:
templates_dict[year][region] = utils.rebin_hist(
templates_dict[year][region],
"bbFatJetParticleNetMass",
list(range(50, 240, 10)) + [250],
)


def get_templates(
templates_dir: str,
years: List[str],
sig_separate: bool,
scale: float = None,
combine_lasttwo: bool = False,
):
"""Loads templates, combines bg and sig templates if separate, sums across all years"""
templates_dict: Dict[str, Dict[str, Hist]] = {}

Expand Down Expand Up @@ -538,6 +566,9 @@ def get_templates(templates_dir: str, years: List[str], sig_separate: bool, scal
for key in templates_dict[year]:
templates_dict[year][key] = templates_dict[year][key] * scale

if combine_lasttwo:
_combine_last_two_bins(templates_dict)

templates_summed: Dict[str, Hist] = sum_templates(templates_dict) # sum across years
return templates_dict, templates_summed

Expand Down
8 changes: 4 additions & 4 deletions src/HHbbVV/postprocessing/PlotFits.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@
"source": [
"MAIN_DIR = \"../../../\"\n",
"\n",
"plot_dir = \"../../../plots/PostFit/23Jul6\"\n",
"plot_dir = \"../../../plots/PostFit/23Sep1Nonres\"\n",
"_ = os.system(f\"mkdir -p {plot_dir}\")"
]
},
Expand All @@ -50,10 +50,10 @@
"metadata": {},
"outputs": [],
"source": [
"cards_dir = \"f_tests/23May14HP\"\n",
"cards_dir = \"23May16RemFeats/rem_feats_3/\"\n",
"file = uproot.open(\n",
" # f\"/uscms/home/rkansal/hhcombine/cards/{cards_dir}/FitShapes.root\"\n",
" f\"/uscms/home/rkansal/eos/bbVV/cards/{cards_dir}/FitShapes.root\"\n",
" f\"/uscms/home/rkansal/hhcombine/cards/{cards_dir}/FitShapes.root\"\n",
" # f\"/uscms/home/rkansal/eos/bbVV/cards/{cards_dir}/FitShapes.root\"\n",
")"
]
},
Expand Down
35 changes: 20 additions & 15 deletions src/HHbbVV/postprocessing/PlotFitsRes.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,18 +40,18 @@
"source": [
"MAIN_DIR = \"../../../\"\n",
"\n",
"plot_dir = f\"{MAIN_DIR}/plots/PostFit/23Jul12\"\n",
"plot_dir = f\"{MAIN_DIR}/plots/PostFit/23Sep1\"\n",
"_ = os.system(f\"mkdir -p {plot_dir}\")"
]
},
{
"cell_type": "code",
"execution_count": 9,
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"# cards_dir = \"f_tests/23May2/nTF1_2_nTF2_1\"\n",
"cards_dir = \"biastests/23Jul11ResScaled100/NMSSM_XToYHTo2W2BTo4Q2B_MX-3000_MY-250\"\n",
"cards_dir = \"f_tests/23May2/nTF1_2_nTF2_1\"\n",
"# cards_dir = \"biastests/23Jul11ResScaled100/NMSSM_XToYHTo2W2BTo4Q2B_MX-3000_MY-250\"\n",
"file = uproot.open(\n",
" f\"/uscms/home/rkansal/hhcombine/cards/{cards_dir}/FitShapes.root\"\n",
" # f\"/uscms/home/rkansal/eos/bbVV/cards/{cards_dir}/FitShapes.root\"\n",
Expand All @@ -60,7 +60,7 @@
},
{
"cell_type": "code",
"execution_count": 11,
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -86,7 +86,7 @@
},
{
"cell_type": "code",
"execution_count": 12,
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -108,7 +108,7 @@
},
{
"cell_type": "code",
"execution_count": 13,
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
Expand Down Expand Up @@ -156,33 +156,38 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 44,
"metadata": {},
"outputs": [],
"source": [
"pass_ylim = 1000\n",
"fail_ylim = 170000\n",
"pass_ylims = [25, 4]\n",
"fail_ylims = [17000, 1700]\n",
"figs = []\n",
"\n",
"for shape, shape_label in shapes.items():\n",
" for region, region_label in selection_regions.items():\n",
" pass_region = region.startswith(\"pass\")\n",
" # if region == \"pass\" or region == \"fail\":\n",
" # continue\n",
" for i, shape_var in enumerate(shape_vars):\n",
" plot_params = {\n",
" \"hists\": hists[shape][region].project(0, i + 1),\n",
" \"sig_keys\": [f\"X[{mx}]->H(bb)Y[{my}](VV)\"],\n",
" \"bg_keys\": [\"QCD\", \"V+Jets\", \"TT\", \"ST\"],\n",
" \"sig_scale_dict\": None,\n",
" \"show\": True,\n",
" \"sig_scale_dict\": {f\"X[{mx}]->H(bb)Y[{my}](VV)\": 10},\n",
" \"show\": False,\n",
" \"year\": \"all\",\n",
" \"ylim\": pass_ylim if pass_region else fail_ylim,\n",
" \"ylim\": pass_ylims[i] if pass_region else fail_ylims[i],\n",
" \"plot_data\": region != \"pass\",\n",
" \"title\": f\"{shape_label} {region_label} Region\",\n",
" \"name\": f\"{plot_dir}/{shape}_{region}_{shape_var.var}.pdf\",\n",
" \"divide_bin_width\": True,\n",
" }\n",
"\n",
" plotting.ratioHistPlot(**plot_params)\n",
"\n",
" break\n",
" break"
" # break\n",
" # break"
]
},
{
Expand Down
Loading

0 comments on commit 653681e

Please sign in to comment.