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.

1261 lines
1.4 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
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
1 year ago
1 year ago
1 year ago
1 year 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. "51"
  31. ]
  32. },
  33. "execution_count": 2,
  34. "metadata": {},
  35. "output_type": "execute_result"
  36. }
  37. ],
  38. "source": [
  39. "file = uproot.open(\"tracking_losses_ntuple_Dst0ToD0EE.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\")\n",
  40. "\n",
  41. "\n",
  42. "#selektiere nur elektronen von D*0->D0ee und nur solche mit einem momentum von unter 5 GeV \n",
  43. "allcolumns = file.arrays()\n",
  44. "found = allcolumns[(allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p < 5e3)] #D: 2591\n",
  45. "lost = allcolumns[(allcolumns.isElectron) & (allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p < 5e3)] #D: 1908\n",
  46. "\n",
  47. "#ak.num(lost, axis=0)\n",
  48. "ak.num(allcolumns[(allcolumns.fromPairProd) & (allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal)],axis=0)\n",
  49. "#found[0]"
  50. ]
  51. },
  52. {
  53. "cell_type": "code",
  54. "execution_count": 3,
  55. "metadata": {},
  56. "outputs": [
  57. {
  58. "name": "stdout",
  59. "output_type": "stream",
  60. "text": [
  61. "eff all = 0.5759057568348522 +/- 0.007367987865301135\n"
  62. ]
  63. }
  64. ],
  65. "source": [
  66. "def t_eff(found, lost, axis = 0):\n",
  67. " sel = ak.num(found, axis=axis)\n",
  68. " des = ak.num(lost, axis=axis)\n",
  69. " return sel/(sel + des)\n",
  70. "\n",
  71. "def eff_err(found, lost):\n",
  72. " n_f = ak.num(found, axis=0)\n",
  73. " n_all = ak.num(found, axis=0) + ak.num(lost,axis=0)\n",
  74. " return 1/n_all * np.sqrt(np.abs(n_f*(1-n_f/n_all)))\n",
  75. "\n",
  76. "\n",
  77. "print(\"eff all = \", t_eff(found, lost), \"+/-\", eff_err(found, lost))"
  78. ]
  79. },
  80. {
  81. "cell_type": "code",
  82. "execution_count": 4,
  83. "metadata": {},
  84. "outputs": [
  85. {
  86. "data": {
  87. "text/plain": [
  88. "0.7960893854748603"
  89. ]
  90. },
  91. "execution_count": 4,
  92. "metadata": {},
  93. "output_type": "execute_result"
  94. }
  95. ],
  96. "source": [
  97. "\n",
  98. "#statistics\n",
  99. "\n",
  100. "nobrem_found = found[found[\"brem_photons_pe_length\"]==0]\n",
  101. "nobrem_lost = lost[lost[\"brem_photons_pe_length\"]==0]\n",
  102. "\n",
  103. "\"\"\"\n",
  104. "die effizienz mit der wir elektronen finden, die keine bremsstrahlung gemacht haben, ist nicht besonders gut, aber trotzdem besser als\n",
  105. "für alle elektronen.\n",
  106. "Auch hier handelt es sich um eine recht geringe sample size (~350)\n",
  107. "\"\"\"\n",
  108. "\n",
  109. "t_eff(nobrem_found, nobrem_lost)\n",
  110. "#ak.num(nobrem_lost, axis=0)"
  111. ]
  112. },
  113. {
  114. "cell_type": "code",
  115. "execution_count": 5,
  116. "metadata": {},
  117. "outputs": [
  118. {
  119. "data": {
  120. "text/html": [
  121. "<pre>{energy: 2.04e+03,\n",
  122. " photon_length: 4,\n",
  123. " brem_photons_pe: [519, 484, 173, 466],\n",
  124. " brem_vtx_z: [9.84e+03, 9.84e+03, 9.85e+03, 9.85e+03]}\n",
  125. "------------------------------------------------------\n",
  126. "type: {\n",
  127. " energy: float64,\n",
  128. " photon_length: int64,\n",
  129. " brem_photons_pe: var * float64,\n",
  130. " brem_vtx_z: var * float64\n",
  131. "}</pre>"
  132. ],
  133. "text/plain": [
  134. "<Record {energy: 2.04e+03, ...} type='{energy: float64, photon_length: int6...'>"
  135. ]
  136. },
  137. "execution_count": 5,
  138. "metadata": {},
  139. "output_type": "execute_result"
  140. }
  141. ],
  142. "source": [
  143. "#try excluding all photons that originate from a vtx @ z>9500mm\n",
  144. "#ignore all brem vertices @ z>9500mm \n",
  145. "\n",
  146. "#found\n",
  147. "\n",
  148. "brem_e_f = found[\"brem_photons_pe\"]\n",
  149. "brem_z_f = found[\"brem_vtx_z\"]\n",
  150. "e_f = found[\"energy\"]\n",
  151. "length_f = found[\"brem_vtx_z_length\"]\n",
  152. "\n",
  153. "brem_f = ak.ArrayBuilder()\n",
  154. "\n",
  155. "for itr in range(ak.num(found,axis=0)):\n",
  156. " brem_f.begin_record()\n",
  157. " #[:,\"energy\"] energy\n",
  158. " brem_f.field(\"energy\").append(e_f[itr])\n",
  159. " #[:,\"photon_length\"] number of vertices\n",
  160. " brem_f.field(\"photon_length\").integer(length_f[itr])\n",
  161. " #[:,\"brem_photons_pe\",:] photon energy \n",
  162. " brem_f.field(\"brem_photons_pe\").append(brem_e_f[itr])\n",
  163. " #[:,\"brem_vtx_z\",:] brem vtx z\n",
  164. " brem_f.field(\"brem_vtx_z\").append(brem_z_f[itr])\n",
  165. " brem_f.end_record()\n",
  166. "\n",
  167. "brem_f = ak.Array(brem_f)\n",
  168. "\n",
  169. "#lost\n",
  170. "\n",
  171. "brem_e_l = lost[\"brem_photons_pe\"]\n",
  172. "brem_z_l = lost[\"brem_vtx_z\"]\n",
  173. "e_l = lost[\"energy\"]\n",
  174. "length_l = lost[\"brem_vtx_z_length\"]\n",
  175. "\n",
  176. "brem_l = ak.ArrayBuilder()\n",
  177. "\n",
  178. "for itr in range(ak.num(lost,axis=0)):\n",
  179. " brem_l.begin_record()\n",
  180. " #[:,\"energy\"] energy\n",
  181. " brem_l.field(\"energy\").append(e_l[itr])\n",
  182. " #[:,\"photon_length\"] number of vertices\n",
  183. " brem_l.field(\"photon_length\").integer(length_l[itr])\n",
  184. " #[:,\"brem_photons_pe\",:] photon energy \n",
  185. " brem_l.field(\"brem_photons_pe\").append(brem_e_l[itr])\n",
  186. " #[:,\"brem_vtx_z\",:] brem vtx z\n",
  187. " brem_l.field(\"brem_vtx_z\").append(brem_z_l[itr])\n",
  188. " brem_l.end_record()\n",
  189. "\n",
  190. "brem_l = ak.Array(brem_l)\n",
  191. "\n",
  192. "\n",
  193. "\n",
  194. "\n",
  195. "brem_f[0]"
  196. ]
  197. },
  198. {
  199. "cell_type": "code",
  200. "execution_count": 6,
  201. "metadata": {},
  202. "outputs": [],
  203. "source": [
  204. "acc_brem_found = ak.ArrayBuilder()\n",
  205. "\n",
  206. "for itr in range(ak.num(brem_f, axis=0)):\n",
  207. " acc_brem_found.begin_record()\n",
  208. " acc_brem_found.field(\"energy\").real(brem_f[itr,\"energy\"])\n",
  209. " \n",
  210. " acc_brem_found.field(\"brem_photons_pe\")\n",
  211. " acc_brem_found.begin_list()\n",
  212. " for jentry in range(brem_f[itr, \"photon_length\"]):\n",
  213. " if brem_f[itr, \"brem_vtx_z\", jentry]>9500:\n",
  214. " continue\n",
  215. " else:\n",
  216. " acc_brem_found.real(brem_f[itr,\"brem_photons_pe\", jentry])\n",
  217. " \n",
  218. " #acc_brem_found.field(\"brem_vtx_z\").real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  219. " acc_brem_found.end_list()\n",
  220. " \n",
  221. " acc_brem_found.field(\"brem_vtx_z\")\n",
  222. " acc_brem_found.begin_list()\n",
  223. " for jentry in range(brem_f[itr, \"photon_length\"]):\n",
  224. " if brem_f[itr, \"brem_vtx_z\", jentry]>9500:\n",
  225. " continue\n",
  226. " else:\n",
  227. " acc_brem_found.real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  228. " acc_brem_found.end_list()\n",
  229. " \n",
  230. "\n",
  231. " \n",
  232. " acc_brem_found.end_record()\n",
  233. "\n",
  234. "acc_brem_found = ak.Array(acc_brem_found)\n",
  235. "\n",
  236. "\n",
  237. "\n",
  238. "acc_brem_lost = ak.ArrayBuilder()\n",
  239. "\n",
  240. "for itr in range(ak.num(brem_l, axis=0)):\n",
  241. " acc_brem_lost.begin_record()\n",
  242. " acc_brem_lost.field(\"energy\").real(brem_l[itr,\"energy\"])\n",
  243. " \n",
  244. " acc_brem_lost.field(\"brem_photons_pe\")\n",
  245. " acc_brem_lost.begin_list()\n",
  246. " for jentry in range(brem_l[itr, \"photon_length\"]):\n",
  247. " if brem_l[itr, \"brem_vtx_z\", jentry]>9500:\n",
  248. " continue\n",
  249. " else:\n",
  250. " acc_brem_lost.real(brem_l[itr,\"brem_photons_pe\", jentry])\n",
  251. " \n",
  252. " #acc_brem_found.field(\"brem_vtx_z\").real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  253. " acc_brem_lost.end_list()\n",
  254. " \n",
  255. " acc_brem_lost.field(\"brem_vtx_z\")\n",
  256. " acc_brem_lost.begin_list()\n",
  257. " for jentry in range(brem_l[itr, \"photon_length\"]):\n",
  258. " if brem_l[itr, \"brem_vtx_z\", jentry]>9500:\n",
  259. " continue\n",
  260. " else:\n",
  261. " acc_brem_lost.real(brem_l[itr, \"brem_vtx_z\",jentry])\n",
  262. " acc_brem_lost.end_list()\n",
  263. " \n",
  264. " acc_brem_lost.end_record()\n",
  265. "\n",
  266. "acc_brem_lost = ak.Array(acc_brem_lost)\n"
  267. ]
  268. },
  269. {
  270. "cell_type": "code",
  271. "execution_count": 7,
  272. "metadata": {},
  273. "outputs": [
  274. {
  275. "name": "stdout",
  276. "output_type": "stream",
  277. "text": [
  278. "\n",
  279. "cutoff = 50 MeV, sample size: 789\n",
  280. "eff = 0.7427122940430925 +/- 0.0156\n",
  281. "\n",
  282. "cutoff = 100 MeV, sample size: 789\n",
  283. "eff = 0.7427122940430925 +/- 0.0156\n",
  284. "\n",
  285. "cutoff = 150 MeV, sample size: 1025\n",
  286. "eff = 0.7346341463414634 +/- 0.0138\n",
  287. "\n",
  288. "cutoff = 200 MeV, sample size: 1202\n",
  289. "eff = 0.7346089850249584 +/- 0.0127\n",
  290. "\n",
  291. "cutoff = 250 MeV, sample size: 1337\n",
  292. "eff = 0.7270007479431563 +/- 0.0122\n",
  293. "\n",
  294. "cutoff = 300 MeV, sample size: 1474\n",
  295. "eff = 0.7164179104477612 +/- 0.0117\n",
  296. "\n",
  297. "cutoff = 350 MeV, sample size: 1599\n",
  298. "eff = 0.7148217636022514 +/- 0.0113\n",
  299. "\n",
  300. "cutoff = 400 MeV, sample size: 1721\n",
  301. "eff = 0.7094712376525276 +/- 0.0109\n",
  302. "\n",
  303. "cutoff = 450 MeV, sample size: 1854\n",
  304. "eff = 0.7076591154261057 +/- 0.0106\n",
  305. "\n",
  306. "cutoff = 500 MeV, sample size: 1966\n",
  307. "eff = 0.7004069175991862 +/- 0.0103\n",
  308. "\n",
  309. "cutoff = 550 MeV, sample size: 2065\n",
  310. "eff = 0.6968523002421307 +/- 0.0101\n",
  311. "\n",
  312. "cutoff = 600 MeV, sample size: 2152\n",
  313. "eff = 0.6909851301115242 +/- 0.01\n",
  314. "\n",
  315. "cutoff = 650 MeV, sample size: 2252\n",
  316. "eff = 0.6873889875666075 +/- 0.0098\n",
  317. "\n",
  318. "cutoff = 700 MeV, sample size: 2346\n",
  319. "eff = 0.6837169650468883 +/- 0.0096\n",
  320. "\n",
  321. "cutoff = 750 MeV, sample size: 2441\n",
  322. "eff = 0.6780008193363376 +/- 0.0095\n",
  323. "\n",
  324. "cutoff = 800 MeV, sample size: 2534\n",
  325. "eff = 0.6764009471191792 +/- 0.0093\n",
  326. "\n",
  327. "cutoff = 850 MeV, sample size: 2624\n",
  328. "eff = 0.6714939024390244 +/- 0.0092\n",
  329. "\n",
  330. "cutoff = 900 MeV, sample size: 2707\n",
  331. "eff = 0.6693756926486886 +/- 0.009\n",
  332. "\n",
  333. "cutoff = 950 MeV, sample size: 2785\n",
  334. "eff = 0.6660682226211849 +/- 0.0089\n",
  335. "\n",
  336. "cutoff = 1000 MeV, sample size: 2839\n",
  337. "eff = 0.6657273687918281 +/- 0.0089\n",
  338. "\n",
  339. "cutoff energy = 350MeV, sample size: 1599\n",
  340. "eff = 0.7148 +/- 0.0113\n"
  341. ]
  342. }
  343. ],
  344. "source": [
  345. "#finden wir die elektronen die keine bremsstrahlung gemacht haben mit hoher effizienz?\n",
  346. "efficiencies_found = []\n",
  347. "\n",
  348. "\n",
  349. "\n",
  350. "for cutoff_energy in range(50,1050,50):\n",
  351. "\tnobrem_f = acc_brem_found[ak.sum(acc_brem_found[\"brem_photons_pe\"],axis=-1,keepdims=False)<cutoff_energy]\n",
  352. "\tnobrem_l = acc_brem_lost[ak.sum(acc_brem_lost[\"brem_photons_pe\"],axis=-1,keepdims=False)<cutoff_energy]\n",
  353. "\n",
  354. "\tif ak.num(nobrem_f,axis=0)+ak.num(nobrem_l,axis=0)==0:\n",
  355. "\t\tcontinue\n",
  356. "\teff = t_eff(nobrem_f,nobrem_l)\n",
  357. "\tefficiencies_found.append(eff)\n",
  358. "\n",
  359. "\tprint(\"\\ncutoff = \",str(cutoff_energy) ,\"MeV, sample size: \",ak.num(nobrem_f,axis=0)+ak.num(nobrem_l,axis=0))\n",
  360. "\tprint(\"eff = \",eff, \"+/-\", np.round(eff_err(nobrem_f, nobrem_l),4))\n",
  361. "\n",
  362. "\"\"\"\n",
  363. "we see that a cutoff energy of xxxMeV is ideal because the efficiency drops significantly for higher values\n",
  364. "\"\"\"\n",
  365. "cutoff_energy = 350.0 #MeV\n",
  366. "\n",
  367. "\"\"\"\n",
  368. "better statistics: cutoff=xxxMeV - sample size: xxx events and efficiency=xxxx\n",
  369. "\"\"\"\n",
  370. "nobrem_found = acc_brem_found[ak.sum(acc_brem_found[\"brem_photons_pe\"],axis=-1,keepdims=False)<cutoff_energy]\n",
  371. "nobrem_lost = acc_brem_lost[ak.sum(acc_brem_lost[\"brem_photons_pe\"],axis=-1,keepdims=False)<cutoff_energy]\n",
  372. "\n",
  373. "print(\"\\ncutoff energy = 350MeV, sample size:\",ak.num(nobrem_found,axis=0)+ak.num(nobrem_lost,axis=0))\n",
  374. "print(\"eff = \",np.round(t_eff(nobrem_found, nobrem_lost),4), \"+/-\", np.round(eff_err(nobrem_found, nobrem_lost),4))"
  375. ]
  376. },
  377. {
  378. "cell_type": "code",
  379. "execution_count": 8,
  380. "metadata": {},
  381. "outputs": [
  382. {
  383. "data": {
  384. "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHICAYAAACyBMv/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABJbElEQVR4nO3deVxU9f4/8NcAw4CyuCCLirgvuIuaW5m7pmaZ5ZaSolcFLbU0/WZudUO7N691L5qaSmkulcttsZJMy61E1Ewx98J0EEUETYUB3r8//M1cx2FnZg6c83o+HjxqznzOmc/5zDi8+JzP53N0IiIgIiIiUgkXpStAREREZE8MN0RERKQqDDdERESkKgw3REREpCoMN0RERKQqDDdERESkKgw3REREpCoMN0RERKQqDDdERESkKgw3Cjt+/DiSkpKUrgYREZFqMNwobOXKlfj++++VrgYREZFqMNwo4Mcff8T48eORkpJi2ZaYmIjBgwfj6tWrCtaMiIio/GO4KaEdO3ZAp9NZftzc3FC1alU89thjePfdd5GZmZnvvmFhYQgICECLFi2wa9cu/Oc//0GvXr3Qt29fVKtWrVycg5bs2bPHqp0e/Pnpp59syh8/fhwRERGoV68ePD094enpiQYNGmDChAk4fPhwsV//6aefhqenJ27evJlvmZEjR0Kv1zs0HM+fPx86nQ7Xr1+32zEPHDiA+fPnF3hu5V1OTg78/f3xr3/9K98yjmjb8kDN7/+hQ4fQp08feHt7w8vLC926dcP+/futyhT3u+X27duYOnUqqlevDg8PD7Rq1QqbNm0qcTk1c1O6AuXVkSNHAABbt25FUFAQsrOzcfXqVezcuROvvPIKYmNjsXv3blSqVMlm34oVK2LhwoXIzc1FdHQ0XFxc8PHHH2PYsGHl5hy06K233kK3bt2stjVr1szq8YoVKzB58mQ0atQIL730Epo2bQqdTodTp05h48aNaNeuHc6dO4d69eoV+XUjIiKwfft2bNiwAZGRkTbPp6enY9u2bRgwYAACAgJKdnIKOXDgABYsWIAXXnhBtZ+zH3/8EdeuXcPgwYOVrkqZo9b3Pz4+Ho899hjat2+PdevWQUTw9ttvo0ePHti9ezc6duxoVb4o3y0AMHjwYMTHx2PRokVo2LAhNmzYgOHDhyM3NxcjRowodjlVEyqRp59+Wry8vCQ3N9fmudjYWAEgU6ZMyXPf48ePS8uWLaVPnz7yxBNPyLhx46RZs2byxBNPyLVr1xxddYvSnIOaXLlyRdLT0/N9fvfu3QJAPv300wKPs2/fPnFxcZGBAwdKZmZmnmU++eQTuXz5crHql52dLdWrV5ewsLA8n1++fLkAkC+++KJYxy2uefPmCQC7fkb/8Y9/CAC5ePGi3Y5Z1kRGRkrbtm0LLFOatv3rr79KWjXFldf3v7DvjD59+khAQIDVe5ORkSF+fn7SqVMny7aifreIiHz11VcCQDZs2GC1vVevXlK9enXJzs4uVjm1Y7gpoVq1akmXLl3yfb5mzZpSvXr1PJ+7dOmS5RdRVFSUrF27VrKzsyU2NlaysrJKXbclS5bItm3bCi1X0nPYu3ev9OrVS3x8fKRSpUryxBNPyJkzZ0pcrqgGDBggYWFhsnLlSmnRooV4eHhIzZo1Ze7cuZKTk1OsY924cUM++OAD6dGjh7i4uMjRo0fzLVvUL6AnnnhC9Hq9XLlypVh1OXPmjAwfPlyqVasm7u7u0rhxY/nPf/5jVWb27NkCQI4fP26zf/v27SUoKKhEX1rmX6pHjhyRp59+Wry9vcXHx0dGjhwpKSkpeZY9ceKEDBs2THx8fMTf31/GjBkjN2/etDn23r17pXv37uLl5SWenp7SsWNH+fLLL22O9/DP7t27i7R/ceuVkpIi48ePl5o1a4q7u7vlF01cXFy+7XPixAkBIJ988oll2+HDhwWAhIaGWpUdOHCgtGnTxmpbbm6uBAUFSXR0dL6v8eA5FPY+mMslJCTIM888I5UqVZLAwEDL80X5LJmP8csvv8iQIUPEx8dHKleuLNOmTROTySS//fab9OnTR7y8vCQkJEQWL15cYN0ftm3bNgEg3333nc1zy5Yts7x2fu//119/La1atZJ69epZvX9Go1ECAgKka9euRf6s53V8809xA1VxvjO8vLxk6NChNtsHDx4sACzfEcUJN+PGjRMvLy8xmUxW2zds2CAAZP/+/cUqJ1L6z4KjP0ulwXBTAtevXy+0V6N79+6i0+kK/aVrDjf2NGLECNHr9QUGnJKew7x588TFxUXGjh0rX331lXz22WfSvHlzCQ4Ollu3bhW7XHEEBQVJxYoVpUmTJrJu3TrZuXOnDBs2TADIqlWrCt3/r7/+kk2bNsmTTz4p7u7u4unpKc8884x8+umn+fa0iPzvC8jf319cXV3F29tbevfuLXv37rWUyc7OtvwCLo6TJ0+Kr6+vNG/eXD766CPZuXOnvPzyy+Li4iLz58+3lDt79qzodDqZOnWqzf4AZNasWcV6XTPzl1NISIjMmDFDvv32W1myZIlUrFhRWrdubRW2zWUbNWokc+fOlbi4OFmyZIkYDAYZM2aM1XH37Nkjer1ewsLCZPPmzbJ9+3bp3bu36HQ62bRpk4jcD/lTpkwRALJ161Y5ePCgHDx4UNLT04u0f3Hr1adPH6lWrZqsXLlS9uzZI9u3b5e5c+faHO9hQUFB8re//c3yeNGiReLp6SkALL1wJpNJfHx8ZObMmVb77tu3TwAUGuqL+j48WO7VV1+VuLg42b59u4gU/bP0YHu98cYbEhcXJzNnzhQAMnnyZGncuLG89957EhcXJ2PGjBEAsmXLlgLr/yCTyST+/v4ycuRIm+fat29vCYAFvf9nzpwRb29vGTx4sIiI5OTkSPfu3cXf379YfzyYj2n++f7776VGjRoSGBhYYM+LWUm/M9zd3WX06NE224cPHy4A5NtvvxWRon23mHXo0EHatWtns90cwFesWFGsciKl/yw4+rNUGgw3JbBz504BIGvWrMm3TOfOnaVixYpOrNX/ZGdnFxpwSnIOX3zxhQCQt99+26rcmTNnBICsX7++WOWK488//xQAUrduXau/5rKysiQwMFAGDBiQ535ZWVny5ZdfyogRI6RixYri7u4uAwYMkPXr1xc5ZB05ckReeukl2bZtm/z444+yZs0aadKkibi6uso333wjIiLJyckCQIYNG2azf3Z2tphMJsvPg5cB+/TpIzVr1rT5op08ebJ4eHjIjRs3LNu6du0qfn5+VoHj5ZdfLtIvz/yYv5ymTZtmtf3jjz+2ea/MZR9+XyMjI8XDw8PqvDp06CD+/v5WbZydnS3NmjWTmjVrWsrmd1miqPsXp15eXl424bAonn/+ealbt67lcc+ePWX8+PFSuXJl+fDDD0VEZP/+/QJAdu7cabXv1KlTpXnz5oW+RlHfB3O5uXPn2hyjqJ8l8zHeeecdq3KtWrWyBA0zk8kk1apVs4SMopo+fbp4enpa/VtNTEwUAPLvf//bsq2gy1KbN28WALJ06VKZO3euuLi42LRvcWRnZ8ugQYPEy8tLEhIS8i1nj++MVq1aScOGDa3+MDSZTFK3bl2rS0ZF+W4xa9CggfTp08fmta5cuSIA5K233ipWOZHSfxac8VkqKYabEoiOjrZ0IecnICDApovaHq5du1ZgV+vDP3q9XpKTk+1yDq1bt5Z69epJZmam1S9rk8kknp6esnDhwmKVKw5zV/dHH31k81yXLl2kQ4cONtt/+eUXqVKliri6ukrv3r1l9erVkpaWVuzXzktaWprUrFlTWrRoISIFh5uWLVtavSf/+Mc/RETk7t274ubmJlOmTLFppx07dggA2bFjh+U4H330kQCQzz77TETuf1kEBATIo48+WuLzMH85HT582Gq7yWQSNzc3iYiIsCn722+/WZV9//33BYDlc3b79m3R6XQSGRlp83qLFy8WAHLq1CkRyfuXW3H2L069unfvLpUqVZI33nhDDh48WORLwGvXrhUAcuHCBbl79654eHjIli1bZPDgwfL888+LiMiCBQvEYDD
  385. "text/plain": [
  386. "<Figure size 640x480 with 1 Axes>"
  387. ]
  388. },
  389. "metadata": {},
  390. "output_type": "display_data"
  391. }
  392. ],
  393. "source": [
  394. "x_ = np.arange(50,1050,step=50)\n",
  395. "plt.scatter(x_,efficiencies_found)\n",
  396. "plt.xlabel(\"cutoff energy [MeV]\")\n",
  397. "plt.ylabel(r\"$\\epsilon$\")\n",
  398. "plt.title(r'$D^\\ast\\rightarrow Dee$, $p<5$GeV, photons w/ brem_vtx_z$<9500$mm')\n",
  399. "plt.ylim([0.5,0.95])\n",
  400. "plt.yticks(np.arange(0.5,0.96,step=0.01),minor=True)\n",
  401. "plt.xticks(np.arange(0,1100,step=50),minor=True)\n",
  402. "plt.grid()\n",
  403. "plt.show()"
  404. ]
  405. },
  406. {
  407. "cell_type": "code",
  408. "execution_count": 9,
  409. "metadata": {},
  410. "outputs": [
  411. {
  412. "name": "stdout",
  413. "output_type": "stream",
  414. "text": [
  415. "eff = 0.4993 +/- 0.0093\n"
  416. ]
  417. },
  418. {
  419. "data": {
  420. "text/html": [
  421. "<pre>[{energy: 4.19e+03, brem_photons_pe: [1e+03, ..., 158], brem_vtx_z: [...]},\n",
  422. " {energy: 4.12e+03, brem_photons_pe: [278, 468], brem_vtx_z: [...]},\n",
  423. " {energy: 4.53e+03, brem_photons_pe: [445], brem_vtx_z: [9.33e+03]},\n",
  424. " {energy: 4.22e+03, brem_photons_pe: [426, ..., 207], brem_vtx_z: [...]},\n",
  425. " {energy: 2.54e+03, brem_photons_pe: [362, ...], brem_vtx_z: [...]},\n",
  426. " {energy: 3.59e+03, brem_photons_pe: [715], brem_vtx_z: [2.59e+03]},\n",
  427. " {energy: 3.78e+03, brem_photons_pe: [198, ..., 180], brem_vtx_z: [...]},\n",
  428. " {energy: 3.07e+03, brem_photons_pe: [1.26e+03, ...], brem_vtx_z: [...]},\n",
  429. " {energy: 3.89e+03, brem_photons_pe: [486], brem_vtx_z: [5.2e+03]},\n",
  430. " {energy: 4.66e+03, brem_photons_pe: [906, ..., 538], brem_vtx_z: [...]},\n",
  431. " ...,\n",
  432. " {energy: 3.08e+03, brem_photons_pe: [575, 337], brem_vtx_z: [...]},\n",
  433. " {energy: 2.39e+03, brem_photons_pe: [163, 971], brem_vtx_z: [...]},\n",
  434. " {energy: 2.57e+03, brem_photons_pe: [184, 798], brem_vtx_z: [...]},\n",
  435. " {energy: 2.19e+03, brem_photons_pe: [161, 337], brem_vtx_z: [...]},\n",
  436. " {energy: 3.62e+03, brem_photons_pe: [1.19e+03], brem_vtx_z: [3.98e+03]},\n",
  437. " {energy: 3.96e+03, brem_photons_pe: [135, ...], brem_vtx_z: [...]},\n",
  438. " {energy: 3.48e+03, brem_photons_pe: [3.17e+03], brem_vtx_z: [8.58e+03]},\n",
  439. " {energy: 3.19e+03, brem_photons_pe: [240, ...], brem_vtx_z: [...]},\n",
  440. " {energy: 3.6e+03, brem_photons_pe: [142, ...], brem_vtx_z: [...]}]\n",
  441. "---------------------------------------------------------------------------\n",
  442. "type: 1452 * {\n",
  443. " energy: float64,\n",
  444. " brem_photons_pe: var * float64,\n",
  445. " brem_vtx_z: var * float64\n",
  446. "}</pre>"
  447. ],
  448. "text/plain": [
  449. "<Array [{energy: 4.19e+03, ...}, ..., {...}] type='1452 * {energy: float64,...'>"
  450. ]
  451. },
  452. "execution_count": 9,
  453. "metadata": {},
  454. "output_type": "execute_result"
  455. }
  456. ],
  457. "source": [
  458. "#wie viel energie relativ zur anfangsenergie verlieren die elektronen durch bremstrahlung und hat das einen einfluss darauf ob wir sie finden oder nicht?\n",
  459. "\n",
  460. "cutoff_energy=350\n",
  461. "\n",
  462. "brem_found = acc_brem_found[ak.sum(acc_brem_found[\"brem_photons_pe\"],axis=-1,keepdims=False)>=cutoff_energy]\n",
  463. "energy_found = ak.to_numpy(brem_found[\"energy\"])\n",
  464. "eph_found = ak.to_numpy(ak.sum(brem_found[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  465. "residual_found = energy_found - eph_found\n",
  466. "energyloss_found = eph_found/energy_found\n",
  467. "\n",
  468. "brem_lost = acc_brem_lost[ak.sum(acc_brem_lost[\"brem_photons_pe\"],axis=-1,keepdims=False)>=cutoff_energy]\n",
  469. "energy_lost = ak.to_numpy(brem_lost[\"energy\"])\n",
  470. "eph_lost = ak.to_numpy(ak.sum(brem_lost[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  471. "residual_lost = energy_lost - eph_lost\n",
  472. "energyloss_lost = eph_lost/energy_lost\n",
  473. "\n",
  474. "print(\"eff = \", np.round(t_eff(brem_found,brem_lost),4), \"+/-\", np.round(eff_err(brem_found, brem_lost),4))\n",
  475. "brem_lost"
  476. ]
  477. },
  478. {
  479. "cell_type": "code",
  480. "execution_count": null,
  481. "metadata": {},
  482. "outputs": [],
  483. "source": []
  484. },
  485. {
  486. "cell_type": "code",
  487. "execution_count": 10,
  488. "metadata": {},
  489. "outputs": [
  490. {
  491. "name": "stdout",
  492. "output_type": "stream",
  493. "text": [
  494. "mean energyloss relative to initial energy (found): 0.3135\n",
  495. "mean energyloss relative to initial energy (lost): 0.4443\n"
  496. ]
  497. }
  498. ],
  499. "source": [
  500. "mean_energyloss_found = ak.mean(energyloss_found)\n",
  501. "mean_energyloss_lost = ak.mean(energyloss_lost)\n",
  502. "print(\"mean energyloss relative to initial energy (found): \", np.round(mean_energyloss_found,4))\n",
  503. "print(\"mean energyloss relative to initial energy (lost): \", np.round(mean_energyloss_lost,4))"
  504. ]
  505. },
  506. {
  507. "cell_type": "code",
  508. "execution_count": 11,
  509. "metadata": {},
  510. "outputs": [
  511. {
  512. "data": {
  513. "image/png": "iVBORw0KGgoAAAANSUhEUgAAAioAAAHMCAYAAAD7xYOxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABLmUlEQVR4nO3deVwU9f8H8Ncut9wqIIpHKKLgVV6h5pV3Hh1WHhl5lWkeaWpWJpZ91fpZaWmmqXR5lFdm5lGJt3lnhiaoJeaBigh4ICzv3x/G5rIL7A67zACv5+OxD52Zz3zmPfPZ483M5zOjExEBERERkQbp1Q6AiIiIKD9MVIiIiEizmKgQERGRZjFRISIiIs1iokJERESaxUSFiIiINIuJChEREWkWExUiIiLSLCYqREREpFlMVOzo6NGjOHv2rNphEBERlRpMVOxowYIF+OWXX9QOg4iIqNRgolJE27dvx9ChQ5GcnGycFx8fj8cffxyXLl1SMTIiIqKSj4kKgA0bNkCn0xlfzs7OqFChAlq3bo3Zs2cjMzMz33UbN26MoKAgNGjQAD///DM+/vhjdOzYEV26dEFAQECJ2IeyJC4uzuQ43fvau3evWfmjR49i8ODBqFmzJjw8PODh4YGwsDC88MILOHDggM3bf+yxx+Dh4YHU1NR8y/Tv3x8uLi4OTXRjYmKg0+lw5coVu9W5e/duxMTEFLhvJZ3BYEBgYCA++OCDfMs44tiWBKW5/fft24fOnTvD29sbXl5eaNeuHXbt2mVSxtbvloyMDIwZMwaVK1eGu7s7GjVqhOXLlysuV5o5qx2AFhw6dAgAsHr1agQHByM7OxuXLl3C5s2b8corryA2NhZbt26Fn5+f2bqenp546623kJOTg+nTp0Ov1+Prr79Gnz59Ssw+lEX/+9//0K5dO5N59erVM5n+9NNP8dJLLyE8PByjR49GZGQkdDodjh8/jmXLlqFp06ZITExEzZo1rd7u4MGDsXbtWixduhTDhw83W379+nWsWbMG3bt3R1BQkLKdU8nu3bsxdepUPPfcc6X2fbZ9+3ZcvnwZjz/+uNqhaE5pbf/9+/ejdevWaNasGb788kuICN599108/PDD2Lp1K6KiokzKW/PdAgCPP/449u/fjxkzZqB27dpYunQp+vbti5ycHPTr18/mcqWakDz22GPi5eUlOTk5ZstiY2MFgIwcOdLiukePHpWGDRtK586dpVu3bjJkyBCpV6+edOvWTS5fvuzo0I2Ksg+lyfnz5+X69ev5Lt+6dasAkG+//bbAenbu3Cl6vV569OghmZmZFst888038s8//9gUX3Z2tlSuXFkaN25scfknn3wiAOT777+3qV5bTZkyRQDY9T363nvvCQA5c+aM3erUmuHDh0uTJk0KLFOUY3vjxg2loamupLZ/Yd8ZnTt3lqCgIJO2SUtLk4oVK0qLFi2M86z9bhER+eGHHwSALF261GR+x44dpXLlypKdnW1TudKOiYqIVKtWTVq1apXv8pCQEKlcubLFZUlJScYflREjRsiSJUskOztbYmNj5c6dO0WO7f3335c1a9YUWk7pPuzYsUM6duwoPj4+4ufnJ926dZOTJ08qLmet7t27S+PGjWXBggXSoEEDcXd3l5CQEHnzzTfFYDDYVFdKSop89tln8vDDD4ter5fDhw/nW9baL5Nu3bqJi4uLnD9/3qZYTp48KX379pWAgABxdXWVOnXqyMcff2xSZtKkSQJAjh49arZ+s2bNJDg4WNEXUO4P5KFDh+Sxxx4Tb29v8fHxkf79+0tycrLFsseOHZM+ffqIj4+PBAYGysCBAyU1NdWs7h07dkj79u3Fy8tLPDw8JCoqStavX29WX97X1q1brVrf1riSk5Nl6NChEhISIq6ursYfjS1btuR7fI4dOyYA5JtvvjHOO3DggACQiIgIk7I9evSQBx54wGReTk6OBAcHy/Tp0/Pdxr37UFg75JY7ePCgPPHEE+Ln5yeVKlUyLrfmvZRbx2+//Sa9e/cWHx8f8ff3l5dfflmysrLkxIkT0rlzZ/Hy8pLq1avLzJkzC4w9rzVr1ggA+emnn8yWzZs3z7jt/Nr/xx9/lEaNGknNmjVN2u/ChQsSFBQkbdq0sfq9bqn+3JetyZEt3xleXl7y9NNPm81//PHHBYDxO8KWRGXIkCHi5eUlWVlZJvOXLl0qAGTXrl02lRMp+nvB0e+loijzicqVK1cKPdvQvn170el0hf6A5iYq9tSvXz9xcXEpMFlRug9TpkwRvV4vgwYNkh9++EFWrlwp9evXl6pVq0p6errN5WwRHBwsnp6eUrduXfnyyy9l8+bN0qdPHwEgCxcuLHT9GzduyPLly6Vnz57i6uoqHh4e8sQTT8i3336b7xkQkf++TAIDA8XJyUm8vb2lU6dOsmPHDmOZ7Oxs44+pLf744w/x9fWV+vXryxdffCGbN2+WcePGiV6vl5iYGGO5hIQE0el0MmbMGLP1Acirr75q03Zz5X7RVK9eXcaPHy+bNm2S999/Xzw9PeX+++83SZxzy4aHh8ubb74pW7Zskffff1/c3Nxk4MCBJvXGxcWJi4uLNG7cWFasWCFr166VTp06iU6nk+XLl4vI3YR95MiRAkBWr14te/bskT179sj169etWt/WuDp37iwBAQGyYMECiYuLk7Vr18qbb75pVl9ewcHB8vzzzxunZ8yYIR4eHgLAeHYsKytLfHx8ZMKECSbr7ty5UwAUmqBb2w73lps4caJs2bJF1q5dKyLWv5fuPV5vv/22bNmyRSZMmCAA5KWXXpI6derInDlzZMuWLTJw4EABIKtWrSow/ntlZWVJYGCg9O/f32xZs2bNjMlcQe1/8uRJ8fb2lscff1xERAwGg7Rv314CAwNt+kMgt87c1y+//CJVqlSRSpUqFXhGJJfS7wxXV1d59tlnzeb37dtXAMimTZtExLrvllwPPvigNG3a1Gx+bjL96aef2lROpOjvBUe/l4qizCcqmzdvFgCyePHifMu0bNlSPD09izGq/2RnZxearCjZh++//14AyLvvvmtS7uTJkwJAvvrqK5vK2eLcuXMCQEJDQ03+yrpz545UqlRJunfvbnG9O3fuyPr166Vfv37i6ekprq6u0r17d/nqq6+sTpgOHToko0ePljVr1sj27dtl8eLFUrduXXFycpKNGzeKiMjFixcFgPTp08ds/ezsbMnKyjK+7r3U1rlzZwkJCTH70nzppZfE3d1dUlJSjPPatGkjFStWNEkexo0bZ9UPYX5yv2hefvllk/lff/21WVvlls3brsOHDxd3d3eT/XrwwQclMDDQ5BhnZ2dLvXr1JCQkxFg2v1P/1q5vS1xeXl5miZ41nnnmGQkNDTVOd+jQQYYOHSr+/v7y+eefi4jIrl27BIBs3rzZZN0xY8ZI/fr1C92Gte2QW+7NN980q8Pa91JuHbNmzTIp16hRI2PSkCsrK0sCAgKMCYO1xo4dKx4eHiaf1fj4eAEgH330kXFeQZd+VqxYIQDkww8/lDfffFP0er3Z8bVFdna29OrVS7y8vOTgwYP5lrPHd0ajRo2kdu3aJn/kZWVlSWhoqMllGWu+W3KFhYVJ586dzbZ1/vx5ASD/+9//bConUvT3QnG8l5Qq84nK9OnTjadp8xMUFGR2GtgeLl++XODpzLwvFxcXuXjxol324f7775eaNWtKZmamyQ9vVlaWeHh4yFtvvWVTOVvknk7+4osvzJa1atVKHnzwQbP5v/32m5QvX16cnJykU6dOsmjRIrl27ZrN27bk2rVrEhISIg0aNBCRghOVhg0bmrTJe++9JyIit27dEmdnZxk5cqTZcdqwYYMAkA0bNhjr+eKLLwSArFy5UkTufvCDgoLkoYceUrwfuV80Bw4cMJmflZUlzs7OMnjwYLOyJ06cMCk7f/58AWB8n2VkZIhOp5Phw4ebbW/mzJkCQI4fPy4iln+obFnflrjat28vfn5+8vbbb8uePXusvsy6ZMkSASCnT5+WW7duibu7u6xatUoef/xxeeaZZ0REZOrUqeLm5iY3b940WbdatWomZzPyY20
  514. "text/plain": [
  515. "<Figure size 640x480 with 1 Axes>"
  516. ]
  517. },
  518. "metadata": {},
  519. "output_type": "display_data"
  520. }
  521. ],
  522. "source": [
  523. "plt.hist(energyloss_lost, bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=\"lost\")\n",
  524. "plt.hist(energyloss_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  525. "plt.xticks(np.arange(0,1.1,0.1), minor=True,)\n",
  526. "plt.yticks(np.arange(0,5.5,0.5), minor=True)\n",
  527. "plt.xlabel(r\"$E_\\gamma/E_0$\")\n",
  528. "plt.ylabel(\"counts (normed)\")\n",
  529. "plt.title(r'$D^\\ast\\rightarrow Dee$, $p<5$GeV, photons w/ brem_vtx_z$<9500$mm')\n",
  530. "plt.legend()\n",
  531. "plt.grid(alpha=0.5)\n",
  532. "\n",
  533. "\"\"\"\n",
  534. "found: elektronen verlieren durchschnittlich 0.44 ihrer anfangsenergie durch bremsstrahlung\n",
  535. "lost: elektronen verlieren durchschnittlich 0.59 ihrer anfangsenergie durch bremsstrahlung\n",
  536. "\n",
  537. "-> lost electrons lose slightly more energy than found electrons. This is however nowhere near as extreme as for the B decay\n",
  538. "\"\"\"\n",
  539. "\n",
  540. "#check if pt near 1 have small eta; perhaps exclude pt with eta outside range; or 2d plot\n",
  541. "\n",
  542. "plt.show()"
  543. ]
  544. },
  545. {
  546. "cell_type": "code",
  547. "execution_count": 12,
  548. "metadata": {},
  549. "outputs": [
  550. {
  551. "data": {
  552. "image/png": "iVBORw0KGgoAAAANSUhEUgAABiUAAAJOCAYAAADYhwTwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAACkVElEQVR4nOzdd3hUVeL/8c8lCWkkIaEk9CYgUgRBQ3ElobNEVBSUIIKVFRBRUOwJFhBchRVErOACAXRV1kYElGBBiiAWcLGhgBKKCaFIS3J+f/jN/BySkLk3mcmQvF/PM8+T3Dn9ljlnzpx7LWOMEQAAAAAAAAAAgJdVKe8CAAAAAAAAAACAyoFJCQAAAAAAAAAA4BNMSgAAAAAAAAAAAJ9gUgIAAAAAAAAAAPgEkxIAAAAAAAAAAMAnmJQAAAAAAAAAAAA+waQEAAAAAAAAAADwCSYlAAAAAAAAAACATzApAQAAAAAAAAAAfIJJCQAAAD/x1VdfaefOneVdDAAAAAAAvIZJCQAAAD/x/PPP68MPPyzvYgAAAAAA4DVMSgAAAJSjjz76SDfffLP27dvn2rZt2zYNGjRIe/fuLceSAQAAAABQ9piUAAAAKKX33ntPlmW5XoGBgapRo4YuueQS/etf/9KJEyeKjduxY0fFxsaqXbt2+uCDDzR79mz17t1b/fr1U61atc6KOlQmGRkZbu3019e6desKhf/qq6904403qlmzZgoNDVVoaKiaN2+uUaNG6fPPP7ed/xVXXKHQ0FAdPHiw2DDDhg1TUFCQVye1UlNTZVmWDhw4UGZprl27VqmpqWes29kuLy9PtWvX1owZM4oN4422PRtU5P2/YcMG9e3bVxEREapWrZoSExP16aefuoWxe205cuSIxo8fr7p16yokJETt27fXkiVLHIcDAADwpcDyLgAAAMDZbvPmzZKkN954Q3Xq1FFubq727t2rFStWaOLEiZo/f75Wr16t6tWrF4obHh6uhx9+WPn5+Zo6daqqVKmiRYsW6Zprrjlr6lAZTZkyRYmJiW7b2rRp4/b/c889p7Fjx6ply5a6/fbb1bp1a1mWpW+//VaLFy/WhRdeqB9++EHNmjXzON8bb7xRy5YtU1pamkaPHl3o/ZycHL355ptKSkpSbGyss8qVk7Vr12ry5MkaOXJkhT3OPvroI+3fv1+DBg0q76L4nYq6/zdu3KhLLrlEF110kRYsWCBjjKZPn66ePXtq9erV6tKli1t4T64tkjRo0CBt3LhRjz/+uFq0aKG0tDQNHTpU+fn5Sk5Oth0OAADAl5iUAAAAKKXNmzerWrVquvzyy2VZlmv7lVdeqa5du2rkyJF66KGH9PTTTxeK+/XXX2v48OGKi4vT3//+d9WtW1ePPfaYFixYoFdeeUU1a9b0+zpUJHv27FF4eLgiIyPPGK558+bq3Llzse9/+umnGj16tAYMGKD//Oc/qlq1quu9Hj16aMyYMXrttdcUGhpqq3z9+/dX3bp19fLLLxc5KbF48WIdO3ZMN954o6104Rv/+c9/1KlTJzVq1Mgr6f/xxx8KCwvzStooWknXjAcffFDVq1dXenq6a9/06tVLTZs21cSJEwutmCjp2iL9ubJt5cqVrgkGSUpMTNQvv/yiu+66S1dffbUCAgI8DgcAAOBr3L4JAACglDZt2qT27du7fZlfYMSIEapfv75ef/31IuNGR0fr0UcfVXp6upo0aaJu3bppy5YtGjJkiKKiokpdthkzZmjZsmUlhnNah08++UR9+vRRVFSUoqOjNWDAAH3//feOw3nq0ksvVadOnfTCCy/o/PPPV2hoqBo0aKCUlBTl5+fbSis7O1svvfSSevXqpfr16+unn35yXK4CU6ZMUUBAgJ577jm3CYm/Gjx4sOrWreu27fvvv1dycrJq166t4OBgtWrVSs8884zr/YCAAI0YMUKbNm3S119/XSjNefPmqU6dOurfv7/tMhfcNuiLL77QoEGDFBkZqaioKF177bXav39/kXH27t2roUOHKioqSrGxsbrhhhuUk5NTKNwnn3yinj17KiIiQmFhYerataveffddt7zvuusuSVKTJk1ct63JyMjwKP7pddi6desZy7V//37dcsstatCggYKDg1WrVi1169ZNq1atKrZ9tm7dKsuy9Nprr7m2bdq0SZZlqXXr1m5hBw4cqI4dO7ptM8bozTff1JVXXllsHn+1a9euM+6Hgrpu3rxZV111laKjo91W3ZR0LP01ja+++kqDBw9WVFSUYmJidOeddyo3N1fbt29Xv379FBERocaNG2v69Okelb3AsmXLZFmWPvjgg0LvPfvss668i9v/6enp6tChg8455xy3/ZeZmam4uDglJCQoLy/Po7IUd2sky7L0888/26qXnWvGp59+qoSEBLfJooiICF1yySVau3at9uzZYytvSXrzzTdVrVo1DR482G379ddfr99++03r16+3FU4q/bHg7WMJAABULExKAAAAlMLvv/+unTt3qkOHDsWGadGihfbs2VPkl+X169dXUlKS27aCL56DgoJKXb7PP/9cQ4YMOePEhNM6pKamqnv37mrQoIEWL16sF198Ubt27VLPnj115MgR2+Hs2LRpk/73v/9pxowZuuuuu/TWW2/p4osv1sMPP6yXX365xPh//PGHli5dqssuu0xxcXG67bbbVL16dS1dulTnnXdeifHHjBmjwMBARUZGqm/fvvrkk09c7+Xl5Wn16tXq1KmT6tSp43Gdtm3bpgsvvFDffPONnnzySb3zzjsaMGCAxo0bp8mTJ7vC3XDDDbIsq1A9t23bpg0bNmjEiBGl+vXzFVdcoXPOOUf/+c9/lJqaqmXLlqlv3746depUobBXXnmlWrRooddff1333HOP0tLSdMcdd7iFWbNmjXr06KGcnBy99NJLWrx4sSIiInTppZdq6dKlkqSbbrpJt912m6Q/byH22Wef6bPPPtMFF1zgUXy75Ro+fLiWLVumhx56SCtWrNCLL76oXr166ffffy+2XVq3bq06deq4TVysWrVKoaGh2rZtm3777TdJUm5urtasWaNevXq5xS/4AtrTSQlP98OgQYN0zjnn6LXXXtPcuXMleX4sFRgyZIjOP/98vf7667r55ps1Y8YM3XHHHbr88ss1YMAAvfnmm+rRo4cmTZqkN954w6PyS1JSUpJq166tefPmFXpv/vz5uuCCC9SuXbti93/Xrl316quvat++fbrhhhskSfn5+Ro2bJiMMVq8eLHHx3pBmgWvDz/8UPXq1VNcXJxiYmJKjO/0mnHy5EkFBwcX2l6w7fTJxTNdWwp88803atWqlQID3W980K5dO9f7dsL9VWmPBW8dSwAAoIIxAAAAcGzFihVGknn55ZeLDdOtWzcTHh7uw1L9f7m5uSY5OdkEBQWZN998s8gwTurw9ttvG0lm+vTpbuG+++47I8ksXLjQVjg7du/ebSSZpk2bmoMHD7q2nzx50sTFxZmkpKQi4508edK88847Jjk52YSHh5uqVauapKQks3DhQnP48GGP8t68ebO5/fbbzZtvvmk++ugj8/LLL5tWrVqZgIAAk56ebowxJjMz00gy11xzTaH4ubm55tSpU65Xfn6+672+ffua+vXrm5ycHLc4Y8eONSEhISYrK8u1rXv37qZmzZrm5MmTrm0TJkwwksx3333nUV1Ol5KSYiSZO+64w237okWLCu2rgrCn79fRo0ebkJAQt3p17tzZ1K5d262Nc3NzTZs2bUz9+vVdYZ944gkjyezYscMtTU/j2ylXtWrVzPjx4+00jzHGmGuvvdY0bdrU9X+vXr3MzTffbKKjo80rr7xijDHm008/NZLMihUr3OKOHz/etG3btsQ8PN0PBeEeeuihQml4eiwVpPHkk0+6hWvfvr2RZN544w3XtlOnTplatWqZQYMGlViHv7rzzjtNaGio27m6bds2I8nMmjXLta24/W+MMUuXLjWSzMyZM81DDz1kqlSpUqh97cjNzTWXXXaZqVatmtm0aVOx4crimtG+fXvTokULk5eX59p26tQp07RpUyPJpKWlGWM8u7YUaN6
  553. "text/plain": [
  554. "<Figure size 2000x600 with 3 Axes>"
  555. ]
  556. },
  557. "metadata": {},
  558. "output_type": "display_data"
  559. }
  560. ],
  561. "source": [
  562. "#energyloss in abh von der energie der elektronen\n",
  563. "\n",
  564. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20,6))\n",
  565. "\n",
  566. "a0=ax0.hist2d(energyloss_found, energy_found, bins=(np.linspace(0,1,50), np.linspace(0,6e3,50)), cmap=plt.cm.jet, cmin=1, vmax=8)\n",
  567. "ax0.set_ylim(0,6e3)\n",
  568. "ax0.set_xlim(0,1)\n",
  569. "ax0.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  570. "ax0.set_ylabel(r\"$E_0$\")\n",
  571. "ax0.set_title(\"found energyloss wrt electron energy\")\n",
  572. "\n",
  573. "a1=ax1.hist2d(energyloss_lost, energy_lost, bins=(np.linspace(0,1,50), np.linspace(0,6e3,50)), cmap=plt.cm.jet, cmin=1, vmax=8)\n",
  574. "ax1.set_ylim(0,6e3)\n",
  575. "ax1.set_xlim(0,1)\n",
  576. "ax1.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  577. "ax1.set_ylabel(r\"$E_0$\")\n",
  578. "ax1.set_title(\"lost energyloss wrt electron energy\")\n",
  579. "\n",
  580. "fig.colorbar(a1[3],ax=ax1)\n",
  581. "fig.suptitle(r\"$D^\\ast\\rightarrow D ee$, $p<5$GeV, photons w/ brem_vtx_z$<9500$mm\")\n",
  582. "\n",
  583. "\n",
  584. "plt.show()"
  585. ]
  586. },
  587. {
  588. "cell_type": "code",
  589. "execution_count": null,
  590. "metadata": {},
  591. "outputs": [],
  592. "source": []
  593. },
  594. {
  595. "cell_type": "code",
  596. "execution_count": 13,
  597. "metadata": {},
  598. "outputs": [],
  599. "source": [
  600. "#ist die shape der teilspur im scifi anders? (koenntest du zum beispiel durch vergleich der verteilungen der fit parameter studieren,\n",
  601. "#in meiner thesis findest du das fitmodell -- ist einfach ein polynom dritten grades)\n",
  602. "z_ref=8520 #mm\n",
  603. "\n",
  604. "def scifi_track(z, a, b, c, d):\n",
  605. " return a + b*(z-z_ref) + c*(z-z_ref)**2 + d*(z-z_ref)**3\n",
  606. "\n",
  607. "def z_mag(xv, zv, tx, a, b):\n",
  608. " \"\"\" optical centre of the magnet is defined as the intersection between the trajectory tangents before and after the magnet\n",
  609. "\n",
  610. " Args:\n",
  611. " xv (double): velo x track\n",
  612. " zv (double): velo z track\n",
  613. " tx (double): velo x slope\n",
  614. " a (double): ax parameter of track fit\n",
  615. " b (double): bx parameter of track fit\n",
  616. "\n",
  617. " Returns:\n",
  618. " double: z_mag\n",
  619. " \"\"\"\n",
  620. " return (xv-tx*zv-a+b*z_ref)/(b-tx)"
  621. ]
  622. },
  623. {
  624. "cell_type": "code",
  625. "execution_count": 14,
  626. "metadata": {},
  627. "outputs": [],
  628. "source": [
  629. "scifi_found = found[found[\"scifi_hit_pos_x_length\"]>3]\n",
  630. "scifi_lost = lost[lost[\"scifi_hit_pos_x_length\"]>3]\n",
  631. "\n",
  632. "scifi_x_found = scifi_found[\"scifi_hit_pos_x\"]\n",
  633. "scifi_z_found = scifi_found[\"scifi_hit_pos_z\"]\n",
  634. "\n",
  635. "tx_found = scifi_found[\"velo_track_tx\"]\n",
  636. "\n",
  637. "scifi_x_lost = scifi_lost[\"scifi_hit_pos_x\"]\n",
  638. "scifi_z_lost = scifi_lost[\"scifi_hit_pos_z\"]\n",
  639. "\n",
  640. "tx_lost = scifi_lost[\"velo_track_tx\"]\n",
  641. "\n",
  642. "xv_found = scifi_found[\"velo_track_x\"]\n",
  643. "zv_found = scifi_found[\"velo_track_z\"]\n",
  644. "\n",
  645. "xv_lost = scifi_lost[\"velo_track_x\"]\n",
  646. "zv_lost = scifi_lost[\"velo_track_z\"]\n",
  647. "\n",
  648. "\n",
  649. "\n",
  650. "#ak.num(scifi_found[\"energy\"], axis=0)\n",
  651. "#scifi_found.snapshot()"
  652. ]
  653. },
  654. {
  655. "cell_type": "code",
  656. "execution_count": 15,
  657. "metadata": {},
  658. "outputs": [],
  659. "source": [
  660. "scifi_fitpars_found = ak.ArrayBuilder()\n",
  661. "\n",
  662. "for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
  663. " popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]))\n",
  664. " scifi_fitpars_found.begin_list()\n",
  665. " scifi_fitpars_found.real(popt[0])\n",
  666. " scifi_fitpars_found.real(popt[1])\n",
  667. " scifi_fitpars_found.real(popt[2])\n",
  668. " scifi_fitpars_found.real(popt[3])\n",
  669. " scifi_fitpars_found.end_list()\n",
  670. "\n",
  671. "scifi_fitpars_lost = ak.ArrayBuilder()\n",
  672. "\n",
  673. "for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
  674. " popt, pcov = curve_fit(scifi_track,ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]))\n",
  675. " scifi_fitpars_lost.begin_list()\n",
  676. " scifi_fitpars_lost.real(popt[0])\n",
  677. " scifi_fitpars_lost.real(popt[1])\n",
  678. " scifi_fitpars_lost.real(popt[2])\n",
  679. " scifi_fitpars_lost.real(popt[3])\n",
  680. " scifi_fitpars_lost.end_list()\n",
  681. "\n",
  682. "\n",
  683. "scifi_fitpars_lost = scifi_fitpars_lost.to_numpy()\n",
  684. "scifi_fitpars_found = scifi_fitpars_found.to_numpy()\n",
  685. "\n",
  686. "\n",
  687. "\n",
  688. "dtx_found = scifi_fitpars_found[:,1] - tx_found\n",
  689. "dtx_lost = scifi_fitpars_lost[:,1] - tx_lost\n"
  690. ]
  691. },
  692. {
  693. "cell_type": "code",
  694. "execution_count": 16,
  695. "metadata": {},
  696. "outputs": [
  697. {
  698. "data": {
  699. "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
  700. "text/plain": [
  701. "<Figure size 1500x1000 with 4 Axes>"
  702. ]
  703. },
  704. "metadata": {},
  705. "output_type": "display_data"
  706. }
  707. ],
  708. "source": [
  709. "fig, ((ax0, ax1), (ax2, ax3)) = plt.subplots(nrows=2, ncols=2, figsize=(15,10))\n",
  710. "\n",
  711. "ax0.hist(scifi_fitpars_found[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$a_x$ found\")\n",
  712. "ax0.hist(scifi_fitpars_lost[:,0], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$a_x$ lost\")\n",
  713. "ax0.set_xlabel(\"a\")\n",
  714. "ax0.set_ylabel(\"normed\")\n",
  715. "ax0.set_title(\"fitparameter a der scifi track\")\n",
  716. "ax0.legend()\n",
  717. "\n",
  718. "ax1.hist(scifi_fitpars_found[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$b_x$ found\")\n",
  719. "ax1.hist(scifi_fitpars_lost[:,1], bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$b_x$ lost\")\n",
  720. "ax1.set_xlabel(\"b\")\n",
  721. "ax1.set_ylabel(\"normed\")\n",
  722. "ax1.set_title(\"fitparameter b der scifi track\")\n",
  723. "ax1.legend()\n",
  724. "\n",
  725. "ax2.hist(scifi_fitpars_found[:,2], bins=200, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$c_x$ found\")\n",
  726. "ax2.hist(scifi_fitpars_lost[:,2], bins=1000, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$c_x$ lost\")\n",
  727. "ax2.set_xlim([-4e-5,4e-5])\n",
  728. "ax2.set_xticks(np.arange(-4e-5,4.5e-5,1e-5),minor=False)\n",
  729. "ax2.set_xlabel(\"c\")\n",
  730. "ax2.set_ylabel(\"normed\")\n",
  731. "ax2.set_title(\"fitparameter c der scifi track\")\n",
  732. "ax2.legend()\n",
  733. "\n",
  734. "ax3.hist(scifi_fitpars_found[:,3], bins=200, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=r\"$d_x$ found\")\n",
  735. "ax3.hist(scifi_fitpars_lost[:,3], bins=1000, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=r\"$d_x$ lost\")\n",
  736. "ax3.set(xlim=(-5e-8,5e-8))\n",
  737. "#ax3.text(-4e-8,3e8,\"d negligible <1e-7\")\n",
  738. "ax3.set_xlabel(\"d\")\n",
  739. "ax3.set_ylabel(\"normed\")\n",
  740. "ax3.set_title(\"fitparameter d der scifi track\")\n",
  741. "ax3.legend()\n",
  742. "\n",
  743. "\"\"\"\n",
  744. "a_x: virtual hit on the reference plane\n",
  745. "\"\"\"\n",
  746. "\n",
  747. "plt.show()"
  748. ]
  749. },
  750. {
  751. "cell_type": "code",
  752. "execution_count": 17,
  753. "metadata": {},
  754. "outputs": [
  755. {
  756. "name": "stdout",
  757. "output_type": "stream",
  758. "text": [
  759. "found\n",
  760. "a = 18.04180550457298\n",
  761. "b = 0.0057651944487250376\n",
  762. "c = 9.480277789440454e-08\n",
  763. "d = -4.452015498411874e-11\n",
  764. "lost\n",
  765. "a = -35.48371609662243\n",
  766. "b = -0.010382326385451382\n",
  767. "c = -6.208301288913938e-07\n",
  768. "d = 9.580933791481267e-11\n"
  769. ]
  770. }
  771. ],
  772. "source": [
  773. "print(\"found\")\n",
  774. "print(\"a = \", str(np.mean(scifi_fitpars_found[:,0])))\n",
  775. "print(\"b = \", str(np.mean(scifi_fitpars_found[:,1])))\n",
  776. "print(\"c = \", str(np.mean(scifi_fitpars_found[:,2])))\n",
  777. "print(\"d = \", str(np.mean(scifi_fitpars_found[:,3])))\n",
  778. "\n",
  779. "print(\"lost\")\n",
  780. "print(\"a = \", str(np.mean(scifi_fitpars_lost[:,0])))\n",
  781. "print(\"b = \", str(np.mean(scifi_fitpars_lost[:,1])))\n",
  782. "print(\"c = \", str(np.mean(scifi_fitpars_lost[:,2])))\n",
  783. "print(\"d = \", str(np.mean(scifi_fitpars_lost[:,3])))"
  784. ]
  785. },
  786. {
  787. "cell_type": "code",
  788. "execution_count": 18,
  789. "metadata": {},
  790. "outputs": [
  791. {
  792. "data": {
  793. "text/plain": [
  794. "-1.5438992626615335e-08"
  795. ]
  796. },
  797. "execution_count": 18,
  798. "metadata": {},
  799. "output_type": "execute_result"
  800. }
  801. ],
  802. "source": [
  803. "np.min(scifi_fitpars_found[:,3])"
  804. ]
  805. },
  806. {
  807. "cell_type": "code",
  808. "execution_count": 19,
  809. "metadata": {},
  810. "outputs": [
  811. {
  812. "data": {
  813. "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
  814. "text/plain": [
  815. "<Figure size 1500x600 with 2 Axes>"
  816. ]
  817. },
  818. "metadata": {},
  819. "output_type": "display_data"
  820. }
  821. ],
  822. "source": [
  823. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(15,6))\n",
  824. "\n",
  825. "for i in range(0,ak.num(scifi_found[\"energy\"], axis=0)):\n",
  826. " z_coord = np.linspace(scifi_z_found[i,0],12000,300)\n",
  827. " fit = scifi_track(z_coord, *scifi_fitpars_found[i])\n",
  828. " ax0.plot(z_coord, fit, \"-\", lw=0.5)\n",
  829. " ax0.errorbar(ak.to_numpy(scifi_z_found[i,:]),ak.to_numpy(scifi_x_found[i,:]),fmt=\".\",ms=2)\n",
  830. "\n",
  831. "#ax0.legend()\n",
  832. "ax0.set_xlabel(\"z [mm]\")\n",
  833. "ax0.set_ylabel(\"x [mm]\")\n",
  834. "ax0.set_title(\"found tracks of scifi hits\")\n",
  835. "ax0.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
  836. "ax0.grid()\n",
  837. "\n",
  838. "for i in range(0,ak.num(scifi_lost[\"energy\"], axis=0)):\n",
  839. " z_coord = np.linspace(scifi_z_lost[i,0],12000,300)\n",
  840. " fit = scifi_track(z_coord, *scifi_fitpars_lost[i])\n",
  841. " ax1.plot(z_coord, fit, \"-\", lw=0.5)\n",
  842. " ax1.errorbar(ak.to_numpy(scifi_z_lost[i,:]),ak.to_numpy(scifi_x_lost[i,:]),fmt=\".\",ms=2)\n",
  843. "\n",
  844. "#ax1.legend()\n",
  845. "ax1.set_xlabel(\"z [mm]\")\n",
  846. "ax1.set_ylabel(\"x [mm]\")\n",
  847. "ax1.set_title(\"lost tracks of scifi hits\")\n",
  848. "ax1.set(xlim=(7e3,12000), ylim=(-4000,4000))\n",
  849. "ax1.grid()\n",
  850. "\n",
  851. "plt.show()"
  852. ]
  853. },
  854. {
  855. "cell_type": "code",
  856. "execution_count": 20,
  857. "metadata": {},
  858. "outputs": [
  859. {
  860. "name": "stdout",
  861. "output_type": "stream",
  862. "text": [
  863. "found \n",
  864. "zmag = 5318.452650242005\n",
  865. "lost \n",
  866. "zmag = 5425.137423441005\n"
  867. ]
  868. }
  869. ],
  870. "source": [
  871. "#vergleich der zmag werte\n",
  872. "zmag_found = z_mag(xv_found, zv_found, tx_found, scifi_fitpars_found[:,0], scifi_fitpars_found[:,1])\n",
  873. "zmag_lost = z_mag(xv_lost, zv_lost, tx_lost, scifi_fitpars_lost[:,0], scifi_fitpars_lost[:,1])\n",
  874. "zmag_lost = zmag_lost[~np.isnan(zmag_lost)]\n",
  875. "zmag_found = zmag_found[~np.isnan(zmag_found)]\n",
  876. "\n",
  877. "print(\"found \\nzmag = \", str(np.mean(zmag_found)))\n",
  878. "print(\"lost \\nzmag = \", str(np.mean(zmag_lost)))"
  879. ]
  880. },
  881. {
  882. "cell_type": "code",
  883. "execution_count": 21,
  884. "metadata": {},
  885. "outputs": [
  886. {
  887. "data": {
  888. "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
  889. "text/plain": [
  890. "<Figure size 640x480 with 1 Axes>"
  891. ]
  892. },
  893. "metadata": {},
  894. "output_type": "display_data"
  895. }
  896. ],
  897. "source": [
  898. "plt.hist(zmag_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  899. "plt.hist(zmag_lost, bins=3000, density=True, alpha=0.5, histtype=\"bar\",color=\"darkorange\", label=\"lost\")\n",
  900. "plt.xlabel(\"$z_{mag}$ [mm]\")\n",
  901. "plt.ylabel(\"normed\")\n",
  902. "plt.title(\"magnet kick position $z_{mag}$ calculated w fitparameters\")\n",
  903. "plt.legend()\n",
  904. "plt.xticks(np.arange(5100,5805,5), minor=True)\n",
  905. "plt.yticks(np.arange(0,0.011,0.001), minor=True)\n",
  906. "plt.xlim(5100,5800)\n",
  907. "\n",
  908. "\"\"\"\n",
  909. "wir können definitiv einen unterschied zwischen den z_mag werten, bzw deren verteilungen für lost and found, erkennen\n",
  910. "\"\"\"\n",
  911. "\n",
  912. "plt.show()"
  913. ]
  914. },
  915. {
  916. "cell_type": "code",
  917. "execution_count": 22,
  918. "metadata": {},
  919. "outputs": [
  920. {
  921. "name": "stdout",
  922. "output_type": "stream",
  923. "text": [
  924. "2591\n",
  925. "1908\n",
  926. "4499\n",
  927. "5767\n"
  928. ]
  929. }
  930. ],
  931. "source": [
  932. "both = ak.concatenate([found,lost],axis=0)\n",
  933. "print(ak.num(found,axis=0))\n",
  934. "print(ak.num(lost,axis=0))\n",
  935. "print(ak.num(both,axis=0))\n",
  936. "\n",
  937. "tracks = allcolumns[(allcolumns.fromSignal) & (allcolumns.isElectron) & (~allcolumns.lost)]\n",
  938. "alltracks = allcolumns[(allcolumns.fromSignal) & (allcolumns.isElectron)]\n",
  939. "print(ak.num(tracks,axis=0))"
  940. ]
  941. },
  942. {
  943. "cell_type": "code",
  944. "execution_count": 23,
  945. "metadata": {},
  946. "outputs": [
  947. {
  948. "name": "stdout",
  949. "output_type": "stream",
  950. "text": [
  951. "#events w/ shared track electrons from found and lost: 35\n",
  952. "event_count: [359, 359]\n",
  953. "velo idx: [25, 25]\n",
  954. "mcp_index: [2926, 2936]\n",
  955. "\n",
  956. "velo x: [5.89, 5.89]\n",
  957. "velo y: [9.81, 9.81]\n",
  958. "\n",
  959. "velo tx: [0.00824, 0.00824]\n",
  960. "velo ty: [0.0135, 0.0135]\n",
  961. "percentage of e with shared tracks: 0.7488\n"
  962. ]
  963. }
  964. ],
  965. "source": [
  966. "#versuche teilchen von denselben events mit shared tracks zu finden. \n",
  967. "#idee: alle teilchen eines events sind durch event_count auffindbar. \n",
  968. "a_f_itr = tracks[\"event_count\"].to_numpy()\n",
  969. "f_itr = np.unique(a_f_itr)\n",
  970. "\n",
  971. "shared = ak.ArrayBuilder()\n",
  972. "count = 0\n",
  973. "\n",
  974. "for itr in f_itr:\n",
  975. " temp = alltracks[alltracks[\"event_count\"]==itr]\n",
  976. " if ak.num(temp,axis=0)>1:\n",
  977. " #iterate over cols in temp and append all with duplicate velo_track_idx, such that possibles is array with possible shared tracks particles\n",
  978. " #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",
  979. " #concatenate to the array of other shared track particles\n",
  980. " _jitr = temp[\"velo_track_idx\"].to_numpy()\n",
  981. " jitr = np.unique(_jitr)\n",
  982. " for jentry in jitr:\n",
  983. " jtem = temp[temp[\"velo_track_idx\"]==jentry]\n",
  984. " if ak.num(jtem,axis=0)>1:\n",
  985. " shared.append(jtem)\n",
  986. " else:\n",
  987. " continue\n",
  988. " else:\n",
  989. " continue\n",
  990. "shared = ak.Array(shared)\n",
  991. "\n",
  992. "idx=0\n",
  993. "print(\"#events w/ shared track electrons from found and lost: \",ak.num(shared, axis=0))\n",
  994. "\n",
  995. "print(\"event_count: \", shared[idx,:,\"event_count\"])\n",
  996. "print(\"velo idx: \" ,shared[idx,:,\"velo_track_idx\"])\n",
  997. "print(\"mcp_index: \", shared[idx,:,\"mcp_idx\"])\n",
  998. "\n",
  999. "print(\"\\nvelo x: \" ,shared[idx,:,\"velo_track_x\"])\n",
  1000. "print(\"velo y: \" ,shared[idx,:,\"velo_track_y\"])\n",
  1001. "\n",
  1002. "print(\"\\nvelo tx: \" ,shared[idx,:,\"velo_track_tx\"])\n",
  1003. "print(\"velo ty: \" ,shared[idx,:,\"velo_track_ty\"])\n",
  1004. "\n",
  1005. "\n",
  1006. "\n",
  1007. "print(\"percentage of e with shared tracks: \", np.round((ak.num(shared,axis=0)*2)/(ak.num(alltracks,axis=0))*100,4)) #error for percentage"
  1008. ]
  1009. },
  1010. {
  1011. "cell_type": "code",
  1012. "execution_count": 24,
  1013. "metadata": {},
  1014. "outputs": [
  1015. {
  1016. "name": "stdout",
  1017. "output_type": "stream",
  1018. "text": [
  1019. "#velo_track_idx in all events: 217\n",
  1020. "velo idx: [0, 0, 0]\n",
  1021. "mcp_index: [1066, 1251, 666]\n",
  1022. "event_count: [5735, 7049, 7378]\n",
  1023. "\n",
  1024. "velo x: [-2.21, 9.45, -18.5]\n",
  1025. "velo y: [-21.6, -33.4, 17.4]\n",
  1026. "\n",
  1027. "velo tx: [-0.0022, 0.0153, -0.0224]\n",
  1028. "velo ty: [-0.0263, -0.0469, 0.0212]\n"
  1029. ]
  1030. }
  1031. ],
  1032. "source": [
  1033. "#electrons with same velo_track_idx from all events\n",
  1034. "temp_ = found[\"velo_track_idx\"].to_numpy()\n",
  1035. "temp = np.unique(temp_)\n",
  1036. "count=0\n",
  1037. "psb=ak.ArrayBuilder()\n",
  1038. "for jentry in temp:\n",
  1039. " jtem = found[found[\"velo_track_idx\"]==jentry]\n",
  1040. " if ak.num(jtem,axis=0)>1:\n",
  1041. " psb.append(jtem)\n",
  1042. " else:\n",
  1043. " continue\n",
  1044. "\n",
  1045. "psb = ak.Array(psb)\n",
  1046. "\n",
  1047. "print(\"#velo_track_idx in all events: \",ak.num(psb,axis=0))\n",
  1048. "idx = 0\n",
  1049. "print(\"velo idx: \" ,psb[idx,:,\"velo_track_idx\"])\n",
  1050. "print(\"mcp_index: \", psb[idx,:,\"mcp_idx\"])\n",
  1051. "print(\"event_count: \", psb[idx,:,\"event_count\"])\n",
  1052. "print(\"\\nvelo x: \" ,psb[idx,:,\"velo_track_x\"])\n",
  1053. "print(\"velo y: \" ,psb[idx,:,\"velo_track_y\"])\n",
  1054. "\n",
  1055. "print(\"\\nvelo tx: \" ,psb[idx,:,\"velo_track_tx\"])\n",
  1056. "print(\"velo ty: \" ,psb[idx,:,\"velo_track_ty\"])\n",
  1057. "\n"
  1058. ]
  1059. },
  1060. {
  1061. "cell_type": "code",
  1062. "execution_count": null,
  1063. "metadata": {},
  1064. "outputs": [],
  1065. "source": []
  1066. },
  1067. {
  1068. "cell_type": "code",
  1069. "execution_count": 25,
  1070. "metadata": {},
  1071. "outputs": [
  1072. {
  1073. "name": "stdout",
  1074. "output_type": "stream",
  1075. "text": [
  1076. "#events w/ shared track electrons from Photon Conversions: 951\n",
  1077. "shared_idx: 0\n",
  1078. "event_count: [11, 11]\n",
  1079. "velo idx: [27, 27]\n",
  1080. "mcp_index: [1211, 1215]\n",
  1081. "\n",
  1082. "velo x: [-20.4, -20.4]\n",
  1083. "velo y: [-24.4, -24.4]\n",
  1084. "\n",
  1085. "velo tx: [-0.0234, -0.0234]\n",
  1086. "velo ty: [-0.028, -0.028]\n",
  1087. "percentage of e with shared tracks: 5.59\n"
  1088. ]
  1089. }
  1090. ],
  1091. "source": [
  1092. "#generell: wie viele elektronen von photon conversions teilen sich einen track?\n",
  1093. "conv = allcolumns[(allcolumns.fromPairProd) & (allcolumns.isElectron) & (~allcolumns.lost)]\n",
  1094. "\n",
  1095. "conv_it = conv[\"event_count\"].to_numpy()\n",
  1096. "conv_itr = np.unique(conv_it)\n",
  1097. "\n",
  1098. "cshared = ak.ArrayBuilder()\n",
  1099. "count = 0\n",
  1100. "\n",
  1101. "for itr in conv_itr:\n",
  1102. " temp = conv[conv[\"event_count\"]==itr]\n",
  1103. " if ak.num(temp,axis=0)>1:\n",
  1104. " #iterate over cols in temp and append all with duplicate velo_track_idx, such that possibles is array with possible shared tracks particles\n",
  1105. " #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",
  1106. " #concatenate to the array of other shared track particles\n",
  1107. " _jitr = temp[\"velo_track_idx\"].to_numpy()\n",
  1108. " jitr = np.unique(_jitr)\n",
  1109. " for jentry in jitr:\n",
  1110. " jtem = temp[temp[\"velo_track_idx\"]==jentry]\n",
  1111. " if ak.num(jtem,axis=0)>1:\n",
  1112. " cshared.append(jtem)\n",
  1113. " else:\n",
  1114. " continue\n",
  1115. " else:\n",
  1116. " continue\n",
  1117. "cshared = ak.Array(cshared)\n",
  1118. "\n",
  1119. "\n",
  1120. "print(\"#events w/ shared track electrons from Photon Conversions: \",ak.num(cshared, axis=0))\n",
  1121. "idx = 0\n",
  1122. "print(\"shared_idx: \", idx)\n",
  1123. "print(\"event_count: \", cshared[idx,:,\"event_count\"])\n",
  1124. "print(\"velo idx: \" ,cshared[idx,:,\"velo_track_idx\"])\n",
  1125. "print(\"mcp_index: \", cshared[idx,:,\"mcp_idx\"])\n",
  1126. "\n",
  1127. "print(\"\\nvelo x: \" ,cshared[idx,:,\"velo_track_x\"])\n",
  1128. "print(\"velo y: \" ,cshared[idx,:,\"velo_track_y\"])\n",
  1129. "\n",
  1130. "print(\"\\nvelo tx: \" ,cshared[idx,:,\"velo_track_tx\"])\n",
  1131. "print(\"velo ty: \" ,cshared[idx,:,\"velo_track_ty\"])\n",
  1132. "\n",
  1133. "print(\"percentage of e with shared tracks: \", np.round((ak.num(cshared,axis=0)*2)/(ak.num(conv,axis=0))*100,4))\n",
  1134. "\n",
  1135. "#constrain z pos of brem vtx; only interested in conversions in velo z<770\n"
  1136. ]
  1137. },
  1138. {
  1139. "cell_type": "code",
  1140. "execution_count": 26,
  1141. "metadata": {},
  1142. "outputs": [
  1143. {
  1144. "data": {
  1145. "text/plain": [
  1146. "34025"
  1147. ]
  1148. },
  1149. "execution_count": 26,
  1150. "metadata": {},
  1151. "output_type": "execute_result"
  1152. }
  1153. ],
  1154. "source": [
  1155. "ak.num(conv,axis=0)"
  1156. ]
  1157. },
  1158. {
  1159. "cell_type": "code",
  1160. "execution_count": null,
  1161. "metadata": {},
  1162. "outputs": [],
  1163. "source": []
  1164. },
  1165. {
  1166. "cell_type": "code",
  1167. "execution_count": null,
  1168. "metadata": {},
  1169. "outputs": [],
  1170. "source": []
  1171. },
  1172. {
  1173. "cell_type": "code",
  1174. "execution_count": null,
  1175. "metadata": {},
  1176. "outputs": [],
  1177. "source": []
  1178. },
  1179. {
  1180. "cell_type": "code",
  1181. "execution_count": null,
  1182. "metadata": {},
  1183. "outputs": [],
  1184. "source": []
  1185. },
  1186. {
  1187. "cell_type": "code",
  1188. "execution_count": null,
  1189. "metadata": {},
  1190. "outputs": [],
  1191. "source": []
  1192. },
  1193. {
  1194. "cell_type": "code",
  1195. "execution_count": 27,
  1196. "metadata": {},
  1197. "outputs": [
  1198. {
  1199. "data": {
  1200. "text/plain": [
  1201. "'\\nwas meinst du? velo track sharen? wie überprüfe ich das? versuche ich einfach e teilchen mit identischen velo tracks und slopes zu finden?\\nbenutze ich mother key um elektronen von gleichen events zu finden? meinst du wie viele elektronen vom selben event teilen sich einen track,\\noder wie viele teilchen von unterschiedlichen events haben dieselben tracks?\\n'"
  1202. ]
  1203. },
  1204. "execution_count": 27,
  1205. "metadata": {},
  1206. "output_type": "execute_result"
  1207. }
  1208. ],
  1209. "source": [
  1210. "\"\"\"\n",
  1211. "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",
  1212. "benutze ich mother key um elektronen von gleichen events zu finden? meinst du wie viele elektronen vom selben event teilen sich einen track,\n",
  1213. "oder wie viele teilchen von unterschiedlichen events haben dieselben tracks?\n",
  1214. "\"\"\"\n"
  1215. ]
  1216. },
  1217. {
  1218. "cell_type": "code",
  1219. "execution_count": null,
  1220. "metadata": {},
  1221. "outputs": [],
  1222. "source": []
  1223. },
  1224. {
  1225. "cell_type": "code",
  1226. "execution_count": null,
  1227. "metadata": {},
  1228. "outputs": [],
  1229. "source": []
  1230. },
  1231. {
  1232. "cell_type": "code",
  1233. "execution_count": null,
  1234. "metadata": {},
  1235. "outputs": [],
  1236. "source": []
  1237. }
  1238. ],
  1239. "metadata": {
  1240. "kernelspec": {
  1241. "display_name": "env1",
  1242. "language": "python",
  1243. "name": "python3"
  1244. },
  1245. "language_info": {
  1246. "codemirror_mode": {
  1247. "name": "ipython",
  1248. "version": 3
  1249. },
  1250. "file_extension": ".py",
  1251. "mimetype": "text/x-python",
  1252. "name": "python",
  1253. "nbconvert_exporter": "python",
  1254. "pygments_lexer": "ipython3",
  1255. "version": "3.11.5"
  1256. },
  1257. "orig_nbformat": 4
  1258. },
  1259. "nbformat": 4,
  1260. "nbformat_minor": 2
  1261. }