MakePostfitFigures.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6604e7db",
   "metadata": {},
   "outputs": [],
   "source": [
    "import uproot\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import mplhep as hep\n",
    "import hist\n",
    "\n",
    "plt.style.use(hep.style.CMS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6eca2491",
   "metadata": {},
   "outputs": [],
   "source": [
    "def add_uncertainty(hist, ax, newaxis, ratio=False):\n",
    "    opts = {\n",
    "        \"step\": \"post\",\n",
    "        \"label\": \"Uncertainty\",\n",
    "        \"hatch\": \"///\",\n",
    "        \"facecolor\": \"none\",\n",
    "        \"edgecolor\": (0, 0, 0, 0.5),\n",
    "        \"linewidth\": 0,\n",
    "        \"zorder\": 10.0,\n",
    "    }\n",
    "\n",
    "    if ratio:\n",
    "        down = np.ones(len(hist.counts())) - hist.errors() / hist.counts()\n",
    "        up = np.ones(len(hist.counts())) + hist.errors() / hist.counts()\n",
    "    else:\n",
    "        down = hist.counts() - hist.errors()\n",
    "        up = hist.counts() + hist.errors()\n",
    "    up = np.nan_to_num(up)[: len(newaxis) - 1]\n",
    "    down = np.nan_to_num(down)[: len(newaxis) - 1]\n",
    "    ax.fill_between(x=newaxis, y1=np.r_[down, down[-1]], y2=np.r_[up, up[-1]], **opts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ba901658",
   "metadata": {},
   "outputs": [],
   "source": [
    "fitDiagnostics = uproot.open(\"Fitting/fitDiagnostics.prefit-asimov.root\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3cccb2b8",
   "metadata": {},
   "outputs": [],
   "source": [
    "shapetype = \"shapes_prefit\"\n",
    "fitDiagnostics[shapetype].keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ced38b9e",
   "metadata": {},
   "outputs": [],
   "source": [
    "processes = [\n",
    "    \"ttgamma\",\n",
    "    \"NonPrompt\",\n",
    "    \"MisID\",\n",
    "    \"WG\",\n",
    "    \"ZG\",\n",
    "    \"other\",\n",
    "]\n",
    "\n",
    "processes.reverse()\n",
    "\n",
    "colors = {\n",
    "    \"ttgamma\": \"\",  # FIXME\n",
    "}\n",
    "\n",
    "labels = {\n",
    "    \"ttgamma\": r\"$t\\bar{t}\\gamma$\",\n",
    "    \"NonPrompt\": r\"nonprompt $\\gamma$\",\n",
    "    \"MisID\": r\"MisID $\\gamma$\",\n",
    "    \"WG\": r\"$W\\gamma$\",\n",
    "    \"ZG\": r\"$Z\\gamma$\",\n",
    "    \"other\": r\"other\",\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c439442a",
   "metadata": {},
   "outputs": [],
   "source": [
    "groupingCategory = {\n",
    "    \"NonPrompt\": [3j, 4j],\n",
    "    \"MisID\": [2j],\n",
    "    \"Prompt\": [1j],\n",
    "}\n",
    "\n",
    "groupingFineCategory = {\n",
    "    \"Genuine Photon\": [1j],\n",
    "    \"Misidentified Electron\": [2j],\n",
    "    \"Hadronic Photon\": [3j],\n",
    "    \"Hadronic Fake\": [4j],\n",
    "}\n",
    "\n",
    "groupingMCDatasets = {\n",
    "    \"ZG\": [\n",
    "        \"ZGamma_01J_5f_lowMass\",\n",
    "    ],\n",
    "    \"WG\": [\n",
    "        \"WGamma\",\n",
    "    ],\n",
    "    \"other\": [\n",
    "        \"TTbarPowheg_Dilepton\",\n",
    "        \"TTbarPowheg_Semilept\",\n",
    "        \"TTbarPowheg_Hadronic\",\n",
    "        \"W2jets\",\n",
    "        \"W3jets\",\n",
    "        \"W4jets\",\n",
    "        \"DYjetsM50\",\n",
    "        \"ST_s_channel\",\n",
    "        \"ST_tW_channel\",\n",
    "        \"ST_tbarW_channel\",\n",
    "        \"ST_tbar_channel\",\n",
    "        \"ST_t_channel\",\n",
    "        \"TTWtoLNu\",\n",
    "        \"TTWtoQQ\",\n",
    "        \"TTZtoLL\",\n",
    "        \"GJets_HT200To400\",\n",
    "        \"GJets_HT400To600\",\n",
    "        \"GJets_HT600ToInf\",\n",
    "        \"ZZ\",\n",
    "        \"WZ\",\n",
    "        \"WW\",\n",
    "        \"TGJets\",\n",
    "    ],\n",
    "    \"ttgamma\": [\n",
    "        \"TTGamma_Dilepton\",\n",
    "        \"TTGamma_SingleLept\",\n",
    "        \"TTGamma_Hadronic\",\n",
    "    ],\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "96abb647",
   "metadata": {},
   "outputs": [],
   "source": [
    "# axes hardcoded from save_to_root.py\n",
    "inputshapes = {\n",
    "    \"m3\": \"RootFiles/M3_Output.root\",\n",
    "    \"misidEle\": \"RootFiles/MisID_Output_electron.root\",\n",
    "    \"Vgamma\": \"RootFiles/MisID_Output_muon.root\",\n",
    "    \"chIso\": \"RootFiles/Isolation_Output.root\",\n",
    "}\n",
    "newaxes = {\n",
    "    channel: uproot.open(fn)[\"data_obs\"].axes[0].edges()\n",
    "    for channel, fn in inputshapes.items()\n",
    "}\n",
    "axeslabels = {\n",
    "    \"m3\": r\"$M_3$ [GeV]\",\n",
    "    \"misidEle\": r\"$m_{e\\gamma}$ [GeV]\",\n",
    "    \"Vgamma\": r\"$m_{\\mu\\gamma}$ [GeV]\",\n",
    "    \"chIso\": \"Charged Hadron Isolation [GeV]\",\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "db2a0b85",
   "metadata": {},
   "outputs": [],
   "source": [
    "channel = \"Vgamma\"\n",
    "dc_x = fitDiagnostics[f\"{shapetype}/{channel}\"]\n",
    "h = dc_x[\"data\"]\n",
    "\n",
    "fig, (ax, rax) = plt.subplots(\n",
    "    2,\n",
    "    1,\n",
    "    figsize=(10, 10),\n",
    "    gridspec_kw={\"height_ratios\": (3, 1), \"hspace\": 0.05},\n",
    "    sharex=True,\n",
    ")\n",
    "newaxis = newaxes[channel]\n",
    "xlabel = axeslabels[channel]\n",
    "\n",
    "hep.cms.label(\n",
    "    \"Preliminary\",\n",
    "    data=True,\n",
    "    lumi=35.9,\n",
    "    loc=0,\n",
    "    ax=ax,\n",
    ")\n",
    "\n",
    "hep.histplot(\n",
    "    [dc_x[x].counts()[: len(newaxis) - 1] for x in processes],\n",
    "    bins=newaxis,\n",
    "    # w2=[ dc_sr[x].errors() for x in processes ],  # not needed\n",
    "    histtype=\"fill\",\n",
    "    stack=True,\n",
    "    label=[labels[x] for x in processes],\n",
    "    # color=['red', 'green'],\n",
    "    ax=ax,\n",
    ")\n",
    "\n",
    "# this is for the data points\n",
    "wq = hep.histplot(\n",
    "    dc_x[\"data\"].values()[1][: len(newaxis) - 1],\n",
    "    bins=newaxis,\n",
    "    w2=dc_x[\"data\"].values()[1][: len(newaxis) - 1],\n",
    "    histtype=\"errorbar\",\n",
    "    stack=False,\n",
    "    label=\"Observation\",\n",
    "    color=\"black\",\n",
    "    ax=ax,\n",
    ")\n",
    "\n",
    "ratio_val = dc_x[\"data\"].values()[1] / dc_x[\"total\"].counts()\n",
    "ratio_err_hep = np.sqrt(dc_x[\"data\"].values()[1]) / dc_x[\"total\"].counts()\n",
    "\n",
    "# this is for the ratio plot\n",
    "hep.histplot(\n",
    "    ratio_val[: len(newaxis) - 1],\n",
    "    bins=newaxis,\n",
    "    yerr=ratio_err_hep[: len(newaxis) - 1],\n",
    "    histtype=\"errorbar\",\n",
    "    color=\"black\",\n",
    "    ax=rax,\n",
    ")\n",
    "\n",
    "rax.set_ylim(0, 1.99)\n",
    "rax.set_xlabel(xlabel)\n",
    "rax.set_ylabel(r\"Data/Pred.\")\n",
    "ax.set_ylabel(r\"Events\")\n",
    "\n",
    "add_uncertainty(dc_x[\"total\"], ax, newaxis)\n",
    "add_uncertainty(dc_x[\"total\"], rax, newaxis, ratio=True)\n",
    "\n",
    "# ax.legend()\n",
    "ax.legend(prop={\"size\": 15}, bbox_to_anchor=(1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8b7d66dd-353c-463a-928a-f4c9e4c4cf63",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}