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.

863 lines
383 KiB

  1. {
  2. "cells": [
  3. {
  4. "cell_type": "code",
  5. "execution_count": 125,
  6. "metadata": {},
  7. "outputs": [],
  8. "source": [
  9. "import uproot\t\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. "from mpl_toolkits.axes_grid1 import ImageGrid\n",
  20. "%matplotlib inline"
  21. ]
  22. },
  23. {
  24. "cell_type": "code",
  25. "execution_count": 126,
  26. "metadata": {},
  27. "outputs": [
  28. {
  29. "data": {
  30. "text/plain": [
  31. "10522"
  32. ]
  33. },
  34. "execution_count": 126,
  35. "metadata": {},
  36. "output_type": "execute_result"
  37. }
  38. ],
  39. "source": [
  40. "file = uproot.open(\"tracking_losses_ntuple_Bd2KstEE.root:PrDebugTrackingLosses.PrDebugTrackingTool/Tuple;1\")\n",
  41. "\n",
  42. "#selektiere nur elektronen von B->K*ee und nur solche mit einem momentum von ueber 5 GeV \n",
  43. "allcolumns = file.arrays()\n",
  44. "found = allcolumns[(allcolumns.isElectron) & (~allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p > 5e3)] #B: 9056\n",
  45. "lost = allcolumns[(allcolumns.isElectron) & (allcolumns.lost) & (allcolumns.fromSignal) & (allcolumns.p > 5e3)] #B: 1466\n",
  46. "\n",
  47. "ak.num(found, axis=0) + ak.num(lost, axis=0)\n",
  48. "#ak.count(found, axis=None)"
  49. ]
  50. },
  51. {
  52. "cell_type": "code",
  53. "execution_count": 127,
  54. "metadata": {},
  55. "outputs": [
  56. {
  57. "name": "stdout",
  58. "output_type": "stream",
  59. "text": [
  60. "eff all = 0.8606728758791105 +/- 0.003375885792719708\n"
  61. ]
  62. }
  63. ],
  64. "source": [
  65. "def t_eff(found, lost, axis = 0):\n",
  66. " sel = ak.num(found, axis=axis)\n",
  67. " des = ak.num(lost, axis=axis)\n",
  68. " return sel/(sel + des)\n",
  69. "\n",
  70. "def eff_err(found, lost):\n",
  71. " n_f = ak.num(found, axis=0)\n",
  72. " n_all = ak.num(found, axis=0) + ak.num(lost,axis=0)\n",
  73. " return 1/n_all * np.sqrt(np.abs(n_f*(1-n_f/n_all)))\n",
  74. "\n",
  75. "\n",
  76. "print(\"eff all = \", t_eff(found, lost), \"+/-\", eff_err(found, lost))"
  77. ]
  78. },
  79. {
  80. "cell_type": "code",
  81. "execution_count": 128,
  82. "metadata": {},
  83. "outputs": [],
  84. "source": [
  85. "#try excluding all photons that originate from a vtx @ z>9500mm\n",
  86. "#ignore all brem vertices @ z>9500mm \n",
  87. "\n",
  88. "#found\n",
  89. "\n",
  90. "brem_e_f = found[\"brem_photons_pe\"]\n",
  91. "brem_z_f = found[\"brem_vtx_z\"]\n",
  92. "e_f = found[\"energy\"]\n",
  93. "length_f = found[\"brem_vtx_z_length\"]\n",
  94. "\n",
  95. "brem_f = ak.ArrayBuilder()\n",
  96. "\n",
  97. "for itr in range(ak.num(found,axis=0)):\n",
  98. " brem_f.begin_record()\n",
  99. " #[:,\"energy\"] energy\n",
  100. " brem_f.field(\"energy\").append(e_f[itr])\n",
  101. " #[:,\"photon_length\"] number of vertices\n",
  102. " brem_f.field(\"photon_length\").integer(length_f[itr])\n",
  103. " #[:,\"brem_photons_pe\",:] photon energy \n",
  104. " brem_f.field(\"brem_photons_pe\").append(brem_e_f[itr])\n",
  105. " #[:,\"brem_vtx_z\",:] brem vtx z\n",
  106. " brem_f.field(\"brem_vtx_z\").append(brem_z_f[itr])\n",
  107. " brem_f.end_record()\n",
  108. "\n",
  109. "brem_f = ak.Array(brem_f)\n",
  110. "\n",
  111. "#lost\n",
  112. "\n",
  113. "brem_e_l = lost[\"brem_photons_pe\"]\n",
  114. "brem_z_l = lost[\"brem_vtx_z\"]\n",
  115. "e_l = lost[\"energy\"]\n",
  116. "length_l = lost[\"brem_vtx_z_length\"]\n",
  117. "\n",
  118. "brem_l = ak.ArrayBuilder()\n",
  119. "\n",
  120. "for itr in range(ak.num(lost,axis=0)):\n",
  121. " brem_l.begin_record()\n",
  122. " #[:,\"energy\"] energy\n",
  123. " brem_l.field(\"energy\").append(e_l[itr])\n",
  124. " #[:,\"photon_length\"] number of vertices\n",
  125. " brem_l.field(\"photon_length\").integer(length_l[itr])\n",
  126. " #[:,\"brem_photons_pe\",:] photon energy \n",
  127. " brem_l.field(\"brem_photons_pe\").append(brem_e_l[itr])\n",
  128. " #[:,\"brem_vtx_z\",:] brem vtx z\n",
  129. " brem_l.field(\"brem_vtx_z\").append(brem_z_l[itr])\n",
  130. " brem_l.end_record()\n",
  131. "\n",
  132. "brem_l = ak.Array(brem_l)"
  133. ]
  134. },
  135. {
  136. "cell_type": "code",
  137. "execution_count": 129,
  138. "metadata": {},
  139. "outputs": [],
  140. "source": [
  141. "cut_brem_found = ak.ArrayBuilder()\n",
  142. "\n",
  143. "for itr in range(ak.num(brem_f, axis=0)):\n",
  144. " cut_brem_found.begin_record()\n",
  145. " cut_brem_found.field(\"energy\").real(brem_f[itr,\"energy\"])\n",
  146. " \n",
  147. " cut_brem_found.field(\"brem_photons_pe\")\n",
  148. " cut_brem_found.begin_list()\n",
  149. " for jentry in range(brem_f[itr, \"photon_length\"]):\n",
  150. " if brem_f[itr, \"brem_vtx_z\", jentry]>9500:\n",
  151. " continue\n",
  152. " else:\n",
  153. " cut_brem_found.real(brem_f[itr,\"brem_photons_pe\", jentry])\n",
  154. " \n",
  155. " #cut_brem_found.field(\"brem_vtx_z\").real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  156. " cut_brem_found.end_list()\n",
  157. " \n",
  158. " cut_brem_found.field(\"brem_vtx_z\")\n",
  159. " cut_brem_found.begin_list()\n",
  160. " for jentry in range(brem_f[itr, \"photon_length\"]):\n",
  161. " if brem_f[itr, \"brem_vtx_z\", jentry]>9500:\n",
  162. " continue\n",
  163. " else:\n",
  164. " cut_brem_found.real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  165. " cut_brem_found.end_list()\n",
  166. " \n",
  167. "\n",
  168. " \n",
  169. " cut_brem_found.end_record()\n",
  170. "\n",
  171. "cut_brem_found = ak.Array(cut_brem_found)\n",
  172. "\n",
  173. "\n",
  174. "\n",
  175. "cut_brem_lost = ak.ArrayBuilder()\n",
  176. "\n",
  177. "for itr in range(ak.num(brem_l, axis=0)):\n",
  178. " cut_brem_lost.begin_record()\n",
  179. " cut_brem_lost.field(\"energy\").real(brem_l[itr,\"energy\"])\n",
  180. " \n",
  181. " \n",
  182. " cut_brem_lost.field(\"brem_photons_pe\")\n",
  183. " cut_brem_lost.begin_list()\n",
  184. " for jentry in range(brem_l[itr, \"photon_length\"]):\n",
  185. " if brem_l[itr, \"brem_vtx_z\", jentry]>9500:\n",
  186. " continue\n",
  187. " else:\n",
  188. " cut_brem_lost.real(brem_l[itr,\"brem_photons_pe\", jentry])\n",
  189. " \n",
  190. " #cut_brem_found.field(\"brem_vtx_z\").real(brem_f[itr, \"brem_vtx_z\",jentry])\n",
  191. " cut_brem_lost.end_list()\n",
  192. " \n",
  193. " cut_brem_lost.field(\"brem_vtx_z\")\n",
  194. " cut_brem_lost.begin_list()\n",
  195. " for jentry in range(brem_l[itr, \"photon_length\"]):\n",
  196. " if brem_l[itr, \"brem_vtx_z\", jentry]>9500:\n",
  197. " continue\n",
  198. " else:\n",
  199. " cut_brem_lost.real(brem_l[itr, \"brem_vtx_z\",jentry])\n",
  200. " cut_brem_lost.end_list()\n",
  201. " \n",
  202. " cut_brem_lost.end_record()\n",
  203. "\n",
  204. "cut_brem_lost = ak.Array(cut_brem_lost)\n"
  205. ]
  206. },
  207. {
  208. "cell_type": "code",
  209. "execution_count": 130,
  210. "metadata": {},
  211. "outputs": [
  212. {
  213. "data": {
  214. "text/html": [
  215. "<pre>{energy: 9.36e+03,\n",
  216. " brem_photons_pe: [2.47e+03, 170, 224, 388, 3.23e+03, 809, 172, 224],\n",
  217. " brem_vtx_z: [400, 501, 638, 667, 677, 709, 8.58e+03, 9.28e+03]}\n",
  218. "---------------------------------------------------------------------\n",
  219. "type: {\n",
  220. " energy: float64,\n",
  221. " brem_photons_pe: var * float64,\n",
  222. " brem_vtx_z: var * float64\n",
  223. "}</pre>"
  224. ],
  225. "text/plain": [
  226. "<Record {energy: 9.36e+03, ...} type='{energy: float64, brem_photons_pe: va...'>"
  227. ]
  228. },
  229. "execution_count": 130,
  230. "metadata": {},
  231. "output_type": "execute_result"
  232. }
  233. ],
  234. "source": [
  235. "#data in cut_brem_found and cut_brem_lost\n",
  236. "\n",
  237. "cut_length_found = ak.num(cut_brem_found[\"brem_photons_pe\"],axis=-1)\n",
  238. "cut_length_lost = ak.num(cut_brem_lost[\"brem_photons_pe\"], axis=-1)\n",
  239. "\n",
  240. "cut_brem_found[1]\n"
  241. ]
  242. },
  243. {
  244. "cell_type": "code",
  245. "execution_count": 131,
  246. "metadata": {},
  247. "outputs": [
  248. {
  249. "data": {
  250. "text/plain": [
  251. "8"
  252. ]
  253. },
  254. "execution_count": 131,
  255. "metadata": {},
  256. "output_type": "execute_result"
  257. }
  258. ],
  259. "source": [
  260. "cut_length_found[1]"
  261. ]
  262. },
  263. {
  264. "cell_type": "markdown",
  265. "metadata": {},
  266. "source": [
  267. "### Split in Upstream and Downstream Events and analyse separately"
  268. ]
  269. },
  270. {
  271. "cell_type": "code",
  272. "execution_count": 132,
  273. "metadata": {},
  274. "outputs": [],
  275. "source": [
  276. "#try to find a split between energy lost before and after the magnet (z~5000mm)\n",
  277. "\n",
  278. "upstream_found = ak.ArrayBuilder()\n",
  279. "downstream_found = ak.ArrayBuilder()\n",
  280. "\n",
  281. "for itr in range(ak.num(cut_brem_found, axis=0)):\n",
  282. " upstream_found.begin_record()\n",
  283. " upstream_found.field(\"energy\").real(cut_brem_found[itr,\"energy\"])\n",
  284. " \n",
  285. " downstream_found.begin_record()\n",
  286. " downstream_found.field(\"energy\").real(cut_brem_found[itr,\"energy\"])\n",
  287. " \n",
  288. " upstream_found.field(\"brem_photons_pe\")\n",
  289. " downstream_found.field(\"brem_photons_pe\")\n",
  290. " upstream_found.begin_list()\n",
  291. " downstream_found.begin_list()\n",
  292. " for jentry in range(cut_length_found[itr]):\n",
  293. " if (cut_brem_found[itr, \"brem_vtx_z\", jentry]>5000):\n",
  294. " if cut_brem_found[itr, \"brem_vtx_z\", jentry]<=9500:\n",
  295. " downstream_found.real(cut_brem_found[itr,\"brem_photons_pe\",jentry])\n",
  296. " else:\n",
  297. " continue\n",
  298. " else:\n",
  299. " upstream_found.real(cut_brem_found[itr,\"brem_photons_pe\", jentry]) \n",
  300. " upstream_found.end_list()\n",
  301. " downstream_found.end_list()\n",
  302. " \n",
  303. " upstream_found.field(\"brem_vtx_z\")\n",
  304. " downstream_found.field(\"brem_vtx_z\")\n",
  305. " upstream_found.begin_list()\n",
  306. " downstream_found.begin_list()\n",
  307. " for jentry in range(cut_length_found[itr]):\n",
  308. " if cut_brem_found[itr, \"brem_vtx_z\", jentry]>5000:\n",
  309. " if cut_brem_found[itr,\"brem_vtx_z\",jentry]<=9500:\n",
  310. " downstream_found.real(cut_brem_found[itr,\"brem_vtx_z\",jentry])\n",
  311. " else:\n",
  312. " continue\n",
  313. " else:\n",
  314. " upstream_found.real(cut_brem_found[itr, \"brem_vtx_z\",jentry])\n",
  315. " upstream_found.end_list()\n",
  316. " downstream_found.end_list()\n",
  317. " upstream_found.end_record()\n",
  318. " downstream_found.end_record()\n",
  319. " \n",
  320. "\n",
  321. "upstream_found = ak.Array(upstream_found)\n",
  322. "downstream_found = ak.Array(downstream_found)\n",
  323. "\n",
  324. "\n",
  325. "upstream_lost = ak.ArrayBuilder()\n",
  326. "downstream_lost = ak.ArrayBuilder()\n",
  327. "\n",
  328. "for itr in range(ak.num(cut_brem_lost, axis=0)):\n",
  329. " upstream_lost.begin_record()\n",
  330. " upstream_lost.field(\"energy\").real(cut_brem_lost[itr,\"energy\"])\n",
  331. " \n",
  332. " downstream_lost.begin_record()\n",
  333. " downstream_lost.field(\"energy\").real(cut_brem_lost[itr,\"energy\"])\n",
  334. " \n",
  335. " upstream_lost.field(\"brem_photons_pe\")\n",
  336. " downstream_lost.field(\"brem_photons_pe\")\n",
  337. " upstream_lost.begin_list()\n",
  338. " downstream_lost.begin_list()\n",
  339. " for jentry in range(cut_length_lost[itr]):\n",
  340. " if (cut_brem_lost[itr, \"brem_vtx_z\", jentry]>5000):\n",
  341. " if cut_brem_lost[itr, \"brem_vtx_z\", jentry]<=9500:\n",
  342. " downstream_lost.real(cut_brem_lost[itr,\"brem_photons_pe\",jentry])\n",
  343. " else:\n",
  344. " continue\n",
  345. " else:\n",
  346. " upstream_lost.real(cut_brem_lost[itr,\"brem_photons_pe\", jentry]) \n",
  347. " upstream_lost.end_list()\n",
  348. " downstream_lost.end_list()\n",
  349. " \n",
  350. " upstream_lost.field(\"brem_vtx_z\")\n",
  351. " downstream_lost.field(\"brem_vtx_z\")\n",
  352. " upstream_lost.begin_list()\n",
  353. " downstream_lost.begin_list()\n",
  354. " for jentry in range(cut_length_lost[itr]):\n",
  355. " if cut_brem_lost[itr, \"brem_vtx_z\", jentry]>5000:\n",
  356. " if cut_brem_lost[itr,\"brem_vtx_z\",jentry]<=9500:\n",
  357. " downstream_lost.real(cut_brem_lost[itr,\"brem_vtx_z\",jentry])\n",
  358. " else:\n",
  359. " continue\n",
  360. " else:\n",
  361. " upstream_lost.real(cut_brem_lost[itr, \"brem_vtx_z\",jentry])\n",
  362. " upstream_lost.end_list()\n",
  363. " downstream_lost.end_list()\n",
  364. " upstream_lost.end_record()\n",
  365. " downstream_lost.end_record()\n",
  366. " \n",
  367. "\n",
  368. "upstream_lost = ak.Array(upstream_lost)\n",
  369. "downstream_lost = ak.Array(downstream_lost)\n"
  370. ]
  371. },
  372. {
  373. "cell_type": "code",
  374. "execution_count": 133,
  375. "metadata": {},
  376. "outputs": [
  377. {
  378. "data": {
  379. "text/html": [
  380. "<pre>{energy: 4.62e+04,\n",
  381. " brem_photons_pe: [3.26e+03, 4.45e+03, 178, 1.45e+04, 1.1e+03, 3.79e+03],\n",
  382. " brem_vtx_z: [162, 187, 387, 487, 1.34e+03, 2.32e+03]}\n",
  383. "-------------------------------------------------------------------------\n",
  384. "type: {\n",
  385. " energy: float64,\n",
  386. " brem_photons_pe: var * float64,\n",
  387. " brem_vtx_z: var * float64\n",
  388. "}</pre>"
  389. ],
  390. "text/plain": [
  391. "<Record {energy: 4.62e+04, ...} type='{energy: float64, brem_photons_pe: va...'>"
  392. ]
  393. },
  394. "execution_count": 133,
  395. "metadata": {},
  396. "output_type": "execute_result"
  397. }
  398. ],
  399. "source": [
  400. "upstream_found[0]"
  401. ]
  402. },
  403. {
  404. "cell_type": "code",
  405. "execution_count": null,
  406. "metadata": {},
  407. "outputs": [],
  408. "source": []
  409. },
  410. {
  411. "cell_type": "code",
  412. "execution_count": 134,
  413. "metadata": {},
  414. "outputs": [
  415. {
  416. "name": "stdout",
  417. "output_type": "stream",
  418. "text": [
  419. "\n",
  420. "upstream: cutoff energy = 200MeV, sample size: 1279\n",
  421. "eff = 0.921 +/- 0.008\n"
  422. ]
  423. }
  424. ],
  425. "source": [
  426. "#plot efficiency against cutoff energy \n",
  427. "up_efficiencies = []\n",
  428. "\n",
  429. "\n",
  430. "\n",
  431. "for cutoff_energy in range(0,1000,50):\n",
  432. "\tup_nobrem_f = upstream_found[ak.all(upstream_found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  433. "\tup_nobrem_l = upstream_lost[ak.all(upstream_lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  434. "\teff = t_eff(up_nobrem_f,up_nobrem_l)\n",
  435. "\tup_efficiencies.append(eff)\n",
  436. "\n",
  437. "\n",
  438. "\t#print(\"\\ncutoff = \",str(cutoff_energy),\"MeV, sample size: \",ak.num(up_nobrem_f,axis=0)+ak.num(up_nobrem_l,axis=0))\n",
  439. "\t#print(\"eff = \",np.round(eff,4), \"+/-\", np.round(eff_err(up_nobrem_f, up_nobrem_l),4))\n",
  440. "\n",
  441. "\"\"\"\n",
  442. "we see that a cutoff energy of xxxMeV is ideal because the efficiency drops significantly for higher values\n",
  443. "\"\"\"\n",
  444. "cutoff_energy = 200.0 #MeV\n",
  445. "\n",
  446. "\"\"\"\n",
  447. "better statistics: cutoff=xxxMeV - sample size: xxx events and efficiency=xxxx\n",
  448. "\"\"\"\n",
  449. "up_nobrem_found = upstream_found[ak.all(upstream_found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  450. "up_nobrem_lost = upstream_lost[ak.all(upstream_lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  451. "\n",
  452. "print(\"\\nupstream: cutoff energy = 200MeV, sample size:\",ak.num(up_nobrem_found,axis=0)+ak.num(up_nobrem_lost,axis=0))\n",
  453. "print(\"eff = \",np.round(t_eff(up_nobrem_found, up_nobrem_lost),4), \"+/-\", np.round(eff_err(up_nobrem_found, up_nobrem_lost),3))\n"
  454. ]
  455. },
  456. {
  457. "cell_type": "code",
  458. "execution_count": 135,
  459. "metadata": {},
  460. "outputs": [
  461. {
  462. "name": "stdout",
  463. "output_type": "stream",
  464. "text": [
  465. "\n",
  466. "downstream: cutoff energy = 200MeV, sample size: 4634\n",
  467. "eff = 0.8869 +/- 0.005\n"
  468. ]
  469. }
  470. ],
  471. "source": [
  472. "down_efficiencies = []\n",
  473. "for cutoff_energy in range(0,1000,50):\n",
  474. "\tdown_nobrem_f = downstream_found[ak.all(downstream_found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  475. "\tdown_nobrem_l = downstream_lost[ak.all(downstream_lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  476. "\teff = t_eff(down_nobrem_f,down_nobrem_l)\n",
  477. "\tdown_efficiencies.append(eff)\n",
  478. "\n",
  479. "\n",
  480. "\t#print(\"\\ncutoff = \",str(cutoff_energy),\"MeV, sample size: \",ak.num(down_nobrem_f,axis=0)+ak.num(down_nobrem_l,axis=0))\n",
  481. "\t#print(\"eff = \",np.round(eff,4), \"+/-\", np.round(eff_err(down_nobrem_f, down_nobrem_l),4))\n",
  482. "\n",
  483. "\"\"\"\n",
  484. "we see that a cutoff energy of xxxMeV is ideal because the efficiency drops significantly for higher values\n",
  485. "\"\"\"\n",
  486. "cutoff_energy = 200.0 #MeV\n",
  487. "\n",
  488. "\"\"\"\n",
  489. "better statistics: cutoff=xxxMeV - sample size: xxx events and efficiency=xxxx\n",
  490. "\"\"\"\n",
  491. "down_nobrem_found = downstream_found[ak.all(downstream_found[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  492. "down_nobrem_lost = downstream_lost[ak.all(downstream_lost[\"brem_photons_pe\"]<cutoff_energy,axis=1)]\n",
  493. "\n",
  494. "print(\"\\ndownstream: cutoff energy = 200MeV, sample size:\",ak.num(down_nobrem_found,axis=0)+ak.num(down_nobrem_lost,axis=0))\n",
  495. "print(\"eff = \",np.round(t_eff(down_nobrem_found, down_nobrem_lost),4), \"+/-\", np.round(eff_err(down_nobrem_found, down_nobrem_lost),3))\n"
  496. ]
  497. },
  498. {
  499. "cell_type": "code",
  500. "execution_count": 136,
  501. "metadata": {},
  502. "outputs": [
  503. {
  504. "data": {
  505. "image/png": "iVBORw0KGgoAAAANSUhEUgAABcMAAAIhCAYAAACCO6JIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB0KElEQVR4nO3deXxU1f3/8fcQQhZIwhJJwpqAFsI3soUtAVSQRQSUtrbBhc1AxVARYykgIgQXClaKGykoiywVqqUoNlUiLgWBRmOoxSC4oGkxMRKQsJhkSO7vD36ZOk5Agpm5ycnr+Xjw+Dp3ztxzTz/JfE/ec+Zch2VZlgAAAAAAAAAAMFgDuy8AAAAAAAAAAABvIwwHAAAAAAAAABiPMBwAAAAAAAAAYDzCcAAAAAAAAACA8QjDAQAAAAAAAADGIwwHAAAAAAAAABiPMBwAAAAAAAAAYDzCcAAAAAAAAACA8QjDAQAAAAAAAADGIwwHAMNlZGRowYIFdl8GAAAAUKctWLBADofD7su4JI888oi2bt1q92UAgO0IwwHAcBkZGUpLS7P7MgAAAADYhDAcAM4hDAcAuFiWpW+//dbuywAAAABgk2+//VaWZdl9GQDgFYThAOBDEydOVHR0tMfx73/l0uFw6Ne//rVWrFihn/zkJwoICFCXLl20adMmt9edOXNGv/nNbxQTE6PAwEA1b95cvXr10vPPP+/q7+mnn3ads/Lf559/7tbPH//4R8XGxiogIEDPPfecJOnjjz/WLbfcopYtWyogIECxsbGuc1UqKSnRvffeq+7duyssLEzNmzdXQkKCXnrpJY8xVva1Zs0aderUSUFBQerVq5f27t0ry7L06KOPKiYmRk2aNNHgwYP1ySefXPL/zgAAAMCP8be//U3du3dXQECAYmJi9Pvf/96jTUlJiebMmaOYmBg1atRIrVu31rRp0/TNN9+42sycOVNhYWEqLy93HbvrrrvkcDj06KOPuo4VFRWpQYMGevLJJyVJb731lhwOh55//nnNnTtXrVq1UmhoqIYMGaKDBw+6XUdOTo5GjRrlmre3atVKI0eO1H//+19J5+bhp0+f1nPPPef6e+Caa66RJK1du1YOh0Pbt2/X7bffrssuu0zBwcEqLS2VJG3evFkJCQlq3LixmjRpouHDhysnJ8et//fee09jx45VdHS0goKCFB0drZtvvllffPGFW7vKvt544w1NmTJFLVq0UGhoqMaPH6/Tp0+roKBAv/zlL9W0aVNFRUXpN7/5jZxOZzUrBwAX1tDuCwAAVO3ll1/Wm2++qYULF6px48Zavny5br75ZjVs2FA33XSTJCk1NVXr16/XQw89pB49euj06dPav3+/ioqKJEnz5s3T6dOn9eKLL2rPnj2uc0dFRbn+e+vWrdq5c6ceeOABRUZGqmXLlsrNzVViYqLatWunxx57TJGRkXrttdc0ffp0HT16VPPnz5cklZaW6tixY/rNb36j1q1bq6ysTK+//rp+9rOfac2aNRo/frzbmF555RXl5OTod7/7nRwOh2bNmqWRI0dqwoQJ+uyzz/TUU0/pxIkTSk1N1c9//nPt27evzu7LCAAAgLppx44duvHGG5WQkKBNmzapvLxcS5Ys0VdffeVqY1mWxowZox07dmjOnDkaOHCgPvjgA82fP1979uzRnj17FBAQoCFDhuj3v/+9srKylJCQIEl6/fXXFRQUpMzMTM2cOdPVp2VZGjJkiNu13Hffferfv7+effZZFRcXa9asWRo9erQOHDggPz8/nT59WkOHDlVMTIyefvppRUREqKCgQG+++aZOnjwpSdqzZ48GDx6sQYMGad68eZKk0NBQt35uv/12jRw5UuvXr9fp06fl7++vRx55RPfff78mTZqk+++/X2VlZXr00Uc1cOBAZWVlqUuXLpKkzz//XJ06ddLYsWPVvHlz5efnKz09Xb1791Zubq7Cw8Pd+po8ebJ+9rOfadOmTcrJydF9992ns2fP6uDBg/rZz36mX/3qV3r99de1ePFitWrVSqmpqTVYXQD1ngUA8JkJEyZY7du39zg+f/5867tvyZKsoKAgq6CgwHXs7NmzVufOna3LL7/cdSwuLs4aM2bMBfucNm2adb63e0lWWFiYdezYMbfjw4cPt9q0aWOdOHHC7fivf/1rKzAw0KP9d6/R6XRaycnJVo8ePTz6ioyMtE6dOuU6tnXrVkuS1b17d6uiosJ1fNmyZZYk64MPPrjg2AAAAICa1rdvX6tVq1bWt99+6zpWXFxsNW/e3DWvfvXVVy1J1pIlS9xeu3nzZkuStXLlSsuyLOv06dNWo0aNrIULF1qWZVn//e9/LUnWrFmzrKCgIKukpMSyLMuaMmWK1apVK9d53nzzTUuSdf3117ud/89//rMlydqzZ49lWZb13nvvWZKsrVu3XnBMjRs3tiZMmOBxfM2aNZYka/z48W7H8/LyrIYNG1p33XWX2/GTJ09akZGR1i9/+cvz9nX27Fnr1KlTVuPGja3HH3/co6/vn3PMmDGWJGvp0qVux7t372717NnzguMCgOpimxQAqKWuvfZaRUREuB77+fkpKSlJn3zyiesrj3369NHf//53zZ49W2+99dYl7fc9ePBgNWvWzPW4pKREO3bs0E9/+lMFBwfr7Nmzrn/XX3+9SkpKtHfvXlf7F154Qf3791eTJk3UsGFD+fv7a9WqVTpw4IBHX4MGDVLjxo1dj2NjYyVJI0aMcFsBXnn8+1+tBAAAALzp9OnTevfdd/Wzn/1MgYGBruMhISEaPXq06/Ebb7wh6dy2hN/1i1/8Qo0bN9aOHTskScHBwUpISNDrr78uScrMzFTTpk01c+ZMlZWVadeuXZLOrRb//qpwSbrhhhvcHnft2lXS/+bJl19+uZo1a6ZZs2bpj3/8o3Jzcy9p3D//+c/dHr/22ms6e/asxo8f7/b3QGBgoK6++mq99dZbrranTp3SrFmzdPnll6thw4Zq2LChmjRpotOnT1f5N8GoUaPcHlfO/UeOHOlxnL8HANQ0wnAAqKUiIyPPe6xyG5QnnnhCs2bN0tatWzVo0CA1b95cY8aM0ccff3zR/Xx3y5TKc589e1ZPPvmk/P393f5df/31kqSjR49KkrZs2aJf/vKXat26tTZs2KA9e/bo3Xff1e23366SkhKPvpo3b+72uFGjRhc8XtU5AAAAAG85fvy4KioqLjgXl87NmRs2bKjLLrvMrY3D4VBkZKRrvi5JQ4YM0d69e3X69Gm9/vrrGjx4sFq0aKH4+Hi9/vrrOnz4sA4fPlxlGN6iRQu3xwEBAZLkWgQTFhamt99+W927d9d9992n//u//1OrVq00f/78au23/f2/CSq3hOndu7fH3wSbN292/T0gSbfccoueeuopTZ48Wa+99pqysrL07rvv6rLLLqtysU51/ibg7wEANY09wwHAhwIDA103o/mu704mKxUUFJz3WOWkuHHjxkpLS1NaWpq++uor1yrx0aNH66OPPrqoa/r+ntzNmjWTn5+fxo0bp2nTplX5mpiYGEnShg0bFBMTo82bN7udp6oxAgAAALVds2bN5HA4LjgXl87Nx8+ePauvv/7aLRC3LEsFBQXq3bu369i1116refPm6R//+Id27Njhuv/Otddeq+3bt7vm1tdee+0lXfOVV16pTZs2ybIsffDBB1q7dq0WLlyooKAgzZ49+6LO8f2/CSr3+X7xxRfVvn37877uxIkTeuWVVzR//ny3virvLQQAtQ0rwwHAh6Kjo1VYWOh2852ysjK99tprHm137Njh1q68vFybN29Wx44d1aZNG4/2ERERmjhxom6++WYdPHhQZ86ckeS5euSHBAcHa9CgQcrJyVHXrl3Vq1cvj3+VYbzD4VCjRo3cJs8FBQV66aWXLqovAAAAoDZp3Lix+vTpoy1btritSj558qS2bdvmelwZXG/YsMHt9X/5y190+vRpt2C7T58+Cg0N1bJly1RQUKChQ4dKOrdiPCcnR3/+85/VpUsXtWrV6kddu8PhULdu3fSHP/xBTZs21fvvv+96LiAgoFpbKg4fPlwNGzb
  506. "text/plain": [
  507. "<Figure size 1800x600 with 2 Axes>"
  508. ]
  509. },
  510. "metadata": {},
  511. "output_type": "display_data"
  512. }
  513. ],
  514. "source": [
  515. "#plot efficiencies wrt cutoff energy\n",
  516. "fig, ax = plt.subplots(nrows=1,ncols=2,figsize=(18,6))\n",
  517. "x_ = np.linspace(0,1000,20)\n",
  518. "ax[0].scatter(x_,up_efficiencies)\n",
  519. "ax[0].set(xlabel=\"cutoff energy [MeV]\",ylabel=r\"$\\epsilon$\",title=\"upstream\", ylim=[0.85,0.95])\n",
  520. "ax[0].set_yticks(np.arange(0.85,0.96,step=0.01),minor=False)\n",
  521. "ax[0].grid()\n",
  522. "\n",
  523. "ax[1].scatter(x_,down_efficiencies)\n",
  524. "ax[1].set(xlabel=\"cutoff energy [MeV]\",ylabel=r\"$\\epsilon$\",title=\"downstream\", ylim=[0.85,0.95])\n",
  525. "ax[1].set_yticks(np.arange(0.85,0.96,step=0.01),minor=False)\n",
  526. "ax[1].grid()\n",
  527. "\n",
  528. "plt.show()"
  529. ]
  530. },
  531. {
  532. "cell_type": "code",
  533. "execution_count": null,
  534. "metadata": {},
  535. "outputs": [],
  536. "source": []
  537. },
  538. {
  539. "cell_type": "code",
  540. "execution_count": 137,
  541. "metadata": {},
  542. "outputs": [
  543. {
  544. "name": "stdout",
  545. "output_type": "stream",
  546. "text": [
  547. "upstream eff = 0.852 +/- 0.004\n",
  548. "downstream eff = 0.84 +/- 0.005\n"
  549. ]
  550. }
  551. ],
  552. "source": [
  553. "cutoff_energy=200\n",
  554. "#possibly: instead of checking if any photons exceed the cutoff, use the sum of all photon energies to separate nobrem and brem\n",
  555. "\n",
  556. "upstream_brem_found = upstream_found[ak.any(upstream_found[\"brem_photons_pe\"]>=cutoff_energy,axis=1)]\n",
  557. "up_energy_found = ak.to_numpy(upstream_brem_found[\"energy\"])\n",
  558. "up_eph_found = ak.to_numpy(ak.sum(upstream_brem_found[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  559. "up_residual_found = up_energy_found - up_eph_found\n",
  560. "up_energyloss_found = up_eph_found/up_energy_found\n",
  561. "\n",
  562. "\n",
  563. "upstream_brem_lost = upstream_lost[ak.any(upstream_lost[\"brem_photons_pe\"]>=cutoff_energy,axis=1)]\n",
  564. "up_energy_lost = ak.to_numpy(upstream_brem_lost[\"energy\"])\n",
  565. "up_eph_lost = ak.to_numpy(ak.sum(upstream_brem_lost[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  566. "up_residual_lost = up_energy_lost - up_eph_lost\n",
  567. "up_energyloss_lost = up_eph_lost/up_energy_lost\n",
  568. "\n",
  569. "\n",
  570. "print(\"upstream eff = \", np.round(t_eff(upstream_brem_found,upstream_brem_lost),3), \"+/-\", np.round(eff_err(upstream_brem_found, upstream_brem_lost),3))\n",
  571. "\n",
  572. "\n",
  573. "downstream_brem_found = downstream_found[ak.any(downstream_found[\"brem_photons_pe\"]>=cutoff_energy,axis=1)]\n",
  574. "down_energy_found = ak.to_numpy(downstream_brem_found[\"energy\"])\n",
  575. "down_eph_found = ak.to_numpy(ak.sum(downstream_brem_found[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  576. "down_residual_found = down_energy_found - down_eph_found\n",
  577. "down_energyloss_found = down_eph_found/down_energy_found\n",
  578. "\n",
  579. "\n",
  580. "downstream_brem_lost = downstream_lost[ak.any(downstream_lost[\"brem_photons_pe\"]>=cutoff_energy,axis=1)]\n",
  581. "down_energy_lost = ak.to_numpy(downstream_brem_lost[\"energy\"])\n",
  582. "down_eph_lost = ak.to_numpy(ak.sum(downstream_brem_lost[\"brem_photons_pe\"], axis=-1, keepdims=False))\n",
  583. "down_residual_lost = down_energy_lost - down_eph_lost\n",
  584. "down_energyloss_lost = down_eph_lost/down_energy_lost\n",
  585. "\n",
  586. "\n",
  587. "print(\"downstream eff = \", np.round(t_eff(downstream_brem_found,downstream_brem_lost),3), \"+/-\", np.round(eff_err(downstream_brem_found, downstream_brem_lost),3))"
  588. ]
  589. },
  590. {
  591. "cell_type": "code",
  592. "execution_count": 138,
  593. "metadata": {},
  594. "outputs": [
  595. {
  596. "name": "stdout",
  597. "output_type": "stream",
  598. "text": [
  599. "upstream:\n",
  600. "mean energyloss relative to initial energy (found): 0.3207102540525612\n",
  601. "mean energyloss relative to initial energy (lost): 0.5602258293743071\n",
  602. "downstream:\n",
  603. "mean energyloss relative to initial energy (found): 0.17552539358035377\n",
  604. "mean energyloss relative to initial energy (lost): 0.2870828762276071\n"
  605. ]
  606. }
  607. ],
  608. "source": [
  609. "print(\"upstream:\\nmean energyloss relative to initial energy (found): \",ak.mean(up_energyloss_found))\n",
  610. "print(\"mean energyloss relative to initial energy (lost): \", ak.mean(up_energyloss_lost))\n",
  611. "\n",
  612. "print(\"downstream:\\nmean energyloss relative to initial energy (found): \",ak.mean(down_energyloss_found))\n",
  613. "print(\"mean energyloss relative to initial energy (lost): \", ak.mean(down_energyloss_lost))"
  614. ]
  615. },
  616. {
  617. "cell_type": "code",
  618. "execution_count": 139,
  619. "metadata": {},
  620. "outputs": [
  621. {
  622. "data": {
  623. "image/png": "iVBORw0KGgoAAAANSUhEUgAABbYAAAJNCAYAAADtbSO1AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB4hklEQVR4nOzdZ3hUVff38d+kJxBCDb036ShNikpHELCBCooIWAERsQAWCDbao2IBUURAkaJSROVGUAlSBQFFhVtQUKOA9BB6yn5e+M/cDkkgJ0w7J9/PdeXS2XPKmllDsmbNnn1cxhgjAAAAAAAAAABsIiTQAQAAAAAAAAAAYAWNbQAAAAAAAACArdDYBgAAAAAAAADYCo1tAAAAAAAAAICt0NgGAAAAAAAAANgKjW0AAAAAAAAAgK3Q2AYAAAAAAAAA2AqNbQAAAAAAAACArdDYBgAAAAAAAADYCo1tAACQZ9u2bdMff/wR6DAAAAAAAPkMjW0AAJBnb731lr766qtAhwEAAAAAyGdobAMAAEu+/vpr3XPPPTpw4IB7bPv27brpppv0999/BzAyAAAAAEB+QWMbAIAAWrp0qVwul8dPoUKF1KhRI33wwQd+Ofe7777rMX706FF17txZERERev3117Ps16hRI5UsWVL169fXl19+qddff10dOnTQtddeqxIlSvg0ZjtJTEzMktvMnw0bNmS7z7Zt2zRgwABVrVpV0dHRio6OVvXq1XXffffp22+/tRzDjTfeqOjoaB07dizHbW6//XaFh4cHxYcSM2fOlMvl0m+//ebX8yYkJMjlcunQoUNeO+a6deuUkJBwwefe7tLT0xUfH6+XX345x2188dzagZPzv3HjRnXq1EmxsbEqWLCg2rRpo7Vr12bZzsrvwBMnTmjo0KEqU6aMoqKi1LBhQ82bNy/b81vZFgAAOFtYoAMAACA/27JliyTp448/Vnx8vIwx+uOPPzR69Gj16tVLl112merXr+/Tczdq1Mg9tm3bNt144406efKkvvrqK7Vq1SrLfgUKFNAzzzyjjIwMjR07ViEhIXr//fd12223+SROu3vhhRfUpk0bj7G6detm2e7NN9/U4MGDVbNmTT300EOqU6eOXC6XduzYoblz56pJkyb65ZdfVLVq1Vyfe8CAAVq8eLHmzJmjgQMHZrk/OTlZixYtUteuXVWyZEnrDw45WrduncaMGaO77rpLhQsXDnQ4PvH111/r4MGDuummmwIdStBxav43bdqkq6++Wk2bNtV7770nY4wmTJigdu3aaeXKlWrevHmWfXLzO/Cmm27Spk2bNG7cONWoUUNz5sxRr169lJGRod69e+d5WwAA4Gw0tgEACKAtW7YoLi5O3bt3d481b95caWlpuuOOO7R161afNrajo6N12WWXSZLmzZunAQMGqH79+lqwYIHKlCmT7X4//PCD+vTpo1KlSqlLly4qU6aMnn/+eb333nuaNWuWihcv7pN4/enIkSPKyMjwymOpXr26rrzyygtus3btWg0cOFDXXXedPvroI0VERLjva9u2rQYNGqQPP/xQ0dHRls7duXNnlSlTRu+88062je25c+fq9OnTGjBggKXjApL00UcfqXHjxqpYsaJPjn/q1CnFxMT45NjI2b59+1SgQAEVKlQoy31PP/20ChcurGXLlrlz0759e1WpUkWPPvpotjO3L/Y7cOnSpVqxYoW7QS1Jbdq00e+//67HHntMt956q0JDQy1vCwAAnI+lSAAACKDNmzerYcOGWcb//PNPSVKtWrUsH/Pll1/W4sWLLZ370UcfVa9evXT77bdr1apVOTa1JalIkSJ67rnntGzZMlWuXFktW7bUd999p1tuuUVxcXE57rdmzRp17NhRcXFxKlKkiK677jrt2rUrz9vlVrdu3dS4cWNNmzZNDRo0UHR0tMqXL6/Ro0crIyMj2322bdum0qVLq3Pnznr33XeVkpKS5/PnxgsvvKDQ0FC9+eabHk3tf+vZs2eWvOzatUu9e/dWfHy8IiMjVatWLU2ePNl9f2hoqPr27avNmzfrhx9+yHLMGTNmuB9nXqxZs0bt2rVTbGysYmJi1KJFC3322Wce22QuRfHTTz+pV69eiouLU8mSJdW/f38lJydf8PirV6+Wy+XS3Llzs9z37rvvyuVyadOmTTnun3nurVu36qabblKhQoUUFxenO+64QwcPHsyy/d9//52rGC/2uBMSEvTYY49JkipXruxefiExMdEnz93Bgwd17733qnz58oqMjFSJEiXUsmVLffHFFzk+Nz/99JNcLpc+/PBD99jmzZvlcrlUp04dj227d+/u8c0OSTLGaNGiRbr55ptzPMe/JSUlXTAHmY91y5Yt6tGjh4oUKeLx7YSLvdbPP862bdvUs2dPxcXFqWjRoho2bJjS0tL0888/69prr1VsbKwqVaqkCRMm5Cp+SVq8eLFcLpe+/PLLLPe98cYb7vNeKP9nzpzR5ZdfrmrVqnnkcP/+/SpVqpRat26t9PT0XMWT0zIfeVnO5+jRo5o+fbrat2+vcuXKaffu3dlut3btWrVu3drjA4fY2FhdffXVWrdunfbt22fpvJK0aNEiFSxYUD179vQY79evn/bu3atvvvnG8rbeeB348rUEAAC8xAAAgIA4dOiQkWSGDBliUlNTTWpqqvn777/Nu+++a2JjY83dd9+dp+P27t3bhIeHm0WLFl303Lfeeqtp27atiYyMNNOmTbN8rkGDBpkZM2ZcdLvRo0ebkJAQ079/f/PZZ5+Zjz76yNSrV8+UL1/epKSkWN7OitKlS5sCBQqYWrVqmffee88sX77c3HbbbUZSjo/59OnT5v333zfdunUzERERJioqyvTo0cMsWLDAnDlzJlfnXblypZFk4uPjTWhoqImNjTUdO3Y0q1ev9tguLS3NREdHm+bNm1t6XD/99JOJi4sz9erVM++++65Zvny5eeSRR0xISIhJSEhwb7dr1y7jcrnM0KFDs+wvyYwYMcLSeTMlJiaa8PBw06hRIzN//nyzePFi07FjR+Nyucy8efPc240ePdpIMjVr1jSjRo0yK1asMC+99JKJjIw0/fr18zjmjBkzjCSzZ88e99jll19uWrZsmeX8TZo0MU2aNLlgjJnnrlixonnsscfM559/bl566SVToEABc/nll5tz585ZjjE3jzspKck8+OCDRpJZuHChWb9+vVm/fr1JTk72yXPXqVMnU6JECfPWW2+ZxMREs3jxYjNq1CiPY2WndOnS5t5773XfHjdunImOjjaSzF9//WWMMSY1NdUUKlTIPP744x77rlmzxkgyO3fu9GoOKlasaIYPH25WrFhhFi9ebIzJ/Wv9/Ofs2WefNStWrDCPP/64kWQGDx5sLrvsMvPqq6+aFStWmH79+hlJZsGCBRd8DJlSU1NNfHy8uf3227Pc17RpU3PFFVcYYy6e/507d5rY2Fhz0003GWOMSU9PN23btjXx8fFm7969uYrFGOM+bubPV199ZcqWLWtKlSrlPteFnDx50sybN890797dREREmOjoaHPzzTebDz/80Jw9ezbbfSIiIsydd96ZZbxXr15Gkvn888/dY7n9HXjllVdm+2/5xx9/NJLMm2++aXlbb7wOfPlaAgAA3kFjGwCAAFm+fLmRlOUnLCzMPPfcc3k+blpa2kWb2/8+d1RUlNmwYUOez3cxn3zyiZFkJkyY4DG+c+dOI8nMnj3b0nZW/Pnnn0aSqVKlijl27Jh7/Ny5c6ZUqVKma9euFz3G0aNHzTvvvGM6duxowsLCTFxcnLnrrrvM559/btLS0nLcb8uWLeahhx4yixYtMl9//bV55513TK1atUxoaKhZtmyZe7v9+/cbSea2227Lcoy0tDT3hx6pqakmIyPDfV+nTp1MuXLlsjSwBg8ebKKiosyRI0fcY9dcc40pXry4u4lojDGPPPJIrhqTObnyyitNfHy8xwcOaWlppm7duqZcuXLuWDObQ+fndeDAgSYqKsrjMWXX2M4c27p1q3ts48aNRpKZNWvWBWPMPPfDDz/sMf7+++97vKasxJjbxz1x4sQsj8X
  624. "text/plain": [
  625. "<Figure size 1800x600 with 2 Axes>"
  626. ]
  627. },
  628. "metadata": {},
  629. "output_type": "display_data"
  630. }
  631. ],
  632. "source": [
  633. "#in abhängigkeit von der energie der elektronen\n",
  634. "fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(18,6))\n",
  635. "\n",
  636. "\n",
  637. "ax[0].hist(up_energyloss_lost, bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=\"lost\")\n",
  638. "ax[0].hist(up_energyloss_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  639. "ax[0].set_xticks(np.arange(0,1.1,0.1), minor=True,)\n",
  640. "ax[0].set_yticks(np.arange(0,11,1), minor=True)\n",
  641. "ax[0].set_xlabel(r\"$E_\\gamma/E_0$\")\n",
  642. "ax[0].set_ylabel(\"counts (normed)\")\n",
  643. "ax[0].set_title(\"Upstream\")\n",
  644. "ax[0].legend()\n",
  645. "ax[0].grid()\n",
  646. "\n",
  647. "ax[1].hist(down_energyloss_lost, bins=100, density=True, alpha=0.5, histtype='bar', color=\"darkorange\", label=\"lost\")\n",
  648. "ax[1].hist(down_energyloss_found, bins=100, density=True, alpha=0.5, histtype='bar', color=\"blue\", label=\"found\")\n",
  649. "ax[1].set_xticks(np.arange(0,1.1,0.1), minor=True,)\n",
  650. "ax[1].set_yticks(np.arange(0,11,1), minor=True)\n",
  651. "ax[1].set_xlabel(r\"$E_\\gamma/E_0$\")\n",
  652. "ax[1].set_ylabel(\"counts (normed)\")\n",
  653. "ax[1].set_title(\"Downstream\")\n",
  654. "ax[1].legend()\n",
  655. "ax[1].grid()\n",
  656. "\n",
  657. "\"\"\"\n",
  658. "most electrons lose little energy relative to their initial energy downstream\n",
  659. "\"\"\"\n",
  660. "fig.suptitle(r\"$B\\rightarrow K^\\ast ee$, $p>5$GeV, only photons w/ brem_vtx_z$<9500$mm\")\n",
  661. "plt.show()"
  662. ]
  663. },
  664. {
  665. "cell_type": "code",
  666. "execution_count": 140,
  667. "metadata": {},
  668. "outputs": [
  669. {
  670. "data": {
  671. "image/png": "iVBORw0KGgoAAAANSUhEUgAABj8AAAJOCAYAAADoCxXRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADUXUlEQVR4nOzde5zN1eL/8fc2d2NmzJAZ466QMBQll5P7LZciFDWhkpPKETqlGzq5ROQcktJFudcJJ3QmRCR3kiN9qZNCDMqYQRgzs35/+M0+trnszzaf2fbMfj0fj/0on732un8++/OZtddaDmOMEQAAAAAAAAAAQDFR4lpnAAAAAAAAAAAAwE4MfgAAAAAAAAAAgGKFwQ8AAAAAAAAAAFCsMPgBAAAAAAAAAACKFQY/AAAAAAAAAABAscLgBwAAAAAAAAAAKFYY/AAAAAAAAAAAAMUKgx8AAAAAAAAAAKBYYfADAAAAAAAAAAAUKwx+AAAA+Ljdu3fr4MGD1zobAAAAAAAUGQx+AAAA+Li3335ba9asudbZAAAAAACgyGDwAwAAwAetX79eAwcO1PHjx53H9u7dqx49eujYsWPXMGcAAAAAAPg+Bj8AAECx99lnn8nhcLi8IiMj1bBhQ3300UdeSfvDDz90OZ6SkqJOnTopODhY06dPz/G5hg0bKjY2VgkJCfriiy80ffp0tWvXTh07dtR1111XqHkuSr788sscbZv92rx5c66f2b17tx5++GFdf/31CgsLU1hYmGrUqKFBgwZp+/btHuehe/fuCgsL06lTp/IMc//99ysoKMjjgavRo0fL4XDot99+y/X9unXrqmXLlh7F6YmNGzdq9OjR+ZbNH7hrB0/5Q71mZmaqXLlyev311/MMY3e9FhXFuf23bt2qDh06KCIiQqVKlVKrVq309ddf5wjnybX7zJkzGjp0qOLj4xUaGqoGDRpo4cKFuabvSVgAAFD8BV7rDAAAABS2nTt3SpL+9a9/qVy5cjLG6ODBgxo1apT69OmjG2+8UQkJCYWadsOGDZ3Hdu/ere7du+vs2bNas2aNmjdvnuNz4eHhevnll5WVlaXx48erRIkSmjdvnu67775CyWdRN27cOLVq1crlWN26dXOEe+utt/TEE0+oVq1a+stf/qI6derI4XDo+++/14IFC3Trrbfqxx9/1PXXX2857YcfflhLly7V/PnzNXjw4Bzvp6amasmSJerSpYtiY2M9L9w1tHHjRo0ZM0b9+/dX6dKlr3V2ig1/qNf169frxIkT6tGjx7XOis8pru2/bds23XHHHbrttts0Z84cGWM0ceJEtWnTRmvXrlWTJk1yfMbKtbtHjx7atm2bJkyYoJo1a2r+/Pnq06ePsrKy1Ldv36sOCwAAij8GPwAAQLG3c+dORUVFqVu3bs5jTZo0UUZGhh544AF98803hTr4ERYWphtvvFGStHDhQj388MNKSEjQJ598ovj4+Fw/95///EeJiYmKi4vTnXfeqfj4eI0dO1Zz5szRBx98oLJlyxZKfr3p5MmTysrKsqUsNWrU0O23355vmK+//lqDBw9W586d9c9//lPBwcHO91q3bq3HH39cH3/8scLCwjxKu1OnToqPj9d7772X6+DHggULdO7cOT388MMexVvU/PHHHypZsuS1zgZ8xD//+U81atRIVapUKbQ06HPed/ToUYWHhysyMjLHey+++KJKly6tpKQkZ7u0bdtW1atX14gRI3KdAeLu2v3ZZ59p1apVzkEMSWrVqpV++eUXPf3007r33nsVEBDgcVgAAOAfWPYKAAAUezt27FCDBg1yHD98+LAkqXbt2h7H+frrr2vp0qUepT1ixAj16dNH999/v9atW5fnwIckRUdH65VXXlFSUpKqVaumZs2aadeuXerdu7eioqLy/NyGDRvUvn17RUVFKTo6Wp07d9YPP/xw1eGs6tq1qxo1aqRZs2apfv36CgsLU6VKlTRq1ChlZWXl+pndu3erfPny6tSpkz788EOdPn36qtO3Yty4cQoICNBbb73lMvBxuV69euVolx9++EF9+/ZVuXLlFBISotq1a+uNN95wvh8QEKB+/fppx44d+s9//pMjzvfff99ZzsKWvYzQN998ox49eigyMlJRUVF64IEHdOLECZewJ06c0KOPPqpKlSopJCRE1113nZo1a6bVq1c743r66aclSdWqVXMuSfPll18609m5c6d69uyp6Ohol9ky7upMkn788UcNGDBANWrUUMmSJVWhQgV17do1Rx1mp7V792716tVLUVFRiomJ0bBhw5SRkaF9+/apY8eOioiIUNWqVTVx4kTb6yrbsWPH1KdPH0VFRSk2NlYPPfSQUlNTXcJs2LBBbdq0UUREhEqWLKmmTZtqxYoVLunmVa9W47g8/999953bPLlr69x89913cjgc+vjjj53HduzYIYfDoTp16riE7datm8vsNmOMlixZonvuuSfP+C936NAht21gR58rrL4kSUuXLpXD4dAXX3yR470333zTmW5+7X/+/HndfPPNuuGGG1zaMDk5WXFxcWrZsqUyMzMt5SevJaUcDod+/vlny+WSLi3T+O6776pt27aqWLGifvrpp1zDff3112rZsqXLgFRERITuuOMObdy4UUePHvUoXUlasmSJSpUqpV69erkcHzBggI4cOaItW7ZcVdiC9oXC7EsAAMA+DH4AAIBi7ffff9fBgwdVv359ZWRkKCMjQ8ePH9ecOXM0duxYPfLII7rttts8jnf79u3q3bt3vgMg2WlXrlxZ7du31/Tp0zVr1iy9/fbbef7xPVvFihXVpUsXl2PZf2QPCgrK9TOjR49WixYtVKlSJS1YsEDvvPOODh06pDZt2ujMmTMeh/PEjh079H//9396/fXX9fTTT+vTTz9V8+bN9fLLL+u9997L9TO33367PvjgAwUFBWngwIEqV66cevXqpcWLF+vChQsepf/4448rMDBQkZGR6tChgzZs2ODyfmZmptauXatGjRqpfPnyluPdu3evbr31Vu3Zs0eTJ0/W8uXL1blzZw0ZMkRjxoxxhnvooYfkcDhylHXv3r3aunWr+vXr59VfHHfv3l033HCD/vnPf2r06NFaunSpOnTooIsXLzrDJCYmaunSpXrppZe0cuVKvfPOO2rbtq1+//13SdIjjzyiJ598UpK0ePFibdq0SZs2bdItt9zijKNHjx664YYb9PHHH2vmzJnOMlupsyNHjqhMmTKaMGGCkpKS9MYbbygwMFCNGzfWvn37cpSpd+/eql+/vj755BMNHDhQr7/+up566indfffd6ty5s5YsWaLWrVvrmWee0eLFi22tq2z33HOPatasqU8++UTPPvus5s+fr6eeesr5/rp169S6dWulpqbq3Xff1YIFCxQREaGuXbtq0aJFlurVShye5Ely39a5qVOnjsqXL+8yQLJ69WqFhYVp7969OnLkiCQpIyND69atU9u2bZ3hsv/QbXXww5M2KEify1YYfalLly4qV66c3n///RzvzZ49W7fccosSEhLybf/Q0FB99NFHOn78uB566CFJUlZWlu6//34ZY7RgwQLL15HseLNfa9asUYUKFRQXF6eYmBi3n//jjz+0aNEi3XXXXYqLi9OTTz6p0qVLa9GiRbrpppty/Ux6erpCQkJyHM8+ltvgsLtr9549e1S7dm0FBrouWpE9W3PPnj1XFTZbQftCYV2XAACATQwAAEAxtnLlSiMpxyswMNC88sorVx1vRkaG6du3rwkKCjJLlixxm3ZoaKjZvHnzVafnzrJly4wkM3HiRJfj+/fvN5LM3LlzPQrnicOHDxtJpnr16ubUqVPO4+np6SYuLs506dLFbRwpKSnmvffeM+3btzeBgYEmKirK9O/f33z++ecmIyMjz8/t3LnT/OUvfzFLliwx69evN++9956pXbu2CQgIMElJSc5wycnJRpK57777csSRkZFhLl686HxlZWU53+vQoYOpWLGiSU1NdfnME088YUJDQ83Jkyedx1q0aGHKli1r0tPTnceGDx9uJJn9+/e7rYPcjBo1ykgyJ06cyPX9OnX
  672. "text/plain": [
  673. "<Figure size 2000x600 with 3 Axes>"
  674. ]
  675. },
  676. "metadata": {},
  677. "output_type": "display_data"
  678. }
  679. ],
  680. "source": [
  681. "#energyloss in abh von der energie der elektronen\n",
  682. "#upstream\n",
  683. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20,6))\n",
  684. "\n",
  685. "a0=ax0.hist2d(up_energyloss_found, up_energy_found, bins=(np.linspace(0,1,80), np.linspace(0,1.5e5,80)), cmap=plt.cm.jet, cmin=1, vmax=15)\n",
  686. "ax0.set_ylim(0,1.5e5)\n",
  687. "ax0.set_xlim(0,1)\n",
  688. "ax0.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  689. "ax0.set_ylabel(r\"$E_0$\")\n",
  690. "ax0.set_title(\"found energyloss wrt electron energy\")\n",
  691. "\n",
  692. "a1=ax1.hist2d(up_energyloss_lost, up_energy_lost, bins=(np.linspace(0,1,50), np.linspace(0,1.5e5,50)), cmap=plt.cm.jet, cmin=1, vmax=15)\n",
  693. "ax1.set_ylim(0,1.5e5)\n",
  694. "ax1.set_xlim(0,1)\n",
  695. "ax1.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  696. "ax1.set_ylabel(r\"$E_0$\")\n",
  697. "ax1.set_title(\"lost energyloss wrt electron energy\")\n",
  698. "\n",
  699. "fig.colorbar(a1[3],ax=ax1)\n",
  700. "fig.suptitle(r\"$B\\rightarrow K^\\ast ee$, $p>5$GeV, Upstream photons w/ brem_vtx_z$<9500$mm\")\n",
  701. "\n",
  702. "\n",
  703. "plt.show()"
  704. ]
  705. },
  706. {
  707. "cell_type": "code",
  708. "execution_count": 141,
  709. "metadata": {},
  710. "outputs": [
  711. {
  712. "data": {
  713. "image/png": "iVBORw0KGgoAAAANSUhEUgAABj8AAAJOCAYAAADoCxXRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADHH0lEQVR4nOzdeXwURf7/8feQm5AMCUdCkFO5BAIKyrlyX3IoCChoBFRkRUUEXMWDw5VDkGMXRJRVUG5dgRVxIyAIIoccIov4BV1RQAighHAIhCT1+8NfZhlyTE/ohEnm9Xw85qF0f7qquqtmpjs1VeUwxhgBAAAAAAAAAAAUEcWudwEAAAAAAAAAAADsROcHAAAAAAAAAAAoUuj8AAAAAAAAAAAARQqdHwAAAAAAAAAAoEih8wMAAAAAAAAAABQpdH4AAAAAAAAAAIAihc4PAAAAAAAAAABQpND5AQAAAAAAAAAAihQ6PwAAAAAAAAAAQJFC5wcAAEAB27Nnjw4dOnS9iwEAAAAAQJFF5wcAAEABe+utt7Ru3brrXQwAAAAAAIosOj8AAAAKwMaNGzVw4ECdOHHCtW3fvn3q0aOHjh8/fh1LBgAAAABA0UPnBwAAKHQ++eQTORwOt1dkZKQaNGig999/v0Dyfu+999y2Jycnq1OnTgoODtbMmTOzHNegQQPFxMQoPj5en332mWbOnKl27dqpY8eOKlOmTL6WuTD5/PPPs9Rt5mvr1q3ZHrNnzx49/PDDuvHGGxUWFqawsDBVq1ZNgwYN0o4dO7wuQ/fu3RUWFqbTp0/nGHP//fcrKCjI646refPmuZ1TaGioYmNj1apVK02YMMGtc6ww2rx5s8aMGZPrtfMHY8aMkcPh0K+//mpbmv5wbdPT01W2bFlNmzYtx5j8uLaFQVGu/6+++kodOnRQRESESpQooVatWunLL790i/H2u+HcuXMaOnSo4uLiFBoaqvr162vJkiV5jgMAAIVT4PUuAAAAgLd27dolSfrXv/6lsmXLyhijQ4cOafTo0erTp49q1qyp+Pj4fM27QYMGrm179uxR9+7ddf78ea1bt07NmzfPclx4eLhefvllZWRkaMKECSpWrJgWLlyo++67L1/KWdiNHz9erVq1cttWp06dLHFvvvmmnnjiCdWoUUNPPfWUateuLYfDoe+++06LFy/Wbbfdph9++EE33nij5bwffvhhrVixQosWLdLgwYOz7E9JSdHy5cvVpUsXxcTEeH9ykubOnauaNWvq8uXLOnHihDZt2qRXX31Vr732mpYuXaq2bdvmKd3rbfPmzRo7dqz69++vkiVLXu/iFCn+cG03btyokydPqkePHte7KD6nqNb/9u3bdccdd+j222/X/PnzZYzRpEmT1KZNG61fv15NmjRxi7f63dCjRw9t375dEydOVPXq1bVo0SL16dNHGRkZ6tu3r9dxAACgcKLzAwAAFDq7du2S0+lUt27dXNuaNGmitLQ0PfDAA/r666/ztfMjLCxMNWvWlCQtWbJEDz/8sOLj4/Xhhx8qLi4u2+P+85//KCEhQbGxsbrzzjsVFxencePGaf78+Xr33XdVunTpfClvQTp16pQyMjJsOZdq1aqpcePGucZ8+eWXGjx4sDp37qx//vOfCg4Odu1r3bq1Hn/8cX3wwQcKCwvzKu9OnTopLi5O77zzTradH4sXL9aFCxf08MMPe5XulerUqaOGDRu6/n3PPffo6aefVvPmzdWjRw99//33ee5YKUx+//13FS9e/HoXAz7in//8pxo2bKhKlSrlS/q0t4J37NgxhYeHKzIyMtv9L730kkqWLKnExERX3bRt21ZVq1bViBEjsowAsfLd8Mknn2jNmjWujgxJatWqlX7++Wc988wzuvfeexUQEGA5DgAAFF5MewUAAAqdnTt3qn79+lm2HzlyRJJUq1Ytr9OcNm2aVqxY4VXeI0aMUJ8+fXT//fdrw4YNOXZ8SFJUVJReeeUVJSYmqkqVKmrWrJl2796t3r17y+l05njcpk2b1L59ezmdTkVFRalz5876/vvv8xxnVdeuXdWwYUPNmTNH9erVU1hYmCpUqKDRo0crIyMj22P27NmjcuXKqVOnTnrvvfd09uzZPOdvxfjx4xUQEKA333zTrePjSr169cpSL99//7369u2rsmXLKiQkRLVq1dLrr7/u2h8QEKB+/fpp586d+s9//pMlzblz57rO004VK1bUlClTdPbsWb355ptu+zZt2qQ2bdooIiJCxYsXV9OmTbVq1SrX/m+//VYOh0MffPCBa9vOnTvlcDhUu3Ztt7S6devmNnIpcxqhb7/9Vn369JHT6VRMTIweeughpaSkuOJOnjypRx99VBUqVFBISIjKlCmjZs2aae3ata50nnnmGUlSlSpVXFPSfP7552757Nq1Sz179lRUVJRrRI6nOsn0ww8/aMCAAapWrZqKFy+u8uXLq2vXrlnqKTOvPXv2qFevXnI6nYqOjtawYcOUlpam/fv3q2PHjoqIiFDlypU1adIkj/WTmebXX3+tHj16KDIyUk6nUw888IBOnjyZ7THHjx/P9Zpm8lS/nq6tp+OvPodrrevsXEsblCRjjJYvX6577rknxzyudPjw4VzrIbf2Jllrc/nRjq60YsUKORwOffbZZ1n2vfHGG668c6r/xMRE3XLLLbrpppvc6i8pKUmxsbFq2bKl0tPTLZUlpymlHA6HfvrpJ6/OKzk5WW+//bbatm2rG264QT/++GOOsV9++aVatmzp1ikVERGhO+64Q5s3b9axY8e8yluSli9frhIlSqhXr15u2wcMGKCjR49q27ZtXsVJ194W8rstAQCA7NH5AQAACpXffvtNhw4dUr169ZSWlqa0tDSdOHFC8+fP17hx4/TII4/o9ttv9zrdHTt2qHfv3rl2gGTmXbFiRbVv314zZ87UnDlz9NZbb+X4x/dMN9xwg7p06eK2LfOP7EFBQdkeM2bMGLVo0UIVKlTQ4sWL9Y9//EOHDx9WmzZtdO7cOa/jvLFz50793//9n6ZNm6ZnnnlGH330kZo3b66XX35Z77zzTrbHNG7cWO+++66CgoI0cOBAlS1bVr169dKyZct06dIlr/J//PHHFRgYqMjISHXo0EGbNm1y25+enq7169erYcOGKleunOV09+3bp9tuu0179+7VlClT9PHHH6tz584aMmSIxo4d64p76KGH5HA4spzrvn379NVXX6lfv3758ovgO++8UwEBAdq4caNr24YNG9S6dWulpKTo7bff1uLFixUREaGuXbtq6dKlkqTatWurXLlybn+cXrt2rcLCwrRv3z4dPXpUkpSWlqYNGzZkO63WPffco+rVq+vDDz/Uc889p0WLFunpp5927U9ISNCKFSs0atQorV69Wv/4xz/Utm1b/fbbb5KkRx55RE8++aQkadmyZdqyZYu2bNmiW2+91S2fHj166KabbtIHH3yg2bNnW64TSTp69KhKlSqliRMnKjExUa+//roCAwPVqFEj7d+/P8s59e7dW/Xq1dOHH36ogQMHatq0aXr66ad19913q3Pnzlq+fLlat26tZ599VsuWLbNUR927d9dNN92kf/7znxozZoxWrFihDh066PLly15fU8la/eZ2ba0c7225PNV1dq61DWb+odtq54fVeri6vUnWPwcy5Uc7kqQuXbqobNmymjt3bpZ98+bN06233qr4+Pgc679p06Z6//33deLECT300EOSpIyMDN1///0yxmjx4sWWP6cy08x8rVu3TuXLl1dsbKyio6M9Hv/7779r6dKluuuuuxQbG6snn3xSJUuW1NKlS3XzzTfneFxqaqpCQkKybM/cdnXHpqfvBknau3evatWqpcBA94kuMkeE7t2716u4K11rW8ivtgQAAHJgAAAACpHVq1cbSVlegYGB5pVXXslzumlpaaZv374mKCjILF++3GPeoaGhZuvWrXnOz5OVK1caSWbSpElu2w8cOGAkmQULFngV540jR44YSaZq1arm9OnTru2pqakmNjbWdOnSxWMaycnJ5p133jHt27c3gYGBxul0mv79+5tPP/3UpKW
  714. "text/plain": [
  715. "<Figure size 2000x600 with 3 Axes>"
  716. ]
  717. },
  718. "metadata": {},
  719. "output_type": "display_data"
  720. }
  721. ],
  722. "source": [
  723. "#downstream\n",
  724. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20,6))\n",
  725. "\n",
  726. "a0=ax0.hist2d(down_energyloss_found, down_energy_found, bins=(np.linspace(0,1,80), np.linspace(0,1.5e5,80)), cmap=plt.cm.jet, cmin=1, vmax=15)\n",
  727. "ax0.set_ylim(0,1.5e5)\n",
  728. "ax0.set_xlim(0,1)\n",
  729. "ax0.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  730. "ax0.set_ylabel(r\"$E_0$\")\n",
  731. "ax0.set_title(\"found energyloss wrt electron energy\")\n",
  732. "\n",
  733. "a1=ax1.hist2d(down_energyloss_lost, down_energy_lost, bins=(np.linspace(0,1,50), np.linspace(0,1.5e5,50)), cmap=plt.cm.jet, cmin=1, vmax=15)\n",
  734. "ax1.set_ylim(0,1.5e5)\n",
  735. "ax1.set_xlim(0,1)\n",
  736. "ax1.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  737. "ax1.set_ylabel(r\"$E_0$\")\n",
  738. "ax1.set_title(\"lost energyloss wrt electron energy\")\n",
  739. "\n",
  740. "fig.colorbar(a1[3],ax=ax1)\n",
  741. "fig.suptitle(r\"$B\\rightarrow K^\\ast ee$, $p>5$GeV, Downstream photons w/ brem_vtx_z$<9500$mm\")\n",
  742. "\n",
  743. "\n",
  744. "plt.show()"
  745. ]
  746. },
  747. {
  748. "cell_type": "code",
  749. "execution_count": 142,
  750. "metadata": {},
  751. "outputs": [
  752. {
  753. "data": {
  754. "image/png": "iVBORw0KGgoAAAANSUhEUgAABk4AAAJOCAYAAADxgPt3AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADUZklEQVR4nOzdd3gU1f7H8c+mh5AsoSWEDlKlKUhVA9KlWQAFjTSRKyAiYJdmARUEFGzXq4J0kQ6KEKr0LoJegSsICgFUSAAJpJzfH5D9saTthF3T3q/n2Qd25jtnzpyZzM7Zs+ccmzHGCAAAAAAAAAAAAPLK7gwAAAAAAAAAAADkFDScAAAAAAAAAAAAXEPDCQAAAAAAAAAAwDU0nAAAAAAAAAAAAFxDwwkAAAAAAAAAAMA1NJwAAAAAAAAAAABcQ8MJAAAAAAAAAADANTScAAAAAAAAAAAAXEPDCQAAAAAAAAAAwDU0nAAAAAAAAAAAAFxDwwkAAMhRevbsqalTp2Z3NuCiffv26dixY9mdDQAAAAAA3IaGEwAAAGTZv//9b61Zsya7swEAAAAAgNvQcAIAAPKUuXPn6tZbb1VgYKBsNpv27t2b3VnK0Ndffy2bzeb0CgkJUd26dfXll1/+Y/v/4osvnJafPXtWbdu2lZ+fn6ZMmeK0bsOGDerbt69Onz7tWPbjjz/qgQce0KlTpzye59xk3bp1qc5vymvr1q2p4vft26c+ffqoYsWKCgwMVGBgoCpVqqR+/fpp586dWcrD/fffr8DAQJ07dy7dmEceeUS+vr6Wz9+oUaNks9n0xx9/pLm+Ro0aatq0qaU0XbV582aNGjUqw+PKDzI7B1mRH8o2KSlJxYsX18SJE9ON8UTZ5gZ5+fxv375drVu3VnBwsAoWLKhmzZpp06ZNTjFW79sXLlzQ4MGDFRERoYCAANWpU0dz5sxJc/9WYgEAQPai4QQAAGS79u3bq1ChQipUqJBmzZql/v37O96/+eabLqdz5swZRUVFqWLFilqxYoW2bNmiypUrezDnN2/37t2SpMWLF2vLli3avHmzPvnkE128eFHdunXTvn37/pH9161b17Fs3759qlevnvbs2aM1a9Zo4MCBTtvUrVtXYWFhqlWrllavXq0pU6aoZcuWatOmjYoVK+bR/OZWY8aM0ZYtW5xeNWrUcIr5+OOPVbduXW3btk1PP/20li1bpuXLl2vw4ME6cOCA7rjjDv3vf/+zvO8+ffooPj5es2bNSnN9bGysFi5cqPbt2yssLCxLx5cdNm/erNGjR+fJL3ezW34o2w0bNujMmTN64IEHsjsrOU5ePf87duzQ3XffrUuXLmn69OmaPn264uPj1bx5c23ZsiVVvCv3bUl64IEHNG3aNI0cOVLffPON7rjjDnXr1i3Ne66VWAAAkL18sjsDAAAAy5Ytc/y/Z8+eatq0qXr27Gk5nYMHDyohIUGPPvqoIiMj0437+++/VaBAgaxk1e12794tu92ujh07OpY1atRIiYmJevTRR7Vnzx7VqlXLo/sPDAxU1apVJUlz5sxRnz59VKtWLc2fP18RERGptgkKCtKrr76q5ORkjR07Vl5eXpo5c6Yefvhhj+UzO/z1119KTk5W0aJFbzqtSpUqqWHDhumu37Rpk/r376927drpq6++kp+fn2PdPffcowEDBmjevHkKDAy0vO+2bdsqIiJCn332mfr3759q/ezZs3Xp0iX16dPHctq5SU76u0f2++qrr1SvXj2VLVvWI+lzvWWPkydPKigoSCEhIanWDR8+XIUKFdKKFSsc56ZFixaqUKGChg0blqrnSWb3belqr81Vq1Zp1qxZ6tatmySpWbNm+vXXX/Xss8/qoYcekre3t+VYAACQ/ehxAgAAPGLjxo1q1aqV7Ha7QkND1a5dOx06dMhj++vZs6fuvPNOSdJDDz0km82mpk2bOoZa2b17tzp37qzQ0FBVrFjRKZ/NmzdXcHCwChQooMaNG2v58uVOaaeksW/fPnXp0kV2u12FCxfWkCFDlJiYqJ9//llt2rRRcHCwypUrp7ffftvlfO/atUt16tRJtfy3336TJFWrVi0LpSFNnDhRixYtsrT/YcOGqVu3bnrkkUe0fv36NBtNJOmHH37Q7bffrt27d+vee+9V79699cYbb6hdu3aZDmnj6nXh7uunQ4cOqlevnj755BPVrl1bgYGBKl26tEaOHKnk5OQ0t9m3b59KlCihtm3b6osvvtD58+ezvP/MjBkzRt7e3vr444+dGk2u16VLl1Tn5NChQ+revbuKFy8uf39/VatWTe+//75TjLe3t3r06KFdu3bphx9+SJXu559/7jhOT0v5W9qzZ48eeOABhYSEyG6369FHH9WZM2cccWfOnNETTzyh0qVLy9/fX8WKFVOTJk0UHR3tSOfZZ5+VJJUvX94xjM66deuc9pPe370r5Xb48GH16tVLlSpVUoECBVSyZEl16NAhVRl66v7galmlOHXqlLp16ya73a6wsDD17t1bsbGxqeIyu+dlVrau3DOvz/+BAwcyzFdm5zotBw4ckM1m07x58xzLdu3aJZvNpltvvdUptmPHjk496iTJGKOFCxfqwQcfTHcf1zt+/HiG58Ad15snP2ckadGiRbLZbFq9enWqdR9++KFj3+md/xUrVui2227TLbfc4nT+YmJiFB4erqZNmyopKcnl/KQ3FJbNZtPRo0ctHdvZs2f16aefqkWLFipVqpR++eWXNOM2bdqkpk2bOjVoBQcH6+6779bmzZt18uRJS/uVpIULF6pgwYLq0qWL0/JevXrpxIkT2rZtm+XYm70WPH0tAQCQbxgAAAA3GzlypPHy8jK9e/c2y5cvN1999ZWpWbOmKV26tDl//rxH9nn48GHz/vvvG0lmzJgxZsuWLebAgQNm5MiRRpIpW7asef75582qVavMokWLjDHGrFu3zvj6+pq6deuauXPnmkWLFplWrVoZm81m5syZ43Q8kkyVKlXMa6+9ZlatWmWee+45I8kMHDjQVK1a1bz33ntm1apVplevXkaSmT9/fqZ5/uOPP4wkM2jQIJOQkGASEhLMqVOnzBdffGGCg4PN448/nuXy6N69u/H19TULFy7MdP8PPfSQueeee4y/v7/55JNPMk37+PHjZunSpcYYYwYMGGA+//xzk5iYaKZOnWquXLmS7nauXheeuH5KlChhgoKCTLVq1cz06dPNypUrzcMPP2wkpXvMly5dMjNnzjQdOnQwfn5+JiAgwHTu3NnMnz/fxMfHu7TftWvXGkmmePHixtvb2wQHB5tWrVqZ7777zhGTmJhoAgMDTaNGjSwd04EDB4zdbjc1a9Y0X3zxhVm5cqUZOnSo8fLyMqNGjXKKPXTokLHZbGbw4MGp0pBkXnjhBUv7TpHyt3HmzJk01996660mMjIyVXzZsmXNs88+a7799lszYcIEExQUZG677TbH9dO6dWtTrFgx8+9//9usW7fOLFq0yIwYMcLxd3n8+HHz1FNPGUlmwYIFZsuWLWbLli0mNjY21X5u/Lt3tdzWr19vhg4dar766iuzfv16s3DhQnPfffeZwMBA89///jfVMbn7/uBqWV2//xEjRphVq1aZCRMmGH9/f9OrVy+nNF2552VUtq7eM63kK7NznZ4SJUqYJ554wvH+zTffNIGBgUaS+f33340xxiQkJJiQkBDz3HPPOW27ceNGI8kcPHjQrefgZq43T11HKRISEkzx4sXNI488kmpd/fr1ze23326Myfj8Hzx40AQHB5sHHnjAGGNMUlKSueeee0zx4sXNiRMnXM6LMcaRbsprzZo1pmTJkiY8PNzxd5yRixcvmjlz5piOHTsaPz8/ExgYaB588EEzb948c/ny5TS38fPzM4899liq5d26dTOSzLfffmuMce2+naJhw4bmjjvuSLV8//79RpL5+OOPLcfe7LXg6WsJAID8goYTAADgVkuXLjWSzNtvv+20/ODBg0aSmTFjRqpt2rRpY4KCgtJ8vfHGGy7vO+XLjnnz5jmWpXyBMGLEiFTxDRs2NMWLF3f6Mj4xMdHUqFHDlCpVyiQ
  755. "text/plain": [
  756. "<Figure size 2000x600 with 3 Axes>"
  757. ]
  758. },
  759. "metadata": {},
  760. "output_type": "display_data"
  761. }
  762. ],
  763. "source": [
  764. "#plot residual energy against energyloss and try to find a good split (eg energyloss before and after the magnet)\n",
  765. "#upstream\n",
  766. "nbins=60\n",
  767. "\n",
  768. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20,6))\n",
  769. "\n",
  770. "a0=ax0.hist2d(up_energyloss_found, up_residual_found, bins=(np.linspace(0,1,nbins), np.linspace(0,1.5e5,nbins)), cmap=plt.cm.jet, cmin=1, vmax=20)\n",
  771. "ax0.set_ylim(0,1.5e5)\n",
  772. "ax0.set_xlim(0,1)\n",
  773. "ax0.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  774. "ax0.set_ylabel(r\"$E_0-E_\\gamma$\")\n",
  775. "ax0.set_title(\"found energyloss wrt residual electron energy\")\n",
  776. "\n",
  777. "a1=ax1.hist2d(up_energyloss_lost, up_residual_lost, bins=(np.linspace(0,1,nbins), np.linspace(0,1.5e5,nbins)), cmap=plt.cm.jet, cmin=1, vmax=20) \n",
  778. "ax1.set_ylim(0,1.5e5)\n",
  779. "ax1.set_xlim(0,1)\n",
  780. "ax1.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  781. "ax1.set_ylabel(r\"$E_0-E_\\gamma$\")\n",
  782. "ax1.set_title(\"lost energyloss wrt residual electron energy\")\n",
  783. "\n",
  784. "fig.colorbar(a1[3],ax=ax1)\n",
  785. "fig.suptitle(r\"$e^\\pm$ from $B\\rightarrow K^\\ast ee$, $p>5$GeV, Upstream photons w/ brem_vtx_z$<9500$mm\")\n",
  786. "\n",
  787. "\"\"\"\n",
  788. "\"\"\"\n",
  789. "plt.show()"
  790. ]
  791. },
  792. {
  793. "cell_type": "code",
  794. "execution_count": 143,
  795. "metadata": {},
  796. "outputs": [
  797. {
  798. "data": {
  799. "image/png": "iVBORw0KGgoAAAANSUhEUgAABkEAAAJOCAYAAAAAi6D6AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAADE+klEQVR4nOzdd3hUxf7H8c+mh5AsoYQQOkivCggIGpAuTRFU0EgTuQIiAnYlgAoqRbxgu14FpAgqTYp0BBGQLoJewSsICgHEEIoEUub3h7/sZUnbk+zCJnm/nmcf2LPfM2fOnJPdmZ2dGZsxxggAAAAAAAAAACCf8bnRGQAAAAAAAAAAAPAEOkEAAAAAAAAAAEC+RCcIAAAAAAAAAADIl+gEAQAAAAAAAAAA+RKdIAAAAAAAAAAAIF+iEwQAAAAAAAAAAORLdIIAAAAAAAAAAIB8iU4QAAAAAAAAAACQL9EJAgAAAAAAAAAA8iU6QQAAAAAAAAAAQL5EJwgAAPCoPn36aMaMGTc6G3DRvn37dPTo0RudDQAAAAAA3IJOEAAAADj861//0vr16290NgAAAAAAcAs6QQAAgFebP3++atWqpeDgYNlsNu3du/dGZylLK1askM1mc3qEhYWpQYMG+vTTT6/b8T/++GOn7fHx8erQoYMCAgI0bdo0p9c2bdqkAQMG6NSpU45tP/zwg7p166aTJ096PM95yVdffZXu+qY9tm3bli5+37596t+/vypXrqzg4GAFBwerSpUqGjhwoHbu3JmjPNxzzz0KDg7W2bNnM4158MEH5e/vb/n6zZgxw+mcgoKCFBkZqZYtW2r8+PFO90hetGXLFo0ePTrLsisIRo8eLZvNpj/++MMt6RWEck1JSVFERITefPPNTGPcXa55RX6+/tu3b1e7du0UGhqqwoULq2XLlvrmm2/SxVn5bLhw4YKGDRumqKgoBQUFqX79+po3b16Gx7cSCwAAMkcnCAAAcLtOnTqpSJEiKlKkiObOnatBgwY5nr/22msup3P69GnFxMSocuXKWrlypbZu3aqqVat6MOe5t3v3bknSkiVLtHXrVm3ZskUffPCBLl68qJ49e2rfvn3X5fgNGjRwbNu3b58aNmyoPXv2aP369RoyZIjTPg0aNFDJkiVVt25drVu3TtOmTVObNm3Uvn17lShRwqP5zavGjRunrVu3Oj1q167tFPP++++rQYMG+vbbb/XEE09o2bJlWr58uYYNG6YDBw6oUaNG+u9//2v52P3791diYqLmzp2b4esJCQlatGiROnXqpJIlS+bo/KZPn66tW7dqzZo1evvtt1W/fn29/vrrqlGjhtauXZujNL3Bli1bNGbMmHz5Ze2NVBDKddOmTTp9+rS6det2o7PidfLr9d+xY4fuuOMOXbp0SbNmzdKsWbOUmJioVq1aaevWrRnu48pnQ7du3TRz5kzFxsbqyy+/VKNGjdSzZ88M39OtxAIAgMz53egMAACA/GfZsmWO//fp00ctWrRQnz59LKdz8OBBJSUl6aGHHlJ0dHSmcX/99ZcKFSqUk6y63e7du2W329WlSxfHtqZNmyo5OVkPPfSQ9uzZo7p163r0+MHBwapevbokad68eerfv7/q1q2rBQsWKCoqKt0+ISEhGjt2rFJTUzV+/Hj5+Phozpw5euCBBzyWzxvhzz//VGpqqooXL57rtKpUqaImTZpk+vo333yjQYMGqWPHjvr8888VEBDgeO3OO+/U4MGD9dlnnyk4ONjysTt06KCoqCh99NFHGjRoULrXP/nkE126dEn9+/e3nHaa2rVrq2HDho7n9957r5588kk1b95c3bp106FDh3LcwZJXeNP7Cm68zz//XA0bNlT58uU9dgzuuevvxIkTCgkJUVhYWLrXXnrpJRUpUkQrV650XJfWrVurUqVKGjlyZIYjQrL7bFixYoXWrFmjuXPnqmfPnpKkli1b6tdff9VTTz2l+++/X76+vpZjAQBA1hgJAgAAXLJ582a1bdtWdrtd4eHh6tixow4dOuSx4/Xp00fNmzeXJN1///2y2Wxq0aKFY7qR3bt3q3v37goPD1flypWd8tmqVSuFhoaqUKFCuu2227R8+XKntNPS2Ldvn3r06CG73a6iRYtq+PDhSk5O1k8//aT27dsrNDRUFSpU0BtvvOFyvnft2qX69eun2/7bb79JkmrUqJGD0pDefPNNLV682NLxR44cqZ49e+rBBx/Uxo0bM+wAkaTvv/9et9xyi3bv3q277rpL/fr106uvvqqOHTtmO62Lq/eFu++fzp07q2HDhvrggw9Ur149BQcHq2zZsoqNjVVqamqG++zbt0+lSpVShw4d9PHHH+v8+fM5Pn52xo0bJ19fX73//vtOHSBX69GjR7prcujQIfXq1UsREREKDAxUjRo19PbbbzvF+Pr6qnfv3tq1a5e+//77dOlOnz7dcZ7uVK5cOU2aNEnnz5/X+++/79juyt/cgQMHZLPZ9Nlnnzm27dq1SzabTbVq1XKK7dKli2MkU9rf6oEDB9SzZ0/Z7XaVLFlS/fr1U0JCgtN+p0+f1qOPPqqyZcsqMDBQJUqUULNmzRwjV0aPHq2nnnpKklSxYkXHVDVfffVVtu8rrlyXn3/+WX379lWVKlVUqFAhlS5dWp07d053jTz1/pOW7p49e9StWzeFhYXJbrfroYce0unTpzPc5+TJk9mWa3bXN6tydTWNq/PvjmudkZzeg5JkjNGiRYt07733Zpr+1Y4dO5btNXDHPefJz7LFixfLZrNp3bp16V579913HcfN6vonJibq5ptv1k033eR0DePi4hQZGakWLVooJSXFpfxkNtWUzWbTkSNHXD4v6e/pIT/88EO1bt1aZcqU0S+//JJh3DfffKMWLVo4dUyFhobqjjvu0JYtW3TixAlLx5WkRYsWqXDhwurRo4fT9r59++r48eP69ttvcxSb23vBk/cSAABewQAAAGQjNjbW+Pj4mH79+pnly5ebzz//3NSpU8eULVvWnD9/3iPH/Pnnn83bb79tJJlx48aZrVu3mgMHDpjY2FgjyZQvX94888wzZs2aNWbx4sXGGGO++uor4+/vbxo0aGDmz59vFi9ebNq2bWtsNpuZN2+e0/lIMtWqVTMvv/yyWbNmjXn66aeNJDNkyBBTvXp1889//tOsWbPG9O3b10gyCxYsyDbPf/zxh5Fkhg4dapKSkkxSUpI5efKk+fjjj01oaKh55JFHclwevXr1Mv7+/mbRokXZHv/+++83d955pwkMDDQffPBBtmkfO3bMLF261BhjzODBg8306dNNcnKymTFjhrly5Uqm+7l6X3ji/ilVqpQJCQkxNWrUMLNmzTKrV682DzzwgJGU6TlfunTJzJkzx3Tu3NkEBASYoKAg0717d7NgwQKTmJjo0nE3bNhgJJmIiAjj6+trQkNDTdu2bc3XX3/tiElOTjbBwcGmadOmls7pwIEDxm63mzp16piPP/7YrF692owYMcL4+PiY0aNHO8UeOnTI2Gw2M2zYsHRpSDLPPvuspWOnmT59upFkduzYkeHrFy5cML6+vqZVq1bGGNf/5oz5+5o9+uijjuevvfaaCQ4ONpLM77//bowxJikpyYSFhZmnn37aGOP8tzpq1CizZs0aM3nyZBMYGGj69u3rlH67du1MiRIlzL/+9S/z1VdfmcWLF5tRo0Y58nHs2DHz+OOPG0lm4cKFZuvWrWbr1q0mISEhy/cVV6/Lxo0bzYgRI8znn39uNm7caBYtWmTuvvtuExwcbP7zn/844jz1/nP1OTz11FNm1apVZvLkySYkJMTcfPPNTn/LrparK9c3q3K1co+481pnJif3oDHGbN682UgyBw8edPs1yM0956l7Ka0cIiIizIMPPpjutVtvvdXccsstxpjsr//BgwdNaGio6datmzHGmJSUFHPnnXeaiIgIc/z4cZfyYoxxpJv2WL9+vSldurSJjIx0HCsrFy9eNPPmzTNdunQxAQEBJjg42Nx7773ms88+M5cvX85wn4CAAPPwww+n296zZ08jyaxatcqxzZXPBmOMadKkiWnUqFG6NPf
  800. "text/plain": [
  801. "<Figure size 2000x600 with 3 Axes>"
  802. ]
  803. },
  804. "metadata": {},
  805. "output_type": "display_data"
  806. }
  807. ],
  808. "source": [
  809. "#downstream\n",
  810. "fig, ((ax0, ax1)) = plt.subplots(nrows=1, ncols=2, figsize=(20,6))\n",
  811. "\n",
  812. "a0=ax0.hist2d(down_energyloss_found, down_residual_found, bins=(np.linspace(0,1,60), np.linspace(0,1.5e5,60)), cmap=plt.cm.jet, cmin=1, vmax=25)\n",
  813. "ax0.set_ylim(0,1.5e5)\n",
  814. "ax0.set_xlim(0,1)\n",
  815. "ax0.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  816. "ax0.set_ylabel(r\"$E_0-E_\\gamma$\")\n",
  817. "ax0.set_title(\"found energyloss wrt residual electron energy\")\n",
  818. "\n",
  819. "a1=ax1.hist2d(down_energyloss_lost, down_residual_lost, bins=(np.linspace(0,1,60), np.linspace(0,1.5e5,60)), cmap=plt.cm.jet, cmin=1, vmax=20) \n",
  820. "ax1.set_ylim(0,1.5e5)\n",
  821. "ax1.set_xlim(0,1)\n",
  822. "ax1.set_xlabel(r\"energyloss $E_\\gamma/E_0$\")\n",
  823. "ax1.set_ylabel(r\"$E_0-E_\\gamma$\")\n",
  824. "ax1.set_title(\"lost energyloss wrt residual electron energy\")\n",
  825. "\n",
  826. "fig.colorbar(a0[3],ax=ax1)\n",
  827. "fig.suptitle(r\"$e^\\pm$ from $B\\rightarrow K^\\ast ee$, $p>5$GeV, Downstream photons w/ brem_vtx_z$<9500$mm\")\n",
  828. "\n",
  829. "\"\"\"\n",
  830. "\"\"\"\n",
  831. "plt.show()"
  832. ]
  833. },
  834. {
  835. "cell_type": "code",
  836. "execution_count": null,
  837. "metadata": {},
  838. "outputs": [],
  839. "source": []
  840. }
  841. ],
  842. "metadata": {
  843. "kernelspec": {
  844. "display_name": "env1",
  845. "language": "python",
  846. "name": "python3"
  847. },
  848. "language_info": {
  849. "codemirror_mode": {
  850. "name": "ipython",
  851. "version": 3
  852. },
  853. "file_extension": ".py",
  854. "mimetype": "text/x-python",
  855. "name": "python",
  856. "nbconvert_exporter": "python",
  857. "pygments_lexer": "ipython3",
  858. "version": "3.11.5"
  859. }
  860. },
  861. "nbformat": 4,
  862. "nbformat_minor": 2
  863. }