diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index aada7b4d2..99e8a17d6 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -315,7 +315,7 @@ jobs:
- name: bsm rabbit setup
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run_python.sh scripts/rabbit/setupRabbit.py
- -i $HIST_FILE --lumiScale $LUMI_SCALE --excludeProcGroups QCD WtoNMu_10 WtoNMu_50 --breitwignerWMassWeights
+ -i $HIST_FILE --lumiScale $LUMI_SCALE --addBSM WtoNMu_5 --breitwignerWMassWeights
--postfix bsm -o $WREMNANTS_OUTDIR
- name: bsm rabbit fit
@@ -772,13 +772,13 @@ jobs:
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py --config utilities/styles/styles.py
$WREMNANTS_OUTDIR/ZMassWLike_eta_pt_charge/fitresults_uncorr.hdf5 --result uncorr --title CMS --titlePos 1 --subtitle Preliminary
- -o $WEB_DIR/$PLOT_DIR --legCols 6 --yscale 1.2 --prefit --legCol 6 --legSize large -m Basemodel
+ -o $WEB_DIR/$PLOT_DIR --legCols 6 --yscale 1.2 --prefit --legCol 6 --legSize large -m BaseMapping
- name: wlike postfit plot
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py --config utilities/styles/styles.py
$WREMNANTS_OUTDIR/ZMassWLike_eta_pt_charge/fitresults_uncorr.hdf5 --result uncorr --title CMS --titlePos 1 --subtitle Preliminary
- -o $WEB_DIR/$PLOT_DIR --legCols 6 --yscale 1.2 --rrange 1.03 0.97 --binSeparationLines -2.2 -2.1 1.2 1.3 -m Basemodel
+ -o $WEB_DIR/$PLOT_DIR --legCols 6 --yscale 1.2 --rrange 1.03 0.97 --binSeparationLines -2.2 -2.1 1.2 1.3 -m BaseMapping
- name: wlike rabbit impacts
run: >-
@@ -873,13 +873,13 @@ jobs:
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_fit.py
'$WREMNANTS_OUTDIR/ZMassDilepton_etaAbsEta_mll/ZMassDilepton.hdf5'
- -t -1 --computeHistErrors --saveHists --saveHistsPerProcess --computeVariations --doImpacts -m Basemodel -m Project ch0 mll
+ -t -1 --computeHistErrors --saveHists --saveHistsPerProcess --computeVariations --doImpacts -m BaseMapping -m Project ch0 mll
-o '$WREMNANTS_OUTDIR/ZMassDilepton_etaAbsEta_mll/'
- name: dilepton postfit mll-etaAbsEta plot
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py --config utilities/styles/styles.py
- $WREMNANTS_OUTDIR/ZMassDilepton_etaAbsEta_mll/fitresults.hdf5 -o $WEB_DIR/$PLOT_DIR -m Basemodel
+ $WREMNANTS_OUTDIR/ZMassDilepton_etaAbsEta_mll/fitresults.hdf5 -o $WEB_DIR/$PLOT_DIR -m BaseMapping
--yscale 1.4 --binSeparationLines -0.6 0.0 --titlePos 1 --title CMS --subtitle Preliminary
- name: dilepton postfit mll plot
@@ -1145,12 +1145,12 @@ jobs:
- name: dilepton prefit plot
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py $WREMNANTS_OUTDIR/ZMassDilepton_ptll_yll/fitresults.hdf5
- --config utilities/styles/styles.py -o $WEB_DIR/$PLOT_DIR --yscale 1.2 --prefit --title CMS --subtitle Preliminary -m Basemodel
+ --config utilities/styles/styles.py -o $WEB_DIR/$PLOT_DIR --yscale 1.2 --prefit --title CMS --subtitle Preliminary -m BaseMapping
- name: dilepton postfit plot
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py $WREMNANTS_OUTDIR/ZMassDilepton_ptll_yll/fitresults.hdf5
- --config utilities/styles/styles.py -o $WEB_DIR/$PLOT_DIR --yscale 1.2 --title CMS --subtitle Preliminary -m Basemodel
+ --config utilities/styles/styles.py -o $WEB_DIR/$PLOT_DIR --yscale 1.2 --title CMS --subtitle Preliminary -m BaseMapping
- name: dilepton plotting ptll
run: >-
@@ -1190,7 +1190,7 @@ jobs:
run: >-
scripts/ci/run_with_singularity.sh scripts/ci/setup_and_run.sh rabbit_plot_hists.py --config utilities/styles/styles.py
$WREMNANTS_OUTDIR/ZMassDilepton_ptll/fitresults_from_ZMassWLike_eta_pt_charge.hdf5 --result asimov_uncorr --postfix uncorr_fromWLike
- -o $WEB_DIR/$PLOT_DIR --legCols 1 --yscale 1.2 --rrange 1.03 0.97 --title CMS --subtitle Preliminary -m Basemodel
+ -o $WEB_DIR/$PLOT_DIR --legCols 1 --yscale 1.2 --rrange 1.03 0.97 --title CMS --subtitle Preliminary -m BaseMapping
gen:
# The type of runner that the job will run on
diff --git a/rabbit b/rabbit
index 79953d0c7..99e39bffa 160000
--- a/rabbit
+++ b/rabbit
@@ -1 +1 @@
-Subproject commit 79953d0c7581383db2d614109f2ddc4dfcddfa46
+Subproject commit 99e39bffa108c5aa25e770f452823727a6a95ba2
diff --git a/scripts/corrections/make_theory_corr.py b/scripts/corrections/make_theory_corr.py
index d7bcf0a8a..9b3d8248f 100644
--- a/scripts/corrections/make_theory_corr.py
+++ b/scripts/corrections/make_theory_corr.py
@@ -250,6 +250,7 @@ def main():
eventgen_procs = [
"WtoNMu_MN-5-V-0p001",
"WtoNMu_MN-10-V-0p001",
+ "WtoNMu_MN-30-V-0p001",
"WtoNMu_MN-50-V-0p001",
]
diff --git a/scripts/histmakers/mw_with_mu_eta_pt.py b/scripts/histmakers/mw_with_mu_eta_pt.py
index f2fc8bb11..fac876061 100644
--- a/scripts/histmakers/mw_with_mu_eta_pt.py
+++ b/scripts/histmakers/mw_with_mu_eta_pt.py
@@ -76,6 +76,11 @@
default=common.get_default_mtcut(analysis_label),
help="Value for the transverse mass cut in the event selection",
)
+parser.add_argument(
+ "--mtGenCut",
+ action="store_true",
+ help="Apply mt cut at generator level",
+)
parser.add_argument(
"--vetoGenPartPt",
type=float,
@@ -675,11 +680,7 @@ def build_graph(df, dataset):
logger.info(f"build graph for dataset: {dataset.name}")
results = []
isW = dataset.name in common.wprocs
- isBSM = dataset.name in [
- "WtoNMu_MN-5-V-0p001",
- "WtoNMu_MN-10-V-0p001",
- "WtoNMu_MN-50-V-0p001",
- ]
+ isBSM = dataset.name.startswith("WtoNMu")
isWmunu = isBSM or dataset.name in [
"WplusmunuPostVFP",
"WminusmunuPostVFP",
@@ -785,7 +786,7 @@ def build_graph(df, dataset):
)
cols = [*cols, "run4axis"]
- if isBSM:
+ if isWmunu:
# to compute inclusive cross section
unfolding_tools.add_xnorm_histograms(
results,
@@ -809,7 +810,7 @@ def build_graph(df, dataset):
"pt_min": template_minpt,
"pt_max": template_maxpt,
"abseta_max": template_maxeta,
- "mtw_min": None,
+ "mtw_min": args.mtCut if args.mtGenCut else None,
}
if hasattr(dataset, "out_of_acceptance"):
df = unfolding_tools.select_fiducial_space(
diff --git a/scripts/histmakers/mz_dilepton.py b/scripts/histmakers/mz_dilepton.py
index d1f83c551..36884d92e 100644
--- a/scripts/histmakers/mz_dilepton.py
+++ b/scripts/histmakers/mz_dilepton.py
@@ -273,6 +273,7 @@
"nonTrigMuons_charge0": hist.axis.Regular(
2, -2.0, 2.0, underflow=False, overflow=False, name="nonTrigMuons_charge0"
),
+ "ptll_resolution": hist.axis.Regular(1000, -1, 1, name="ptll_resolution"),
}
auxiliary_gen_axes = [
@@ -1068,6 +1069,25 @@ def build_graph(df, dataset):
)
# test plots
+ if args.validationHists:
+ # resolution plot
+ df = df.Define("ptll_relResolution", "(ptll - postfsrPTV)/postfsrPTV")
+ df = df.Define("ptll_resolution", "(ptll - postfsrPTV)")
+ results.append(
+ df.HistoBoost(
+ f"nominal_relResolution",
+ [all_axes["ptll_resolution"], all_axes["ptll"], axis_absYll],
+ ["ptll_resolution", "postfsrPTV", "absYll", "nominal_weight"],
+ )
+ )
+ results.append(
+ df.HistoBoost(
+ f"nominal_resolution",
+ [all_axes["ptll_resolution"], all_axes["ptll"], axis_absYll],
+ ["ptll_resolution", "postfsrPTV", "absYll", "nominal_weight"],
+ )
+ )
+
if args.validationHists and args.useDileptonTriggerSelection:
df_plusTrig = df.Filter("trigMuons_passTrigger0")
df_minusTrig = df.Filter("nonTrigMuons_passTrigger0")
diff --git a/scripts/plotting/inclusive_xsec_summary.py b/scripts/plotting/inclusive_xsec_summary.py
index 717a46b4d..2659656a7 100644
--- a/scripts/plotting/inclusive_xsec_summary.py
+++ b/scripts/plotting/inclusive_xsec_summary.py
@@ -47,14 +47,12 @@
pdf_results = {}
comp_result = {}
-pdf_lumis = {}
+# pdf_lumis = {}
for pdf_file in args.pdfFiles:
pdf_name = pdf_file.split("/")[-2].split("_")[-1]
pdf_result, pdf_meta = rabbit.io_tools.get_fitresult(pdf_file, meta=True)
- pdf_lumis[pdf_name] = pdf_meta["meta_info_input"]["channel_info"]["ch0"]["lumi"]
-
pdf_model = pdf_result["physics_models"]
if "CompositeModel" in pdf_model.keys():
@@ -79,9 +77,15 @@
(r"$\mathrm{W}^{+}$", "Project ch1_masked qGen", "ch1_masked", {"qGen": 1}),
(r"$\mathrm{W}$", "Project ch1_masked", "ch1_masked", None),
(r"$\mathrm{Z}$", "Project ch0_masked", "ch0_masked", None),
+ # (
+ # r"$\mathrm{W}^{-}/\mathrm{W}^{+}$",
+ # "Ratio ch1_masked ch1_masked qGen:0,ptGen:sum,absEtaGen:sum qGen:1,ptGen:sum,absEtaGen:sum",
+ # "ch1_masked",
+ # None,
+ # ),
(
r"$\mathrm{W}^{+}/\mathrm{W}^{-}$",
- "Ratio ch1_masked ch1_masked qGen:0,ptGen:sum,absEtaGen:sum qGen:1,ptGen:sum,absEtaGen:sum",
+ "Ratio ch1_masked ch1_masked qGen:1,ptGen:sum,absEtaGen:sum qGen:0,ptGen:sum,absEtaGen:sum",
"ch1_masked",
None,
),
@@ -130,18 +134,13 @@
h1 = h1[{"yield": hist.sum}]
hi = hi[{"yield": hist.sum}]
- if model.startswith("Ratio"):
- scale = 1
- else:
- scale = 1 / (lumi * 1000)
-
- prefit = hp.value * scale
- prefit_error = hp.variance**0.5 * scale
+ prefit = hp.value
+ prefit_error = hp.variance**0.5
- value = h1.value * scale
- error = h1.variance**0.5 * scale
+ value = h1.value
+ error = h1.variance**0.5
- impacts = hi.values() * scale
+ impacts = hi.values()
labels = np.array(hi.axes["impacts"])
mask = np.isin(labels, grouping)
@@ -190,22 +189,26 @@
df["prefit_error"] = prefit_error
for pdf_name, pdf_res in pdf_results.items():
- hr = pdf_res[model.replace("_masked", "")]["channels"][
+ channel_models = pdf_res[model.replace("_masked", "")]["channels"][
channel.replace("_masked", "")
- ]["hist_prefit_inclusive"].get()
+ ]
+ hr = channel_models["hist_prefit_inclusive"].get()
+ hr_impacts = channel_models[
+ "hist_prefit_inclusive_global_impacts_grouped"
+ ].get()
if selection is not None:
hr = hr[selection]
+ hr_impacts = hr_impacts[selection]
if getattr(hr, "axes", False) and "yield" in hr.axes.name:
hr = hr[{"yield": hist.sum}]
+ hr_impacts = hr_impacts[{"yield": hist.sum}]
- if model.startswith("Ratio"):
- scale = 1
- else:
- scale = 1 / (pdf_lumis[pdf_name] * 1000)
-
- df[pdf_name] = hr.value * scale
- df[f"{pdf_name}_error"] = hr.variance**0.5 * scale
+ df[pdf_name] = hr.value
+ df[f"{pdf_name}_error"] = hr.variance**0.5
+ df[f"{pdf_name}_pdf"] = hr_impacts[
+ {"impacts": f"pdf{pdf_name.replace('aN3LO','an3lo')}"}
+ ]
# Convert 'labels' column to categorical with the custom order
df["label"] = pd.Categorical(df["label"], categories=custom_order, ordered=True)
@@ -305,6 +308,17 @@
marker="o",
label=pdf_name if i == 0 else None,
)
+ # # only plot PDF uncertainties
+ # pdf_error_pdf = df_g[f"{pdf_name}_pdf"].values[0] / norm
+ # ax.errorbar(
+ # [pdf_value],
+ # [i + 1 - (j + 1) / (nPDFs + 1)],
+ # xerr=pdf_error_pdf,
+ # color=pdf_colors[pdf_name],
+ # capsize=5,
+ # capthick=2,
+ # marker="o",
+ # )
# round to two significant digits in total uncertainty
sig_digi = 2 - int(math.floor(math.log10(abs(total)))) - 1
@@ -387,7 +401,7 @@
ax.set_xlim([lo, hi])
ax.set_ylim([0, len(norms) + 2])
-ax.set_xlabel("1./Measurement", fontsize=20)
+ax.set_xlabel("Prediction / Measurement")
# Disable ticks on the top and right axes
ax.tick_params(top=False)
@@ -442,10 +456,10 @@ def plot_cov_ellipse(cov, pos, nstd=2, **kwargs):
("WZ", xsec_keys[2], xsec_keys[3], "pb"),
("R", xsec_keys[4], xsec_keys[5], None),
):
- ckey1 = (
+ ckey0 = (
channel0[1].replace("_masked", "") + " " + channel0[2].replace("_masked", "")
)
- ckey2 = (
+ ckey1 = (
channel1[1].replace("_masked", "") + " " + channel1[2].replace("_masked", "")
)
@@ -456,10 +470,6 @@ def plot_cov_ellipse(cov, pos, nstd=2, **kwargs):
for pdf_name, result in comp_result.items():
ibin = 0
- if name == "R":
- scale = 1
- else:
- scale = 1 / (1000 * pdf_lumis[pdf_name])
for k, r in result["channels"].items():
fittype = "postfit" if f"hist_postfit_inclusive" in r.keys() else "prefit"
@@ -467,41 +477,42 @@ def plot_cov_ellipse(cov, pos, nstd=2, **kwargs):
if getattr(hi, "axes", False) and "yield" in hi.axes.name:
hi = hi[{"yield": hist.sum}]
- if k == ckey1:
+ if k == ckey0:
sel = channel0[-1]
if sel is not None:
- x = hi[sel].value * scale
+ x = hi[sel].value
ix = ibin + [i for i in sel.values()][0]
else:
- x = hi.value * scale
+ x = hi.value
ix = ibin
- if k == ckey2:
+ if k == ckey1:
sel = channel1[-1]
if sel is not None:
- y = hi[channel0[-1]].value * scale
+ y = hi[sel].value
iy = ibin + [i for i in sel.values()][0]
else:
- y = hi.value * scale
+ y = hi.value
iy = ibin
ibin += hi.size if hasattr(hi, "size") else 1
cov = result[f"hist_{fittype}_inclusive_cov"].get().values()
- cov = cov[np.ix_([ix, iy], [ix, iy])] * scale**2
+ cov = cov[np.ix_([ix, iy], [ix, iy])]
# for pos, cov in zip(points, covs):
if fittype == "postfit":
- icol = "grey"
ell = plot_cov_ellipse(
cov,
np.array([x, y]),
nstd=2,
edgecolor="none",
- facecolor=icol,
+ facecolor="grey",
label="Measurement",
)
+ ax.add_patch(ell)
+ ax.plot(x, y, color="black", marker="P")
else:
icol = pdf_colors[pdf_name]
ell = plot_cov_ellipse(
@@ -513,8 +524,8 @@ def plot_cov_ellipse(cov, pos, nstd=2, **kwargs):
linewidth=2,
label=pdf_name,
)
- ax.add_patch(ell)
- ax.plot(x, y, color=icol, marker="o", alpha=0) # measurement center
+ ax.add_patch(ell)
+ ax.plot(x, y, color=icol, marker="o", alpha=0)
xlim = ax.get_xlim()
ylim = ax.get_ylim()
@@ -550,7 +561,7 @@ def plot_cov_ellipse(cov, pos, nstd=2, **kwargs):
padding_loc="auto",
)
- plot_tools.add_cms_decor(ax, args.cmsDecor, data=True, lumi=lumi, loc=args.logoPos)
+ plot_tools.add_cms_decor(ax, args.cmsDecor, data=True, loc=args.logoPos)
outname = f"summary_2D_{name}"
if args.postfix:
diff --git a/scripts/plotting/makeDataMCStackPlot.py b/scripts/plotting/makeDataMCStackPlot.py
index d4b976c56..81cab0124 100644
--- a/scripts/plotting/makeDataMCStackPlot.py
+++ b/scripts/plotting/makeDataMCStackPlot.py
@@ -82,6 +82,7 @@
"QCD",
"WtoNMu_5",
"WtoNMu_10",
+ "WtoNMu_30",
"WtoNMu_50",
],
)
diff --git a/scripts/plotting/response_matrix.py b/scripts/plotting/response_matrix.py
index 4c147c893..83dd3e2e3 100644
--- a/scripts/plotting/response_matrix.py
+++ b/scripts/plotting/response_matrix.py
@@ -51,6 +51,8 @@
args = parser.parse_args()
+cmap = plt.get_cmap("tab10")
+
logger = logging.setup_logger(__file__, args.verbose, args.noColorLogger)
outdir = output_tools.make_plot_dir(args.outpath, args.outfolder, eoscp=args.eoscp)
@@ -66,14 +68,18 @@
) # set gen axes empty to not integrate over when loading
if "Wmunu" in groups.groups:
- groups.copyGroup(
- "Wmunu", "Wmunu_qGen0", member_filter=lambda x: x.name.startswith("Wminus")
- )
- groups.copyGroup(
- "Wmunu", "Wmunu_qGen1", member_filter=lambda x: x.name.startswith("Wplus")
- )
+ if "minus" in args.channels:
+ groups.copyGroup(
+ "Wmunu", "Wmunu_qGen0", member_filter=lambda x: x.name.startswith("Wminus")
+ )
+
+ if "plus" in args.channels:
+ groups.copyGroup(
+ "Wmunu", "Wmunu_qGen1", member_filter=lambda x: x.name.startswith("Wplus")
+ )
- groups.deleteGroup("Wmunu")
+ if "all" not in args.channels:
+ groups.deleteGroup("Wmunu")
if "Zmumu" in groups.groups:
groups.groups["Zmumu"].deleteMembers(
@@ -94,12 +100,12 @@
"eta": r"$\mathrm{Reco}\ \eta$",
"abs(eta)": r"$\mathrm{Reco}\ |\eta|$",
"absEtaGen": r"$\mathrm{Gen}\ |\eta|$",
- "ptll": r"$\mathrm{Reco}\ p_\mathrm{T}(\ell\ell)\ [\mathrm{GeV}]$",
- "ptW": r"$\mathrm{Reco}\ p_\mathrm{T}(\ell\nu)\ [\mathrm{GeV}]$",
- "ptVGen": r"$\mathrm{Gen}\ p_\mathrm{T}(V)\ [\mathrm{GeV}]$",
- "yll": r"$\mathrm{Reco}\ Y(\mathrm{V})$",
- "abs(yll)": r"$\mathrm{Reco}\ |Y(\mathrm{V})|$",
- "absYVGen": r"$\mathrm{Gen}\ |Y(\mathrm{V})|$",
+ "ptll": r"$\mathrm{Reco}\ p_\mathrm{T}^{\ell\ell}\ [\mathrm{GeV}]$",
+ "ptW": r"$\mathrm{Reco}\ p_\mathrm{T}^{\ell\nu}\ [\mathrm{GeV}]$",
+ "ptVGen": r"$\mathrm{Gen}\ p_\mathrm{T}^Z\ [\mathrm{GeV}]$",
+ "yll": r"$\mathrm{Reco}\ y^{\ell\ell}$",
+ "abs(yll)": r"$\mathrm{Reco}\ |y^{\ell\ell}|$",
+ "absYVGen": r"$\mathrm{Gen}\ |y^\mathrm{Z}|$",
"cosThetaStarll": r"$\mathrm{Reco}\ \cos{\theta^{\star}_{\ell\ell}}$",
"phiStarll": r"$\mathrm{Reco}\ \phi^{\star}_{\ell\ell}$",
"helicitySig": r"Helicity",
@@ -187,7 +193,7 @@ def plot_resolution(
ax.set_ylabel(ylabel)
xedges = None
- for sel2, idx2 in selections_slices:
+ for j, (sel2, idx2) in enumerate(selections_slices):
if (
(
idx2 not in [hist.underflow, hist.overflow]
@@ -208,11 +214,19 @@ def plot_resolution(
else:
edges = h2d.axes[sel2].edges
var2 = translate_label[sel2].replace(r"[\mathrm{GeV}]", "")
+ unit = "GeV" if r"[\mathrm{GeV}]" in translate_label[sel2] else ""
if idx2 == hist.overflow:
- label = f"{var2} > {edges[-1]}"
+ lo = edges[-1]
+ if int(lo) == lo:
+ lo = int(lo)
+ label = f"{var2} > {lo} {unit}"
elif idx2 + 1 < len(edges):
lo, hi = edges[idx2], edges[idx2 + 1]
- label = f"{lo} < {var2} < {hi}"
+ if int(lo) == lo:
+ lo = int(lo)
+ if int(hi) == hi:
+ hi = int(hi)
+ label = f"{lo} < {var2} < {hi} {unit}"
h1d = h2d[{sel2: idx2}]
if len(axes_reco) > 1:
h1d = hh.unrolledHist(h1d, binwnorm=None, obs=axes_reco)
@@ -224,7 +238,9 @@ def plot_resolution(
if normalize:
values /= abs(values).sum()
- ax.stairs(values, edges=h1d.axes[0].edges, label=label)
+ ax.stairs(
+ values, edges=h1d.axes[0].edges, label=label, color=cmap(j % cmap.N)
+ )
if xedges is None:
xedges = h1d.axes[0].edges
@@ -243,14 +259,19 @@ def plot_resolution(
if sel is not None:
lo, hi = histo.axes[sel].edges[idx], histo.axes[sel].edges[idx + 1]
- var = translate_label[sel].replace(r"[\mathrm{GeV}]", "")
+ var = translate_label[sel].replace(r"\ [\mathrm{GeV}]", "")
+ unit = "GeV" if r"[\mathrm{GeV}]" in translate_label[sel] else ""
+ if int(hi) == hi:
+ hi = int(hi)
+ if int(lo) == lo:
+ lo = int(lo)
if sel.startswith("abs") and lo == 0:
- title = f"{var} < {hi}"
+ title = f"{var} < {hi} {unit}"
else:
- title = f"{lo} < {var} < {hi}"
+ title = f"{lo} < {var} < {hi} {unit}"
plt.text(
- 0.06,
+ 0.04,
0.94,
title,
horizontalalignment="left",
@@ -292,6 +313,93 @@ def plot_resolution(
else {"charge": -1.0j if channel == "minus" else 1.0j}
)
+ if "resolution" in args.histName or "relResolution" in args.histName:
+
+ if "relResolution" in args.histname:
+ relative = True
+ else:
+ relative = False
+
+ fig = plt.figure(figsize=(10, 6))
+ ax = fig.add_subplot()
+
+ vreco = translate_label["ptll"]
+ vgen = translate_label["ptVGen"]
+ binlabel = translate_label["abs(yll)"]
+
+ ax.set_xlabel(vgen)
+ vreco = (
+ vreco.replace(r"[\mathrm{GeV}]", "")
+ .replace("$", "")
+ .replace("(V)", "(Z)")
+ .replace("Reco", "")
+ .replace(r"\ ", "")
+ )
+ vgen = (
+ vgen.replace(r"[\mathrm{GeV}]", "")
+ .replace("$", "")
+ .replace("(V)", "(Z)")
+ .replace("Gen", "")
+ .replace(r"\ ", "")
+ )
+ if relative:
+ ylabel = rf"$\sigma(({vreco} - {vgen}) / {vgen})$"
+ else:
+ ylabel = rf"$\sigma({vreco} - {vgen})$"
+ ax.set_ylabel(ylabel)
+
+ x = histo.axes["ptll"].edges
+ x_centers = histo.axes["ptll"].centers
+ r_centers = histo.axes["ptll_resolution"].centers
+
+ for idx, ibin in enumerate(histo.axes["absYll"]):
+ values = histo[{"absYll": idx}].values()
+
+ mean_res = np.zeros(len(x_centers))
+ std_res = np.zeros(len(x_centers))
+
+ for i in range(len(x_centers)):
+ weights = values[:, i]
+
+ if weights.sum() == 0:
+ mean_res[i] = np.nan
+ std_res[i] = np.nan
+ continue
+
+ # weighted mean
+ mean = np.average(r_centers, weights=weights)
+
+ # weighted variance
+ var = np.average((r_centers - mean) ** 2, weights=weights)
+
+ mean_res[i] = mean
+ std_res[i] = np.sqrt(var)
+
+ ax.stairs(std_res, x, label=f"{binlabel}[{idx}]")
+
+ plot_tools.add_cms_decor(
+ ax, args.cmsDecor, data=True, lumi=None, loc=args.logoPos
+ )
+ plot_tools.addLegend(
+ ax, ncols=args.legCols, loc=args.legPos, text_size=args.legSize
+ )
+
+ ax.set_ylim(0, 0.9)
+ ax.set_xlim(0, 100)
+
+ outfile = f"resolution_1d_{g_name}_pt"
+
+ plot_tools.save_pdf_and_png(outdir, outfile)
+
+ output_tools.write_index_and_log(
+ outdir,
+ outfile,
+ analysis_meta_info={args.infile: groups.getMetaInfo()},
+ args=args,
+ )
+
+ continue
+
# plot slices of resolution
if all(x in histo.axes.name for x in ["pt", "ptGen", "absEtaGen"]):
plot_resolution(
@@ -299,7 +407,13 @@ def plot_resolution(
axes_reco="pt",
axis_gen="ptGen",
selections_global=(("absEtaGen", 0), ("absEtaGen", 17)),
- selections_slices=(("ptGen", 0), ("ptGen", 6), ("ptGen", 13)),
+ selections_slices=(
+ ("ptGen", 0),
+ ("ptGen", 3),
+ ("ptGen", 6),
+ ("ptGen", 10),
+ ("ptGen", 13),
+ ),
)
if all(x in histo.axes.name for x in ["ptGen", "eta", "absEtaGen"]):
plot_resolution(
@@ -309,7 +423,9 @@ def plot_resolution(
selections_global=(("ptGen", 0), ("ptGen", 13)),
selections_slices=(
("absEtaGen", 0),
- ("absEtaGen", 8),
+ ("absEtaGen", 3),
+ ("absEtaGen", 7),
+ ("absEtaGen", 10),
("absEtaGen", 17),
),
)
@@ -321,11 +437,16 @@ def plot_resolution(
axis_gen="ptVGen",
selections_global=(
("absYVGen", 0),
- ("absYVGen", 3),
+ ("absYVGen", 4),
),
selections_slices=(
("ptVGen", 0),
+ ("ptVGen", 4),
("ptVGen", 8),
+ ("ptVGen", 12),
+ ("ptVGen", 14),
+ ("ptVGen", 16),
+ ("ptVGen", 18),
("ptVGen", hist.overflow),
),
)
@@ -341,6 +462,7 @@ def plot_resolution(
selections_slices=(
("absYVGen", 0),
("absYVGen", 2),
+ ("absYVGen", 4),
("absYVGen", hist.overflow),
),
)
diff --git a/scripts/rabbit/setupRabbit.py b/scripts/rabbit/setupRabbit.py
index 28f8e6257..0e24dd3a7 100644
--- a/scripts/rabbit/setupRabbit.py
+++ b/scripts/rabbit/setupRabbit.py
@@ -166,7 +166,7 @@ def make_parser(parser=None):
type=str,
nargs="*",
help="Don't run over processes belonging to these groups (only accepts exact group names)",
- default=["QCD", "WtoNMu_5", "WtoNMu_10", "WtoNMu_50"],
+ default=["QCD"],
)
parser.add_argument(
"--filterProcGroups",
@@ -175,6 +175,13 @@ def make_parser(parser=None):
help="Only run over processes belonging to these groups",
default=[],
)
+ parser.add_argument(
+ "--addBSM",
+ type=str,
+ default=None,
+ help="Add BSM model",
+ choices=["WtoNMu_0", "WtoNMu_5", "WtoNMu_10", "WtoNMu_30", "WtoNMu_50"],
+ )
parser.add_argument(
"-x",
"--excludeNuisances",
@@ -224,10 +231,14 @@ def make_parser(parser=None):
)
parser.add_argument(
"--rebin",
- type=int,
+ type=parsing.str_to_list_or_int,
nargs="*",
default=[],
- help="Rebin axis by this value (default, 1, does nothing)",
+ help="""
+ Rebin axis by this value (default, 1, does nothing);
+ use integer 'n' for merging 'n' bins;
+ use comma separated list with new edges, use a leading space in case the bin edge is negative (e.g. " -2.4")
+ """,
)
parser.add_argument(
"--absval",
@@ -927,6 +938,7 @@ def setup(
xnorm=any(
inputBaseName.startswith(x) for x in ["gen", "xnorm", "prefsr", "postfsr"]
),
+ bsm_model=args.addBSM,
)
datagroups.fit_axes = fitvar
@@ -994,6 +1006,15 @@ def setup(
else:
base_group = "Zee" if datagroups.flavor == "ee" else "Zmumu"
+ if args.addBSM == "WtoNMu_0":
+ # an empty group should have been created, add SM W members as proxy for BSM
+ datagroups.groups["WtoNMu_0"].addMembers(datagroups.groups["Wmunu"].members)
+
+ # normalize total cross section to BSM default
+ total_xsec = sum([m.xsec for m in datagroups.groups["WtoNMu_0"].members])
+ for member in datagroups.groups["WtoNMu_0"].members:
+ member.xsec = common.xsec_WtoNMu * member.xsec / total_xsec
+
if args.addTauToSignal:
# add tau signal processes to signal group
datagroups.groups[base_group].addMembers(
diff --git a/utilities/common.py b/utilities/common.py
index ed8bb6180..6aa4e42ba 100644
--- a/utilities/common.py
+++ b/utilities/common.py
@@ -58,6 +58,7 @@
"WminusCharmToMuNu",
"WtoNMu_MN-5-V-0p001",
"WtoNMu_MN-10-V-0p001",
+ "WtoNMu_MN-30-V-0p001",
"WtoNMu_MN-50-V-0p001",
]
zprocs = [
@@ -353,6 +354,7 @@ def get_dilepton_ptV_binning(fine=False):
[
0,
1,
+ 1.5,
2,
2.5,
3,
diff --git a/utilities/differential.py b/utilities/differential.py
index 075b23aae..a79398051 100644
--- a/utilities/differential.py
+++ b/utilities/differential.py
@@ -104,7 +104,7 @@ def get_dilepton_axes(
gen_level,
add_out_of_acceptance_axis=False,
flow_y=False,
- rebin_pt=False,
+ rebin_pt=None,
):
"""
construct axes, columns, and selections for differential Z dilepton measurement from correponding reco edges. Currently supported: pT(Z), |yZ|
@@ -129,13 +129,8 @@ def get_dilepton_axes(
if v == "ptVGen":
edges = reco_edges["ptll"]
- if rebin_pt:
- # use 2 ptll bin for each ptVGen bin, last bin is overflow
- if len(edges) % 2:
- # in case it's an odd number of edges, last two bins are overflow
- edges = edges[:-1]
- # 1 gen bin for 2 reco bins
- edges = edges[::2]
+ if rebin_pt is not None:
+ edges = rebin_pt(edges)
axes.append(
hist.axis.Variable(
diff --git a/utilities/parsing.py b/utilities/parsing.py
index aba18ad6c..1a7e1eef0 100644
--- a/utilities/parsing.py
+++ b/utilities/parsing.py
@@ -29,6 +29,17 @@ def str_to_complex_or_int(value):
raise argparse.ArgumentTypeError(f"Invalid integer: '{value}'")
+def str_to_list_or_int(value):
+ # .strip() handles the "leading space" trick if you use it
+ value = value.strip()
+ if "," in value:
+ return [float(x) for x in value.split(",")]
+ try:
+ return int(value)
+ except ValueError:
+ raise argparse.ArgumentTypeError(f"Invalid integer: '{value}'")
+
+
def set_parser_attribute(parser, argument, attribute, newValue):
# change an argument of the parser, must be called before parse_arguments
logger = logging.child_logger(__name__)
diff --git a/utilities/styles/styles.py b/utilities/styles/styles.py
index 87a9ef761..e2a3c70e0 100644
--- a/utilities/styles/styles.py
+++ b/utilities/styles/styles.py
@@ -12,12 +12,14 @@ def translate_html_to_latex(n):
# transform html style formatting into latex style
if "", r"\mathit{")
+ n.replace("", r"\mathit{")
+ .replace("", "}")
.replace("", "_{")
.replace("", "^{")
- .replace("", "}")
.replace("", "}")
.replace("", "}")
+ .replace(r"μ", r"\mu")
+ .replace(r"ε", r"\epsilon")
.replace(" ", r"\ ")
)
return n
@@ -52,9 +54,10 @@ def translate_html_to_latex(n):
"Fake_e": "#964A8B",
"Fake_mu": "#964A8B",
"Prompt": "#E42536",
- "WtoNMu_5": "#409C3D",
- "WtoNMu_10": "#2BC74D",
- "WtoNMu_50": "#00FF80",
+ "WtoNMu_5": "#006400",
+ "WtoNMu_10": "#009933",
+ "WtoNMu_30": "#00CC99",
+ "WtoNMu_50": "#00FFFF",
"BSM": "#409C3D",
}
@@ -116,9 +119,10 @@ def translate_html_to_latex(n):
"Fake_e": "Nonprompt (e)",
"Fake_mu": r"Nonprompt ($\mu$)",
"Prompt": "Prompt",
- "WtoNMu_5": r"W$^{\pm}\to\mathrm{N}\mu (5GeV)$",
- "WtoNMu_10": r"W$^{\pm}\to\mathrm{N}\mu (10GeV)$",
- "WtoNMu_50": r"W$^{\pm}\to\mathrm{N}\mu (50GeV)$",
+ "WtoNMu_5": r"W$^{\pm}\to\mu\mathrm{N} (5GeV)$",
+ "WtoNMu_10": r"W$^{\pm}\to\mu\mathrm{N} (10GeV)$",
+ "WtoNMu_30": r"W$^{\pm}\to\mu\mathrm{N} (30GeV)$",
+ "WtoNMu_50": r"W$^{\pm}\to\mu\mathrm{N} (50GeV)$",
}
axis_labels = {
@@ -132,7 +136,7 @@ def translate_html_to_latex(n):
"muonJetPt": {"label": r"$\mathit{p}_{T}^\mathrm{jet[\mu]}$", "unit": "GeV"},
"recoil_perp": {"label": r"$\it{u}_{\mathrm{T}}^{\perp}$", "unit": "GeV"},
"recoil_para": {"label": r"$\it{u}_{\mathrm{T}}^{\parallel}$", "unit": "GeV"},
- "qGen": r"$|\mathit{q}^{\mu}|$",
+ "qGen": r"$\mathit{q}^{\mu}$",
"eta": r"$\mathit{\eta}^{\mu}$",
"etaGen": r"$\mathit{\eta}^{\mu}$",
"abseta": r"$|\mathit{\eta}^{\mu}|$",
@@ -376,6 +380,7 @@ def translate_html_to_latex(n):
"binByBinStatZtautau",
"binByBinStatWtoNMu_5",
"binByBinStatWtoNMu_10",
+ "binByBinStatWtoNMu_30",
"binByBinStatWtoNMu_50",
],
"unfolding": [
@@ -406,6 +411,10 @@ def translate_html_to_latex(n):
"pdfCT18Z",
"pTModeling",
"theory_ew",
+ "massShift",
+ "widthW",
+ "widthZ",
+ "sin2thetaZ",
],
"unfolding_min": [
"Total",
@@ -433,6 +442,30 @@ def translate_html_to_latex(n):
"bcQuarkMass",
"theory_ew",
],
+ "xsecs": [
+ "Total",
+ "stat",
+ "binByBinStat",
+ "binByBinStatW",
+ "binByBinStatZ",
+ # "expNoLumi",
+ "luminosity",
+ "angularCoeffs",
+ "pdfCT18Z",
+ "pTModeling",
+ "theory_ew",
+ "massShift",
+ "widthW",
+ "widthZ",
+ "sin2thetaZ",
+ "muon_eff_syst",
+ "muon_eff_stat",
+ "prefire",
+ "muonCalibration",
+ "Fake",
+ "recoil",
+ "CMS_background",
+ ],
}
text_dict = {
@@ -465,9 +498,10 @@ def translate_html_to_latex(n):
}
impact_labels = {
- "WtoNMu_5": "μW→μN(5GeV)",
- "WtoNMu_10": "μW→μN(10GeV)",
- "WtoNMu_50": "μW→μN(50GeV)",
+ "WtoNMu_5": "μW→Nμ(5GeV)",
+ "WtoNMu_10": "μW→Nμ(10GeV)",
+ "WtoNMu_30": "μW→Nμ(30GeV)",
+ "WtoNMu_50": "μW→Nμ(50GeV)",
"massShiftZ100MeV": "mZ",
"massShiftW100MeV": "mW",
"widthZ": "ΓmZ",
@@ -483,7 +517,7 @@ def translate_html_to_latex(n):
"QCDscaleWPtHelicityMiNNLO": "μR μF scale (W)",
"QCDscaleZPtChargeHelicityMiNNLO": "μR μF scale (Z)",
"QCDscaleWPtChargeHelicityMiNNLO": "μR μF scale (W)",
- "binByBinStat": "Bin-by-bin stat.",
+ "binByBinStat": "Simulation stat.",
"binByBinStatW": "Bin-by-bin stat. (W)",
"binByBinStatZ": "Bin-by-bin stat. (Z)",
"binByBinStatWmunu": "Bin-by-bin stat. (W→μν)",
@@ -496,9 +530,10 @@ def translate_html_to_latex(n):
"binByBinStatTop": "Bin-by-bin stat. (top)",
"binByBinStatWtoNMu_5": "Bin-by-bin stat. (BSM)",
"binByBinStatWtoNMu_10": "Bin-by-bin stat. (BSM)",
+ "binByBinStatWtoNMu_30": "Bin-by-bin stat. (BSM)",
"binByBinStatWtoNMu_50": "Bin-by-bin stat. (BSM)",
"recoil": "recoil",
- "CMS_background": "Bkg.",
+ "CMS_background": "Other bkg.",
"FakeHighMT": "FakeHighMT",
"FakeLowMT": "FakeLowMT",
"rFake": "fakerate",
diff --git a/wremnants/datasets/datagroup.py b/wremnants/datasets/datagroup.py
index aed47966c..c21725857 100644
--- a/wremnants/datasets/datagroup.py
+++ b/wremnants/datasets/datagroup.py
@@ -69,7 +69,7 @@ def addMember(self, member, member_operation=None):
self.memberOp = [None] * len(self.members)
self.memberOp.append(deepcopy(member_operation))
- self.members.append(member)
+ self.members.append(deepcopy(member))
def deleteMembers(self, members=None):
if members is None:
diff --git a/wremnants/datasets/datagroups.py b/wremnants/datasets/datagroups.py
index 59a7eb3f5..2b9875900 100644
--- a/wremnants/datasets/datagroups.py
+++ b/wremnants/datasets/datagroups.py
@@ -895,9 +895,10 @@ def getHistForUnfolding(
if self.groups[group_name].memberOp is not None:
base_member_op = self.groups[group_name].memberOp[base_member_idx]
- return base_member_op(nominal_hist)
- else:
- return nominal_hist
+ if base_member_op is not None:
+ return base_member_op(nominal_hist)
+
+ return nominal_hist
def getPOINames(self, gen_bin_indices, axes_names, base_name, flow=True):
poi_names = []
diff --git a/wremnants/datasets/datagroups2016.py b/wremnants/datasets/datagroups2016.py
index 8104eeeab..e4c8ed9b9 100644
--- a/wremnants/datasets/datagroups2016.py
+++ b/wremnants/datasets/datagroups2016.py
@@ -4,7 +4,12 @@
def make_datagroups_2016(
- dg, combine=False, pseudodata_pdfset=None, excludeGroups=None, filterGroups=None
+ dg,
+ combine=False,
+ pseudodata_pdfset=None,
+ excludeGroups=None,
+ filterGroups=None,
+ bsm_model=None,
):
# reset datagroups
dg.groups = {}
@@ -65,18 +70,12 @@ def make_datagroups_2016(
"QCD",
members=dg.get_members_from_results(startswith=["QCD"]),
)
- dg.addGroup(
- "WtoNMu_5",
- members=dg.get_members_from_results(startswith=["WtoNMu_MN-5-"]),
- )
- dg.addGroup(
- "WtoNMu_10",
- members=dg.get_members_from_results(startswith=["WtoNMu_MN-10-"]),
- )
- dg.addGroup(
- "WtoNMu_50",
- members=dg.get_members_from_results(startswith=["WtoNMu_MN-50-"]),
- )
+ if bsm_model is not None:
+ model, mass = bsm_model.split("_")
+ dg.addGroup(
+ bsm_model,
+ members=dg.get_members_from_results(startswith=[f"{model}_MN-{mass}-"]),
+ )
else:
dg.addGroup(
"Other",
diff --git a/wremnants/datasets/datagroupsLowPU.py b/wremnants/datasets/datagroupsLowPU.py
index e06cb67b5..440af8072 100644
--- a/wremnants/datasets/datagroupsLowPU.py
+++ b/wremnants/datasets/datagroupsLowPU.py
@@ -3,7 +3,9 @@
logger = logging.child_logger(__name__)
-def make_datagroups_lowPU(dg, combine=False, excludeGroups=None, filterGroups=None):
+def make_datagroups_lowPU(
+ dg, combine=False, excludeGroups=None, filterGroups=None, bsm_model=None
+):
# reset datagroups
dg.groups = {}
diff --git a/wremnants/datasets/datasetDict_v9.py b/wremnants/datasets/datasetDict_v9.py
index 7b295ec29..c4d45cee7 100644
--- a/wremnants/datasets/datasetDict_v9.py
+++ b/wremnants/datasets/datasetDict_v9.py
@@ -231,6 +231,13 @@
"xsec": common.xsec_WtoNMu,
"group": "WtoNMu_10",
},
+ "WtoNMu_MN-30-V-0p001": {
+ "filepaths": [
+ "{BASE_PATH}/WtoNMu_MN-30-V-0p001_TuneCP5_13TeV_madgraph-pythia8/"
+ ],
+ "xsec": common.xsec_WtoNMu,
+ "group": "WtoNMu_30",
+ },
"WtoNMu_MN-50-V-0p001": {
"filepaths": [
"{BASE_PATH}/WtoNMu_MN-50-V-0p001_TuneCP5_13TeV_madgraph-pythia8/"
diff --git a/wremnants/helicity_utils.py b/wremnants/helicity_utils.py
index 7b823ffd0..7384f3bdc 100644
--- a/wremnants/helicity_utils.py
+++ b/wremnants/helicity_utils.py
@@ -11,6 +11,7 @@
from utilities.io_tools import input_tools
from wremnants.correctionsTensor_helper import makeCorrectionsTensor
from wremnants.theory_tools import helicity_xsec_to_angular_coeffs
+from wums import boostHistHelpers as hh
from wums import logging
logger = logging.child_logger(__name__)
@@ -32,7 +33,7 @@
def make_helicity_weight_helper(
is_z=False,
filename=f"{common.data_dir}/angularCoefficients/w_z_helicity_xsecs_theoryAgnosticBinning_scetlib_dyturboCorr_maxFiles_m1.hdf5",
- rebi_ptVgen=False,
+ rebin_ptVgen_edges=None,
):
with h5py.File(filename, "r") as ff:
@@ -41,8 +42,10 @@ def make_helicity_weight_helper(
hist_helicity_xsec_scales = out["Z"] if is_z else out["W"]
# rebinning must happen *before* calculating Ai's
- if rebi_ptVgen:
- hist_helicity_xsec_scales = hist_helicity_xsec_scales[{"ptVgen": hist.rebin(2)}]
+ if rebin_ptVgen_edges is not None:
+ hist_helicity_xsec_scales = hh.rebinHist(
+ hist_helicity_xsec_scales, "ptVgen", rebin_ptVgen_edges
+ )
corrh = helicity_xsec_to_angular_coeffs(hist_helicity_xsec_scales)
diff --git a/wremnants/theory_tools.py b/wremnants/theory_tools.py
index 143064723..a4ceb573a 100644
--- a/wremnants/theory_tools.py
+++ b/wremnants/theory_tools.py
@@ -657,7 +657,7 @@ def define_postfsr_vars(df, mode=None):
)
df = df.Define(
"postfsrOtherLep_mass",
- "isEvenEvent ? wrem::get_pdgid_mass(GenPart_pdgId[postfsrLep][postfsrAntiLep_idx]) : wrem::get_pdgid_mass(GenPart_pdgId[postfsrAntiLep][postfsrLep_idx])",
+ "isEvenEvent ? wrem::get_pdgid_mass(GenPart_pdgId[postfsrAntiLep][postfsrAntiLep_idx]) : wrem::get_pdgid_mass(GenPart_pdgId[postfsrLep][postfsrLep_idx])",
)
df = df.Define(
diff --git a/wremnants/unfolding_tools.py b/wremnants/unfolding_tools.py
index f8df87bf2..2923acc1f 100644
--- a/wremnants/unfolding_tools.py
+++ b/wremnants/unfolding_tools.py
@@ -346,16 +346,16 @@ def __init__(
self.poi_as_noi = poi_as_noi
self.unfolding_levels = unfolding_levels
- # self.qcdScaleByHelicity_helper =
+ def rebin_pt(edges):
+ # use 2 ptll bin for each ptVGen bin, except first and last
+ # first gen bin same size as reco bin, then 1 gen bin for 2 reco bins
+ new_edges = np.array([*edges[:2], *edges[3::2]])
+ if len(new_edges) % 2:
+ # in case it's an odd number of edges, last two bins are overflow
+ edges = edges[:-1]
+ return new_edges
- if self.add_helicity_axis:
- # helper to derive helicity xsec shape from event by event reweighting
- self.weightsByHelicity_helper_unfolding = helicity_utils.make_helicity_weight_helper(
- is_z=True,
- filename=f"{common.data_dir}/angularCoefficients/w_z_helicity_xsecs_maxFiles_m1_alphaSunfoldingBinning_helicity.hdf5",
- # filename=f"{common.data_dir}/angularCoefficients/w_z_helicity_xsecs_maxFiles_m1_nnpdf31_alphaSunfoldingBinning_helicity.hdf5",
- rebi_ptVgen=True,
- )
+ self.weightsByHelicity_helper_unfolding = None
self.unfolding_axes = {}
self.unfolding_cols = {}
@@ -368,13 +368,22 @@ def __init__(
level,
flow_y=self.poi_as_noi,
add_out_of_acceptance_axis=self.poi_as_noi,
- rebin_pt=True,
+ rebin_pt=rebin_pt,
)
self.unfolding_axes[level] = a
self.unfolding_cols[level] = c
self.unfolding_selections[level] = s
if self.add_helicity_axis:
+ if self.weightsByHelicity_helper_unfolding is None:
+ edges = [ax for ax in a if ax.name == "ptVGen"][0].edges
+ # helper to derive helicity xsec shape from event by event reweighting
+ self.weightsByHelicity_helper_unfolding = helicity_utils.make_helicity_weight_helper(
+ is_z=True,
+ filename=f"{common.data_dir}/angularCoefficients/w_z_helicity_xsecs_maxFiles_m1_alphaSunfoldingBinning_helicity.hdf5",
+ rebin_ptVgen_edges=edges,
+ )
+
for ax in a:
if ax.name == "acceptance":
continue
@@ -383,6 +392,7 @@ def __init__(
wbh_axis = self.weightsByHelicity_helper_unfolding.hist.axes[
ax.name.replace("Gen", "gen")
]
+
if any(ax.edges != wbh_axis.edges):
raise RuntimeError(
f"""