Projektpraktikum/notebooks/trackinglosses_photons.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import uproot\t\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from mpl_toolkits import mplot3d\n",
    "import awkward as ak\n",
    "from scipy.optimize import curve_fit\n",
    "import mplhep\n",
    "mplhep.style.use([\"LHCbTex2\"])\n",
    "\n",
    "plt.rcParams[\"savefig.dpi\"] = 600\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "file = uproot.open(\n",
    "    \"/work/cetin/Projektpraktikum/trackinglosses_B_photon_cuts.root\")\n",
    "\n",
    "# selektiere nur elektronen von B->K*ee\n",
    "allcolumns = []\n",
    "for i in range(11):\n",
    "    allcolumns.append(file[\"Tree\" + str(i)].arrays())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<pre>{oneCut_event_id: 1,\n",
       " oneCut_lost: False,\n",
       " oneCut_rad_length_frac: 0.148,\n",
       " oneCut_energy: 1.28e+04,\n",
       " noneCut_brem_photons_pe: 1,\n",
       " oneCut_brem_photons_pe: [7.42e+03],\n",
       " noneCut_brem_vtx_x: 1,\n",
       " oneCut_brem_vtx_x: [-3.61],\n",
       " noneCut_brem_vtx_z: 1,\n",
       " oneCut_brem_vtx_z: [35.6],\n",
       " oneCut_photon_length: 1}\n",
       "------------------------------------------\n",
       "type: {\n",
       "    oneCut_event_id: int64,\n",
       "    oneCut_lost: bool,\n",
       "    oneCut_rad_length_frac: float64,\n",
       "    oneCut_energy: float64,\n",
       "    noneCut_brem_photons_pe: int32,\n",
       "    oneCut_brem_photons_pe: var * float64,\n",
       "    noneCut_brem_vtx_x: int32,\n",
       "    oneCut_brem_vtx_x: var * float64,\n",
       "    noneCut_brem_vtx_z: int32,\n",
       "    oneCut_brem_vtx_z: var * float64,\n",
       "    oneCut_photon_length: int64\n",
       "}</pre>"
      ],
      "text/plain": [
       "<Record {oneCut_event_id: 1, ...} type='{oneCut_event_id: int64, oneCut_los...'>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "allcolumns[1][1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cutdict():\n",
    "    basedict = {\n",
    "        \"0\": {},\n",
    "        \"1\": {},\n",
    "        \"2\": {},\n",
    "        \"3\": {},\n",
    "        \"4\": {},\n",
    "        \"5\": {},\n",
    "        \"6\": {},\n",
    "        \"7\": {},\n",
    "        \"8\": {},\n",
    "        \"9\": {},\n",
    "        \"10\": {},\n",
    "    }\n",
    "\n",
    "    basedict[\"0\"] = \"no\"\n",
    "    basedict[\"1\"] = \"one\"\n",
    "    basedict[\"2\"] = \"two\"\n",
    "    basedict[\"3\"] = \"three\"\n",
    "    basedict[\"4\"] = \"four\"\n",
    "    basedict[\"5\"] = \"five\"\n",
    "    basedict[\"6\"] = \"six\"\n",
    "    basedict[\"7\"] = \"seven\"\n",
    "    basedict[\"8\"] = \"eight\"\n",
    "    basedict[\"9\"] = \"nine\"\n",
    "    basedict[\"10\"] = \"ten\"\n",
    "\n",
    "    return basedict\n",
    "\n",
    "\n",
    "Cuts = cutdict()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# electrons = []\n",
    "# for jcut in range(11):\n",
    "\n",
    "jcut = 4  # cut 0.2*E\n",
    "\n",
    "energy_emissions = ak.ArrayBuilder()\n",
    "\n",
    "for jelec in range(ak.num(allcolumns[jcut], axis=0)):\n",
    "    energy_emissions.begin_record()\n",
    "    energy_emissions.field(\"lost\").boolean(\n",
    "        allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"lost\"])\n",
    "    energy_emissions.field(\"rad_length_frac\").real(\n",
    "        allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"rad_length_frac\"])\n",
    "    energy_emissions.field(\"energy\").real(\n",
    "        allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"energy\"])\n",
    "\n",
    "    tmp_velo = 0\n",
    "    tmp_richut = 0\n",
    "    tmp_neither = 0\n",
    "    tmp_velo_length = 0\n",
    "    tmp_richut_length = 0\n",
    "    tmp_neither_length = 0\n",
    "\n",
    "    for jphoton in range(\n",
    "            ak.num(\n",
    "                allcolumns[jcut][jelec][Cuts[str(jcut)] + \"Cut_\" +\n",
    "                                        \"brem_photons_pe\"],\n",
    "                axis=0,\n",
    "            )):\n",
    "        if (allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"brem_vtx_z\",\n",
    "                             jphoton] <= 770):\n",
    "            tmp_velo += allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" +\n",
    "                                         \"brem_photons_pe\", jphoton]\n",
    "            tmp_velo_length += 1\n",
    "        elif (allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"brem_vtx_z\",\n",
    "                               jphoton]\n",
    "              > 770) and (allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" +\n",
    "                                           \"brem_vtx_z\", jphoton] <= 2700):\n",
    "            tmp_richut += allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" +\n",
    "                                           \"brem_photons_pe\", jphoton]\n",
    "            tmp_richut_length += 1\n",
    "        else:\n",
    "            tmp_neither += allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" +\n",
    "                                            \"brem_photons_pe\", jphoton]\n",
    "            tmp_neither_length += 1\n",
    "\n",
    "    energy_emissions.field(\"velo_length\").integer(tmp_velo_length)\n",
    "    energy_emissions.field(\"velo\").real(tmp_velo)\n",
    "\n",
    "    energy_emissions.field(\"rich_length\").integer(tmp_richut_length)\n",
    "    energy_emissions.field(\"rich\").real(tmp_richut)\n",
    "\n",
    "    energy_emissions.field(\"neither_length\").integer(tmp_neither_length)\n",
    "    energy_emissions.field(\"downstream\").real(tmp_neither)\n",
    "\n",
    "    energy_emissions.field(\"photon_length\").integer(tmp_richut_length +\n",
    "                                                    tmp_velo_length)\n",
    "\n",
    "    if ((tmp_velo == 0) and (tmp_richut == 0)\n",
    "            or (allcolumns[jcut][jelec, Cuts[str(jcut)] + \"Cut_\" + \"energy\"] -\n",
    "                tmp_velo < 3000)):\n",
    "        energy_emissions.field(\"quality\").integer(0)\n",
    "    else:\n",
    "        energy_emissions.field(\"quality\").integer(1)\n",
    "\n",
    "    energy_emissions.end_record()\n",
    "\n",
    "energy_emissions = ak.Array(energy_emissions)\n",
    "# electrons.append(energy_emissions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABMEAAAOWCAYAAADm+IXMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAACq5UlEQVR4nOz9zW8b6Z3vf3/KafyQHw5GKsmDAxgWkGaxfTazGKsoB8jyWGQn24lJqXeziEW2J4sbUGLR8h/QNtndBu5FJiblXgT35kik3duMSbmXAdJk2b2YVZolNWBDwAGaKjM4mGAwUd0LHzKiRD2SYhWl9wsgIpJXXfUtPXSsj77XVYbv+74AAAAAAACAc+xS0AUAAAAAAAAAZ40QDAAAAAAAAOceIRgAAAAAAADOPUIwAAAAAAAAnHuEYAAAAAAAADj3CMEAAAAAAABw7hGCAQAAAAAA4NwjBAMAAAAAAMC5RwgGAAAAAACAc48QDAAAAAAAAOceIRgAAAAAAADOPUIwAAAAAAAAnHuEYAAAAAAAADj3CMEAAAAAAABw7hGCAQAAAAAA4NwjBAMAAAAAAMC5RwgGAAAAAACAc48QDAAAAAAAAOceIRgAAAAAAADOPUIwAAAAAAAAnHuEYAAAAAAAADj3CMEAAAAAAABw7hGCASfkeZ7K5XLQZQAAAAAAgBMgBAN6cBxHhmH0fExMTMg0zaBLBAAAAAAAJ0AIBvTw4MGDA9+zLEvxeHyI1QAAAAAAgH4Zvu/7QRcBhInruopGoyqVSj07vizLkmVZwy8MAAAAAACcGiEYsEcmk1Gz2VSpVAq6FAAAAAAAMCDvBV0AECae56lYLBKAAQAAAABwzrAnGLBLey+wVColwzAUi8WUz+fleV6whQEAAAAAgL4QgqEvrusqkUioXC73NU8+n1cikdDExIQMw1A0GlUqlVK1Wh1QpcevYzfHcZTNZjUxMaFsNjvUWgAAAAAAwOAQguFUPM9TKpVSNBpVtVpVs9k81TzVarUrYCqVSmo0GsrlcnIcR4lEQolEYmidWNvb26rX6yoUCkomk13vtYM6AAAAAAAwetgYHyfieZ4ePHiwr2OqUCgonU6faK5qtdoJldLptAqFwr4xsVhMjuPIsizV6/Wed2s8S72uN5fLaWlpaah1AAAAAACA/tAJhmPL5/OdUKpf7U4ySbIsq2cAJqmzQb3rup3xw2SapnK5XFcAx7JIAAAAAABGDyEYjsVxHMXjcTUaDVUqlQNDq+NKpVKdJY6HhUqWZXWWJVarVRWLxb7Oe1q2bSuXy3WeD3uvMgAAAAAA0B9CMByLbduybbvzfGZm5tRzua7bFSLNzc0dOn5+fr7zca/AzHEcGYZxqsdJusvS6XSnG4y7RQIAAAAAMFreC7oAnJ7rupLedUv1q1qtKh6PH3t8P3tz7e6oisfjR861e4N6z/NULpe7XjNNc98m9sd10o3uZ2ZmVK1Wh743GQAAAAAA6A8h2IhyXVexWEySVK/X+wrCUqmUyuXy0DZ8372kcXd32WEsy+qEfqurq12hl2VZnb3DhqWfTjgAAAAAADB8LIccUa7ryvM8eZ6nWCzWCYhOqh2ASdLXX389yBJ72rup/o0bN4513O6wrF1vEFzXlWmadIIBAAAAADBiCMFGVDwe73Q/nTYI2x2AJZPJoXRT7d1Q/rgdbHvHDeIOlSflOI5c19Xy8vLQzw0AAAAAAPpDCDbCdgdXJw3CggjApP3dZsftqLp8+XLX81qtNqiSOrLZrBKJxIGdZgsLC0omk0NZMgoAAAAAAAaLEGzEnSYICyoAk7SvttN2gjUajYHV1JbP51WtVpVKpRSLxVStVjsb8Uej0a7uOwAAAAAAMFoIwc6BkwRhQQZg0v4Q7LQ8zxvIPLvV6/XO3Sodx1EqldLs7Ky+/vprVSqVrrtaAgAAAACA0cLdIc+JdqCVSqU6Qdjeu0YGHYBJpw+v9i6bbDab/Rezh23bqlQqA58XAAAAAAAEj06wc+SwjrBMJhN4ADZIZ9EJBgAAAAAAzi9CsHNmbxCWSCSUSqVULBb3vR+E426EP6x5AAAAAADAxUAIdg7tDrpc1w1VB9jk5GSo5gEAAAAAABcDIdg5lUwmFY/Hu15bWVkJqJq/OW0H197lj3SCAQAAAACAkyAEO6dSqZSq1WrXa7FYLPC9tGZmZrqeH7eevRvhR6PRQZUEAAAAAAAuAEKwc+igu0C6rht4EBaLxbqetzfuP0qj0eh6vrfLDQAAAAAA4DCEYOdMrwAsTEHY3k6w44Zgu+s1TVOWZQ2yLAAAAAAAcM4Rgp0jB3WA7X0eZBBm23bXfl5ff/31sY6r1Wqdj/cGaQAAAAAAAEchBDsnMpnMkXeBDEsQNjc31/nYcZxjHbN7XDabHXhNAAAAAADgfCMEOwcymYyKxaKkgwOwtjAEYZlMpvPx3s37e9k9xrIs9gMDAAAAAAAnRgg24k4SgLUNIgjrJzizbbsryGp3sB1k9zXRBQYAAAAAAE6DEGyEnSYAa+s3CNu7of1JQ7FCodD5+MGDBweO8zyvc43xeFzpdPpE5wEAAAAAAJAIwUZWsVg8dQDWtjcIm52dPdZxnuft68haXV090bkty+qc23Ec5fP5nuPaNZmmeaprBAAAAAAAkCTD930/6CJwOolEYiDhULlc1sLCgur1uizL6jnG8zwtLCzI87xD9/GKx+MyTVPLy8uybfvIc1erVaVSKXmep2QyqeXlZVmWpVqtpmw2K8dxZNu21tfXu+4qCQAAAAAAcBKEYAiFfD6v1dVVua4rz/NkmqZmZmaUyWSUTCaDLg8AAAAAAIw4QjAAAAAAAACce+wJBgAAAAAAgHPvvaALwMn9t//23/SXv/xFP/jBD/Tf//t/D7ocAAAAAACAQ/3v//2/9de//lU//OEP9X/+z/8JpAaWQ46gS5cuiS8bAAAAAAAYNYZhaGdnJ5Bz0wk2gnYHYFevXg2wktPb2trSzs6OLl26pCtXrgRdzomNev0S1xAWXEM4cA3hMOrXMOr1S1xDWHAN4cA1hMOoX8Oo1y9xDWFxHq7hzZs3khRoUw8h2Ai6dOlS55v/9evXQZdzKlNTU3rz5o2uXLkyktcw6vVLXENYcA3hwDWEw6hfw6jXL3ENYcE1hAPXEA6jfg2jXr/ENYTFebiGH/zgB50sIygjF4Jtbm7KcRw1m015nqfvv/9eknT58mWZpinLsjQzM6OxsbGAKwUAAAAAAEBYjEQI9uzZM1UqFa2trcnzvGMdY1mW4vG4UqmUbt68ebYFAgAAAAAAINRCHYI9efJEuVxOrutK6l43ahhGz2PaYxqNhlzXVbFYlGVZunfvnn7xi1+cfdEAAAAAAAAIneAWYh7ixYsXunbtmjKZTM8ArP2816PXGNd1lU6n9T/+x//QN998M7TrAAAAAAAAQDiErhNseXlZ+Xx+X6AVj8dl27YuX74sy7JkmqYkaXJyUpLUbDYlSZ7nyXVdff/993JdV9VqtbOE8ttvv5Vt28pms/rkk0+Gdk0AAAAAAAAIVmhCsFarpdnZWdXrdUnvQq9EIqF4PK7p6em+5n779q1WV1dVLBblOI5yuZyq1aqq1Sob6AMAAAAAAFwAoVgO+erVK0UiEdXrdaXTaTUaDT1//lx3797tOwCTpPHxcaXTadVqNX377be6efOmarWaZmZm9Oc//3kAVwAAAAAAAIAwC7wT7OXLl4rFYrJtW7VaTZFI5EzPZ1mWKpWKqtWqUqmUYrGYarXaSHWE/d3f/Z3evn2rv/u7vwu6lAtrcXFRrVZrpL5v9uIawuE8XMN5cB6+DlxD8Ea9/vPiPHwduIZwOA/XcB6M+tdh1OuXuAYMThiyDMPfu/nWEG1sbCgWiykej2ttbW3o5/c8T7FYTNeuXdPvf//7oZ//tKampvTmzRtdvXpVr1+/DrqcUzkP1wDgHX6egfODn2fgfOFnGjg/zsPPcxiuIdBOsEQiobm5OT1+/DiQ85umqXq9rpmZGX3++ef61a9+FUgdp7W1taWpqalDxywuLmpxcXFIFQEAAAAAgIvm0aNHevTo0aF
      "text/plain": [
       "<Figure size 1200x900 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "velo_found = ak.to_numpy(\n",
    "    energy_emissions[(~energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"velo\"])\n",
    "rich_found = ak.to_numpy(\n",
    "    energy_emissions[(~energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"rich\"])\n",
    "energy_found = ak.to_numpy(\n",
    "    energy_emissions[(~energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"energy\"])\n",
    "\n",
    "velo_lost = ak.to_numpy(\n",
    "    energy_emissions[(energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"velo\"])\n",
    "rich_lost = ak.to_numpy(\n",
    "    energy_emissions[(energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"rich\"])\n",
    "energy_lost = ak.to_numpy(\n",
    "    energy_emissions[(energy_emissions.lost)\n",
    "                     & (energy_emissions.quality == 1)][\"energy\"])\n",
    "\n",
    "diff_found = velo_found - rich_found\n",
    "diff_lost = velo_lost - rich_lost\n",
    "\n",
    "xlim = 20000\n",
    "nbins = 60\n",
    "\n",
    "plt.hist(\n",
    "    diff_lost,\n",
    "    bins=nbins,\n",
    "    density=True,\n",
    "    alpha=0.5,\n",
    "    histtype=\"bar\",\n",
    "    color=\"#F05342\",\n",
    "    label=\"lost\",\n",
    "    range=[-xlim, xlim],\n",
    ")\n",
    "plt.hist(\n",
    "    diff_found,\n",
    "    bins=nbins,\n",
    "    density=True,\n",
    "    alpha=0.5,\n",
    "    histtype=\"bar\",\n",
    "    color=\"#2A9D8F\",  # \"#107E7D\",\n",
    "    label=\"found\",\n",
    "    range=[-xlim, xlim],\n",
    ")\n",
    "# plt.xlim(-20000, 20000)\n",
    "# plt.yscale(\"log\")\n",
    "# plt.title(\"emitted energy difference\")\n",
    "plt.xlabel(r\"$\\Delta E_{VELO} - \\Delta E_{RICH1+UT}$ [MeV]\")\n",
    "plt.ylabel(\"Number of Tracks (normalised)\")\n",
    "plt.legend(loc=\"best\")\n",
    "mplhep.lhcb.text(\"Simulation\")\n",
    "# plt.show()\n",
    "plt.savefig(\n",
    "    \"/work/cetin/Projektpraktikum/thesis/emitted_energy_difference.pdf\",\n",
    "    format=\"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,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "nbins = 6\n",
    "quality_cut = electrons[jcut].quality != -1\n",
    "\n",
    "### all split in velo and rich\n",
    "\n",
    "fig, axs = plt.subplots(3, 3, figsize=(15, 6))\n",
    "ax = axs.ravel()\n",
    "for jcut, ax in enumerate(ax):\n",
    "ax.hist(\n",
    "ak.to_numpy(electrons[jcut][quality_cut][\"velo_length\"]),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"darkorange\",\n",
    "histtype=\"bar\",\n",
    "label=\"velo\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.hist(\n",
    "ak.to_numpy(electrons[jcut][quality_cut][\"rich_length\"]),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"blue\",\n",
    "histtype=\"bar\",\n",
    "label=\"rich\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.set_xlim(0, nbins)\n",
    "ax.set_ylim(0, 1)\n",
    "ax.set_title(\"Photon Cut: \" + str(np.round(jcut \\* 0.05, 2)) + f\"$E_0$\")\n",
    "ax.set_xlabel(\"number of photons\")\n",
    "ax.set_ylabel(\"a.u.\")\n",
    "plt.suptitle(\"number of photons in velo and rich\")\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "\n",
    "### found\n",
    "\n",
    "fig, axs = plt.subplots(3, 3, figsize=(15, 6))\n",
    "ax = axs.ravel()\n",
    "for jcut, ax in enumerate(ax):\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut]~(electrons[jcut].lost) & quality_cut][\"velo_length\"]\n",
    "),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"darkorange\",\n",
    "histtype=\"bar\",\n",
    "label=\"velo\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut]~(electrons[jcut].lost) & quality_cut][\"rich_length\"]\n",
    "),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"blue\",\n",
    "histtype=\"bar\",\n",
    "label=\"rich\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.set_xlim(0, nbins)\n",
    "ax.set_ylim(0, 1)\n",
    "ax.set_title(\"Photon Cut: \" + str(np.round(jcut \\* 0.05, 2)) + f\"$E_0$\")\n",
    "ax.set_xlabel(\"number of photons\")\n",
    "ax.set_ylabel(\"a.u.\")\n",
    "plt.suptitle(\"number of photons of found electrons\")\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "\n",
    "### lost\n",
    "\n",
    "fig, axs = plt.subplots(3, 3, figsize=(15, 6))\n",
    "ax = axs.ravel()\n",
    "for jcut, ax in enumerate(ax):\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut](electrons[jcut].lost) & quality_cut][\"velo_length\"]\n",
    "),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"darkorange\",\n",
    "histtype=\"bar\",\n",
    "label=\"velo\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut](electrons[jcut].lost) & quality_cut][\"rich_length\"]\n",
    "),\n",
    "bins=nbins,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"blue\",\n",
    "histtype=\"bar\",\n",
    "label=\"rich\",\n",
    "range=[0, nbins],\n",
    ")\n",
    "ax.set_xlim(0, nbins)\n",
    "ax.set_ylim(0, 1)\n",
    "ax.set_title(\"Photon Cut: \" + str(np.round(jcut \\* 0.05, 2)) + f\"$E_0$\")\n",
    "ax.set_xlabel(\"number of photons\")\n",
    "ax.set_ylabel(\"a.u.\")\n",
    "plt.suptitle(\"number of photons of lost electrons\")\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "quality_cut = electrons[jcut].quality != -1\n",
    "\n",
    "### all split in lost and found\n",
    "\n",
    "fig, axs = plt.subplots(3, 3, figsize=(15, 6))\n",
    "ax = axs.ravel()\n",
    "for jcut, ax in enumerate(ax):\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut](electrons[jcut].lost) & (quality_cut)][\"photon_length\"]\n",
    "),\n",
    "bins=10,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"darkorange\",\n",
    "histtype=\"bar\",\n",
    "label=\"lost\",\n",
    "range=[0, 10],\n",
    ")\n",
    "ax.hist(\n",
    "ak.to_numpy(\n",
    "electrons[jcut](~electrons[jcut].lost) & (quality_cut)][\"photon_length\"]\n",
    "),\n",
    "bins=10,\n",
    "density=True,\n",
    "alpha=0.5,\n",
    "color=\"blue\",\n",
    "histtype=\"bar\",\n",
    "label=\"found\",\n",
    "range=[0, 10],\n",
    ")\n",
    "ax.set_xlim(0, 10) # ax.set_ylim(0,1) # ax.set_yscale('log')\n",
    "ax.set_title(\"Photon Cut: \" + str(np.round(jcut \\* 0.05, 2)) + f\"$E_0$\")\n",
    "ax.set_xlabel(\"number of photons\")\n",
    "ax.set_ylabel(\"a.u.\")\n",
    "plt.suptitle(\"number of photons in lost and found\")\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "tuner",
   "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}