You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

717 lines
1.3 MiB

1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
11 months ago
1 year ago
1 year ago
1 year ago
11 months ago
1 year ago
11 months ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
11 months ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
11 months ago
1 year ago
1 year ago
1 year ago
1 year ago
1 year ago
11 months ago
1 year ago
11 months ago
1 year ago
11 months ago
1 year ago
1 year ago
11 months ago
1 year ago
11 months ago
1 year ago
1 year ago
11 months ago
1 year ago
  1. {
  2. "cells": [
  3. {
  4. "cell_type": "code",
  5. "execution_count": 1,
  6. "metadata": {},
  7. "outputs": [],
  8. "source": [
  9. "import uproot\n",
  10. "import numpy as np\n",
  11. "import sys\n",
  12. "import os\n",
  13. "import matplotlib\n",
  14. "import matplotlib.pyplot as plt\n",
  15. "from mpl_toolkits import mplot3d\n",
  16. "import itertools\n",
  17. "import awkward as ak\n",
  18. "from scipy.optimize import curve_fit\n",
  19. "%matplotlib inline"
  20. ]
  21. },
  22. {
  23. "cell_type": "code",
  24. "execution_count": 2,
  25. "metadata": {},
  26. "outputs": [
  27. {
  28. "data": {
  29. "text/plain": [
  30. "9056"
  31. ]
  32. },
  33. "execution_count": 2,
  34. "metadata": {},
  35. "output_type": "execute_result"
  36. }
  37. ],
  38. "source": [
  39. "file = uproot.open(\"tracking_losses_ntuple_Bd2KstEE.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\")\n",
  40. "#file = uproot.open(\"tracking_losses_ntuple_Dst0ToD0EE.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\")\n",
  41. "\n",
  42. "\n",
  43. "#selektiere nur elektronen von B->K*ee und nur solche mit einem momentum von ueber 5 GeV \n",
  44. "allcolumns = file.arrays()\n",
  45. "found = allcolumns[(allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p > 5e3)] #B: 9056\n",
  46. "lost = allcolumns[(allcolumns.isElectron) & (allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p > 5e3)] #B: 1466\n",
  47. "\n",
  48. "ak.num(found, axis=0)\n",
  49. "#ak.count(found, axis=None)\n"
  50. ]
  51. },
  52. {
  53. "cell_type": "code",
  54. "execution_count": 3,
  55. "metadata": {},
  56. "outputs": [
  57. {
  58. "data": {
  59. "text/plain": [
  60. "0.8606728758791105"
  61. ]
  62. },
  63. "execution_count": 3,
  64. "metadata": {},
  65. "output_type": "execute_result"
  66. }
  67. ],
  68. "source": [
  69. "def t_eff(found, lost, axis = 0):\n",
  70. " sel = ak.num(found, axis=axis)\n",
  71. " des = ak.num(lost, axis=axis)\n",
  72. " return sel/(sel + des)\n",
  73. "\n",
  74. "t_eff(found, lost)"
  75. ]
  76. },
  77. {
  78. "cell_type": "code",
  79. "execution_count": 4,
  80. "metadata": {},
  81. "outputs": [
  82. {
  83. "name": "stdout",
  84. "output_type": "stream",
  85. "text": [
  86. "sample size: 32\n",
  87. "eff (cutoff = 0 ) = 0.96875\n",
  88. "sample size: 32\n",
  89. "eff (cutoff = 50 ) = 0.96875\n",
  90. "sample size: 32\n",
  91. "eff (cutoff = 100 ) = 0.96875\n",
  92. "sample size: 43\n",
  93. "eff (cutoff = 150 ) = 0.9767441860465116\n",
  94. "sample size: 65\n",
  95. "eff (cutoff = 200 ) = 0.9692307692307692\n",
  96. "sample size: 97\n",
  97. "eff (cutoff = 250 ) = 0.9587628865979382\n",
  98. "sample size: 129\n",
  99. "eff (cutoff = 300 ) = 0.9457364341085271\n",
  100. "sample size: 150\n",
  101. "eff (cutoff = 350 ) = 0.9533333333333334\n",
  102. "sample size: 169\n",
  103. "eff (cutoff = 400 ) = 0.9408284023668639\n",
  104. "sample size: 197\n",
  105. "eff (cutoff = 450 ) = 0.9390862944162437\n",
  106. "sample size: 227\n",
  107. "eff (cutoff = 500 ) = 0.920704845814978\n",
  108. "sample size: 257\n",
  109. "eff (cutoff = 550 ) = 0.9260700389105059\n",
  110. "sample size: 297\n",
  111. "eff (cutoff = 600 ) = 0.9326599326599326\n",
  112. "sample size: 334\n",
  113. "eff (cutoff = 650 ) = 0.9281437125748503\n",
  114. "sample size: 366\n",
  115. "eff (cutoff = 700 ) = 0.9289617486338798\n",
  116. "sample size: 400\n",
  117. "eff (cutoff = 750 ) = 0.925\n",
  118. "sample size: 436\n",
  119. "eff (cutoff = 800 ) = 0.9151376146788991\n",
  120. "sample size: 468\n",
  121. "eff (cutoff = 850 ) = 0.9102564102564102\n",
  122. "sample size: 500\n",
  123. "eff (cutoff = 900 ) = 0.912\n",
  124. "sample size: 533\n",
  125. "eff (cutoff = 950 ) = 0.9136960600375235\n",
  126. "sample size: 562\n",
  127. "eff (cutoff = 1000 ) = 0.9163701067615658\n",
  128. "\n",
  129. "sample size: 150\n"
  130. ]
  131. },
  132. {
  133. "data": {
  134. "text/plain": [
  135. "0.9533333333333334"
  136. ]
  137. },
  138. "execution_count": 4,
  139. "metadata": {},
  140. "output_type": "execute_result"
  141. }
  142. ],
  143. "source": [
  144. "#finden wir die elektronen die keine bremsstrahlung gemacht haben mit hoher effizienz?\n",
  145. "#von energie der photonen abmachen\n",
  146. "#scan ab welcher energie der photonen die effizienz abfällt\n",
  147. "\n",
  148. "#abhängigkeit vom ort der emission untersuchen <- noch nicht gemacht\n",
  149. "\n",
  150. "\n",
  151. "\n",
  152. "#idea: we make an event cut st all events that contain a photon of energy > cutoff_energy are not included\n",
  153. "\"\"\"\n",
  154. "ph_e = found[\"brem_photons_pe\"]\n",
  155. "event_cut = ak.all(ph_e<cutoff_energy,axis=1)\n",
  156. "ph_e = ph_e[event_cut]\n",
  157. "\"\"\"\n",
  158. "\n",
  159. "\n",
  160. "\n",
  161. "for cutoff_energy in range(0,1050,50):\n",
  162. "\tnobrem_f = found[ak.all(found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  163. "\tnobrem_l = lost[ak.all(lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  164. "\tprint(\"sample size: \",ak.num(nobrem_f,axis=0)+ak.num(nobrem_l,axis=0))\n",
  165. "\tprint(\"eff (cutoff = \",str(cutoff_energy),\") = \",str(t_eff(nobrem_f,nobrem_l)))\n",
  166. "\n",
  167. "\"\"\"\n",
  168. "we see that a cutoff energy of 350MeV is ideal because the efficiency drops significantly for higher values\n",
  169. "\"\"\"\n",
  170. "cutoff_energy = 350.0 #MeV\n",
  171. "\n",
  172. "\"\"\"\n",
  173. "better statistics: cutoff=350MeV - sample size: 150 events and efficiency=0.9533\n",
  174. "\"\"\"\n",
  175. "nobrem_found = found[ak.all(found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  176. "nobrem_lost = lost[ak.all(lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  177. "\n",
  178. "print(\"\\nsample size: \",ak.num(nobrem_found,axis=0)+ak.num(nobrem_lost,axis=0))\n",
  179. "t_eff(nobrem_found, nobrem_lost)"
  180. ]
  181. },
  182. {
  183. "cell_type": "code",
  184. "execution_count": null,
  185. "metadata": {},
  186. "outputs": [],
  187. "source": []
  188. },
  189. {
  190. "cell_type": "code",
  191. "execution_count": 5,
  192. "metadata": {},
  193. "outputs": [
  194. {
  195. "name": "stdout",
  196. "output_type": "stream",
  197. "text": [
  198. "31\n"
  199. ]
  200. },
  201. {
  202. "data": {
  203. "text/plain": [
  204. "0.96875"
  205. ]
  206. },
  207. "execution_count": 5,
  208. "metadata": {},
  209. "output_type": "execute_result"
  210. }
  211. ],
  212. "source": [
  213. "#hier wird ohne rücksicht auf energie der photonen getrennt\n",
  214. "\n",
  215. "nobrem_found = found[found[\"brem_photons_pe_length\"]==0]\n",
  216. "nobrem_lost = lost[lost[\"brem_photons_pe_length\"]==0]\n",
  217. "\n",
  218. "\"\"\"\n",
  219. "die effizienz mit der wir elektronen finden, die keine bremsstrahlung gemacht haben, ist gut mit 0.9688.\n",
  220. "allerdings haben wir hier nur sehr wenige teilchen (<100)\n",
  221. "\"\"\"\n",
  222. "print(ak.num(nobrem_found,axis=0))\n",
  223. "t_eff(nobrem_found, nobrem_lost)\n"
  224. ]
  225. },
  226. {
  227. "cell_type": "code",
  228. "execution_count": 6,
  229. "metadata": {},
  230. "outputs": [],
  231. "source": [
  232. "#wie viel energie relativ zur anfangsenergie verlieren die elektronen durch bremstrahlung und hat das einen einfluss darauf ob wir sie finden oder nicht?\n"
  233. ]
  234. },
  235. {
  236. "cell_type": "code",
  237. "execution_count": 7,
  238. "metadata": {},
  239. "outputs": [
  240. {
  241. "data": {
  242. "text/plain": [
  243. "0.8603431839847474"
  244. ]
  245. },
  246. "execution_count": 7,
  247. "metadata": {},
  248. "output_type": "execute_result"
  249. }
  250. ],
  251. "source": [
  252. "#keine rücksicht auf energie der photonen\n",
  253. "brem_found = found[found[\"brem_photons_pe_length\"]!=0]\n",
  254. "energy_found = ak.to_numpy(brem_found[\"energy\"])\n",
  255. "eph_found = ak.to_numpy(ak.sum(brem_found[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  256. "energyloss_found = eph_found/energy_found\n",
  257. "\n",
  258. "\n",
  259. "brem_lost = lost[lost[\"brem_photons_pe_length\"]!=0]\n",
  260. "energy_lost = ak.to_numpy(brem_lost[\"energy\"])\n",
  261. "eph_lost = ak.to_numpy(ak.sum(brem_lost[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  262. "energyloss_lost = eph_lost/energy_lost\n",
  263. "\n",
  264. "t_eff(brem_found,brem_lost)"
  265. ]
  266. },
  267. {
  268. "cell_type": "code",
  269. "execution_count": 8,
  270. "metadata": {},
  271. "outputs": [
  272. {
  273. "name": "stdout",
  274. "output_type": "stream",
  275. "text": [
  276. "mean energyloss relative to initial energy (found): 0.6475128752780828\n",
  277. "mean energyloss relative to initial energy (lost): 0.8241268441538472\n"
  278. ]
  279. }
  280. ],
  281. "source": [
  282. "mean_energyloss_found = ak.mean(energyloss_found)\n",
  283. "mean_energyloss_lost = ak.mean(energyloss_lost)\n",
  284. "print(\"mean energyloss relative to initial energy (found): \", mean_energyloss_found)\n",
  285. "print(\"mean energyloss relative to initial energy (lost): \", mean_energyloss_lost)"
  286. ]
  287. },
  288. {
  289. "cell_type": "code",
  290. "execution_count": 9,
  291. "metadata": {},
  292. "outputs": [
  293. {
  294. "data": {
  295. "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAHOCAYAAABttoiYAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABIp0lEQVR4nO3dTYwbaX7n+V/a3tIIK0iRKV/mMEApaPtgLNClYOrgxgDVbpGusztJydgx4EN3klXIW8FFOudSWYftLKZ1TbjJrAZmgcYCSlK12GObFFyFma2DlYxSX4xBuxkpowbwHFyZUYJ2ZNV0N/eQ9UTzncEgmcxIfj8AISUj4oknHgYj/nzeYqXdbrcFAAAQU7+16AwAAABMg2AGAADEGsEMAACINYIZAAAQawQzAAAg1ghmAABArBHMAACAWCOYAQAAsUYwA8yA67qLzsK58n1flUpFnuctOiuXjud5qlQq8n1/0VnBgi3bdWUaBDNj+L6vvb09raysaGVlRaurq8pms0qn00okEtrb21t0FgOVSkWJRCLIazqdVqPRkHT2pUin08Gy3nwXi0VuTCFls9mgHM2rWCwOXNd1XWWzWSUSCa2uriqZTCqdTqtYLKpSqSiZTE6070aj0fUZJ5NJ1Wq1vvVqtZqSyWSwziwvipVKRbdu3VI+n+eGG0E+nx/6ue/t7SmRSCifz+vk5OScc4ZZ8zxP+Xxe+Xxe2Wx25Hcm7HWlWCwqm80G96Fh1+2w610abYTiOE5bUrtcLgfvVavVtqR2JpNZYM66lUqltqS2bdt9y1KpVNtxnPbp6WnwnjmeUqnUbrVa7Waz2a7X6+eV3dg5PT1t27bdzmQyXa9ms9m3bqFQCM6PVqvVlUYul2tLakf5Cp6engbbdp6PvcrlctuyrInTD8Mc26DjjqrVanWdm5eVZVkjP7dMJtOW1HXOXHaX8bNvNpt935FCodC2LKvvsw17XXEcp+t+U6/X25ZlRV7vMiGYCcmyrIEXmKg3pHkxN8nei2Umk2nncrm+9ZvNZjuXy7VzuVy7VCq1C4XCpbuozFKhUGhXq9Wx65kbUqlUGrpOuVyOfO6YYCKVSg1dx3ym82CC5lleHB3HufQ38Hq93pY08jtmvsOXvSw6XcbP3nGctuM4fe8P+gEc5rpivvO9507vfsKud9lcnLvwBWYi7EG1HSaYmWUAYGpJorBtu+tC2Gq12o7jjP0Fn0ql2oVC4dJE7qenp2OPOwrLstqlUmlkOZkauzAXjqgXl87amWF5sSxrboHprIOZZamNyGQyYz/zZQtm4vbZV6vVscFXq9VqSxr4AzKVSvUdb5jrimVZA+9Bvd/FsOtdNvSZCcH0O0mlUl3vm74KjuPIsqyZ7c+2bZVKpYm3831fnufJsizZtq1araZ0Oq2DgwPlcrmB22SzWdm2rXQ6re3tbT18+FCVSmXaQ1g4y7L0+PFjtVqtmfVt2tvbk+/7KhaLSiaTWl1dHdhfxbRz5/P5sWlWq9W+90w/G9Mva1C7uWVZymQykqTd3d2+5bVaTalUaqLz0nTqTSaTajQawf9NP7FBbf0nJyfBeisrK8pmswPTzefzKhaLSqfTXX25TF5Nnx7Tt6Czj8+47TvzXavV1Gg0Qucnn88rkUiMPecbjUbQj8Gk6bpu0H+pt1+S6f/Q20+hVquFOi8687m6uqrV1dW+7cz323xWveuMOo9MmaXT6aAjdzqd1urqqtLpdPBZmz48q6urQ/uFDeK6rlZXV4Oy6SyHYrEY9Nsb9tmb/Jgyr9Vq8n0/yOMk3+dEIqF0Oh30WensmzJJXzLTJ7Fer+vx48eybXvk8Usa+P0z25l1wlxXXNeV7/tyHKcvPfPew4cPQ68nTX8OzPscmtiio6k4MJF0ZzVg56/vefyiGBRZj2OaLVKpVDuTybRTqVToX+Zha2Wq1Wo7l8u1C4XC0GMPs04Y1Wq1nUql2pZltavVartcLgd9lybtp1Qqldq2bU/VjHZ6etquVqvtQqEQ1ICp59eX+UUmKVLfo2az2dV0ZM6zQb/wOvc1qEp50v139uMxNXXVajX45dx5TppfealUKvhFabbvbNpqNpt9bfXmPO1cz1SN954rYbZvtVpBHjtrGAflp90+qwkoFApd6YVpjjN57CxX02zUmZ5Js/c981mOO/9MvlOpVDuXy3Wd9+bcqFarwTnY+V0ztT7jzqNWq9W1H1Nm5njMvuv1elf5TvKrflifwlar1ZW3YZ99u91f05xKpSauWej97phzd9B3apAo1w5TjoOuU+Z4zTkX5rpizvlhXQXMvsKu125Pfw6cxzk0CYKZEMzJZS4WktqWZc21o6y5kU/CnCymf88k/SXCBBy5XK7ryzmoI2OYdSZhviydwZa5GPTeLMIol8tt27bbuVxu6iYYc7HuvMGZL/Kom1az2WyXSqWul9necZy+L7v5PAelZwLtzrJotVqRO/6aC33vZ9Yb0A9br/OGa45n0HlsvkfmvBt2Qwu7vfkses/53vy022fl2XvuhPmumOCx9wY1qFq/t9O3eS/Md9qc8719KMyNzpwr5jMY9D0Icx6Zm1vvja+3bDvXnbQPlslzp1Kp1HVso4KZ3ht1lO9853XaNM+Oa4I9PT0Ngoso/c46m4F7meMddl0cdF0Z9Vmb89JxnNDrGdOeA+dxDoVFM9MYpjrbcRw1m83gJamvunuWMplMUKUedvirycvjx48lnVXnht12VJWpdHash4eHOjg4kKSg2rez6S3MOpM6OjqSJN2/fz+oJr1//74kRRpqmMvl1Gq1lE6ndffuXeXz+chDFjOZjOr1uiSpXC6H3s5xHOVyORWLxaDaNZVKyfM8ua6r3d3dYEhlZzOJKYtOZvvOZpJSqaTt7e1Ix2T0ng+m+cIcr7G+vt71t2VZQXma4xlU5W3SG1VuUbYfVK3fO8TZtm3t7e11NVUUCoWh+ejcznGcgU2LnucFx+37vnzf7yvDWq02sNlrmN7j7v0MzLHeuXOnLy+TnEe9ZWbyvba21vdeq9UKnf/OPHeW2cOHD4Mm0nHMd6VWq6lUKkVqfu+8/mxubko6+44MOlfMUOpkMqlEIqFWqxXq3OhlWVbQtN9Z9qZJUBp+zR10XQnTXOz7fuj1BuW306TnwDzPobAIZsYwJ1XnF8JxnOCGPe7iZNrko7zK5bIajYZWV1fHBk2mrdRccM0XcBZtlJVKRY1GQwcHB7IsS5VKRcViUfV6PThBw6wTheu6ymQyXRd2c9Po/KJMKpPJqNlsdl3so8yZkkql5DhOkKfOm/uoIMmyrOCYzL+mDb1arXa9Tk9P1W63BwaFqVRKtm0H7deSdHh4OLSPVFQmj5MEfqP6I5hyGpXetNsPU61WZVmWisWiEonERP0mem/OlUpF9+7dk/SbG8/h4WHfdcGsb9aNwnz+vcfceyOJch7NizkPTdm4rtsXAI8TJYAZpNFoqFarBQFSr845WVqt1tTfoXK5rFKpJNd1u/oOmevWqHLova6Ya+iga5R5z7bt0OtdRgQzY5ggIp1Od71vLgjmV9gw5XJZrVYr0qtcLstxHNXr9bEXoN5OyuaXxyxmaTUXcM/zguCo1Wp15SnMOpMyN4Dejo9PnjyRND6QnMQ0k7/Zth1coDpvLINqUgYx23bWaEzClHepVFKtVtP6+vpMO6R35jHKhXDUL8EwAem02/eybVvHx8dBbVgymQzd6b335mxuWKlUKkijXC733QgfPnw4cYfsXmbbcZ9B1PNoHkxH9UajIc/z9PDhw9AdoA0zqMHzvKk68pv9Dup0Py+FQkGtVkvtdlvNZlPr6+vyPE+ZTGbsudB5XTGf+aCJFM17juOEXu8yIpgZwff94FdO702580Ix6xuHdPYLplgs6vHjx6ECAlOD1Bl0mV80k148OpkgKZfLqVAoqFQq9V2ow6wTxaBaMens13AqlZo6UEomkyqXy6pWq6rX65E/x97RA6bcJ/1FaS5Eg5oxJA2tncvlcsHFfnNzc6rPexhzIZxkxmJTJoPybQKURCIxt+2HMTfHer0e3NgmKTNzc+4MHM3Mrnt7ewN/cU/
  296. "text/plain": [
  297. "<Figure size 640x480 with 1 Axes>"
  298. ]
  299. },
  300. "metadata": {},
  301. "output_type": "display_data"
  302. }
  303. ],
  304. "source": [
  305. "#in abhängigkeit von der energie der elektronen\n",
  306. "\n",
  307. "plt.hist(energyloss_lost, bins=200, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=\"lost\")\n",
  308. "plt.hist(energyloss_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  309. "#plt.xticks(np.arange(0,1.1,0.1), minor=True,)\n",
  310. "#plt.yticks(np.arange(0,10,1), minor=True)\n",
  311. "plt.xlabel(r\"$E_\\gamma/E_0$\")\n",
  312. "plt.ylabel(\"counts (normed)\")\n",
  313. "plt.title(r\"$B\\rightarrow K^\\ast ee$, $p>5$GeV, photons w/ brem_vtx_z$<9500$mm\")\n",
  314. "plt.legend(title=\"LHCb Simulation\", title_fontsize=15)\n",
  315. "#plt.grid()\n",
  316. "\n",
  317. "\"\"\"\n",
  318. "found: elektronen verlieren durchschnittlich 0.65 ihrer anfangsenergie durch bremsstrahlung\n",
  319. "lost: elektronen verlieren durchschnittlich 0.82 ihrer anfangsenergie durch bremsstrahlung\n",
  320. "\n",
  321. "-> wir können sofort erkennen, dass verlorene elektronen im schnitt mehr energie durch bremsstrahlung verlieren als gefundene, \n",
  322. "aber auch die rate der gefundenen elektronen steigt für raten nahe 1, wenn auch wesentlich schwächer als für verlorene elektronen.\n",
  323. "die meisten verlorenen elektronen verlieren >0.8 ihrer anfangsenergie.\n",
  324. "\"\"\"\n",
  325. "\n",
  326. "plt.show()"
  327. ]
  328. },
  329. {
  330. "cell_type": "code",
  331. "execution_count": 10,
  332. "metadata": {},
  333. "outputs": [],
  334. "source": [
  335. "#ist die shape der teilspur im scifi anders? (koenntest du zum beispiel durch vergleich der verteilungen der fit parameter studieren,\n",
  336. "#in meiner thesis findest du das fitmodell -- ist einfach ein polynom dritten grades)\n",
  337. "z_ref=8520 #mm\n",
  338. "\n",
  339. "def scifi_track(z, a, b, c, d):\n",
  340. " return a + b*(z-z_ref) + c*(z-z_ref)**2 + d*(z-z_ref)**3\n",
  341. "\n",
  342. "def z_mag(xv, zv, tx, a, b):\n",
  343. " \"\"\" optical centre of the magnet is defined as the intersection between the trajectory tangents before and after the magnet\n",
  344. "\n",
  345. " Args:\n",
  346. " xv (double): velo x track\n",
  347. " zv (double): velo z track\n",
  348. " tx (double): velo x slope\n",
  349. " a (double): ax parameter of track fit\n",
  350. " b (double): bx parameter of track fit\n",
  351. "\n",
  352. " Returns:\n",
  353. " double: z_mag\n",
  354. " \"\"\"\n",
  355. " return (xv-tx*zv-a+b*z_ref)/(b-tx)"
  356. ]
  357. },
  358. {
  359. "cell_type": "code",
  360. "execution_count": 11,
  361. "metadata": {},
  362. "outputs": [],
  363. "source": [
  364. "scifi_found = found[found[\"scifi_hit_pos_x_length\"]>3]\n",
  365. "scifi_lost = lost[lost[\"scifi_hit_pos_x_length\"]>3]\n",
  366. "#should be fulfilled by all candidates\n",
  367. "\n",
  368. "scifi_x_found = scifi_found[\"scifi_hit_pos_x\"]\n",
  369. "scifi_z_found = scifi_found[\"scifi_hit_pos_z\"]\n",
  370. "\n",
  371. "tx_found = scifi_found[\"velo_track_tx\"]\n",
  372. "\n",
  373. "scifi_x_lost = scifi_lost[\"scifi_hit_pos_x\"]\n",
  374. "scifi_z_lost = scifi_lost[\"scifi_hit_pos_z\"]\n",
  375. "\n",
  376. "tx_lost = scifi_lost[\"velo_track_tx\"]\n",
  377. "\n",
  378. "xv_found = scifi_found[\"velo_track_x\"]\n",
  379. "zv_found = scifi_found[\"velo_track_z\"]\n",
  380. "\n",
  381. "xv_lost = scifi_lost[\"velo_track_x\"]\n",
  382. "zv_lost = scifi_lost[\"velo_track_z\"]\n",
  383. "\n",
  384. "\n",
  385. "\n",
  386. "#ak.num(scifi_found[\"energy\"], axis=0)\n",
  387. "#scifi_found.snapshot()"
  388. ]
  389. },
  390. {
  391. "cell_type": "code",
  392. "execution_count": 12,
  393. "metadata": {},
  394. "outputs": [],
  395. "source": [
  396. "#tx_lost.show()"
  397. ]
  398. },
  399. {
  400. "cell_type": "code",
  401. "execution_count": 13,
  402. "metadata": {},
  403. "outputs": [],
  404. "source": [
  405. "scifi_fitpars_found = ak.ArrayBuilder()\n",
  406. "\n",
  407. "for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
  408. " popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]))\n",
  409. " scifi_fitpars_found.begin_list()\n",
  410. " scifi_fitpars_found.real(popt[0])\n",
  411. " scifi_fitpars_found.real(popt[1])\n",
  412. " scifi_fitpars_found.real(popt[2])\n",
  413. " scifi_fitpars_found.real(popt[3])\n",
  414. " scifi_fitpars_found.end_list()\n",
  415. "\n",
  416. "scifi_fitpars_lost = ak.ArrayBuilder()\n",
  417. "\n",
  418. "for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
  419. " popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]))\n",
  420. " scifi_fitpars_lost.begin_list()\n",
  421. " scifi_fitpars_lost.real(popt[0])\n",
  422. " scifi_fitpars_lost.real(popt[1])\n",
  423. " scifi_fitpars_lost.real(popt[2])\n",
  424. " scifi_fitpars_lost.real(popt[3])\n",
  425. " scifi_fitpars_lost.end_list()\n",
  426. "\n",
  427. "\n",
  428. "scifi_fitpars_lost = scifi_fitpars_lost.to_numpy()\n",
  429. "scifi_fitpars_found = scifi_fitpars_found.to_numpy()\n",
  430. "\n",
  431. "\n",
  432. "\n",
  433. "dtx_found = scifi_fitpars_found[:,1] - tx_found\n",
  434. "dtx_lost = scifi_fitpars_lost[:,1] - tx_lost\n"
  435. ]
  436. },
  437. {
  438. "cell_type": "code",
  439. "execution_count": null,
  440. "metadata": {},
  441. "outputs": [],
  442. "source": []
  443. },
  444. {
  445. "cell_type": "code",
  446. "execution_count": 14,
  447. "metadata": {},
  448. "outputs": [
  449. {
  450. "data": {
  451. "image/png": "iVBORw0KGgoAAAANSUhEUgAABNMAAANaCAYAAACjiQrWAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAAC5IUlEQVR4nOz9T4wj+X3nfX6yLVQ5sYVqVrYuPgzQFTQGsE7dwayL5sFUjYrhBuZmJdkF7HgxwDxKcoS89UoZSsyhsy+iSPeVEIIpD7xYPQtUkanFXgawGG13zzwjPDuVDJcvMh77YVT78MAPsO7MUKH05JRGI+6hFNFkJsmMJIMM/nm/AKI7Gf++DAbJb33j92ej1+v1BAAAAAAAAOBKb6QdAAAAAAAAALAsKKYBAAAAAAAAMVFMAwAAAAAAAGKimAYAAAAAAADERDENAAAAAAAAiIliGgAAAAAAABATxTQAAAAAAAAgJoppAAAAAAAAQEwU0wBgSp7nqdFopB3G2nJdV+VyWa1WK+1QAABASnzfV6PRUBAEie6XPC9d5HlYVBTTgBmxbVu1Wk3FYlF37tyRbdtph4SE+b6vYrGoXC4nx3HSDudKvu/rzp07YxPCZbtuW62WbNtWo9HQ6elp2uEAANbAsv1WroNaraZsNqtyuZxYPkCelz7yPCwyimlAwoIgUDab1aNHj7S/v69ms6lqtapWqzX0Tpnv+4nfQVtGy3geDMNQs9lMO4zYgiBQEATqdrtDl13nuh0mjfewUCjo4OBgrscEAKwncrzJzOM87O/vq1AoJLpP8rxB5HnAIIppQMJs21Ymk5FpmtFzpVJJ3W5XmUzm0vrFYpE7LeI8zINpmur1eqpWq5eWXfe6HSat9zBufAAATIMcbzLzOg9bW1szP8YiI88D5otiGpAw13Vjr1ssFuV53gyjWQ6ch/Rd57odhvcQALDqyPGuj/OwGMjzgORRTFtCuVyOJuMLyPM8FYtF+b4fjbFgWZYsy9KdO3e0sbExMIZBq9WKfpTK5XL0IxUEgRqNhnK5nFzXjf7/zp07KhaLl977Wq2mcrks27aVy+UujZPQarVkWVa0rzt37qhcLsfaPozFsiw1Gg35vh+9HsuyoljCcSpGjb0QnhvLspTNZgfWGXUe4mwb5/UNc9U5GycIApXL5ehRq9WGrpd03OFxbdtWuVxWNpu9FHd/bOG1138uW62WisWiisXipTjjXrfDjHsPp7n+4r6uUetvbGyoWCwyYC0AYCrkeMuT44WvrVwu686dO7G3ubgted7gayHPA4boYak4jtOT1Ds7O0s7FIwgqWea5sBz4fvmOM7A8/v7+z1JvW63Gz1XKpV6knqSevl8vre/v99rNpu9QqHQk9QzDOPS9qF2u92T1Gu3271er9drNps9wzB6knqlUqm3v7/fM00ziu+q7bvdbhRPGEun04nWy+fzvVKp1Gu3271utxvF2Ol0on12Op1ePp+P/m42m1E8485DnG2ven3DXPWax+l2u71MJjOwbrVavfSezyLuQqHQ29/fj/52HKdXrVYHYjMMY+AcZjKZXiaTiZaHsfbHFrrOdTvMsPdw2usvzuvq364/znw+32s2m1fGDQBAXOR4i53j9b+eUqnUcxynZ5rmyNznIvK80cjzgMsopi2Rs7MzimlLYNiPVfgjGyfR6vW+/OG+uH4+n+9Jin488vn8wI/N2dlZT9LAj3G4r/7n+vd31fadTudSYtTr9aLkpD/2cN3+H3/TNAcSr17v9Y9k/3U86jzE2Xbc6xsmzmsepVAoxEpQZhF3JpO5tO7F89z/d/8x+r8vrpNkjbpuh7nqWp70+ovzui4mWYVC4dL5BwBgWuR4i53jhcW0i0WWsOBz1Y1T8rzRyPOAy74iTCwIAlUqFUkaOtCj53mqVCoyDENBEMiyrKlmmalUKqpWq9dqqozlZhjGwN/lclmu66rdbqtQKKjZbA4MBHpyciJJA90EwkE77927d2n/cba/uJ/+2DzPGxjsNYw3nEXI9/3oczDMycmJ8vn80GVxtx33+oa5zmu+GE+r1Rr6WZ9H3IZhqFar6a233tL+/r4kRf8Nj3l0dDSwzf7+frROWqa5/q77usLvWcdxLn12AABYJOR4yedKof5B9iVF3Qzb7fbYmMjzro88D+uMYtqEXNeV4zhqtVoqlUqXlvu+r1wup06nE32hZ7NZnZ6eDl0/zvEePXo0ddxYbuG15Pu+pNc/YJlMRq1WS48fPx77gz1sJpzrbD+JcLyDSaYVv+62cWf6mfQ1h+f8qh/vWcXdbDaVy+Vk27Ycx1Gz2Yyuh/CYizzb0STX33Vfl+M48n1frutO9D0LAEBayPFGmza/CQto4bkdhjxvOuR5WEdMQDChfD4/9ku0XC4rn88P3BkJB5ScRLvdvnSXBesnvEMY/tCHRVvf99VsNq99d2ra7ePsv/+/89r2qv1O8prDOK6aEnxWcRuGoefPnyufz0evIRzEdVbHnLWr3ovrvq5yuSzTNFUul5fuXAAA1hs53uyExZpxhTLyvOSR52HVUUybgSAI5LquLMsaeH57e1uSoi/GRqMh27ZHPsIpjGu1mg4ODub7IrCQwh/4XC4nSbIsS1tbWxMnSNNuf5UwaRk1y864abqn2XacSV9zGE+n04m1XtJx+76vTCajdrsdFfLD4nxYaB9V4F/UhOOq9+K6ryuTyUTrXvz+BQBgkZHjxdt2EmG3wnGt88jzkkeeh1W31sW0IAgGpg6+yHXdoVNAXyXsD37x7kf4hdFutyVJpVJJ1Wp15CNskvz48WPdvXs3mt5Zku7evTtyqmYsn6vG6wq1Wi1lMhmVSqVoiuv+ptHhfq66qyZp6u3jCK9h27YvTXE9amrsSbeNY5rX3F8MH/Z+hc/NIm5pcFzGQqEgx3EkvX5N/bFdTOJs2x4Y82TW4l7Lcd6LSV6XYRhRNwDGlwQApI0cb/D4k2w7Ldd1ZRjG2LGryfPiIc8DvrTWxbRMJqNyuTy0su26blTUuq7+sQ7GLY+r0+no7OwsekjS8+fPUx9wEqON+qG5+Hw2m5X05RgAF+9yhT+k4baO40SDdIY/MK1WS41GI2rpKL0eg6DVaikIguiYF48dd/tRSdewpCz8//C/mUwmuk5zuZyKxaJqtZosy1K3240SkmHnIe62o17fMHFf8zAX43FdV77vR9v7vq9arTaTuCXpyZMnA98dQRDIMAwZhqFMJhN9V1mWpWKxKNu2lcvllM1mo++iqxLpuNftMKOu5WmuP0nXel3hf0ulkvL5vBqNBjcdAACJIsdbzBwvjEka/LeW7/vRGGRXbUueNxp5HjBE2tOJLoJ2uz0whfDFv8fRkOmkw6mDh03ZK6lnGMZU8erCFMhYDGdnZ9F7r99OHX52dtZrt9vRFOOGYVyalts0zV4mkxm4jvqnmTZNs1coFHqFQuHSto7j9DKZTM8wjGia6VKpFE2v3Ww2o+nADcO4NPX1Vdt3Op1oqvZMJhNNNe44TvQ6w+mpu91ur1AoROv2H6tarY6NY9R5uGrbq17fMFe95jjbh8c0TbPX7XajffVPF5503Pl8PjrO/v5+r1AoXPoeaDab0bVmmubA9dLpdKL3p//6nPS6Hebiezjt9RfndfUvMwwjukZLpdLANXpxKncAAOIix1uOHC98n/L5fC+fz/dKpVKvVCpd699N5HmjkecBgzZ6vV4vudLc8gpbopXLZTmOE3XFvMrGxoZKpdLA3aVGo6FyuTx0+uWNjQ3l8/nY+8d6qtVqV07hDQAAgOVCjgcAq+EraQewKMICV7FYjLpSTiocK21Uk9mrplwGAAAAAADAYlrrMdP6ua4rz/PUbDbHTkoQRziY4sWx0cK/w1l6AAAAAAAAsFwopunLLp7tdluFQmHkpARxZTIZmaZ5qStnOEvJ+++/P1W8WG1BEETXzlWDpQIAAGA5kOMBwOpY+zHT+gtp/Vqt1pVjpwVBoDt37lwaM016PUtJLpdTt9uNunVms1mVy2Vm4cRYw2al4ZoBAABYbuR4ALA6VqaYZllW1PLrolEDfAZBoN3d3ZF3hlqtlp4
  452. "text/plain": [
  453. "<Figure size 1500x1000 with 4 Axes>"
  454. ]
  455. },
  456. "metadata": {},
  457. "output_type": "display_data"
  458. }
  459. ],
  460. "source": [
  461. "fig, ((ax0, ax1), (ax2, ax3)) = plt.subplots(nrows=2, ncols=2, figsize=(15,10))\n",
  462. "\n",
  463. "ax0.hist(scifi_fitpars_found[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$a_x$ found\")\n",
  464. "ax0.hist(scifi_fitpars_lost[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$a_x$ lost\")\n",
  465. "ax0.set_xlabel(\"a\")\n",
  466. "ax0.set_ylabel(\"normed\")\n",
  467. "ax0.set_title(\"fitparameter a der scifi track\")\n",
  468. "ax0.legend()\n",
  469. "\n",
  470. "ax1.hist(scifi_fitpars_found[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$b_x$ found\")\n",
  471. "ax1.hist(scifi_fitpars_lost[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$b_x$ lost\")\n",
  472. "ax1.set_xlabel(\"b\")\n",
  473. "ax1.set_ylabel(\"normed\")\n",
  474. "ax1.set_title(\"fitparameter b der scifi track\")\n",
  475. "ax1.legend()\n",
  476. "#evtl multiple scattering candidates (lost); findet man einen gewissen endvtx_type (mult scattering)\n",
  477. "#steiler velo winkel (eta)? vertex type? evtl bremsstrahlung?\n",
  478. "\n",
  479. "\n",
  480. "ax2.hist(scifi_fitpars_found[:,2], bins=500, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$c_x$ found\")\n",
  481. "ax2.hist(scifi_fitpars_lost[:,2], bins=500, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$c_x$ lost\")\n",
  482. "ax2.set_xlim([-3e-5,3e-5])\n",
  483. "ax2.set_xticks(np.arange(-3e-5,3.5e-5,1e-5),minor=False)\n",
  484. "ax2.set_xlabel(\"c\")\n",
  485. "ax2.set_ylabel(\"normed\")\n",
  486. "ax2.set_title(\"fitparameter c der scifi track\")\n",
  487. "ax2.legend()\n",
  488. "\n",
  489. "ax3.hist(scifi_fitpars_found[:,3], bins=500, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$d_x$ found\")\n",
  490. "ax3.hist(scifi_fitpars_lost[:,3], bins=500, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$d_x$ lost\")\n",
  491. "ax3.set(xlim=(-5e-8,5e-8))\n",
  492. "ax3.text(-4e-8,3e8,\"d negligible <1e-7\")\n",
  493. "ax3.set_xlabel(\"d\")\n",
  494. "ax3.set_ylabel(\"normed\")\n",
  495. "ax3.set_title(\"fitparameter d der scifi track\")\n",
  496. "ax3.legend()\n",
  497. "\n",
  498. "\"\"\"\n",
  499. "a_x: virtual hit on the reference plane\n",
  500. "\"\"\"\n",
  501. "\n",
  502. "plt.show()"
  503. ]
  504. },
  505. {
  506. "cell_type": "code",
  507. "execution_count": 15,
  508. "metadata": {},
  509. "outputs": [
  510. {
  511. "name": "stdout",
  512. "output_type": "stream",
  513. "text": [
  514. "found\n",
  515. "a = -0.6718207391527037\n",
  516. "b = 0.0013778237292529144\n",
  517. "c = 3.3126998287416195e-08\n",
  518. "d = -1.0330674442255529e-10\n",
  519. "lost\n",
  520. "a = -36.98764338200992\n",
  521. "b = -0.015685137956233643\n",
  522. "c = -8.265859479503501e-07\n",
  523. "d = -1.541510766903436e-11\n"
  524. ]
  525. }
  526. ],
  527. "source": [
  528. "print(\"found\")\n",
  529. "print(\"a = \", str(np.mean(scifi_fitpars_found[:,0])))\n",
  530. "print(\"b = \", str(np.mean(scifi_fitpars_found[:,1])))\n",
  531. "print(\"c = \", str(np.mean(scifi_fitpars_found[:,2])))\n",
  532. "print(\"d = \", str(np.mean(scifi_fitpars_found[:,3])))\n",
  533. "\n",
  534. "print(\"lost\")\n",
  535. "print(\"a = \", str(np.mean(scifi_fitpars_lost[:,0])))\n",
  536. "print(\"b = \", str(np.mean(scifi_fitpars_lost[:,1])))\n",
  537. "print(\"c = \", str(np.mean(scifi_fitpars_lost[:,2])))\n",
  538. "print(\"d = \", str(np.mean(scifi_fitpars_lost[:,3])))"
  539. ]
  540. },
  541. {
  542. "cell_type": "code",
  543. "execution_count": 16,
  544. "metadata": {},
  545. "outputs": [
  546. {
  547. "data": {
  548. "text/plain": [
  549. "-4.6785491318157854e-07"
  550. ]
  551. },
  552. "execution_count": 16,
  553. "metadata": {},
  554. "output_type": "execute_result"
  555. }
  556. ],
  557. "source": [
  558. "np.min(scifi_fitpars_found[:,3])"
  559. ]
  560. },
  561. {
  562. "cell_type": "code",
  563. "execution_count": 17,
  564. "metadata": {},
  565. "outputs": [
  566. {
  567. "data": {
  568. "image/png": "iVBORw0KGgoAAAANSUhEUgAABQYAAAImCAYAAAABqFcVAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd5wd5WHv/8/pvZ89Z3sv0q7qripCgBrdFCOBwXbsxLGIr3OdOM5FdpL7S7+2cHq7RnZyE9vYgGRMB6EVqKCu3VXZ3vvu2dN7n/n9oUhGIEBgQBLM+/XihfacmWeemWdW89UzM88jE0VRRCKRSCQSiUQikUgkEolEIpF8qsivdAUkEolEIpFIJBKJRCKRSCQSycdP6hiUSCQSiUQikUgkEolEIpFIPoWkjkGJRCKRSCQSiUQikUgkEonkU0jqGJRIJBKJRCKRSCQSiUQikUg+haSOQYlEIpFIJBKJRCKRSCQSieRTSOoYlEgkEolEIpFIJBKJRCKRSD6FpI5BiUQikUgkEolEIpFIJBKJ5FNI6hiUSCQSiUQikUgkEolEIpFIPoWkjkGJRHJNaG9vZ8eOHVe6Gu+qtbWVhx9+mF27dl3pqnyihEIhduzYwfDw8JWuikQikUgkkk8gKcN9NKQMJ5FcG6SOQYlEQnt7O9u2bWPLli3U1NTw6KOPXukqXTA8PMyWLVtoaWnhscceu9LVeUe7du1i27Zt7Nixg0AgcKWrA3z87To8PIzNZnvXDtxt27bx6KOPsmXLFmw2G9u2bXvXMnfs2EFVVRUPP/wwoVDo196+RCKRSCSSc87nhJaWFlpaWq50da6YT3uG+yjyG0gZTiK5liivdAUkEsmV1d7ezoYNGwgGg8C5C//Q0NAVrtWvVFdXs3PnTmQy2WWvMzw8jN1ux2q1fnQVe4vNmzcDsGXLlo9tm+/mSrRrKBQiFApdcjuhUIiWlhZ27txJc3MzcC4wbt++ne985zvv2FZbt25laGjosgLxu23/SpwTEolEIpFczc5fjx999NELf/4ovJ9rsJThPv4M91HkN5AynERyLZE6BiWST7nvfve72O32Cz9v3779Ctbmw7FlyxZ27tz5sQeIqymwXIl2bW5uRhTFS363bds2rFbrRf/w2Lp1K1u3bn3Pch0Ox6+9/St1TkgkEolEcjX7KDsEz3s/12Apw338Ge6jym8gZTiJ5FohvUoskXzKtbe3X+kqfKi2bNnyidunD+JqOwatra1XbNvSOSGRSCQSyZXxfq7B0vX6nKvpGFzJ/AbSOSGRfFykjkGJ5FNqx44dbNmyheHh4Qvj+G3ZsuWiABAKhXj44YfZtm0bmzZtYtOmTRd9v2vXLmw2GzKZ7MJFu7W1lS1btiCTyS68knF+4OGWlhZ27dpFa2srLS0tFy3zZue3e/6/yx1XZdeuXRfq8fDDD18UJnbt2nWh/jt27MBms/Hwww9fWPfRRx+9sK8tLS2XHOPkzfU6fzzeLaycX/78fp4f0PrNx/Xhhx+mpqbmssdUea82uZx2fa9y36lO77X/u3bturC989rb299Wn/Prnj93LnffA4HAhfPoUufOpbb/bufEr9MOEolEIpF80r1X5njrMm+9lr7bNfitpAz3wTLctZDfQMpwEslVT5RIJJ9q1dXVYnV19ds+b2trE61Wq9jW1nbhs8cee0wExO3bt1/4bOvWrSJw0XJDQ0MiIG7evPnCz5s3bxYBcePGjeIjjzwitrW1XVj3zeUNDQ2JVqtV3LNnz4XPtm/fLgJic3Pze+7PI488IgLi0NDQhc927twpVldXi4C4detW8ZFHHhGbm5svlHd+nfP27NkjAhfVYWhoSKyurr6oXKvVKlqt1ret99hjj134bOPGjeLOnTsvquPmzZvFRx555KLj+uZj8E4ut01E8Z3b9Z28V53ea/+HhoYutNPGjRvfVv6l2u983d98vC7lzeVu3779kufOu23/UufE5eyzRCKRSCSfBpe6Rl9u5niva+k7XYMvRcpw57yfDHc15zdRlDKcRHKtkJ4YlEgkl/TVr36VZcuWvW1MkebmZrZt28bw8DBw6TFZ3jwuCpybQOSBBx4AYNOmTWzfvp3m5uYLswzv2bPnwrLbtm1j2bJlbNy48cJnjzzyyK+1L5s3b75wZ9lqtbJ9+3ba2tpoa2sDzt0VffN+LFu27G312rJlCw8//DDV1dUXPvvOd75zYbDkS9myZQvbt2+/MKj1eW+983u547Rcbpt8EO9Vp/fa/+rq6vfdTm89T97Lli1beOSRRy557nyQ7X/QdpBIJBKJ5JPucjPHR30tlTLcu7sW8tv5ekgZTiK5ekkdgxKJ5G2Gh4dpb2+/5IDU58PZ+Yv6+3WpjsRAIHBhu+dfF/mwnd/u8uXL3/bdzp07LwRMgJMnTwJcCIvnj8ebOyvhXIelKIpv26dQKHRRB+hbVVdX8+ijj170ivR7BaKPsk3eq07vd/8/KufD/nlWq/XX6gz9IO0gkUgkEskn3fvJHB/HtVTKcO/sWshvIGU4ieRqJ3UMSiSSt3m3MVfOX9h/nYv5Ozlf5pvvan7YLhWCrFYr1dXVF8Y3eev+n//5cgPUY489Rmtr6zuOCXN+ZrVt27ZRU1NzWYMqf9Rt8m51er/7f634IO0gkUgkEskn3fvJHB/ntVTKcO+vPp/U/AZShpNIPmxSx6BEInlHl3q94ny4+CCvEbyX86Ho/BOEH5fh4WFaWloYHh5m586db7vjeL5elxvaHn74YZqbm3n44YcvuU51dTUjIyNs3Ljxwrbfz8DVb/VhtMm71en97v+14tdpB4lEIpFIPukuJ3Nc6Wvppz3DfRrzG1z5804i+aSROgYlEsnbnH/V4VJ3S8+Hmpqamg99u+efFHzzKyEfh02bNmG329/xFYTzx2Pnzp2X/P6tgctqtV5Y9lKvRQ8PD2O1WtmzZ8+F5d48u9671eGjapN3q9P73f9rxQdpB4lEIpFIPuneT+a40tfST3uG+zTmN7jy551E8kkjdQxKJJ9ygUDgbU/oVVdX09zczPDw8NtCw8mTJ7FarRcG+XU4HMDF4eL8n99pQOd3cv51ih07dlxy3fdT3uUue34f3/yaxfl1zx+XN9frraFu27Ztl7zLW11dzWOPPcbw8PDbgsr27dsv/Hnz5s0XxpV5t4D2ftrkfN3fz5OX71anD7L/Hya/3w/8+k+SvvWc+CDtIJFIJBLJJ937yRyXey2VMtxHk+Gu5vwGUoaTSK4VUsegRCK5pPNjd7w5EIVCIbZv384Pf/jDCyHs/N3Ibdu20drayo4dOy5cnFtbWy/cbb2cQGC1Wi/c8W1paaG1tZXh4WG2bdsGnLvYv3mQ4Us5f8f1fKDbtWvXhbq/+f/nnQ9Fu3btYseOHezYsePC9trb2y+sfz6AbNq0iS1btrBt2zZaWlqoqam5cCzeuo2tW7eyceNGduzYcVG9n3rqqYuCy/lZ4d5rbMXLbZMP4t3qdH4WwMvd/3dq63cK+u+3A/nN6715W++0/Xc6Jz5oO0gkEolE8kl3uZnjva6l73QNvhQpw71/12J+O7+ulOEkkquIKJFIPpXa2trErVu3ioAIiFu3bhX37Nlz0TLBYFDcvHmzuHHjRnHr1q3i1q1bxba2treVtX37dtFqtYpWq1V85JFHRFEUxerqavGRRx4R29raxLa2NrG5uVkExOrqanHPnj1iMBi8aPvbt2+/UN5jjz0mVldXi4DY3NwsDg0NXShvaGjoPfetublZtFqt4tatW0VRFMWdO3deKK+6ulp87LHHLlr+scceE61W64VtiKIobt269aL9OV/O+f1obm6+6Hi9+bvq6mpx586dF8o5v4+bN28Wh4aGxI0bN17Y1iOPPCJu3rxZDAaD77lfl9Mml9Oul3I5dXq3/W9raxM3b958UXsGg0ExGAyKjzzyyNs+37Nnz9vOiUvZuXOnaLVaL2yzra3tbefO+fPsUts/763nxOXus0QikUgkn1RvzQxvvXZeTg68nGvppa7B70TKcO8vw12t+e38dqUMJ5FcG2SiKIofTZejRCKRSCQSiUQikUg
  569. "text/plain": [
  570. "<Figure size 1500x600 with 2 Axes>"
  571. ]
  572. },
  573. "metadata": {},
  574. "output_type": "display_data"
  575. }
  576. ],
  577. "source": [
  578. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(15,6))\n",
  579. "\n",
  580. "for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
  581. " z_coord = np.linspace(scifi_z_found[i,0],12000,300)\n",
  582. " fit = scifi_track(z_coord, *scifi_fitpars_found[i])\n",
  583. " ax0.plot(z_coord, fit, \"-\", lw=0.5)\n",
  584. " ax0.errorbar(ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]),fmt=\".\",ms=2)\n",
  585. "\n",
  586. "#ax0.legend()\n",
  587. "ax0.set_xlabel(\"z [mm]\")\n",
  588. "ax0.set_ylabel(\"x [mm]\")\n",
  589. "ax0.set_title(\"found tracks of scifi hits\")\n",
  590. "ax0.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
  591. "ax0.grid()\n",
  592. "\n",
  593. "for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
  594. " z_coord = np.linspace(scifi_z_lost[i,0],12000,300)\n",
  595. " fit = scifi_track(z_coord, *scifi_fitpars_lost[i])\n",
  596. " ax1.plot(z_coord, fit, \"-\", lw=0.5)\n",
  597. " ax1.errorbar(ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]),fmt=\".\",ms=2)\n",
  598. "\n",
  599. "#ax1.legend()\n",
  600. "ax1.set_xlabel(\"z [mm]\")\n",
  601. "ax1.set_ylabel(\"x [mm]\")\n",
  602. "ax1.set_title(\"lost tracks of scifi hits\")\n",
  603. "ax1.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
  604. "ax1.grid()\n",
  605. "\n",
  606. "plt.show()"
  607. ]
  608. },
  609. {
  610. "cell_type": "code",
  611. "execution_count": 18,
  612. "metadata": {},
  613. "outputs": [
  614. {
  615. "name": "stdout",
  616. "output_type": "stream",
  617. "text": [
  618. "found \n",
  619. "zmag = 5215.5640412342\n",
  620. "lost \n",
  621. "zmag = 5450.484726770035\n"
  622. ]
  623. }
  624. ],
  625. "source": [
  626. "#vergleich der zmag werte\n",
  627. "zmag_found = z_mag(xv_found, zv_found, tx_found, scifi_fitpars_found[:,0], scifi_fitpars_found[:,1])\n",
  628. "zmag_lost = z_mag(xv_lost, zv_lost, tx_lost, scifi_fitpars_lost[:,0], scifi_fitpars_lost[:,1])\n",
  629. "zmag_lost = zmag_lost[~np.isnan(zmag_lost)]\n",
  630. "zmag_found = zmag_found[~np.isnan(zmag_found)]\n",
  631. "\n",
  632. "print(\"found \\nzmag = \", str(np.mean(zmag_found)))\n",
  633. "print(\"lost \\nzmag = \", str(np.mean(zmag_lost)))"
  634. ]
  635. },
  636. {
  637. "cell_type": "code",
  638. "execution_count": 19,
  639. "metadata": {},
  640. "outputs": [
  641. {
  642. "data": {
  643. "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAHKCAYAAAATuQ/iAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJ0UlEQVR4nO3dTYwjZ37n+V/NeFRObEGKzPJlbq2gLz6NFMy6GAtImCJHex0nWXUZYA7rIsdIYA5CN6MTGECli9mRrSvhIbO98GK1BxWZvfBxlqyG5AEag1EyLAO7OzvuZZR6gYEN2MoMFcqbTrnt2EMqQmQyGAy+ZDKY/H6AhFSM53niiSfJ4D+ft7gTBEEgAAAAxPpHq64AAABAlhEsAQAAJCBYAgAASECwBAAAkIBgCQAAIAHBEgAAQAKCJQAAgAQES7gVer2eqtWqOp3OWp8DuI1c11Wr1Vp1NYC5ESxh7XU6Hdm2rVarpdPT08S0nudpe3t75hv3LOfIilmvdd62QXqu68q2beXzeeXz+bU9R1qe56lcLiufz6vZbK60LpJk27YODw9VLpe1vb0t27ZXXSWsCYIlZJbnefJ9f2q6Uqmkg4ODVGX6vi/f9zUYDGaqyyznyIqka41r23nbBulZlqXHjx/Ldd1MnyPtZ28a0zTVbrcXLmdRvu8rl8vp8ePHqtVqarfbchxHnU4n9jqXdf3rjnb4DsESMqtcLqfuxTEMI1U6y7IUBIEcx5m5PmnPkRVJ1xrXtou0DdKzLCvz55jls7cObNuWYRgj7VKpVDQYDGI/17ft+udFO3yHYAmZVC6Xr/Wv701G2yLJbXx/9Hq91Glv4/XPg3YYRbC0xnzfV6vVUrFYVKvVkud5KhaL2t7eVrFYjLpPDw8PlcvlJo7RHx4eqlqtRvMcJs1ZqVarqlaryuVyYx+ksC75fF6dTke9Xk/5fF537txRuVweK8t1XZXLZRWLReVyuZF6dTqdqOxqtTrXh9b3fVWr1ej84aTsTqejcrkcW6cwT7VaVbFYVLFYTDzvpHNMShu2T6/Xi/5/e3tb5XI5dkgs/J2Edbl6wx9OE/5ehn93cdea1LZp2mZSfWb9/U/ium7sT1akeY+k/TzNWn6n09H29rbu3LkTvdbr9VQul2dq56T6TfvsJX1u466hWq3q8PAw9fX3ej3duXNn5Hpc11Uul9OdO3eUz+dH6hNeu+d5seWF9fU8L5o/FbZr2JZprn/Wz++090Cn04k+Q61WS9vb26pWq6nyL+u+v8g9eNr7IOn6pt23MivA2hoMBkGlUgkkBYVCIajVakG/3w+63W70WqVSCbrdbjAYDIJSqRRICvr9flRGrVYLht8GYd5utztyLsuyglqtFp1XUiApME0zqFQqI+UP1yWsn+M4UVn9fj8oFArRv9vtdiApqFQqY/UaDAap2iKsd7PZjF4rFApBu90eaS/HcaI6Xm1L0zRHzmcYRmAYxkznmCRsh+H2abfbUZuZphml7ff7gWEYI7+nZrM51o6lUin6nYRpwuNJ1xrXtknp09Rnlt//JP1+P2qj4Z80eW9CmvdI2s+TpMCyrJnLD9tz+HcRfh5LpdLUc6Sp36TPXprP7WAwCAzDGCkvfF9drcsk4fmHywjrOfx+D4LL9+HV1+LEnT98Dw9/nofPP3z9s3x+p7Vxu90OTNOM2q5WqwWWZUX1m5Z/Gff9Re7B0/JOu76k+1aWESytufALZvhNHgSXwc3VN3qYdviNWSgURm7GZ2dnYzel8MMw/GELP6zD5Yfprr7xr34BW5Y1UlYQXH4pSArOzs6CIFg8WCqVSmPnmFSfsE5X6x3e5MM6zXKOOGF5V2/OhUIhkBQFXZZljdUvfH24TQzDGPuimNb2QZDctpPaJk190v7+J2k2myPteXZ2Nva+XqU075E0n6cgiP/yTlN++Lu72k5pg6U09Zv0/kjzuS2VSrG/61mCpUnBn2EYI0FJeL4094i484fv1zTBUhCk//ymaeOwrLhAL03+Re/7i9yD0+RNur40960sYhjulrg6SdE0TUnSzs7O2GvDq53a7bb6/X7075OTE0ka6Vb+/PPPx8433EU+rS6SokmCnufJdV3V6/VoyGd4+CA8/7x831exWJTjOKknuYZ1KhQKI6/XajUFQTB2PfOcY1j4ewiF3dPdbjeqS1y5YbpwCbZpmjo8PBwZ5qjVajPXJ8ks9Qkl/f6TPHr0aOQ8tm1nYrm5lP49kubztEj5i1q0fkmfW8/zouGXRZimKcuyYoe1w+G0sM6+7499nq5b0udXStfG4e/zwYMHY+XP8jua576/yD04bd6k67uJ+9Z1+LVVVwCrZRiGDMNQp9PRJ598EvvmDl/r9XpjX5qzBgthcHVdy4mbzaY8z1Ov11OlUpmpTmm/kOY5R5KwDcMb0SS7u7tROumyDfP5fBRUtNvtpa+0mqU+ixpu/2q1mqlVeWnfI2k+T4uUv6hF65f0uQ3nsC0jeAnnO3U6HZVKJbVaLT169EitVkvNZlOO4+jZs2czzYe7LsOfX2m2No77fc/7O0prkXvwrHnjru8m7lvXgZ6lDed5nvL5vDzPU7vdjo3wS6WSCoWC6vW6er2efN+X4ziq1Woz3xjDG8qyvmCvqlarsixL1Wo19TlmrdM850gS/hU43JZJf0UOp3/x4oUKhUL0e7yuiZJp6rMs4eTWLG3VkPY9kubztEj5i7rO+oXHlrHUPPwjJOxZDAOkQqEQvcebzeZS/lhZ1NXP77xtHFo0f5ryh/97U3lDN3nfWiaCpQ1XLBa1s7Mz9QPZbrdVKBSixxY4jjPXX/7hDWXSyrFZlvjGMQwj+qsn7XBA+FfNpL+Wrt4Y5jlHkvDLJZ/PR3WJa4cwYMnlclG9DMNQt9uN6jO8omYZZqnPMnQ6nbH9cLIg7Xsk7edp3vIXNW/90nxuwzTDQ0iLKJVK6vV66nQ62t3dlWEYqlar8n1fh4eHUc/mqg1/fqX52zi0aP5pFrkHL+P+fRP3retAsLTBwvH/4b/gwy/Aq38dlsvl6K+cWq0295dZOCfDtu2xIZ64vy7m2T3WNM1oqCzNhzC86bZarbEPu23bsT0ns54jSRggVCqVaL7G8NyM0MnJSZRO0kiwWiqVor/C036xpmnbWeqzKM/z1O12Y8vzfV+2bSuXy0VLp4cfzVKtVrW9vT1SxzBPq9XS4eGhtre3x85XrVajZdhJ7ZHmPTLL52me8iXp/v37Ud2Hr2P4XJPMWr/h8tJ8boevIa4us36Ww89VuVyOlqaXSqWoHtf9BZu2vsOf30XeA9Lsv6N5LHIPnjVvnEXvW6tCsLTmJn2A4j5g4f+H/w1vwJ1OR61WS61WK7opua4bPQogvIGXy2UdHh6q1WqN7MMxrS7DDMOI/mLK5/NRmcViUYPBIPowhr0VYUAy7eG14fWG/61UKlGX/fBEwrh2MQwj+gAXi8Xo5pzP55XL5aIbV9pzTDM8adn3fTWbTR0dHUWvtdvt6K/o4XSO4+jo6Ciqz7Nnz8aCA9M0o7/+Jt1kJ7XtpPRp67PIzTycND+pp84wjOh30ev1VCgUVCgUNBgMdHh4KMdxogm/0uWN9+HDh6pWq1HwNdwT4ft+NG8iPP7w4cOJN+w075G0n6d5y5e+64GybTvawyZ8P/V6vcSezrT1i3t/pPncXk3T6/XkeV50Ds/zZvqchGUWCoWRIepSqRQF8bOY1PZXX59270n6/KZt46v3klDa/Ivc9xe5B6fNO+n6pOn3rcxa8Wo8LCDc70JSYBhGtHQ13D9E3y6/7ff7I/ttGIYRLX9tNpvRktxwOWelUhlZ3hnu/xKWOfxjmmZwdnYW9Pv9aNmqaZpBt9uNln6HaYeXhzqOE5VpmubYctwguFyiahjG1OXj7XZ75NxhOwyfu1QqBX/8x38ctUFYn3Cp69VyLMsa2ecl7TmSljEPL6e1LCsolUpBqVQa24MnCC6XC4fLsCuVSlCpVMaW6xYKhej3VqvVglKpFF1Pv99PvNarbTst/bT6zPr7v6pSqUTv5eGl55VKJWrrsL3C816t//DtzLKskfdUpVIZe/8
  644. "text/plain": [
  645. "<Figure size 640x480 with 1 Axes>"
  646. ]
  647. },
  648. "metadata": {},
  649. "output_type": "display_data"
  650. }
  651. ],
  652. "source": [
  653. "plt.hist(zmag_found, bins=5000, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  654. "plt.hist(zmag_lost, bins=400, density=True, alpha=0.5, histtype=\"bar\",color=\"darkorange\", label=\"lost\")\n",
  655. "plt.xlabel(r\"$\\bf{z_{mag}}$ [mm]\")\n",
  656. "plt.ylabel(\"normed\")\n",
  657. "plt.title(\"magnet kick position $z_{mag}$ calculated w fitparameters\")\n",
  658. "plt.legend(title=\"LHCb Simulation\", title_fontsize=15)\n",
  659. "#plt.xticks(np.arange(5100,5800,5), minor=True)\n",
  660. "#plt.yticks(np.arange(0,0.015,0.001), minor=True)\n",
  661. "plt.xlim(5050,5750)\n",
  662. "\n",
  663. "\"\"\"\n",
  664. "wir können einen radikalen unterschied für den z_mag wert erkennen, zwischen den found and lost elektronen.\n",
  665. "\"\"\"\n",
  666. "\n",
  667. "plt.show()"
  668. ]
  669. },
  670. {
  671. "cell_type": "code",
  672. "execution_count": 20,
  673. "metadata": {},
  674. "outputs": [],
  675. "source": [
  676. "\n",
  677. "#file.show()"
  678. ]
  679. },
  680. {
  681. "cell_type": "code",
  682. "execution_count": null,
  683. "metadata": {},
  684. "outputs": [],
  685. "source": []
  686. },
  687. {
  688. "cell_type": "code",
  689. "execution_count": null,
  690. "metadata": {},
  691. "outputs": [],
  692. "source": []
  693. }
  694. ],
  695. "metadata": {
  696. "kernelspec": {
  697. "display_name": "env1",
  698. "language": "python",
  699. "name": "python3"
  700. },
  701. "language_info": {
  702. "codemirror_mode": {
  703. "name": "ipython",
  704. "version": 3
  705. },
  706. "file_extension": ".py",
  707. "mimetype": "text/x-python",
  708. "name": "python",
  709. "nbconvert_exporter": "python",
  710. "pygments_lexer": "ipython3",
  711. "version": "3.10.12"
  712. },
  713. "orig_nbformat": 4
  714. },
  715. "nbformat": 4,
  716. "nbformat_minor": 2
  717. }