Projektpraktikum/D_tasks.ipynb

825 lines
1.3 MiB
Plaintext
Raw Normal View History

2023-09-25 16:36:13 +02:00
{
"cells": [
{
"cell_type": "code",
"execution_count": 60,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [],
"source": [
"import uproot\n",
"import numpy as np\n",
"import sys\n",
"import os\n",
"import matplotlib\n",
"import matplotlib.pyplot as plt\n",
"from mpl_toolkits import mplot3d\n",
"import itertools\n",
"import awkward as ak\n",
"from scipy.optimize import curve_fit\n",
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 61,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"51"
]
},
"execution_count": 61,
"metadata": {},
"output_type": "execute_result"
}
],
2023-09-25 16:36:13 +02:00
"source": [
"file = uproot.open(\"tracking_losses_ntuple_Dst0ToD0EE.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\")\n",
"\n",
"\n",
"#selektiere nur elektronen von D*0->D0ee und nur solche mit einem momentum von unter 5 GeV \n",
"allcolumns = file.arrays()\n",
"found = allcolumns[(allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p < 5e3)] #D: 2591\n",
"lost = allcolumns[(allcolumns.isElectron) & (allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p < 5e3)] #D: 1908\n",
"\n",
"#ak.num(lost, axis=0)\n",
"ak.num(allcolumns[(allcolumns.fromPairProd) & (allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal)],axis=0)"
2023-09-25 16:36:13 +02:00
]
},
{
"cell_type": "code",
"execution_count": 62,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.5759057568348522"
]
},
"execution_count": 62,
2023-09-25 16:36:13 +02:00
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def t_eff(found, lost):\n",
" sel = found[\"energy\"]\n",
" des = lost[\"energy\"]\n",
" return ak.count(sel,axis=None)/(ak.count(sel,axis=None)+ak.count(des,axis=None))\n",
"\n",
"t_eff(found, lost)"
]
},
{
"cell_type": "code",
"execution_count": 63,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.7960893854748603"
]
},
"execution_count": 63,
2023-09-25 16:36:13 +02:00
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#finden wir die elektronen die keine bremsstrahlung gemacht haben mit hoher effizienz?\n",
"nobrem_found = found[found[\"brem_photons_pe_length\"]==0]\n",
"nobrem_lost = lost[lost[\"brem_photons_pe_length\"]==0]\n",
"\n",
"\"\"\"\n",
"die effizienz mit der wir elektronen finden, die keine bremsstrahlung gemacht haben, ist nicht besonders gut, aber trotzdem besser als\n",
"für alle elektronen.\n",
"Auch hier handelt es sich um eine recht geringe sample size (~350)\n",
"\"\"\"\n",
"\n",
"t_eff(nobrem_found, nobrem_lost)\n",
"#ak.num(nobrem_lost, axis=0)"
]
},
{
"cell_type": "code",
"execution_count": 64,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.5568703211784594"
]
},
"execution_count": 64,
2023-09-25 16:36:13 +02:00
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#wie viel energie relativ zur anfangsenergie verlieren die elektronen durch bremstrahlung und hat das einen einfluss darauf ob wir sie finden oder nicht?\n",
"brem_found = found[found[\"brem_photons_pe_length\"]!=0]\n",
"energy_found = ak.to_numpy(brem_found[\"energy\"])\n",
"eph_found = ak.to_numpy(ak.sum(brem_found[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
"energyloss_found = eph_found/energy_found\n",
"\n",
"\n",
"brem_lost = lost[lost[\"brem_photons_pe_length\"]!=0]\n",
"energy_lost = ak.to_numpy(brem_lost[\"energy\"])\n",
"eph_lost = ak.to_numpy(ak.sum(brem_lost[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
"energyloss_lost = eph_lost/energy_lost\n",
"\n",
"t_eff(brem_found,brem_lost)"
]
},
{
"cell_type": "code",
"execution_count": 65,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"mean energyloss relative to initial energy (found): 0.4422334194682742\n",
"mean energyloss relative to initial energy (lost): 0.5885317187224427\n"
]
}
],
"source": [
"mean_energyloss_found = ak.mean(energyloss_found)\n",
"mean_energyloss_lost = ak.mean(energyloss_lost)\n",
"print(\"mean energyloss relative to initial energy (found): \", mean_energyloss_found)\n",
"print(\"mean energyloss relative to initial energy (lost): \", mean_energyloss_lost)"
]
},
{
"cell_type": "code",
"execution_count": 66,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHMCAYAAAAplYnpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABLzElEQVR4nO3de1wU9f4/8NcCy00FFQVR8X5J8BqUomlagUlH7Ry/1flZpqWV3zx5O4ZidoIuX+p0U8skjSQ1y1OoWZrBKRe8lgKWVzTjBCmGVy4Sy8LO7w8P6667izOzl1mG1/Px4BHz2c985j1vxt13M5/Z0QiCIICIiIhIJbyUDoCIiIjImVjcEBERkaqwuCEiIiJVYXFDREREqsLihoiIiFSFxQ0RERGpCosbIiIiUhUWN0RERKQqLG6IiIhIVVjcEBERkaqwuCEiIiJVYXFDRIrYvn07NBqN3Z/09HSnbKe+vh6hoaF4++233b5tIlKGj9IBEFHzlJ+fDwD44osvEBoaavV6ZGSkU7aTm5uL8+fP4y9/+Yvbt01EymBxQ0SKyM/PR1BQEMaPHw+NRuOy7Xz++eeIiYlB165d3b5tIlIGL0sRkSLy8vIwaNAglxYXgiBg8+bNmDRpktu3TUTKYXFDRG538eJFFBcXY8CAAairq7P6EQTBKdvZu3cvSktLLYobd22biJTD4oaI3K5hzst7770HrVZr9XPs2LGbjiEIAoKCgnDhwgW7fT7//HMMGDAAvXv3lrzt8+fP47777kOLFi3Qp08fZGdnO7LLRORGnHNDRG6Xl5cHANi0aRMiIiKsXhczoffUqVNo37492rVrZ7fPpk2b8Pjjj8va9qxZs9ChQwecP38e//73v/Hggw/i559/RkhIyE1jIyJlsbghIrfLz8+Hv78/JkyYAG9vb7v93n//fXz22WfQarX4/vvv0bFjR/zrX/9CZGQk8vPzMXjwYDz11FPYuHEjOnbsiK1bt6JXr14AgB9++AHFxcVW823EbLuqqgpbtmzB6dOnERgYiAkTJmDQoEH44osvrIolIvI8vCxFRG6Xn5+P/v37N1rYAMCRI0dw4MABzJ07F7///juGDh2KJUuWAAAKCgpw8OBBPPbYY7h06RJuvfVWrF692rRuZmYm+vTpg/79+0ve9qlTp9CyZUuLMzsDBgzA0aNH5ewuEbkZixsicqvy8nL88ssvGDRo0E37HjlyBEuWLMHYsWOh1Wrx8MMP4+TJkwCuFTcpKSkYNmwYvLy80LNnT4vJwJmZmVZnbcRuu6qqCkFBQRZtQUFBqKqqErubRKQgXpYiIrfKz8+HIAho0aIF9u/fb/V6p06dTGdMjh49avFtwWVlZaY5NgUFBVi/fr3ptaNHjyIhIQEAcOjQIZw+fdrmJSkx227ZsiUqKiosXquoqEDLli1l7jURuRPP3BCRWzXcrbR8+XLExsZa/WRlZQG4VsicP3/e4huEN2/ejHHjxqGkpAQ+Pj4Wr/30008YOHAggGtnbbp27Yro6GhZ2+7duzeqqqrw22+/mdY9cuQIoqKiXJARInI2jcAvdSAiD/Ttt99i7NixSEtLw9SpU7F27Vq89NJL+PHHH5GTk4MVK1bgm2++AQBUV1ejdevWKC8vR0BAACIjIzFu3Di8+eabsrf/wAMPIDg4GO+88w6+/fZbTJkyBadOnWr07iwi8gy8LEVEHunIkSN4/PHH8emnn2L+/PmIjo5GdnY2goODUVBQYDpLA1y7JNWzZ08EBAQAgKjvybmZ9957D1OnTkVISAg6deqEjRs3srAhaiJ45oaIPNITTzyBmJgYPPXUU0qHQkRNDOfcEJFHOnLkCG655RalwyCiJohnbojIIwUHB+PkyZMICwtTOhQiamI85sxNamoqNBoN5s6da7ePTqeDRqOx+jlx4oT7AiUitygvL2dhQ0SyeMSE4gMHDmDVqlUWEwQbU1hYaPEFW+3bt3dVaERERNTEKH7mpqqqCg8//DBWr16NNm3aiFonNDQUHTp0MP3c7CvciYiIqPlQ/MzNrFmzcN999+Gee+7Byy+/LGqdIUOGoKamBpGRkViyZAnGjBljt69er4derzctG41GXLp0CSEhIdBoNA7HT0RERK4nCAIqKyvRsWNHeHk1fm5G0eLm008/RX5+Pg4cOCCqf3h4OFatWoXo6Gjo9XqsW7cOd999N3Q6HUaNGmVzndTUVKSkpDgzbCIiIlJISUkJOnfu3Ggfxe6WKikpQUxMDLKyskwPsRs9ejQGDx6MpUuXih5n/Pjx0Gg02Lp1q83XbzxzU15eji5duqCoqAitWrUStQ2DwYCdO3dizJgx0Gq1omMjxzDvymDelcG8K4N5dz+5Oa+srET37t1x5coVBAcHN9pXsTM3eXl5KCsrs3j2S319PXJzc/Huu+9Cr9eLmkszbNgwi4fn3cjPzw9+fn5W7W3btrV66q89BoMBgYGBCAkJ4cHvRsy7Mph3ZTDvymDe3U9uzhv6iplSolhxc/fdd+Pw4cMWbY899hhuueUWLFy4UPQk4YKCAoSHh7siRCIiImqCFCtuWrVqhf79+1u0tWjRAiEhIab2pKQknDlzBmvXrgUALF26FN26dUNUVBRqa2uxfv16ZGZmIjMz0+3xExERkWdS/G6pxpSWlqK4uNi0XFtbiwULFuDMmTMICAhAVFQUtm3bhoSEBAWjJCIiIk/iUcWNTqezWM7IyLBYTkxMRGJiovsCIiIiuoHRaERtba3SYTRZBoMBPj4+qKmpQX19vcVrvr6+N73NWwyPKm6IiIg8WW1tLYqKimA0GpUOpckSBAEdOnRASUmJ1eRgLy8vdO/eHb6+vg5tg8UNERGRCIIgoLS0FN7e3oiIiHDKGYbmyGg0oqqqCi1btrTIodFoxNmzZ1FaWoouXbo49EW7LG6IiIhEqKurQ3V1NTp27IjAwEClw2myGi7r+fv7WxWI7du3x9mzZ1FXV+fQrfksO4mIiERomB/i6CUTsq8htzfOxZGKxQ0REZEEfC6h6zgrtyxuiIiISFVY3BAREanY6NGjMXfuXKXDcCtOKCYiInLE3mT3bm+4m7f3XzqdDmPGjMHly5fRunVrRWIQi2duiIiISFVY3BARETUTly9fxqOPPoo2bdogMDAQ48aNw6lTp0yv//rrrxg/fjzatGmDFi1aICoqCtu3b8d//vMfjBkzBgDQpk0baDQaTJs2TaG9uDleliIiImompk2bhlOnTmHr1q0ICgrCwoULkZCQgGPHjkGr1WLWrFmora1Fbm4uWrRogWPHjqFly5aIiIhAZmYmJk2ahMLCQgQFBSEgIEDp3bGLxQ0REVEz0FDU7NmzB8OHDwcAfPzxx4iIiMCWLVvwwAMPoLi4GJMmTcKAAQMAAD169DCt37ZtWwBAaGgo59wQERGR8o4fPw4fHx8MHTrU1BYSEoK+ffvi+PHjAIDZs2fj5ZdfxogRI/DCCy/gp59+Uipch7C4ISIiagYEQbDb3vDleTNmzMAvv/yCKVOm4PDhw4iJicE777zjzjCdgsUNERFRMxAZGYm6ujp8//33praLFy/i5MmT6Nevn6ktIiICM2fOxKZNm/D3v/8dq1evBuC8RyO4A4sbIiKiZqB3796YOHEinnjiCezevRs//vgjHnnkEXTq1AkTJ04EAMydOxfffPMNioqKkJ+fj++++85U+HTt2hUajQZfffUVzp8/j6qqKiV3p1GcUExEROQIhb5UT441a9Zgzpw5+NOf/oTa2lqMGjUK27dvNz2Bu76+HrNmzcJvv/2GoKAg3HvvvXj77bcBAJ06dUJKSgoWLVqExx57DI8++igyMjIU3Bv7WNwQERGpmE6nM/3epk0brF271m7fm82vef755/H88887KzSX4WUpIiIiUhUWN0RERKQqLG6IiIhIVVjcEBERkaqwuCEiIiJVYXFDREREqsLihoiIiFSFxQ0RERGpCosbIiIiUhUWN0RERComCAKefPJJtG3bFhqNBocOHVIsltGjR2PevHku3w4fv0BEROSA5GTP3t6OHTuQkZEBnU6HHj16oF27di6Jy5OwuCE
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.hist(energyloss_lost, bins=150, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=\"lost\")\n",
"plt.hist(energyloss_found, bins=150, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
"plt.xticks(np.arange(0,1.1,0.1), minor=True,)\n",
"plt.yticks(np.arange(0,5,0.5), minor=True)\n",
"plt.xlabel(r\"$E_\\gamma/E_0$\")\n",
"plt.ylabel(\"counts (normed)\")\n",
"plt.title(r'$E_{ph}/E_0$')\n",
"plt.legend()\n",
"plt.grid()\n",
"\n",
"\"\"\"\n",
"found: elektronen verlieren durchschnittlich 0.44 ihrer anfangsenergie durch bremsstrahlung\n",
"lost: elektronen verlieren durchschnittlich 0.59 ihrer anfangsenergie durch bremsstrahlung\n",
"\n",
"-> lost electrons lose slightly more energy than found electrons. This is however nowhere near as extreme as for the B decay\n",
"\"\"\"\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 67,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [],
"source": [
"#ist die shape der teilspur im scifi anders? (koenntest du zum beispiel durch vergleich der verteilungen der fit parameter studieren,\n",
"#in meiner thesis findest du das fitmodell -- ist einfach ein polynom dritten grades)\n",
"z_ref=8520 #mm\n",
"\n",
"def scifi_track(z, a, b, c, d):\n",
" return a + b*(z-z_ref) + c*(z-z_ref)**2 + d*(z-z_ref)**3\n",
"\n",
"def z_mag(xv, zv, tx, a, b):\n",
" \"\"\" optical centre of the magnet is defined as the intersection between the trajectory tangents before and after the magnet\n",
"\n",
" Args:\n",
" xv (double): velo x track\n",
" zv (double): velo z track\n",
" tx (double): velo x slope\n",
" a (double): ax parameter of track fit\n",
" b (double): bx parameter of track fit\n",
"\n",
" Returns:\n",
" double: z_mag\n",
" \"\"\"\n",
" return (xv-tx*zv-a+b*z_ref)/(b-tx)"
]
},
{
"cell_type": "code",
"execution_count": 68,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [],
"source": [
"scifi_found = found[found[\"scifi_hit_pos_x_length\"]>3]\n",
"scifi_lost = lost[lost[\"scifi_hit_pos_x_length\"]>3]\n",
"\n",
"scifi_x_found = scifi_found[\"scifi_hit_pos_x\"]\n",
"scifi_z_found = scifi_found[\"scifi_hit_pos_z\"]\n",
"\n",
"tx_found = scifi_found[\"velo_track_tx\"]\n",
"\n",
"scifi_x_lost = scifi_lost[\"scifi_hit_pos_x\"]\n",
"scifi_z_lost = scifi_lost[\"scifi_hit_pos_z\"]\n",
"\n",
"tx_lost = scifi_lost[\"velo_track_tx\"]\n",
"\n",
"xv_found = scifi_found[\"velo_track_x\"]\n",
"zv_found = scifi_found[\"velo_track_z\"]\n",
"\n",
"xv_lost = scifi_lost[\"velo_track_x\"]\n",
"zv_lost = scifi_lost[\"velo_track_z\"]\n",
"\n",
"\n",
"\n",
"#ak.num(scifi_found[\"energy\"], axis=0)\n",
"#scifi_found.snapshot()"
]
},
{
"cell_type": "code",
"execution_count": 69,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [],
"source": [
"scifi_fitpars_found = ak.ArrayBuilder()\n",
"\n",
"for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
" popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]))\n",
" scifi_fitpars_found.begin_list()\n",
" scifi_fitpars_found.real(popt[0])\n",
" scifi_fitpars_found.real(popt[1])\n",
" scifi_fitpars_found.real(popt[2])\n",
" scifi_fitpars_found.real(popt[3])\n",
" scifi_fitpars_found.end_list()\n",
"\n",
"scifi_fitpars_lost = ak.ArrayBuilder()\n",
"\n",
"for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
" popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]))\n",
" scifi_fitpars_lost.begin_list()\n",
" scifi_fitpars_lost.real(popt[0])\n",
" scifi_fitpars_lost.real(popt[1])\n",
" scifi_fitpars_lost.real(popt[2])\n",
" scifi_fitpars_lost.real(popt[3])\n",
" scifi_fitpars_lost.end_list()\n",
"\n",
"\n",
"scifi_fitpars_lost = scifi_fitpars_lost.to_numpy()\n",
"scifi_fitpars_found = scifi_fitpars_found.to_numpy()\n",
"\n",
"\n",
"\n",
"dtx_found = scifi_fitpars_found[:,1] - tx_found\n",
"dtx_lost = scifi_fitpars_lost[:,1] - tx_lost\n"
]
},
{
"cell_type": "code",
"execution_count": 70,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAABOwAAANVCAYAAADLJ2G2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADqDElEQVR4nOzdeVzU1f7H8ffIjguKC4gLkmWuWUkq5IKZGKlZ6dX0XtxNMzOlm7kmWmmrl8zUFotsUW+praZSiUuQpaJlaouplIGK5q6s398f/ZjrMAMM6wzwej4e86g5c77n+/l+vzPj4TPne47JMAxDAAAAAAAAAJxCNUcHAAAAAAAAAOB/SNgBAAAAAAAAToSEHQAAAAAAAOBESNgBAAAAAAAAToSEHQAAAAAAAOBESNgBAAAAAAAAToSEHQAAAAAAAOBESNgBAAAAAAAAToSEHQAAAAAAAOBESNgBFcDq1avVpk0beXl5yWQyac+ePYqOjpbJZLKot2TJEsXGxjomSCf33nvvKSYmxtFhlFh8fLxMJpPi4+MdHYrdYmNjZTKZdOTIEYvyWbNmqWnTpnJ1dVXt2rUlSWFhYQoLCyu0TWd4r5tMJk2cONGhMQAAKj76eSXniH5eWFiY2rZtW+rtNmvWTCNGjCj1dsuSyWRSdHS0RdmXX36p4OBgVa9eXSaTSR9++GG+fcK8EhISFB0drTNnzpRZzIUZMWKEatSo4bD9AxIJO8DpnTx5UpGRkWrevLk2bNigxMREtWjRQmPGjFFiYqJFXTpy+assCbuKqE+fPkpMTFTDhg3NZR999JGeeuopDRs2TFu2bNEXX3wh6e/38JIlSwptk/c6AKAyoJ9XOujnOVZiYqLGjBljfm4YhgYNGiQ3Nzd9/PHHSkxMVPfu3W32CW1JSEjQ3LlzHZqwA5yBq6MDAFCwn3/+WZmZmfrXv/6l7t27m8u9vb3VuHFjB0Zm7fLly/L09LT6Rbgyu3z5sry8vBwdRrFdunRJ3t7eZbqP+vXrq379+hZl+/btkyRNmjRJDRo0MJe3bt261PefmZkpk8kkV1f+yQMAOBf6ec6tovfzsrOzlZWVJQ8PjzLdT+fOnS2e//nnnzp9+rTuuece9ezZ0+K1vH3C0lDRrxOQH0bYAU5sxIgR6tKliyRp8ODBMplM5tsF894q0axZM/3444/asmWLTCaTTCaTmjVrJul/t1G+8847ioqKkr+/v7y8vNS9e3clJSVZ7HPnzp2677771KxZM3l5ealZs2YaMmSIjh49alEvd0j7pk2bNGrUKNWvX1/e3t5KT0/Xr7/+qpEjR+q6666Tt7e3GjVqpH79+umHH36waCM3rvfee0+PPfaYGjZsqBo1aqhfv346fvy4zp8/r/vvv1/16tVTvXr1NHLkSF24cMGiDcMwtGTJEt14443y8vJSnTp1NHDgQP3222/mOmFhYfrss8909OhR87m5+txlZGToySefVMuWLeXh4aH69etr5MiROnnypMW+mjVrpr59+2rt2rW66aab5Onpqblz5+Z7/eLi4tS/f381btxYnp6euvbaazVu3DilpaXlu83VDh48qDvuuEPe3t6qV6+exo8fr/Pnz9us+8UXX6hnz56qVauWvL29deutt+rLL7+0qJP7ntm9e7cGDhyoOnXqqHnz5vnu/9KlS/r3v/+toKAgeXp6ytfXV8HBwVq5cqVFvR07dqhfv36qW7euPD091bx5c02ePNn8et7bH5o1a6ZZs2ZJkvz8/Cxuo7Dnllh73utvv/22HnnkETVq1EgeHh769ddfdfLkSU2YMEGtW7dWjRo11KBBA912223atm2b1T7S09M1b948tWrVSp6enqpbt6569OihhISEfOMyDEMzZsyQm5ubXnvttQKPAQAA+nkVu5+Xa9u2bercubO8vLzUqFEjzZ49W9nZ2YVul5mZqalTp8rf31/e3t7q0qWLvv32W5t1U1NTNW7cODVu3Fju7u4KCgrS3LlzlZWVZa5z5MgRmUwmPfvss3ryyScVFBQkDw8Pbd68Od8Y3n//fXXq1Ek+Pj7y9vbWNddco1GjRlnUOXPmjB555BFdc8018vDwUIMGDXTnnXfq4MGD5jpX9+Wio6PNyebHHnvM4r1qzy2x0dHRevTRRyVJQUFB5uuZOx1MQdfp5ZdfVrdu3dSgQQNVr15d7dq107PPPqvMzEyr/WzYsEE9e/Y0H3urVq20YMGCfOOSpK+//lr16tVT3759dfHixQLrAqWB4QaAE5s9e7Y6duyoBx98UPPnz1ePHj1Uq1Ytm3XXrVungQMHysfHx3xLYd5f02bMmKGbb75Zr7/+us6ePavo6GiFhYUpKSlJ11xzjaS//7G//vrrdd9998nX11cpKSlaunSpbrnlFu3fv1/16tWzaHPUqFHq06eP3n77bV28eFFubm76888/VbduXT399NOqX7++Tp8+rbfeekudOnVSUlKSrr/+equ4evToodjYWB05ckT//ve/NWTIELm6uqp9+/ZauXKlkpKSNGPGDNWsWVOLFi0ybztu3DjFxsZq0qRJeuaZZ3T69GnNmzdPoaGh2rt3r/z8/LRkyRLdf//9OnTokNatW2ex75ycHPXv31/btm3T1KlTFRoaqqNHj2rOnDkKCwvTzp07LX6x2717tw4cOKBZs2YpKChI1atXz/f6HTp0SCEhIRozZox8fHx05MgRLVy4UF26dNEPP/wgNze3fLc9fvy4unfvLjc3Ny1ZskR+fn569913bc6Z9s4772jYsGHq37+/3nrrLbm5uemVV15R7969tXHjRqtfNu+9917dd999Gj9+fIGdjaioKL399tt68sknddNNN+nixYvat2+fTp06Za6zceNG9evXT61atdLChQvVtGlTHTlyRJs2bcq33XXr1unll1/W8uXLtWHDBvn4+BRpFIE97/Xp06crJCREy5YtU7Vq1dSgQQNzx3zOnDny9/fXhQsXtG7dOoWFhenLL780/5GUlZWliIgIbdu2TZMnT9Ztt92mrKwsffPNN0pOTlZoaKhVTOnp6RoxYoQ+++wzffLJJ7rjjjvsPh4AQNVEP69i9/OkvxNp9913n6ZNm6Z58+bps88+05NPPqm//vpLixcvLnDbsWPHasWKFfr3v/+tXr16ad++fbr33nutfpxNTU1Vx44dVa1aNT3++ONq3ry5EhMT9eSTT+rIkSN68803LeovWrRILVq00PPPP69atWrpuuuus7n/xMREDR48WIMHD1Z0dLQ8PT119OhRffXVV+Y658+fV5cuXXTkyBE99thj6tSpky5cuKCtW7cqJSVFLVu2tGp3zJgxat++ve6991499NBDGjp0aJFG+I0ZM0anT5/WSy+9pLVr15pvn736Loz8rtOhQ4c0dOhQBQUFyd3dXXv37tVTTz2lgwcP6o033jBvv3z5co0dO1bdu3fXsmXL1KBBA/3888/mO0Bs+e9//6thw4Zp1KhReumll+Ti4mL3MQHFZgBwaps3bzYkGe+//75F+Zw5c4y8H+E2bdoY3bt3z7eNm2++2cjJyTGXHzlyxHBzczPGjBmT7/6zsrKMCxcuGNWrVzdefPFFc/mbb75pSDKGDRtW6DFkZWUZGRkZxnXXXWdMmTLFKq5+/fpZ1J88ebIhyZg0aZJF+d133234+vqanycmJhqSjBdeeMGi3u+//254eXkZU6dONZf16dPHCAwMtIpt5cqVhiRjzZo1FuXfffedIclYsmSJuSwwMNBwcXExfvrpp0KPOa+cnBwjMzPTOHr0qCHJ+Oijjwqs/9hjjxkmk8nYs2ePRXmvXr0MScbmzZsNwzCMixcvGr6+vlbnMDs722jfvr3RsWNHc1nue+bxxx+3K+a2bdsad999d4F1mjdvbjRv3ty4fPlyvnVy3yuHDx+2iuXkyZMWdbt3727zPZxXYe/1bt26FdpGVlaWkZmZafTs2dO45557zOUrVqwwJBmvvfZagdtLMh588EHj1KlTRpcuXYxGjRpZXS8AAApCP+9/Klo/r3v37jb7dGPHjjWqVatmHD16NN9tDxw4YEiyOF+GYRj
"text/plain": [
"<Figure size 1500x1000 with 4 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"fig, ((ax0, ax1), (ax2, ax3)) = plt.subplots(nrows=2, ncols=2, figsize=(15,10))\n",
"\n",
"ax0.hist(scifi_fitpars_found[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$a_x$ found\")\n",
"ax0.hist(scifi_fitpars_lost[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$a_x$ lost\")\n",
"ax0.set_xlabel(\"a\")\n",
"ax0.set_ylabel(\"normed\")\n",
"ax0.set_title(\"fitparameter a der scifi track\")\n",
"ax0.legend()\n",
"\n",
"ax1.hist(scifi_fitpars_found[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$b_x$ found\")\n",
"ax1.hist(scifi_fitpars_lost[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$b_x$ lost\")\n",
"ax1.set_xlabel(\"b\")\n",
"ax1.set_ylabel(\"normed\")\n",
"ax1.set_title(\"fitparameter b der scifi track\")\n",
"ax1.legend()\n",
"\n",
"ax2.hist(scifi_fitpars_found[:,2], bins=200, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$c_x$ found\")\n",
"ax2.hist(scifi_fitpars_lost[:,2], bins=1000, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$c_x$ lost\")\n",
"ax2.set_xlim([-4e-5,4e-5])\n",
"ax2.set_xticks(np.arange(-4e-5,4.5e-5,1e-5),minor=False)\n",
"ax2.set_xlabel(\"c\")\n",
"ax2.set_ylabel(\"normed\")\n",
"ax2.set_title(\"fitparameter c der scifi track\")\n",
"ax2.legend()\n",
"\n",
"ax3.hist(scifi_fitpars_found[:,3], bins=200, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$d_x$ found\")\n",
"ax3.hist(scifi_fitpars_lost[:,3], bins=1000, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$d_x$ lost\")\n",
"ax3.set(xlim=(-5e-8,5e-8))\n",
"#ax3.text(-4e-8,3e8,\"d negligible <1e-7\")\n",
"ax3.set_xlabel(\"d\")\n",
"ax3.set_ylabel(\"normed\")\n",
"ax3.set_title(\"fitparameter d der scifi track\")\n",
"ax3.legend()\n",
"\n",
"\"\"\"\n",
"a_x: virtual hit on the reference plane\n",
"\"\"\"\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 71,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"found\n",
"a = 18.04180550457298\n",
"b = 0.0057651944487250376\n",
"c = 9.480277789440454e-08\n",
"d = -4.452015498411874e-11\n",
"lost\n",
"a = -35.48371609662243\n",
"b = -0.010382326385451382\n",
"c = -6.208301288913938e-07\n",
"d = 9.580933791481267e-11\n"
]
}
],
"source": [
"print(\"found\")\n",
"print(\"a = \", str(np.mean(scifi_fitpars_found[:,0])))\n",
"print(\"b = \", str(np.mean(scifi_fitpars_found[:,1])))\n",
"print(\"c = \", str(np.mean(scifi_fitpars_found[:,2])))\n",
"print(\"d = \", str(np.mean(scifi_fitpars_found[:,3])))\n",
"\n",
"print(\"lost\")\n",
"print(\"a = \", str(np.mean(scifi_fitpars_lost[:,0])))\n",
"print(\"b = \", str(np.mean(scifi_fitpars_lost[:,1])))\n",
"print(\"c = \", str(np.mean(scifi_fitpars_lost[:,2])))\n",
"print(\"d = \", str(np.mean(scifi_fitpars_lost[:,3])))"
]
},
{
"cell_type": "code",
"execution_count": 72,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"-1.5438992626615335e-08"
]
},
"execution_count": 72,
2023-09-25 16:36:13 +02:00
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"np.min(scifi_fitpars_found[:,3])"
]
},
{
"cell_type": "code",
"execution_count": 73,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAABQEAAAIhCAYAAAD+cXzUAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdeZhcVYH38e+9VXVr36v3vTvd6U5n6ewLhCSELSQsAirosLg76MygjqOjr7sDijOj46voOCqoiOACKDskgSSQrbMvnXR3Or2v1bXvdbf3Dx/yGhMgOiIE7+d58kedOnXuuXWS9K/PveceQdd1HYPBYDAYDAaDwWAwGAwGg8HwliW+0R0wGAwGg8FgMBgMBoPBYDAYDK8vYxLQYDAYDAaDwWAwGAwGg8FgeIszJgENBoPBYDAYDAaDwWAwGAyGtzhjEtBgMBgMBoPBYDAYDAaDwWB4izMmAQ0Gg8FgMBgMBoPBYDAYDIa3OGMS0GAwGAwGg8FgMBgMBoPBYHiLMyYBDQaDwWAwGAwGg8FgMBgMhrc4YxLQYDAYDAaDwWAwGAwGg8FgeIszJgENBoPBYDAYDAaDwWAwGAyGtzhjEtBg+Bvz0EMP0d7ejt1uRxAEDhw48EZ36ay++MUvIgjCa9Z74IEH+Na3vvX6d+hV1NfXs2HDhje0DwDRaJQbb7yR0tJSBEHg2muvfd2PWV9fz2233XZa2f79+1m1ahVerxdBEPjWt77FCy+8gCAIvPDCC6/a3n333YcgCOzZs+c1j7169WpWr1596nU2m+WLX/ziax7DYDAYDIa/RS//jB0YGHhd2r/nnnu47777zrn+nXfeyaOPPvq69OVcDAwMIAgC//7v//6G9eFlZ8tOr6eXz/2Px+tsvyecaya/7bbbcLlc53R8QRD44he/eOp1V1cXX/ziF1+3v5sGg+H/M7/RHTAYDH894XCYm2++mSuuuIJ77rkHq9VKS0vLG92t/5UHHniAI0eOcMcdd7zRXXnDfeUrX+GRRx7hxz/+MU1NTQQCgdf9mI888ggej+e0sve+971kMhkefPBB/H4/9fX1OBwOduzYwaxZs/5ix77nnntOe53NZvnSl74EcNrkoMFgMBgMhtffPffcQygUOuPi4Cu58847ueGGG/4qFy3f7M6WnV5PFRUV7Nixg6amplNlr/R7wvvf/36uuOKKv+jxd+zYQXV19anXXV1dfOlLX2L16tWv+7kbDH/rjElAg+FvSE9PD7Is83d/93esWrXqje7OX52qqiiKgtVqfaO78ro4cuQITU1NvPvd7/6rHXP+/Pln7ccHPvAB1q1bd1r5smXL/qLH/ktOKBoMBoPBYHjzyuVy2Gy2c7oj7Xz0Stnp9WK1Ws/IZa/0e4LD4Thtwu4v4S+dCQ0Gw7kzlgMbDH8jbrvtNi688EIA3vnOdyIIwml3S/3ud79j+fLlOBwO3G43l156KTt27DijjbNdnTvbMgFBEPjoRz/Kz372M9ra2nA4HMybN4/HH3/8jM8/8cQTdHR0YLVaaWhoOOdlGatXr+aJJ55gcHAQQRBO/YH/v8zh7rvv5qtf/SoNDQ1YrVaef/558vk8n/jEJ+jo6MDr9RIIBFi+fDm//e1vzziGpmn83//7f+no6MBut+Pz+Vi2bBm/+93vXrVv99xzD2azmS984Qunyr73ve8xb948XC4Xbreb1tZWPvOZz7zmeUajUW6//XaqqqqQJInGxkY++9nPUigUTjvXjRs3cuzYsVPfw6sti928eTOrV68mGAxit9upra3l+uuvJ5vNnqpTKBT48pe/TFtbGzabjWAwyJo1a9i+ffupOn+4HPjlZUaKovC9733vtPE41+XAL0ulUvz93/89oVCIYDDIddddx9jY2Gl1/nA58MDAACUlJQB86UtfOnXsl/sWDof54Ac/SE1NDVarlZKSEi644AI2btx4Tv0xGAwGg+Gt6sc//jHz5s3DZrMRCAR429vexrFjx06rc/LkSW688UYqKyuxWq2UlZWxdu3aU4+Vqa+v5+jRo2zZsuXUz+BXu6NLEAQymQw/+clPTtV/+Wf6y3ni2Wef5b3vfS8lJSU4HA4KhQInTpzgPe95D83NzTgcDqqqqrjqqqs4fPjwGceIx+N84hOfoLGxEavVSmlpKVdeeSXHjx9/xX7Jssytt96Ky+U6lVmz2Sz//M//TENDw6nvaNGiRfziF794ze/2yJEjXHPNNfj9fmw2Gx0dHfzkJz859f6rZadXci55cnR09FTukSSJyspKbrjhBiYnJ4EzlwO/2u8J57oc+GUnTpzgyiuvxOVyUVNTwyc+8YlTmfVlf7gc+L777uPtb387AGvWrDn1Hbzct/3797NhwwZKS0uxWq1UVlayfv16RkZGzrlPBoPh/zPuBDQY/kZ87nOfY8mSJXzkIx/hzjvvZM2aNaeWcT7wwAO8+93v5rLLLuMXv/gFhUKBu+++m9WrV7Np06ZToeBP9cQTT9DZ2cmXv/xlXC4Xd999N29729vo7u6msbERgE2bNnHNNdewfPlyHnzwQVRV5e677z4VUl7NPffcwwc/+EH6+vp45JFHzlrn29/+Ni0tLfz7v/87Ho+H5uZmCoUC0WiUf/7nf6aqqopiscjGjRu57rrruPfee7nllltOff62227j/vvv533vex9f/vKXkSSJffv2veIzS3Rd55Of/CTf/va3+eEPf3hqEurBBx/k9ttv5x/+4R/493//d0RR5MSJE3R1db3qOebzedasWUNfXx9f+tKXmDt3Ltu2beOuu+7iwIEDPPHEE6eWdNx+++0kEgl+/vOfA698p9zAwADr169n5cqV/PjHP8bn8zE6OsrTTz9NsVjE4XCgKArr1q1j27Zt3HHHHVx88cUoisLOnTsZGhpixYoVZ7S7fv16duzYwfLly7nhhhv4xCc+8arn9mre//73s379eh544AGGh4f55Cc/yd/93d+xefPms9avqKjg6aef5oorruB973sf73//+wFOTQzefPPN7Nu3j3/7t3+jpaWFeDzOvn37iEQif3YfDQaDwWA4391111185jOf4aabbuKuu+4iEonwxS9+keXLl9PZ2UlzczMAV1555amMVltby/T0NNu3bycejwO/fzzIDTfcgNfrPfW4jldbebFjxw4uvvhi1qxZw+c+9zmAsz5eZP369fzsZz8jk8lgsVgYGxsjGAzyta99jZKSEqLRKD/5yU9YunQp+/fvZ+bMmcDvLyZeeOGFDAwM8KlPfYqlS5eSTqfZunUr4+PjtLa2ntGneDzOddddx7Fjx9iyZQsLFy4E4OMf/zg/+9nP+OpXv8r8+fPJZDIcOXLkNTNEd3c3K1asoLS0lG9/+9sEg0Huv/9+brvtNiYnJ/mXf/mXPzk7nUueHB0dZfHixciyzGc+8xnmzp1LJBLhmWeeIRaLUVZWdka7r/Z7wp9ClmWuvvpq3ve+9/GJT3yCrVu38pWvfAWv18vnP//5s35m/fr13HnnnXzmM5/hu9/9LgsWLACgqamJTCbDpZdeSkNDA9/97ncpKytjYmKC559/nlQq9Sf3z2AwALrBYPib8fzzz+uA/qtf/epUmaqqemVlpT5nzhxdVdVT5alUSi8tLdVXrFhxquzWW2/V6+rqzmj3C1/4gv7H/50AellZmZ5MJk+VTUxM6KIo6nfdddepsqVLl+qVlZV6Lpc7VZZMJvVAIHBGm2ezfv36s/apv79fB/Smpia9WCy+ahuKouiyLOvve9/79Pnz558q37p1qw7on/3sZ1/183V1dfr69ev1bDarX3/99brX69U3btx4Wp2PfvSjus/ne83z+WPf//73dUD/5S9/eVr517/+dR3Qn3322VNlq1at0tvb21+zzV//+tc6oB84cOAV6/z0pz/VAf1//ud/XrWturo6/dZbbz2tDNA/8pGPnFb28t+9559//lXbu/fee3VAv/32208rv/vuu3VAHx8fP1W2atUqfdWqVadeh8NhHdC/8IUvnNGuy+XS77jjjlc9tsFgMBgMb2Uv/4zt7+/XdV3XY7GYbrfb9SuvvPK0ekNDQ7r
"text/plain": [
"<Figure size 1500x600 with 2 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(15,6))\n",
"\n",
"for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
" z_coord = np.linspace(scifi_z_found[i,0],12000,300)\n",
" fit = scifi_track(z_coord, *scifi_fitpars_found[i])\n",
" ax0.plot(z_coord, fit, \"-\", lw=0.5)\n",
" ax0.errorbar(ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]),fmt=\".\",ms=2)\n",
"\n",
"#ax0.legend()\n",
"ax0.set_xlabel(\"z [mm]\")\n",
"ax0.set_ylabel(\"x [mm]\")\n",
"ax0.set_title(\"found tracks of scifi hits\")\n",
"ax0.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
"ax0.grid()\n",
"\n",
"for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
" z_coord = np.linspace(scifi_z_lost[i,0],12000,300)\n",
" fit = scifi_track(z_coord, *scifi_fitpars_lost[i])\n",
" ax1.plot(z_coord, fit, \"-\", lw=0.5)\n",
" ax1.errorbar(ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]),fmt=\".\",ms=2)\n",
"\n",
"#ax1.legend()\n",
"ax1.set_xlabel(\"z [mm]\")\n",
"ax1.set_ylabel(\"x [mm]\")\n",
"ax1.set_title(\"lost tracks of scifi hits\")\n",
"ax1.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
"ax1.grid()\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 74,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"found \n",
"zmag = 5318.452650242005\n",
"lost \n",
"zmag = 5425.137423441005\n"
]
}
],
"source": [
"#vergleich der zmag werte\n",
"zmag_found = z_mag(xv_found, zv_found, tx_found, scifi_fitpars_found[:,0], scifi_fitpars_found[:,1])\n",
"zmag_lost = z_mag(xv_lost, zv_lost, tx_lost, scifi_fitpars_lost[:,0], scifi_fitpars_lost[:,1])\n",
"zmag_lost = zmag_lost[~np.isnan(zmag_lost)]\n",
"zmag_found = zmag_found[~np.isnan(zmag_found)]\n",
"\n",
"print(\"found \\nzmag = \", str(np.mean(zmag_found)))\n",
"print(\"lost \\nzmag = \", str(np.mean(zmag_lost)))"
]
},
{
"cell_type": "code",
"execution_count": 75,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAloAAAHNCAYAAADG0NWMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABOZ0lEQVR4nO3deXwURf7/8feQkwQI4QqJQAigknCIJhgIpwrhUAFBiStGUMTNgiDgsoLKgniAigrItbgIHqvw8wsoKiJBAUEiYrgFlVUUxMRwmXAIuer3B9/Ml0kmdzqT4/V8POZhprq6urqY7vlYVVNtM8YYAQAAoMzVcHUFAAAAqioCLQAAAIsQaAEAAFiEQAsAAMAiBFoAAAAWIdACAACwCIEWAACARQi0AAAALEKgBQAAYBECLQAAAIsQaAEAAFiEQAsVyvbt2zV9+nT98ccfRco/ffp02Ww2nTx5ssB8y5cvl81m088//1ys+hS1/IrC2Xnm16YlbRMUn9VtXZryi3vNFUfO9VOeVq5cqTZt2qhmzZqy2Wzas2dPsa4L0DZljUALFcr27dv11FNPlfkFfuuttyohIUGBgYFlWm5F4+w882vT6tImKJhV15wrnDhxQrGxsWrZsqXWr1+vhIQEXXPNNcW6LkDblDV3V1cAKA8NGzZUw4YNXV0NyxXnPKtLm6D6+OGHH5SRkaF7771XPXr0sKf7+PhUuM/6hQsX5OPj4+pqlKvqeM4SPVqVSk43/L59+3TXXXfJz89P9erV08SJE5WZmanvv/9effv2Ve3atdW8eXO98MILDvv/97//1f3336+rr75aPj4+uuqqq3T77bdr//79eY71wQcfqH379vLy8lKLFi00d+5cp8MAOWnffvut/vKXv8jPz08BAQF64IEHlJqa6pD38OHDuueee9SoUSN5eXkpNDRUCxYscChr0qRJkqSQkBDZbDbZbDZt3ry5WO303XffqUWLFoqMjFRKSook50Mr3333nf7yl78oICBAXl5eatasme677z5dunSp2OU7k9M2u3fv1uDBg1WnTh35+fnp3nvv1YkTJxzybtu2Tbfccotq164tHx8fRUVF6eOPP3bIc+LECT300ENq2rSpvLy81LBhQ3Xp0kUbN26058l9ngW1aX7DTUWpS3H+3Z05dOiQvS65X35+fjLGFFqG1Qr7fBTneipJ+SNGjFDz5s3z7FfU4bii1K+wa66wa/ZKH3/8sTp06CAvLy+FhIRo9uzZRWqHb7/9VjabTe+99549LTExUTabTW3atHHIO2DAAIWHhzstZ8SIEerataskKSYmRjabTT179pRUvOuiONdtUT8DOWXu2rVLd955p/z9/dWyZcsSlVHS+3+O0t6Hi/qZyO+ci3Ifq2ro0aqEhg4dqnvvvVd//etfFR8frxdeeEEZGRnauHGjRo8erb///e9655139Nhjj6lVq1YaPHiwJOm3335T/fr1NWvWLDVs2FCnT5/WG2+8ocjISO3evVvXXnutJGn9+vUaPHiwunfvrpUrVyozM1OzZ8/W77//nm+dhgwZopiYGI0cOVL79+/XlClTJEmvv/66JOngwYOKiopSs2bN9NJLL6lx48b69NNPNW7cOJ08eVLTpk3Tgw8+qNOnT+vVV1/V6tWr7d38YWFhRW6bLVu26I477lD37t31zjvv5Pt/T3v37lXXrl3VoEEDzZgxQ1dffbWSkpK0du1apaeny8vLq1TlX+mOO+7Q0KFDFRcXp2+//VZTp07VwYMHtWPHDnl4eGjLli3q3bu32rdvr6VLl8rLy0sLFy7U7bffrnfffVcxMTGSpNjYWO3atUvPPvusrrnmGv3xxx/atWuXTp06le+xC2pTZ/N5ilqXHIX9u+cnKChICQkJDmkbN27U1KlT9cADD5T7vJ7civL5KOr1VNLyS6so9Svo81GUazbHZ599poEDB6pz585asWKFsrKy9MILLxR4z8jRpk0bBQYGauPGjbrrrrskXf4s1KxZUwcPHtRvv/2moKAgZWZmasuWLYqLi3NaztSpU3XjjTdqzJgxeu6553TTTTepTp06TvMWdN45AUVh121R2/hKgwcP1t133624uDidP3++RGWU9P4vlf4+XJzPRH7nXJL7WKVnUGlMmzbNSDIvvfSSQ3qHDh2MJLN69Wp7WkZGhmnYsKEZPHhwvuVlZmaa9PR0c/XVV5sJEybY0zt27GiaNm1qLl26ZE87e/asqV+/vsn9kcmp0wsvvOCQPnr0aOPt7W2ys7ONMcb06dPHNGnSxKSmpjrke/jhh423t7c5ffq0McaYF1980UgyR44cKUKL/N/xT5w4Yd566y3j6elpxo0bZ7KyshzyLVu2zKHcm2++2dStW9ekpKSUSfkF7Xtl2xpjzH/+8x8jybz99tvGGGM6depkGjVqZM6ePWvPk5mZadq2bWuaNGlib8NatWqZ8ePHF3jM3OdpTP5t6ixvUetS1H/3olqzZo3x8vIyf//734u1n1WK+vm4Un7XkzEl+/wNHz7cBAcH50nPafuCyi9O/fL7fBT1mjXGmMjISBMUFGT+/PNPe1paWpqpV69enro6c++995oWLVrY3/fq1cuMGjXK+Pv7mzfeeMMYY8yXX35pJJkNGzbkW86mTZuMJPPee+85pBfnuijqdetMfm2cU+Y///nPfPctahmluf+X9j5cnM9EfudclPtYVcPQYSV02223ObwPDQ2VzWZTv3797Gnu7u5q1aqVfvnlF3taZmamnnvuOYWFhcnT01Pu7u7y9PTU4cOHdejQIUnS+fPn9c0332jQoEHy9PS071urVi3dfvvt+dZpwIABDu/bt2+vixcvKiUlRRcvXtRnn32mO+64Qz4+PsrMzLS/+vfvr4sXL+qrr74qVZs8++yzGjFihGbNmqW5c+eqRo38P9oXLlzQli1bNHTo0CLP2yhO+bkNGzbM4f3QoUPl7u6uTZs26fz589qxY4fuvPNO1apVy57Hzc1NsbGx+vXXX/X9999Lkm688UYtX75czzzzjL766itlZGQUuQ5FUZy65Cjo372o3nrrLd1111168skn9eKLL5buJMpAUT8fRbmeSlN+aZW0fpKKdc2eP39eO3fu1ODBg+Xt7W0vo3bt2gXeM650yy236KefftKRI0d08eJFbdu2TX379tVNN92k+Ph4SZd7uby8vOzDg1Yr6LrNUdw2HjJkSJ604pZR0vt/ae/DJd0/9zlbfR+riAi0KqF69eo5vPf09JSPj4/DTS4n/eLFi/b3EydO1NSpUzVo0CB9+OGH2rFjh3bu3KnrrrtOf/75pyTpzJkzMsYoICAgz3GdpeWoX7++w/ucoY8///xTp06dUmZmpl599VV5eHg4vPr37y9JpV4+4e2339ZVV12lu+++u9C8Z86cUVZWlpo0aWJJ+bk1btzY4b27u7vq16+vU6dO2dvb2S//goKCJMnepb5y5UoNHz5c//73v9W5c2fVq1dP9913n5KTk4tdJ2eKU5ccBf27F8X8+fN1//3365VXXtGTTz7psG3u3Lm66667dM8996hOnTqKjIxUcnKyxo0bp3r16qlt27b2L5Jjx46pX79+atiwoerWrauHHnpI2dnZki5/kf3jH/9Q48aN1bJlSy1evNg+Pya/dijK56Mo11Npyi+tktZPUrGu2TNnzig7OzvP51zK+9nPT69evSRdDqa2bdumjIwM3XzzzerVq5c+++wz+7YuXbqoZs2aRW6D0ijous1R3DZ2dm0Vt4yS3v9Lex8u6f65z9nq+1hFxBytauTtt9/Wfffdp+eee84h/eTJk6pbt64kyd/fXzabzencipJeCP7+/vZekTFjxjjNExISUqKyc6xfv14xMTHq1q2bPvvsMwUHB+ebt169enJzc9Ovv/5qSfm5JScn66qrrrK/z8zM1KlTp1S/fn35+/urRo0aSkpKyrPfb7/9Jklq0KCB/b9z5szRnDlzdPToUa1du1aTJ09WSkqK1q9fX+T65Kc4dSkLzz77rKZPn67XX39d9913X57t+/bt09dff63Vq1fr9ddfV1RUlG655RbNnTt
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.hist(zmag_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
"plt.hist(zmag_lost, bins=3000, density=True, alpha=0.5, histtype=\"bar\",color=\"darkorange\", label=\"lost\")\n",
"plt.xlabel(\"$z_{mag}$ [mm]\")\n",
"plt.ylabel(\"normed\")\n",
"plt.title(\"magnet kick position $z_{mag}$ calculated w fitparameters\")\n",
"plt.legend()\n",
"plt.xticks(np.arange(5100,5805,5), minor=True)\n",
"plt.yticks(np.arange(0,0.011,0.001), minor=True)\n",
"plt.xlim(5100,5800)\n",
"\n",
"\"\"\"\n",
"wir können definitiv einen unterschied zwischen den z_mag werten, bzw deren verteilungen für lost and found, erkennen\n",
"\"\"\"\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 76,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
2023-09-26 12:18:03 +02:00
"name": "stdout",
"output_type": "stream",
"text": [
"2591\n",
"1908\n",
"4499\n"
]
2023-09-25 16:36:13 +02:00
}
],
"source": [
"both = ak.concatenate([found,lost],axis=0)\n",
2023-09-26 12:18:03 +02:00
"print(ak.num(found,axis=0))\n",
"print(ak.num(lost,axis=0))\n",
"print(ak.num(both,axis=0))"
2023-09-25 16:36:13 +02:00
]
},
{
"cell_type": "code",
"execution_count": 81,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"#events w/ shared track electrons from found and lost: 1\n",
"event_count: [2822, 2822]\n",
"velo idx: [146, 146]\n",
"mcp_index: [5806, 5829]\n",
"\n",
"velo x: [-75.1, -75.1]\n",
"velo y: [5.88, 5.88]\n",
"\n",
"velo tx: [-0.0908, -0.0908]\n",
"velo ty: [0.00703, 0.00703]\n"
]
2023-09-25 16:36:13 +02:00
}
],
"source": [
"#versuche teilchen von denselben events mit shared tracks zu finden. \n",
"#idee: alle teilchen eines events sind durch event_count auffindbar. \n",
"a_f_itr = found[\"event_count\"].to_numpy()\n",
"f_itr = np.unique(a_f_itr)\n",
"\n",
"shared = ak.ArrayBuilder()\n",
2023-09-26 12:18:03 +02:00
"count = 0\n",
"\n",
"for itr in f_itr:\n",
" temp = both[both[\"event_count\"]==itr]\n",
" if ak.num(temp,axis=0)>1:\n",
" #iterate over cols in temp and append all with duplicate velo_track_idx, such that possibles is array with possible shared tracks particles\n",
" #idea: look at one event at a time and at one velo track at a time. if there are at least two e with the same velo_track_idx in the same event,\n",
" #concatenate to the array of other shared track particles\n",
2023-09-26 12:18:03 +02:00
" _jitr = temp[\"velo_track_idx\"].to_numpy()\n",
" jitr = np.unique(_jitr)\n",
" for jentry in jitr:\n",
" jtem = temp[temp[\"velo_track_idx\"]==jentry]\n",
" if ak.num(jtem,axis=0)>1:\n",
" shared.append(jtem)\n",
2023-09-26 12:18:03 +02:00
" else:\n",
" continue\n",
" else:\n",
" continue\n",
"shared = ak.Array(shared)\n",
"\n",
"idx=0\n",
"print(\"#events w/ shared track electrons from found and lost: \",ak.num(shared, axis=0))\n",
"\n",
"print(\"event_count: \", shared[idx,:,\"event_count\"])\n",
"print(\"velo idx: \" ,shared[idx,:,\"velo_track_idx\"])\n",
"print(\"mcp_index: \", shared[idx,:,\"mcp_idx\"])\n",
"\n",
"print(\"\\nvelo x: \" ,shared[idx,:,\"velo_track_x\"])\n",
"print(\"velo y: \" ,shared[idx,:,\"velo_track_y\"])\n",
"\n",
"print(\"\\nvelo tx: \" ,shared[idx,:,\"velo_track_tx\"])\n",
"print(\"velo ty: \" ,shared[idx,:,\"velo_track_ty\"])"
2023-09-26 12:18:03 +02:00
]
},
{
"cell_type": "code",
"execution_count": 82,
2023-09-26 12:18:03 +02:00
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"#velo_track_idx in all events: 217\n",
"velo idx: [0, 0, 0]\n",
"mcp_index: [1066, 1251, 666]\n",
"event_count: [5735, 7049, 7378]\n",
"\n",
"velo x: [-2.21, 9.45, -18.5]\n",
"velo y: [-21.6, -33.4, 17.4]\n",
"\n",
"velo tx: [-0.0022, 0.0153, -0.0224]\n",
"velo ty: [-0.0263, -0.0469, 0.0212]\n"
]
}
],
2023-09-26 12:18:03 +02:00
"source": [
"#electrons with same velo_track_idx from all events\n",
"temp_ = found[\"velo_track_idx\"].to_numpy()\n",
"temp = np.unique(temp_)\n",
"count=0\n",
"psb=ak.ArrayBuilder()\n",
"for jentry in temp:\n",
" jtem = found[found[\"velo_track_idx\"]==jentry]\n",
" if ak.num(jtem,axis=0)>1:\n",
" psb.append(jtem)\n",
" else:\n",
" continue\n",
"\n",
"psb = ak.Array(psb)\n",
"\n",
"print(\"#velo_track_idx in all events: \",ak.num(psb,axis=0))\n",
"idx = 0\n",
"print(\"velo idx: \" ,psb[idx,:,\"velo_track_idx\"])\n",
"print(\"mcp_index: \", psb[idx,:,\"mcp_idx\"])\n",
"print(\"event_count: \", psb[idx,:,\"event_count\"])\n",
"print(\"\\nvelo x: \" ,psb[idx,:,\"velo_track_x\"])\n",
"print(\"velo y: \" ,psb[idx,:,\"velo_track_y\"])\n",
"\n",
"print(\"\\nvelo tx: \" ,psb[idx,:,\"velo_track_tx\"])\n",
"print(\"velo ty: \" ,psb[idx,:,\"velo_track_ty\"])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": 83,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"#events w/ shared track electrons from Photon Conversions: 951\n",
"shared_idx: 0\n",
"event_count: [11, 11]\n",
"velo idx: [27, 27]\n",
"mcp_index: [1211, 1215]\n",
"\n",
"velo x: [-20.4, -20.4]\n",
"velo y: [-24.4, -24.4]\n",
"\n",
"velo tx: [-0.0234, -0.0234]\n",
"velo ty: [-0.028, -0.028]\n"
]
}
],
"source": [
"#generell: wie viele elektronen von photon conversions teilen sich einen track?\n",
"conv = allcolumns[(allcolumns.fromPairProd) & (allcolumns.isElectron) & (~allcolumns.lost)]\n",
"\n",
"conv_it = conv[\"event_count\"].to_numpy()\n",
"conv_itr = np.unique(conv_it)\n",
"\n",
"cshared = ak.ArrayBuilder()\n",
"count = 0\n",
"\n",
"for itr in conv_itr:\n",
" temp = conv[conv[\"event_count\"]==itr]\n",
" if ak.num(temp,axis=0)>1:\n",
" #iterate over cols in temp and append all with duplicate velo_track_idx, such that possibles is array with possible shared tracks particles\n",
" #idea: look at one event at a time and at one velo track at a time. if there are at least two e with the same velo_track_idx in the same event,\n",
" #concatenate to the array of other shared track particles\n",
" _jitr = temp[\"velo_track_idx\"].to_numpy()\n",
" jitr = np.unique(_jitr)\n",
" for jentry in jitr:\n",
" jtem = temp[temp[\"velo_track_idx\"]==jentry]\n",
" if ak.num(jtem,axis=0)>1:\n",
" cshared.append(jtem)\n",
" else:\n",
" continue\n",
" else:\n",
" continue\n",
"cshared = ak.Array(cshared)\n",
"\n",
"\n",
"print(\"#events w/ shared track electrons from Photon Conversions: \",ak.num(cshared, axis=0))\n",
"idx = 0\n",
"print(\"shared_idx: \", idx)\n",
"print(\"event_count: \", cshared[idx,:,\"event_count\"])\n",
"print(\"velo idx: \" ,cshared[idx,:,\"velo_track_idx\"])\n",
"print(\"mcp_index: \", cshared[idx,:,\"mcp_idx\"])\n",
"\n",
"print(\"\\nvelo x: \" ,cshared[idx,:,\"velo_track_x\"])\n",
"print(\"velo y: \" ,cshared[idx,:,\"velo_track_y\"])\n",
"\n",
"print(\"\\nvelo tx: \" ,cshared[idx,:,\"velo_track_tx\"])\n",
"print(\"velo ty: \" ,cshared[idx,:,\"velo_track_ty\"])"
2023-09-26 12:18:03 +02:00
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
2023-09-25 16:36:13 +02:00
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
2023-09-26 12:18:03 +02:00
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
2023-09-25 16:36:13 +02:00
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
2023-09-25 16:36:13 +02:00
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"\"\"\"\n",
"was meinst du? velo track sharen? wie überprüfe ich das? versuche ich einfach e teilchen mit identischen velo tracks und slopes zu finden?\n",
"benutze ich mother key um elektronen von gleichen events zu finden? meinst du wie viele elektronen vom selben event teilen sich einen track,\n",
"oder wie viele teilchen von unterschiedlichen events haben dieselben tracks?\n",
"\"\"\"\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": []
}
],
"metadata": {
"kernelspec": {
"display_name": "env1",
"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.11.5"
},
"orig_nbformat": 4
},
"nbformat": 4,
"nbformat_minor": 2
}