{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Tally a Flux Spectrum\n", "\n", "In this example, we will demonstrate how to get the neutron flux as a function of energy (commonly called a flux spectrum). We will use a pre-built module from the `openmc.examples` package." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "import openmc.examples\n", "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First, we'll generate a pin-cell model:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "model = openmc.examples.pwr_pin_cell()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "By default, the model has no tallies." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "tags": [] }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "model.tallies" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To get the flux spectrum, we just need to create a flux tally with an energy filter. We can take advantage of numpy to get an energy filter specifying equal-lethargy bins. Let's create an energy filter with 500 energy bins." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Create equal-lethargy energies to put in filter\n", "energies = np.logspace(np.log10(1e-5), np.log10(20.0e6), 501)\n", "e_filter = openmc.EnergyFilter(energies)\n", "\n", "# Create tally with energy filter\n", "tally = openmc.Tally()\n", "tally.filters = [e_filter]\n", "tally.scores = ['flux']\n", "\n", "# Set model tallies\n", "model.tallies = [tally]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "OpenMC also has a set of predefined energy group structures that you can take advantage of. Let's see what's available:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "dict_keys(['CASMO-2', 'CASMO-4', 'CASMO-8', 'CASMO-16', 'CASMO-25', 'CASMO-40', 'VITAMIN-J-42', 'SCALE-44', 'CASMO-70', 'XMAS-172', 'VITAMIN-J-175', 'TRIPOLI-315', 'SHEM-361', 'CCFE-709', 'UKAEA-1102', 'ECCO-1968'])" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "openmc.mgxs.GROUP_STRUCTURES.keys()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's also add a flux spectrum tally using the SHEM-361 group structure." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "# Create energy filter using SHEM-361 group structure\n", "energies_shem = openmc.mgxs.GROUP_STRUCTURES['SHEM-361']\n", "shem_filter = openmc.EnergyFilter(openmc.mgxs.GROUP_STRUCTURES['SHEM-361'])\n", "\n", "tally_shem = openmc.Tally()\n", "tally_shem.filters = [shem_filter]\n", "tally_shem.scores = ['flux']\n", "\n", "model.tallies.append(tally_shem)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now let's run the model (making sure to set the number of particles/batches slightly higher than the default values)." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "model.settings.particles = 10000\n", "model.settings.batches = 50\n", "sp_path = model.run(output=False)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Great, the simulation is done. To get our results, we need to load data from the statepoint file. We can get the corresponding tallies from the statepoint file and get the mean values for each energy bin by using the `mean` attribute on the tallies." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "with openmc.StatePoint(sp_path) as sp:\n", " t = sp.tallies[tally.id]\n", " flux500_mean = t.mean.ravel()\n", " flux500_unc = t.std_dev.ravel()\n", " \n", " t_shem = sp.tallies[tally_shem.id]\n", " flux_shem_mean = t_shem.mean.ravel()\n", " flux_shem_unc = t_shem.std_dev.ravel()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can use matplotlib to plot the flux versus the energy. Note that we divide by the energy bin width so that integrating the curve makes sense. This appropriately highlights the fact that most of the spectrum is thermal." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEKCAYAAAAMzhLIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAAun0lEQVR4nO3deXyV9Zn//9cFyI4EiKbKIkHREhpEiT2IS8PYqnWjTAsIOlqkUju1/bVFZ2g7KoOt0unX2mH0a8UFtQ4o6FgUsNZhzKgskWjVmPBFViGURllLlLBevz/O4iEkJyfJOTlL3s/HIw/O/cl97vv6hCRXPvdnM3dHRESkIe1SHYCIiKQ3JQoREYlJiUJERGJSohARkZiUKEREJCYlChERialDqgNIhtzcXB84cGCqw2jQp59+Srdu3VIdRotkQx1A9Ugn2VAHyOx6vP322zvc/aS65VmZKAYOHEhZWVmqw2hQSUkJxcXFqQ6jRbKhDqB6pJNsqANkdj3M7KP6yvXoSUREYlKiEBGRmJQoREQkpqzsoxCRzHfo0CGqqqqora1NdShN0rNnT9asWZPqMGLq3Lkz/fr144QTTojrfCUKEUlLVVVV9OjRg4EDB2JmqQ4nbvv27aNHjx6pDqNB7s7OnTupqqoiPz8/rvfo0ZOIpKXa2lr69OmTUUkiE5gZffr0aVJLTYlCWsX1j5Zy/aOlqQ5DMoySRHI09euqR0+SeGVzqV7xNDtqDkSKbq09DEDFPZ9/yy3vMpruo25mUmBAq4coEo+BAwfSo0cP2rdvT4cOHSLzs3bt2sWECRPYvHkzAwcOZMGCBfTq1Qt35/bbb+e///u/6dq1K0888QTnnntuimvRckoUkhClC++j+7oXABh6sJw8YNPRIfToXP+32NCD5Qw9WM6qJa9SsawDud07kTfqeiia3IpRizTutddeIzc395iyWbNmcckllzB9+nRmzZrFrFmz+NWvfsXLL7/Mhg0bWLduHaWlpXzve9+jtLR5LenDhw/ToUN6/IpOjygkI80r3cKid7cB8ONtC+lvH7G10+lUdCxsvLUQanX0qDnAvtrD9D+wgeoVT5OnRCEZYNGiRZSUlABw4403UlxczK9+9SsWLVrExIkTMTNGjhzJnj172L59O6eccsox77/77rt5+umnOemkk+jfvz8jRozgtttuo7i4mOHDh/Pmm28yceJEhg8fzm233cbhw4c577zzeOihh+jUqVNk9Ync3FzKysq47bbbKCkpYcaMGWzYsIH169ezY8cO/umf/ombb765xfVVopAmm1e6hY+Xv8TIA8u5rZ3RtWN7BrbfQk2vAob+cBkAQxu7SNFk8oomkxe6HkvGM3JXGaUL7yMwblqyqyAZ5l9fqqDyL39L6DULTj2Ru66O/Z1qZlx66aWYGd/97neZOnUqANXV1ZFf/l/4wheorq4GYNu2bfTr1y/y/n79+rFt27ZjEsXq1at5/vnnee+99zh06BDnnnsuI0aMiHz+4MGDlJWVUVtby+DBg1m2bBlnnnkmN9xwAw899BA/+tGPYsb8/vvvs2rVKj799FPOOeccrrzySk499dQmfW3qUme2xKV04X1U3HMh5b+4kEFLxvOjQ48yst0aunZsz9BTetJtwDnBR0fNMCkwACscB0CgYiYV91xI6cL7Ehm+SLO8+eabvPPOO7z88ss8+OCDvP7668edY2ZN6hxevnw5Y8aMoXPnzvTo0YOrr776mM9PmDABgLVr15Kfn8+ZZ54JBFsu9d2/rjFjxtClSxdyc3MZPXo0b731VtyxNUQtColL93Uv0P/ABir9NABWHR3CoiOjKPz6jxiagM7owLhplEbfp3wh8waMU0e3ADT6l3+y9O3bF4CTTz6ZsWPH8tZbb3HxxReTl5cXeaS0fft2Tj755Mj5VVVVkfdXVVVFrhGveFae7dChA0ePHgU4bphr3aSViJFjalFI3Cr9NO7vez8br1zA3T3uYmv+hIT+Ig+Mm8bQn73J/j4FFNhHDFoyXi0LSZlPP/2Uffv2RV7/6U9/4ktf+hIA11xzDU8++SQATz75JGPGjImUz58/H3dn1apV9OzZ87j+iQsuuICXXnqJ2tpaampqWLx4cb33P+uss9i8eTPr168H4Pe//z1f+cpXgOBorLfffhuA559//pj3LVq0iNraWnbu3ElJSQnnnXdei78WalFIvcJzHp7+ToB5pVsYVHuYHp078Ox3zwfg1P0bKS4OJOXeeaOup3rF04zcVQYVMykF9VtIq6uurmbs2LFAcATSpEmTuPzyywGYPn0648eP57HHHuO0005jwYIFAFxxxRX84Q9/4IwzzqBr167MnTv3uOued955XHPNNQwbNoy8vDwKCwvp2bPnced17tyZuXPnMm7cuEhn9i233ALAXXfdxZQpU7jjjjuOW9J82LBhjB49mh07dnDHHXe0uH8ClCiEz0cvjRkebCLXrHiEW3e/CgTnPQyqPUyBfcT+7gWtE1Coo7t04X0EKmaGht0qUUjrGjRoEO+99169n+vTpw/Lli07rtzM+M1vftPoEh633XYbM2bM4LPPPuPiiy+OdGaHR1KFXXLJJfz5z38+7v0XXXQRH374Yb3XHjZsGE899VTM+zeVEoWw6N1tlG7aRemmXQA80/FVCuyjSH9Ej84d2N+9oNmd1c0VGDeNilCfRfXsSzTPQrLG1KlTqayspLa2lhtvvDHtJ+WlfaIws0HAz4Ge7v6tVMeTVcrmQvlz3LlzL/s6Ho5Mjht8tIqOfc9l5OQlKQ4QagaPpbJ8IQU7KzXPQrLGvHnzknLdGTNmJOW6KenMNrPHzexjM/ugTvnlZrbWzNab2XQAd9/o7lNSEWfWK38O/loOBFsNQ0/pydBTetKx79lQmB45OTBuGhuvXECln0ZeaJ6FiLSuVLUongAeACIP0sysPfAg8DWgClhtZi+6e2VKImwrvlDIzIP/AsCzk89PcTD1mxQYQOmWcaD+CpGUSEmLwt1fB3bVKf4ysD7UgjgIPAOMafXg2pDqfbVUbN9L5fbEznhNhsC4aVR0LKT/wQ2akCfSytKpj6IvsDXquAoImFkf4JfAOWb2U3e/t743m9lUYCpAXl7ecaMH0klNTU3S4ivZeoiVfznM+ad2oLh/w7tXlWw9xOU7PgXg1G5HGdK1aTElsw4N2ddzFId3Hib/wAY2rVlIScmIxt/UiFTUIxmyoR5169CzZ8/IPIZMcuTIkYyIu7a2Nv7vGXdPyQcwEPgg6vhbwKNRx/8APNCca48YMcLT2WuvvZa0a4//3Qo/7Z8X+/jfrXB39/9c9ZGP/92KyMdVs9+InLPyjoD/9d//rln3SWYdGvPBLy9wv+tEX7Xg/7T4WqmsRyJlQz3q1qGysjI1gdTxi1/8wgsKCrywsNDPPvtsX7VqlX/lK1/x1atXR87ZtGmTDx061N3dlyxZ4ieeeKKfffbZkY9XX33V3d0Bv+666yLvO3TokOfm5vqVV15Z771vuukmHzZsmBcWFvo3v/lN37dvX+Rzzz77rA8ZMsQLCgp84sSJkfLLLrvMe/bs2eA1w+r7+gJlXs/v1HRqUWwD+kcd9wuVSTPNK93Cz14IdlYH8ntHhr+Gj/MPdCOvR+dUhddsNYPHQkU5AU3GkyRbuXIlixcv5p133qFTp07s2LGDgwcPNvq+iy66qN4Z1926deODDz5g//79dOnShVdffTXmEh/3338/J554IgA/+clPeOCBB5g+fTrr1q3j3nvvZfny5fTq1YuPP/448p7bb7+dzz77jIcffrgZNa5fOi3hsRoYbGb5ZtYRuBZ4McUxZazK7X+LJIl7xhZGZlRHH2dikoDQulBD7wy+rpip/gpJmu3bt5Obm0unTp0AyM3NbfFM5yuuuIIlS4JDz+fPn8/EiRMbPDecJNyd/fv3R9ZteuSRR/j+979Pr169ACJrTUFwkl6i9+xOSYvCzOYDxUCumVUBd7n7Y2Z2K/AK0B543N0rmnWDHetg7pWJCjfhhu/ZA5tyknLtyJyIjh2gI8ENgSo7QyW8fGKw03pI5YlQSXBo7BcKkxJHsoUXEQxUzCRQMZM5VXu0W142e3l6ZCh3wnyhEL4+K+Ypl156KTNnzuTMM8/kq1/9KhMmTIist3TdddfRpUsXILg0eLt2n//d/cYbbzB8+PDI8fPPP8/pp58OwLXXXsvMmTO56qqreP/997npppt44403Goxh8uTJLF26lIKCAu67L/hHUXhW9gUXXMCRI0eYMWNGZHmRZEhJonD3elOouy8FlrZyOFknPCeiriGnnHhswRcK02a+RHNEJ4upe2ezasmrlG4Zp0dRkjDdu3fn7bff5o033uC1115jwoQJzJoVTC7/+Z//SVFREQCbN2/mqquuiryvoUdPEFxiY/PmzcyfP58rrrii0Rjmzp3LkSNH+MEPfsCzzz7L5MmTOXz4MOvWraOkpISqqiouvvhiysvLycnJaXml65FOfRSJkzsY0mBWcUPeLSk5biGvRJn58EogfedEJFpg3DTI733MIoIV616gZvBYJYxs0shf/snUvn17iouLKS4uprCwMLJqbEtcc801kV3pdu7cGSm/7LLLqK6upqioiEcfffSYGK699lr+7d/+jcmTJ9OvXz8CgQAnnHBCZM+KdevWJWSl2PpkZ6KQtiVqEcHu615g6MFyqChXR7e02Nq1a2nXrh2DBw8G4N133+W0007jgw8+aOSdsd10003k5ORQWFh4zBDVV155JfLa3dmwYQNnnHEG7s6LL77IF7/4RQC+8Y1vMH/+fCZPnsyOHTv48MMPGTRoUItiiiWdOrOlha5/tPSYkU1tTXg/i+iObk3Ok5aoqanhxhtvpKCggGHDhlFZWRnXekrhPorwx3PPPXfM5/v168cPf/jDmNdwd2688UYKCwspLCxk+/bt3Hln8Hv7sssuo0+fPhQUFDB69Gh+/etf06dPHyD42GvcuHEsW7aMfv36HZN8mkstiiwxr3QLb67fkeow0kL0bnlqXUhLjBgxghUrVhxXXnei2sCBAyOtjIsuuoi9e/fWe72amprjysKPtepq164dy5cvr/c64eXMf/Ob3xz3uVgd482lFkUWiJ4vAUT2lWjL6mtdzLn/DuaVbklxZCKZRy2KDHf9o6WRlsQ9Yws1PLQOjYwSaTkligwU3qb0isJTIkkikN9bSaIB9Y2Mqt6+NLQRU36qwxNJe0oUGSJ6u9Jwcojuk9DjpkZEjYzy8oXBhLG4jH0nfReSNFRZWs7dI7ORJXGCyzrFT30UGSK8Xemid49d/uqesYVsnnWlWhNxCm+ENKdncMTJ1Z88rFFRaapz587s3Lmzyb/UJDZ3Z+fOnXTuHP8SPmpRpLlwSyI87DV6+OuFZ+QqQTTDpMAACNxN6cKcyBIgmqSXfvr160dVVRWffPJJqkNpktra2ib9Ek6Fzp07069fv7jPV6JIQ3VHMdUnkN+bp78TaKWIslNg3DRe+riaQXtXRIbRRvoutDd3yoVnHWeakpISzjnnnFSHkVB69JSG6j5eCgvk9+aesYUE8nurTyJBegy9IjKMdtXRIeTtKoPFP9LjKJEoShQZZlJgAM9+93w9ckqwun0XWr5c5HNKFCl0/aOlkaGu80q3MOHhlZEJYYH83mye9flS6WpFJN+kwACm/vjuYybpVc++BMrmpjgykdRSH0UKlGw9xEMPrzymYzq6w7o+0RsPSXKFJ+lFD6Mt3bRLHd3SZqlFkQIr/3L4mKQwoU7SiHbhGblceEZua4UmIfU9itISINJWqUWRBqKTRPh1IL83gEY2pVDdYbRT987mpy9+yjx+pD4iaVPUohBpRGDcNLjqtwDce8JjDFoyXq0LaVOUKNJE9NBXSUNFk+Gq31Ldu4iR7dYwde9sBi0Zr5FR0ibo0VMSRa/PBBzzOlogv/dxndUa4ZSGQutFUTb3mAUGtdeFZDsliiSqO5Ipnt3nJgUG6Pl3uotaYDC8BIiShWQzPXoSaabAuGnHzLnQYyjJVmpRJFj0XhH1tSDqjmpqy3tcZ4PojZECFTPZufU5/lrbnuVdRtN91M1qHUpWSPtEYWaDgJ8DPd39W6mOpzHhPSIOHTma4kiktQTGTWNO1R6m7p1Nn79V0geCiwy+PJvy/zmHwp+WpDpEkRZJ6qMnM3vczD42sw/qlF9uZmvNbL2ZTY91DXff6O5Tkhlnc0UvwREPjWjKXt1H3cyqo0OOKy888Gcq7rlQj6UkoyW7j+IJ4PLoAjNrDzwIfB0oACaaWYGZFZrZ4jofJyc5vhZ5c/2OY3aZawqNasoukwID6NE52ECv6FgY6buAYOtCfRiSyZL66MndXzezgXWKvwysd/eNAGb2DDDG3e8FrmruvcxsKjAVIC8vj5KSkuZeqsku+9XLnH/qsV/KPXv21Ps67NT9G9mzZ3/k860ZbyLU1NRkXMz1SWQ9eh0+DMDhw4fZf9IISooXAbCvYilXf/IwgYqZvLdmIVv6XEyPoVck5J5h2fD/kQ11gOypR7RU9FH0BbZGHVcBDa5TYWZ9gF8C55jZT0MJ5TjuPgeYA1BUVOTFrbEP8h+XALB291FycrqHCoOd0zk5ObA76nVYqKy4uJiH1q6E3bvIycmhuDizFv0rKSmhVb7GSZbIelS/3w12wRdyunF29DWLi2HGwwCcfXQNZ3+yhtJP8hI6nDYb/j+yoQ6QPfWIlvad2e6+E7ilte8bPVkunpEr0aOX6vZF1B3pJNkpr0dn2BX6t65BxQCUdrlY269KxklFotgG9I867hcqSyvhyXKlm3Y1KWHEMmZ4Xw2HzWaF3zr232g3BB9DBYBSoPu6F+h/cAMnVsyEipnBRBI6RyTdpGLC3WpgsJnlm1lH4FrgxRTEEbdwsmiJe8YWHpdoxgzvqw2JsknRZJi8pNH9tgPjpjH0Z2+ytePpnxduLFFnt6StpLYozGw+UAzkmlkVcJe7P2ZmtwKvAO2Bx929IplxtLb6Wg31tUa0XEfbtrzL6OB8ixA9jpJ0lexRTxMbKF8KLE3mvdNBIL+3HjVJg5Z1vYKpe2cDUDr0TrqveyGYOCrKtXaUpBWt9SSSItGPHMOPo6LXjmJGTx777Z0NvV2k1ShRJNGY4X21lak0qL7HjoFx0yLbrwJM2fPvmtktKRfz0ZOZnRvHNQ65e3njp7Utgfze6oOQZol+JLXq6BBG6nGUpFhjfRT/S3CUksU4Jx8YmKiARNqU0PyKaGOG94XtwdcjZ646Zt8LdXZLKjSWKFa7+9/FOsHM/ieB8Yi0LfXMnZgUGAAvf34cXso8urNbCUNaU8xE0ViSiPectkJ7TEiyBBPCNJjRE+DY0VH7X2fY7t1Q/HpKY5TsFVdntpmNNbOeUcc5ZvaNpEUl0tYNKq73sVTYA91uBUKjozaW0Hv3e60Tl7RJ8Y56usvd94YP3H0PcFdSIhKR4COpGEt63Hr7LykdeicVHQsjZRpKK8kSb6Ko77y0X1BQJJuF516EHTeU9qkxwQ+RFor3l32Zmf2G4IZDAN8H3k5OSKnVlB3rmuPCM3LZvVt9GJIcx/RdbCwJFpbNbXT9KZFY4m1R/AA4CDwLPAPUEkwWWaclu9bF4+nvBLj9vC5Ju75kuUb6Lo6Z2R1W/lxyY5Ks12iLIrR16WJ3H90K8YhILI0sRR4eShudKKr31ZKX5LAkuzXaonD3I8DR6FFPIpK+AuOmHdPqyNtVpmVApEXi7aOoAcrN7FXg03Chu/+w4beISKsYVAwbS459JBVuedSZd6FNkqQ54u2j+C/gDuB1gp3Y4Q8RSbUbFlFSHHs4bbjvAtAmSdJkcbUo3P3J8Gsz6wX0d/f3kxaViCRU3b4LrRslTRHvzOwSMzvRzHoD7wCPhIbLikiGiO67qOhYyNCD5QQqZqp1IY2Kt4+ip7v/zcy+Azzl7neZmVoUIumu7lDa0OOpoaBVaSVu8SaKDmZ2CjAe+HkS4xGRRIrRb6FVaSVe8SaKmcArwJvuvtrMBgHrkheWiLSG8Kq0pQvv057d0qC4+ijcfaG7D3P3fwwdb3T3byY3NBFpLXXXjVLfhURr8p7ZZvZOMgIRkTQQ1acRqJjJnPvvYF7pltTFI2mhyYmC2Nuiikgmu2ERzNgbmXcxde9sBi0Zr9ZFGxczUZjZefUUL0lSLA3FMMTMfmdmz5nZ91rz3iJtVWDcNLjqt1T3LmJkuzXBkVFaBqTNaqxFMcfM1pnZ3WZWAODu/xLvxc3scTP72Mw+qFN+uZmtNbP1ZjY91jXcfY2730JwxNUF8d5bRFqoaDJ5P1wW2SApPO9CCaPtaWzP7HPM7CzgWuA5MzsEzAeecffNcVz/CeAB4KlwQWg12geBrwFVwGozexFoD9xb5/03ufvHZnYN8D3g9/FUSkQSRyOjxNw9/pPNziaYNMYDf3X3Rv/CN7OBBJcp/1Lo+HxghrtfFjr+KYC7100S9V1ribtf2cDnpgJTAfLy8kY888wzcdWprm//Mbjm4Vm92rF299FI+Vm92vHTQJd6z40+B2Dt7qP1nh9WU1ND9+7dmxVfusiGOoDq0Rz7KpZy9ScPA7C5cwEHBhSz/dTLWnxd/V+k3ujRo99296K65XFvZ2pm7YCTgTygG/BxM2PpC2yNOq4CAjHuWwz8PdAJWNrQee4+B5gDUFRU5MXFxc2L7o/BLpicnByI2okuJyeH4uLz6z03+hwAdu+q//yQkpISmh1fmsiGOoDq0SzFxZQuzMPLFzKythI+rGTPCbktbl3o/yJ9xbNx0UXAROAbQDnBHe5+7O57kxtakLuXACWtcS8RiU9g3DTmDRjH+yseYere2QQqZlK+9gVWdhtN91E3MykwINUhSgLFTBRmthX4iGBymOHuzW1FRNsG9I867hcqE5EMMikwAAJ3U7owJ9i6OFxO4d5y5qwAAnenOjxJoMZGPV3o7he6+wOhTuWuCbjnamCwmeWbWUeCfR4vJuC6IpICgXHT2HjlAub0DO5jNnXvbE3UyzIxE4W7fwTBDmgzqwT+X+j4bDP7v41d3MzmAyuBs8ysysymuPth4FaCa0etARa4e0UL6yEiKTQpMICpP777uIl6ShjZId7O7N8ClxH6y9/d3zOzixt7k7tPbKB8KTE6pkUkMwXGTYP83lSveJqRu8oYuXcNq5a8SumWcRpKm8HiXsLD3bfWKTqS4FhEJBuEJurVndldPfsSKJub6uikGeJNFFvNbBTgZnaCmd1G8LGRiEj9omZ2rzo6hLxdZbD4R0oYGSjeRHEL8H2CcyC2AcNDxyIiMUV3dq86OoQuOyup+NNj6rvIII0tCjjRzPq4+w53v87d89z9ZHe/3t13tlaQIpLZwp3dG69cwNZOp9P/wAatSptBGuvMHgAsNLMTgGXAy8Bb3pR1P0REQiYFBkD7KZHObqL26+akEakOTxrQ2PDYX7n73wFXAO8BNwHvmNk8M7vBzPJaI0gRySINrErb6/XpamGkqXi3Qt3n7i+4+3fd/RzgF8BJRK0KKyLSFOHtV8MJI//IR3j5QvVdpKGmLAo4DBgY9Z5N4RVgRUSaK7yM+eZZ51OwfzMsGU/Fsg7kdu9E3qjroWhyqkNs8+JKFGb2ODAMqADCa2878F9JiktE2pgDA4rZv+MtetQcYF/tYfof2EDFnx7jvSOXaJHBFIu3RTHS3QuSGomItGnbT72MsybdSx4wr3QLPZZNov+BDexbMp45K76mVWlTKN55FCvDW6GKiCTbpMAAhl46hf19ChjZbo3WjkqxeFsUTxFMFn8FDgAGuLsPS1pkItK2FU0mr2gylM3V2lEpFm+ieAz4B4IbFx1t5FwRkcSpJ2FEz79Qwki+eBPFJ+6uPSNEJHVCCaN04X10X/cCQw+WQ0W5EkYriDdR/NnM5gEvEXz0BIC7a9STiLSq8HDa+hKGhtQmR7yJogvBBHFpVJmGx4pIytRNGPtqDweTxuIySjftUgsjgeJKFO6u9CwiaSmcMOaVbuH9FY8wde9sAurDSKh4J9w9Cfx/7r4ndNwLuM/db0pibCIicZsUGACBuyldmKM+jASL99HTsHCSAHD33WZ2TnJCEhFpvob6MOZU7dGkvWaKd8Jdu1ArAgAz600T1okSEWlt0YsOApFJexX3XKiJe00U7y/7+whOuFtIcLLdt4BfJi0qEZEECYybBvm9qV7xND1qDjD0YDlDD5azasmrVCzrwPIuo9XSaES8ndlPmVkZ8Hehor9398rkhSUikkChORh5EJm4V1/SANSfUY+4Hx+FEkOlmU1VkhCRjNVA0gCO6QAPU+JoXj/DLcCcRAfSEDMrBu4muMT5M+5e0lr3FpEsF500INIBHhZOHNXbl7bpiXzNSRQW94nBfSyuAj529y9FlV8O/DvQHnjU3WfFuIwDNUBnoKoZ8YqIxCU8YiqsdOF9ePlCCnevgfLnlCia4OomnPsE8ABRW6aaWXvgQeBrBH/xrzazFwkmjXvrvP8m4A13/9/Q/ty/Aa5rRswiIk0WGDeNCbtG8eNtP2bkR29SPfsSdtQcqPfcbH5EFe+Eu07ANwlthWoWbFS4+8xY73P3181sYJ3iLwPr3X1j6NrPAGPc/V6CrY+G7AY6xYhxKjAVIC8vj5KSklihNWrPnj3HHTd2zej3xDq/pqamxfGlWjbUAVSPdJKudRjS9VDkdZedlezz0+hS5zdn/pGP2LRmISUlI9K2Hi0Rb4tiEbAXeJuoRQGbqS+wNeq4Cgg0dLKZ/T1wGZBDsHVSL3efQ6jvpKioyIuLi5sX3R+XAJCTkwO7d0WKc3JyKC4+v95zo88BYPeu+s8PKSkpodnxpYlsqAOoHukkXetQDFTPngW7oNJPY+OVC44bSltxz4Xsrz3MX7oM4lQ2pmU9WiLeRNHP3S9PaiQNCK1Qq8UHRSRl8kZdz5o/HaSi62im1DPfIrd7J4YeLOf9FY/AOZekIMLkijdRrDCzQncvT8A9twH9o477hcpERNJT0WSGFE1mSAOfzht1PSwu44L9r7Gx4gAVK/4lq5Y8jzdRXAh828w20fKtUFcDg80sn2CCuBaY1IzriIikh6LJVPzpMfof2MDQT4J/T/9tZ1dWvfgwG49ckvGzvuNNFF9vzsXNbD7BR3y5ZlYF3OXuj5nZrcArBEc6Pe7uFc25vohIuqgZPJbK8oXgUNbjEm799AFG2hpYMj7j9/iOdwmPj5pzcXef2ED5UmBpc64pIpKOAuOmcf3eC9m9exdLbv86zAiOvRnZbg1UzMzoSXsxV481s3cau0A854iItAVPfyfA7ed1CR4MKoZBxZQOvZNVR4fQPTRpb17pFiY8vDKjVq9trEUxxMzej/F5A3omMB4RkexwwyIgOPZ/wq5R3Lb9JxRs+TODNo3nx0CPnR2g/ZSMaGE0lii+GMc1jiQiEBGRbDVmeF/+/NlX6br/NXqcAPtqD9P/wAaqVzxNXqYniub2TYiIyOfC27SGzSvdAkvG06PmQGRBwnQW7w53IiKSIJMCA+jRuUOwVTH7Eiibm+qQYlKiEBFJgZrBY6n00+iys5LqFU+nOpyY4koUZlZQT1lxooMREWkrAuOmsfHKBVT6aQ2uSJsu4m1RLDCzf7agLmb2Hxy/JLiIiDRB+BFUuos3UQQIrs+0guASHH8BLkhWUCIibUm691XEmygOAfuBLgR3mtvk7keTFpWISBuRCX0V8SaK1QQTxXnARcBEM1uYtKhERNqI6L6KTTs+TcsZ2/Emiinufqe7H3L37e4+BngxmYGJiLQVkwIDyM/tRoF9xNnLJqXdI6h4E8XHZjYg+gP432QGJiLSluSNup6tnU6PzNhOJ/F2ty8BnODaTp2BfGAtMDRJcYmItC1Fk3nvyCXsS8MZ23G1KNy90N2Hhf4dDHwZWJnc0ERE2pZ0HS7brJnZ7v4OwSGzIiKSYOk2XDau1GVmP4k6bAecS3AuhYiIJFB4p7yC0HDZdFhdNt4WRY+oj04E+yzGJCsoEZG2Kh2X9oh3K9R/TXYgIiISNCkwgIpl6dNXETMSM3uJ4Ginern7NQmPSEREgOAGR/NKtwT3s0ihxlLW/2mVKERE5Bi53Tsx9GA576945JhNj1KhsUSxyd3Tbz65iEiWyxt1PSwu44L9r6U6lEY7s/8QfmFmzyc3lPqZ2UVm9jsze9TMVqQiBhGRVlc0mYqOhamOAmg8UVjU60FNvbiZPW5mH5vZB3XKLzeztWa23symx7qGu7/h7rcAi4EnmxqDiEgmC/dTpFJjicIbeB2vJ4DLowvMrD3wIPB1oIDgSrQFZlZoZovrfJwc9dZJwLxmxCAikpFyu3diZLs11Kx4JKVxNNZHcbaZ/Y1gy6JL6DWhY3f3E2O92d1fN7OBdYq/DKx3940AZvYMMMbd7wWuqu86oUUI97r7vkbiFRHJGunSTxEzUbh7+yTcsy+wNeq4isaXA5kCxJzLbmZTgakAeXl5lJSUtCBE2LNnz3HHjV0z+j2xzq+pqWlxfKmWDXUA1SOdZEMdINH1yKdXuyHsrz3MjN+/SnH/ExJ03aZJnxkdMbj7XXGcMweYA1BUVOTFxcXNu9kflwCQk5MDu3dFinNyciguPr/ec6PPAWD3rvrPDykpKaHZ8aWJbKgDqB7pJBvqAImvR/X73cjbVcb7O16n+B9SM0y2WYsCttA2gvtvh/ULlYmISB15o64HSOnjp1QkitXAYDPLN7OOwLVotzwRkfqlwTDZpCYKM5tPcN+Ks8ysysymuPth4FbgFWANsMDdK5IZh4hIpkvlMNmk9lG4+8QGypcCS5N5bxGRbJHbvRP9D2ygx7JJ0H4KtPLS46l49CQiIk0Q3k974KGNUP5cq99fiUJEJN0VTWZmn19TfmQA1ftqW/32ShQiIhlgzPC+ACnZzEiJQkQkA0wKDKBH59RMfVOiEBHJIEMPlkNZzIUqEk6JQkQkQyzvMjr4opU7tJUoREQyxLKuV7Dq6JBW79BWohARyRCp6tBWohARyRCp6tBWohARkZiUKEREJCYlChGRDNPaQ2SVKEREMkgqhsgqUYiIZJBUDJFVohARySCpGCKrRCEikkHCQ2RbcyMjJQoRkQyT270TAIve3dYq91OiEBHJMHk9OrfqxDslChERiUmJQkREYlKiEBGRmJQoREQy0MBDG7lz5+2tMkNbiUJEJNMUfovNJwxi4KGNrTJDO+0ThZkVmNkCM3vIzL6V6nhERFKuaDIz+/yazScMapXbJTVRmNnjZvaxmX1Qp/xyM1trZuvNbHojl/k68B/u/j3ghqQFKyIi9Ur2QNwngAeAp8IFZtYeeBD4GlAFrDazF4H2wL113n8T8HvgLjO7BuiT5HhFRKSOpCYKd3/dzAbWKf4ysN7dNwKY2TPAGHe/F7iqgUt9P5Rg/itpwYqIZJh9tYep3ldLXpLv0/p76kFfYGvUcRUQaOjkUKL5GdAN+HWM86YCUwHy8vIoKSlpUZB79uw57rixa0a/J9b5NTU1LY4v1bKhDqB6pJNsqAO0Xj2GdD0EwF/3fMqaJN8vFYmiSdx9M6EE0Mh5c4A5AEVFRV5cXNy8G/5xCQA5OTmwe1ekOCcnh+Li8+s9N/ocAHbvqv/8kJKSEpodX5rIhjqA6pFOsqEO0Hr1KAYq7vnX4Osk3y8Vo562Af2jjvuFykREJA2lIlGsBgabWb6ZdQSuBV5MQRwiIhmvNbZFTfbw2PnASuAsM6sysynufhi4FXgFWAMscPeKZMYhIpKNWmtb1GSPeprYQPlSYGky7y0iku2Wdb2CYbtfJT/JI5/Sfma2iIjUr7W2RVWiEBHJUNHbol7/aGnS7qNEISKSBd5cvyNp11aiEBGRmJQoREQkJiUKERGJSYlCRERiUqIQEZGYlChERCQmJQoREYlJiUJERGJK+/0oRESkYR91OJ2/Hq1N6j2UKEREMtiTPW+hdNeuxk9sAT16EhHJYOGFAQG+eMfLSbmHEoWISAabFBjAPWMLAag9dDQp91CiEBHJcJMCA5J6fSUKERGJSYlCRERiUqIQEZGYlChERCQmJQoREYlJiUJERGJSohARkZiUKEREJCZz91THkHBm9gnwUarjiCEX2JHqIFooG+oAqkc6yYY6QGbX4zR3P6luYVYminRnZmXuXpTqOFoiG+oAqkc6yYY6QPbUI5oePYmISExKFCIiEpMSRWrMSXUACZANdQDVI51kQx0ge+oRoT4KERGJSS0KERGJSYlCRERiUqIQEZGYlCjSkJl1M7MyM7sq1bE0h5l9w8weMbNnzezSVMfTFKGv/ZOh+K9LdTzNkclf/7qy4GehnZn90sz+w8xuTHU8zaVEkUBm9riZfWxmH9Qpv9zM1prZejObHsel/hlYkJwoY0tEHdz9D+5+M3ALMCGZ8cajiXX6e+C5UPzXtHqwDWhKHdLt6x+tGd9fKftZaEgT6zAG6AccAqpaO9ZEUaJIrCeAy6MLzKw98CDwdaAAmGhmBWZWaGaL63ycbGZfAyqBj1s7+JAnaGEdot76L6H3pdoTxFkngj/UW0OnHWnFGBvzBPHXISxdvv7RniD+769U/yw05Ani/784C1jh7j8BvtfKcSZMh1QHkE3c/XUzG1in+MvAenffCGBmzwBj3P1e4LjmtJkVA90IfrPtN7Ol7n40mXFHS1AdDJgFvOzu7yQ55EY1pU4E/+rrB7xLGv0h1ZQ6mNka0ujrH62J/xfdSeHPQkOaWIetwMHQOen0h0eTKFEkX18+/wsVgr+IAg2d7O4/BzCzbwM70uEHgybWAfgB8FWgp5md4e6/S2ZwzdRQnWYDD5jZlcBLqQisCRqqQyZ8/aPVWw93vxXS7mehIQ39X/w78B9mdhHweioCSwQlijTl7k+kOobmcvfZBH/hZhx3/xSYnOo4WiKTv/71yfCfhc+AKamOo6XSpmmdxbYB/aOO+4XKMkk21KGubKhTNtQBsqMe2VCHBilRJN9qYLCZ5ZtZR+Ba4MUUx9RU2VCHurKhTtlQB8iOemRDHRqkRJFAZjYfWAmcZWZVZjbF3Q8DtwKvAGuABe5ekco4Y8mGOtSVDXXKhjpAdtQjG+rQVFoUUEREYlKLQkREYlKiEBGRmJQoREQkJiUKERGJSYlCRERiUqIQEZGYlCikTTOzI2b2btRHPMvAJ11UXKfGOOcuM7u3Ttnw0KKAmNlrZlZjZkXJjleym+ZRSJtmZjXu3j3B1+wQmoDVkms0GpeZnQn80d0HRZXNAj5z95mh4xLgNncva0k80rapRSFSDzPbbGb/ambvmFm5mX0xVN4ttHHNW2b2ZzMbEyr/tpm9aGb/Aywzs65mtsDMKs3sBTMrNbMiM7vJzH4bdZ+bzez+OOK51MxWhuJZaGbd3f1DYLeZRa/kOx6Yn9AvhrR5ShTS1nWp8+gpeke4He5+LvAQcFuo7OfA/7j7l4HRwK/NrFvoc+cC33L3rwD/COx29wLgDmBE6JwFwNVmdkLoeDLweKwAzSyX4CZEXw3FUwb8JPTp+QTXFcLMRgK73H1d078MIg3TMuPS1u139+ENfO6/Qv++TXCLVIBLgWvMLJw4OgMDQq9fdfddodcXEtyLAHf/wMzeD72uCbU6rgr1JZzg7uWNxDiS4OY9y4N7QtGR4FpDAM8CK8xsGsGEodaEJJwShUjDDoT+PcLnPysGfNPd10afGHr882mc130U+Bnw/4C5cZxvBJPQxLqfcPetZrYJ+ArwTeD8OGMQiZsePYk0zSvAD0LbvWJm5zRw3nKC/QVYcO/kwvAn3L2U4N4Fk4ivBbAKuMDMzghdr1uoIztsPnA/sNHdq5pWHZHGKVFIW1e3j2JWI+ffDZwAvG9mFaHj+vxf4CQzqwR+AVQAe6M+vwBY7u67GwvQ3T8Bvg3MDz3CWgl8MeqUhcBQ9NhJkkTDY0WSwMzaE+x/qDWz04H/Bs5y94Ohzy8G7nf3ZQ28PyHDdjU8VhJBLQqR5OgKvGlm7wEvAP/o7gfNLMfMPiTYiV5vkgj5W2MT7hpjZq8Bg4BDzb2GCKhFISIijVCLQkREYlKiEBGRmJQoREQkJiUKERGJSYlCRERiUqIQEZGY/n8yGKp7D975FAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots()\n", "ax.step(energies[:-1], flux500_mean/np.diff(energies), where='post', label='500 group')\n", "ax.step(energies_shem[:-1], flux_shem_mean/np.diff(energies_shem), where='post', label='SHEM-361')\n", "ax.set_xscale('log')\n", "ax.set_yscale('log')\n", "ax.set_xlabel('Energy [eV]')\n", "ax.set_ylabel('Flux [n-cm/eV-src]')\n", "ax.grid()\n", "ax.legend()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If we don't divide by the energy bin width, we obtain a plot where the different curves aren't directly comparable." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEKCAYAAAAMzhLIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAABGu0lEQVR4nO2deXhU5dm47+fMJOxLDApKMAHZERCCiFo0uOAultqqWHer1tq931dr61Ktlf6699OqWC1qWdS6FnHXKFUhEARFdpGwKCAxIFtIZub9/XHmzJzZJ8tkJslzX1eu5Jw5y/NOZt7nPOsrxhgURVEUJRFWtgVQFEVRchtVFIqiKEpSVFEoiqIoSVFFoSiKoiRFFYWiKIqSFFUUiqIoSlK82RYgE/Tq1cuUlJRkW4yE7Nu3jy5dumRbjCbRFsYAOo5coi2MAVr3OCorK3caYw6N3t8mFUVJSQlLlizJthgJKS8vp6ysLNtiNIm2MAbQceQSbWEM0LrHISJV8far60lRFEVJiioKRVEUJSmqKBRFUZSkqKJQFEVRkqKKQlEURUmKKgpFUZTGsrkCFvzR/t2GaZPpsYqiKBlncwU8ej7468CTD1e8AP3Gpz5n4wIomZj62BxCFYWiKEpD2VwB5feA/yCYgK0sNi5IPvk3RrHkCKooFEVRGoIz4fsOAgEQy574SyYmP2/jAltJGL/9e/nsVmNdqKJQFEVpCM6ETwCwYEAZlP0idrKPdjOVTLQVir8OLA98MBsCvlZhXaiiUBRFSZfNFbB7M1heW0948hMrCcfNZHlgzLehz2g45mJA7GMqHw1bF6ncVllGFYWiKEoi3FYBRE7+pVfA6EviT/ARbiY/LHnE3i8WeDrAmdPD1kU6bqsso4pCURQlms0VdgzB7R465pLw5B8AehQltgIcN5OvFjDh/U7g+0C17W7auAA6Fdq/IWetClUUiqIobkLBatck76+z/07DCui+ezVs3GdbDduW2crGX09M4NtRCq0gE0oVhaIoiptQsNqxBMSexEdPs3/iZSo5LqpOhYxefqttdTgTv3NOp0LbknCf63ZR+Q7aKbfxYh5ZRhWFoijtF9cEz4Fq+7c7WO0Eot2xiGSBa8Ayfnu/E6Se+NM0XFTBVNsN5VD1fs5ZFqooFEVp+8SriF4yE+b/FAJ+bOtB7N9i2YoiWbDajdsqcGN5Uwep+423lUL5PbaSSLd4r4VRRaFkl+gnulZQfKS0EtyfrZdvjowDQFBJ+FwnBF1NJmArj2TBajdRgWs7+VVgzLT0zu833nY3Vb2fs1lQqiiU7ODOKnECfYht6p/9Rxh3ZZYFVFot0RlLIvbk735aBwgEok50WRQNmKxnf9aHA92u46qae7Hw27aJE9OIQ2VVDQs3VDNhQCGlxQX2TseyyNEsKFUUSssSoSDcAUPsvwM+ePEnsG15ema/osSrdXBnLBkLLItQULpToZ2N5MmzH1IsC46/CTp2T27ZRrmvKqtqmP7SKhZvrOFGzzYCXoMltk2xc+A3+POSTux8awkAh3brwNSxRazZtodbn/sIvwGPwF0XjGTacUfa18/hLChVFErLEeMTToDxw5J/wrI5OfNFUbJIso6r0Y32nFqHiIwlLwyaDF0Ps6ujHTeU5bEt13QeSJzPrgmApwOPDfobt33QJfTywsAw6vGC8REQLzd8OJhKsyniErMWRW77Ddz6/AqG9OkWtiyi+0HlSKxCFYXSMmyuiOMTFvupbsy3oUN3eP9elxIxOZ0uqDSAeJlFcX4f/tlOWFAZqRBSdVyNnljdtQ6Wx1YQ616DNS8FFQnpF8255Xd9dv2+WrZ9+DowJXTIUjOYS+tuYYK1ioWBYSw1gxkrayO24+EPGB54+xMeunycvcPdDyqHYhWqKJTM47RkDriyQiROC4Sh58TGLT55Cz59J+wa0GB36yK606oTB4jzezAG1gVbXDgKIdUTdvTEGl3rsHGBrSTiKZJ0J+Llc0JKwhgIYLEwMCzmsKVmMEv9tkIYK2uZlf9b8vBRj5dL625JqCxeW7md6fNXcfPZwyJjFTn0WVdFoTQP0X5i5wkyojI1ODEkClj3G2//jJ5mK5ZP3iIUt3j3L8nPVTJPQxfdca/ZgBM4Nsl/R6eHpnrCTjSxuuVLpkiSjGP2ok0se+8V7t79eGii9GNxa/2VCSd9hwnWKvLw4ZUAGB8TrFUhJRKPGQs2cPqIPrYLyvke5BCqKJSmE90pE4nMZApNBhYcVZbaleSkC376Tmz6YsBnuwF6D8+5L1ObJl431ES+/UQZbQksiZgaBkchOIrpzOmxAeZopZXos5COIomisqqGXz37Eau27eFGz2LE60cE/Ahz/ZOYGzg1dOz4kgIuGFPEbc+vwBcIx90WBobhEy+CD+PJo8I3HABL4LRhvRnQqwsPvLMhdHzAwMIN1eFYRY6hikJpOhHugURPjgLeDunHG/qNty2HeMFvE8iZIF+7IaYbajDZwOlnhNiKA2KzjpwHhGFTksYo1m7eyZB+vWI7tUbHJhLFLRJZPGk+obuzmBzcQep6vDzjt2WzBH7jzlgCfvXcRzi64gMzmHnHPMiFhZ/iKZnILwKDYlJit31Vy3PLPgud//aaHZEpszmEKgql6UQvyOK2KJwq12RPoIkYd6VtOWxcALVf2cHuYNZJrgT52g0x3VAN+A7AvB8Tcit9MAvGXBqbdZTmA8Ln5eUMmVhmbyz4Y+LYRER/pAN2FtOYy2OL6tL8rFVW1fDA25/w2srtMa/FC1IDXDz+yAgl4fx92/Mr8AcMx+Wt53jPHiiZDP3GUwoxCmBQ724R9nbFxhq++cB7MQooF1BFoTSK7rtXR2aouM17aL5qa/fT4NBzGh/ka6WL2ucMzv94+WxY+i8I1AdfcBWt+Q/C58uT90lKl3ixCXf2lOWxLRuArZXw2TLANKgFRjwLIh7uILXD0Uf0iDlu2nFHMqRPN9578V98r/q3WEvrYfm9CZXWhAGFeCyJcFkFDNzy7Edsqt5nB7dzBFUUSsPZXBHbITPavE9jYpi9aBOP/HcDB3wBunfw8lVtPYjQt0dHBvXuxtSxRZFPYY0N8rWGRe0dv77jwsm2fPEUa+j9F9v1FK8WZutS24occhac+MPG/7+iYxMQ+T8cNBlWzwufY/y2gnKK6pJYnKsXv07l2//h6eqSpEHpop4duXHSIDZV7+PBdza4bSRq9tfFPae0uICCLmuwvqhPWQtRWlzAnVOO5pfPfhTzTj7wzgaOLOySM5aFKgql4WxcgBXwAY1rYBbvSW6r6/WtNQeo2FjDnIpN8c3whloHOVrEFGJzBcw8J9R9lA9mwZXzWk7G6PczqrgsRrGOvsSOT7hrFfZ8bisJjP0+r33ZVhSNkcVJp7UsO07Vb3ysK6rroWGLA2w5z/p/CS1Yx70kWyr468HbGYiPqfnx01YFuP6kARFP9EcWduG251cQMIZ8r8WEAYUJh7Cr59Fpp+A6n+1fPfuR2zYD4InFm1RRKK2YkokELC8ex6JIM17gfFlfX7k9WV12CMcML1+zg+tPPsq2LtwTiUj8J9dEi9o753RK/CVPSmPSQ9M5fuOCYEwnSEsqs2hr68zpkYWRvlrb0om2FqMziTZXwD/PCp/X2ISDjQvCNReBQDjDLVG9RAIrrLKqhqeXbmH99j1s2LmPnXtthXKjZxl53sRpq5OH9w5/1lw4bqWYHk1x+KrH0AbVQjjXjn54yvt8MVv/U07fYyZn/cEm5xWFiAwAfgn0MMZcmG15FKDfeJaPvouxh+xLe9KsrKrhohnv4/OnoyIieXXldl5ftd22LupcE4nBdj+seQnO+ZMd/Ha3CXG7QJwJ0ATgpf8J95JKl4a6rxpyfMnEYN8h5+k4feXbZJbPCQeo/XWw6vmoZnkGljwKe7+IVMjxXI1OllpTEg5KJtqWhCODo3Am/jS2aV7JRDj3L0DwIeSxJXz6xV7q/YaqL/fHvXx0FtPCwDAO7ZrPmCML4ioIN6XFBelnJDXQTVpaXMBTN5zAtx54j4qNNYyVtTzu/S35lb6kcY6WIqOKQkQeAc4FdhhjjnbtPxP4K+AB/mGMmZ7oGsaYDcA1IvLvTMqqpCDq6firHkPByVBJg9+9tCqhkji0az69unbgq9p6DvoD7NwT6/8NGPjlcx/R6+RBTHZPJGC7I178MdR8GmwD4gvvXz3PbuEw5lK7rNYJdgbTO7uPvANIYxwNdV819Pgxl9qTcdfDMhujiC6M/OBfhGINltdOYa16Pyq9NRB+H5O5xNxZao1NGkimcILXC8w8D+Ovox4vv+h6NxW+o9i6qzaty0dnMY2feCbP5FDQeFDvblRsrAkV7Hka6d5tbjJtUcwE7gUec3aIiAe4Dzgd2AIsFpEXsJXGPVHnX22M2ZFhGZVUxHs6bgDT56+iIk5mSbxcdLCD3O6cdAdj4PpyD3NLb+G4lXdHLhRjAvDuX4kbYHW3bnCnd/rr6LlrRXqD6FRou61IswV1uj17WjLQHn2vgadGurwKSmxLy6mNqHws8j32H4x1Q0WTxpN03DbbLmb7T2FZt+kMO/gh7weG8fHsvXTv8A5f1dZzSd2/ucFfh1cCGOPj8JolbPX3TfstOLRrPt5ex7G792n8MjpZIgeYOraIJ5ZsjrB8sPLwZjkdPKOKwhjzjoiURO0eD6wPWgqIyFxgijHmHmzrQ8k14j0dU5rWqdPnr4qoQHVI5AuGsM82XjzDABdVDuWGAX/meu9/KNj0OrHFfVE42TDuxe4DPvDk24HHVGyusHP0AwHbLXLm9NSTebo9e1oq0O5up2ECtvtuzUtEvGc719g/ng625dDnGALzfoIQVha+yn/x8K7j6D7oRGr211HQOZ8Vn+1GgBFH9GDFZ7tZv30PX+6rI89jUe8PcEiXfHp2zrfF+HI/q7ftCdViF3bLp4PXQ/cOXr6o2Uf9O6+we78POCL4A3AglOywxeqEyRN8RkKuo1T06pZPz075XH1i/5wJDieitLiAi8b1Y9YiE7J8FvuGc3NgUJrfuMwgxjTcZ9ygG9iKYp7jehKRC4EzjTHXBrcvA44zxtyU4PxC4G5sC+QfQYUS77jrgOsAevfuXTp37tzmHkqzsXfvXrp27ZptMdKm++7VjF5+KxLwYSwvy0ffxWeeopRjeHBZLe9v88fsP7vEy7eGdkjr3utr/Dyy4iCf7Yv9nArw16K3Obf6H0jwyVeAAEJ14Xjq8gvovG8TPb5ajWAIWHksH30XAD13rWBXz6PTGseRVf+m/6ezEAIEsNjY/1I2FV9I992rQ9f5qsfQtMYTTbz3NtW1yjfV886WesZ713O8tZLXaoeyKX8gRd3zWP2ln8H+tYwOrGR13gg2dxjCxAOvcZPvn1j4sYLvTwALiwAejG25CVjB6wcQHrK+xd/rp3CUby23eh9ntPUJloDPWPzJ903+7p+SRMLG4qiP+IyVtczJ/w15+Ahg8av6qyLaaTgc1gm8IvTpanF2/zwGFngyIGtimvr9Xl/j556KWvyGUAfaTn1HMnL4yGaUMj6TJk2qNMaMi96f88FsY0w1cEMax80AZgCMGzfOlJWVZViyxlNeXk4uyxdLGYwdG3o6HttvPF+lGMP0+at4f1usJXFDVNphGndmzFi7YjXGFQX8YMvJfNgZbjYP4SGAAXx4+Yu5iK921/On/bcgBOzpJ+Cjt7WTvuf9KnSNuONYMtMO6g6bYvvdN3eGR/8N/josTz4DTrmcAQCP3pG2y8hxt3TbsZSOn71HcVE/jusDq48YzYw9f6PHtoWU1w3h49WDGXF4V64/+SjWbNsTqjMZ71nPkNplLDLDeWtffy623uCuvJlYBJhMHpfuu4VtNfC/nnmc7qlEMNQf9PLwvjO53vsiFgYR8Bl4N3A0O013vu55jwDgwwMG8oKWQwDh09qu7A7AUgZzl+8yuxOqKwCcGRIrCYCpngXk40MExAQ42tpIcc/O5HmEQ7rkx6+9yQJN/X6XATs7rKJiwcuhDrTsfB7vUf/JWpwiG4piK9DPtV1EZBq9kos0IItj9qJNPBjH3dRQJeFQWlzAby4YGTduAZB/cBd47RCC3whP+k9i1md9uNHzPOIN2CthGnsC/P57Xdi8/LVw8Ly2lkOXvUO9P0Cex+LkvS/yc98D9oU/eZPHX3oLf153tnMFPT17WdNxNBWz9nCZ7xm+4ztoKyd/HRJ0GTkpwDu+qqV/ry6s2LqbXbX1FO9bwdetBXzT8zYe/Hh2GfwfC8Umj9fqbmGpmRwczQG21hzgVVc7ibBS8HMNHh7ynMV3vC/hxW5Wl2fqmOqxr52PHcgXgXzj4wbvPCS4bTsPLHaa7kz1vBu6vmC4tf4qyqzlnOZZihDgrryZUA9zA6cmbGORaYoP6YwvEAARunfw0mWPh5AXTOCMEb359iWTWkSWlqZbpzyO94Q70PoD2Q1oZ0NRLAYGiUh/bAVxMTAtC3IoGaCyqoZbn18REy244JgjmtSSwB23iO7Jk6hxm70/D0wdgvCQ7yx7kttbxxd7w5lVO7ftAWwz//K8J0HCE+ulvucxPgmvKbC3L1DLRqsDgTzBINQbDz9c0IldH79H5aaaUF/E5Vt2h677r/zfkk996MneGPBgyCN5C+qxspa78maGlIJl/FzvfREirgMnWR+SF3zaduMoCYfPAodwgec9+7Xgfo8JcIjs5UNzFKdRiSd4n7vyZrK2rp/dwiJOG4tUOBN9dEZS34JO9O3REYCtuw6EFMEXNXvo0LkTIw7vHj9+tbkjzCwHfz3iyePQr13VIHlaExMGFPK7N4ZHfK6rOo6mcQ7OppPp9Ng52JZULxHZAtxujHlYRG4CXsH+jD9ijPk4k3IoLcfvXlqFP+qxv7GWRDSlxQU8dPm4mKyoRE+8S81gfl1/GXflzUQIcJX3VV4PjIu7+pjj/3aeyJ3QnWDwiIkozhora7k973GEAAaLX9dfxqsHi+Gr+D2DQqmOYjDGXgLTgrgB2bGylqmeBQA845/IBGuV7ToLKgX7t8HCvQ395IuQ3JHHhseCQJG1M2xhBHf7xcP6LsfQJd+L2ftvjKOUCHBO9084sn8ZK7buBhFOG3oYXx30sXPPQXbtr+OgL8DxAwr5ZOc+dnxVy/EDCunWKS8io8kpfhNI6hpK6bLpNx6ufLFd9OwqLS5gUOmpXFphmOpZgAAdtn7F0GOzI0+ms57iVjQZY+YD8zN5b6UJNLKB3uxFm2LSYE8f3rvZm5s51oW78vaDvZFPvAIcUdCJIaYeqTX203twsidA0PdbjwkuRHOI7MXreiIPKwrbneWe0K/zzKMjdUGfv+EQ2ZtUXrfF48fiKf/JrAiUcIjsjVBsX+vwCQ/Lb8g3trL6lqecN/1jMGIRCE7ejgXhCGcgQuaQ9eB6LQDUd+5DhwM7gsrNCRsLFB9P3mm/Zobzf16SH6ph8Hg6cM23L4d+YxrzbwrRoEK1VOTgoj6ZYurYIn5XKXzDs4A8fFgfvQulaSzdmgFyPpittDCNzOuvrKrhb2+sjdgnwA0nH5URMaMnHydYXNA5n5r9deEn2s1d4dGngn2J8tjcZSyT6z4i31ePB4PBz2/yZ/Jk3hSod+XdSOSSS095zmZy3nousj7ijLolQDjuER3c7dHZy1f7fYwJWi2rO47iNz2mc3r9m3TK83JIt6EM+nwrL+4JK4l8r8Xvxn1FXqU/NNnnGT+TvUtC2UgE5YrODTKu15xVAGXEVPj4GTABLE8HOpzyCzvF13cQSdT6BJqnaE5pNo6TcJzCBOqzFqdQRaFE0oi8/kTtOa4/aUCLZaAkfGp11TN4Sybyf6G+RHMh4EMAL36m+Z6LOM31YI4HuNzMs8vDXWM0Ats6DWZHp1Hg8sN/td/Hsd51/CvvHvKoR+QFOGE6vPwG7DkIXz5HAGFafh6X1t3CB2YwF5YW0feYQgIf/BXjD8dPIpRElGzu7QAW1ok/iFxXfPx3Iif8mk/hvb/ZdRTr34ht2ue2JCf+NPkbrmSchRuqec8/lBst2xo1lpe8LBXeqaJQIkm3otjFg29/EqMkMuFyajRJ+hKZoLLAOC1BBA7pD19GZW2FVtkLT9MCHDlwJC9W/5k/7BvCv+pPCVkRfdmJ19QjTguGVc8HK8Tt+1gY8qnneM8qVspQvjG2CPoVYF31IjWv/T8CO9dTULsFCfiC97WCgYXouhSX3dOxe+QE7x735orwwk9gF965HwJaQyv2dsaEAYX8n2co3667heM9qxh1/LlMbkdZT0ouk25FcZDZizbFZCF5LMmYy6nZCLpY9s+5ii77t4T3i8AJPwy6aYJWgljhCdqTDxO+G7a8PnqSHsCdngUcf8geJu16ljzqcVxAmODvvM52hbg/gDOxW2IoHTaQU742IcIaKvjsv/a1xZV+JQInfB8WPRhuwX38Tfa2vw6DwO4t9oTvXhbUWTxq69LI9cfFinwIyPVW7O2Q0uICbjt3BLc9byAA1nsvcmRhZ4Yee1qLy6KKor2SLGDdgDWGo1NhBbhrytFZL3pKSXD8MWUZJmD76N3tto0reX/MNHulvUUP2stwuhhW8zZ5Uo9XDAaDBAwUnwBbltjtMiwP9B0bWrdBsDjlSA+43yv3hB3Ru8QP1evDq8whthxDz4Hls5HKx2DJI7D0cTj7D6F4hG3BuKMtwXEcf1Pk/7i5WrErzUrN/jpGs4Z/5dmFd4H5z0KfFlyrJEg8F6jS1nHcDG/ebf/eXNGoyzz49icRqbAicPfXc2+93xhc4++0//OoF8WeiA9UR7XbxrYIRk8LT+ZRbDc9AbENAAAMbHo/PPEHfHD4KPB2BPGAxxu2AhxKJhJTDOGweh68fofdq6ryUXsMAF+sCbUwIVAPi+6PcHPF7YG16MHI+/YbbytHy7KV5cs3N/pzoTQfEwYUcoKr8M4K1LN12astLocqivZI3CZ/DWN9jZ/XV0W6nE4b1jv3lQSEF8cxfjuGUHyC60VjT8ShbrHul4KTcclE2zpwvyQexllrg32UJDw1m7CrCROAPsfYVkHp5YCEJ3xnUu43Hg4fnVj2qneDjf2C/7vls+224BHCBDvlhr7eEvXbxP+/H6iObMXeiM+F0ryUFhfQZcgk6vHiMxb1eHnfP7zF5VDXU3ukEQHraN7d6otop2FJ5lJhm51OhThP2+IEqMUKB3oD9fakOeSsqHWZTXgRnag1my3jxxII4OFgtyPptKcq/r23LbPjIxsXhN1a0TGBMZfD1soUg5CgMhBnK8yEG8Mprk6MolNhTOfcmP97M3wulOan+6ATufTjcEHphYfH9OzLOKoo2iMNDFhHU1lVw5ovfRH7Th3WO/fjEg4HqiMVQ9W7ka87gd6Sifbaz04Q2D15dj007qUtIbGSsC9uNx1c/aJ9n+jrgq1IwF4nevdmYlxHVh6MvSy8Qt+yOQR8B7DEghN+ED7fHdTuPdze77jOEsWmmvC5UDJDzf46lmEXlFoCp+6PdXtmGlUU7ZVGVrjOXrSJW5/7yF1OgNfTCrKc3JRMDGYgxfnCWV47ddZ5b656Cd79C+zZZj/pO/tHT4MPZsVeI68T1CWo1BYPbFlsB50dhp6buPBt3JWR2UvblsVfBe+KF9j45mMMOOXyyOs0ZsGpdlT53FqYMKCQfK/F0f7VnOBdzaldvwkMbFEZVFEoaVNZVcNtz68gekXTb43r13qsCbAnwqDrKMJlk2jSXv+mPdluX2k/mTuT6ZUv2kFft5sokZLo0hv2bYdtH0Xu3/N5OBYQb4KOroVwJv5lc8K1Dv3Gs6l4PwOiz4+ORS2fDcvmaq1EK6O0uIDnzs9jwPx78Jh65OXnoE/LthzXYLaSNgs3VOOLaviX5xG7WKy10fWw2H19x8Z++ZIF/vuNt62MVIjHVhLx+PzD9LPPGpqE4FhOSPh3qvM3V9guL814yim6bVuIFajHQ4CAr67FM59UUShpMXvRJmYvivW9f7O1WRMOHboDUd7/eLUD0RlO0cdsW578PmJBjwSKtOeRwXqJOBN3vAnbCTaLJyau0X336tjjt690FdkZ6DM64fmhezZD2rTS/LzvH57VzCd1PSkpmb1oE7c8+1HMfq/VSq0JCE3KEa6nbcviH+u0kjV+eOl/w+4ne2fy+0T0+XZheWHXpuBGsHrbmbiXzAx1cMXTIcLFFDfYvLmC0ctvteVzxyJe/LErk8tnB/GTBau1Ojtn6T9mElcs+SXjzEp2Szdu9KyEzS3XSVYVhZKSl1ZEF6XZE+ydraECOxHdDo/d98Wa2H1OGqtD9AQ6ehoseZRwcRv2E7u75UenHrDbdU0nBTWES5lsrghXhENsT6Z4weaNC7ACPlsGR76tS139q4I4yiHR5KLpsTnNMobgCxhm5f+Wjkt9sPzeFoszqaJQUnLW0YezYN3O0LYlcPmw/NZRXJeIE38Iq+dFtvnzHYw9rmQiePLCE3u8CdTjyqAqPgFGXhS0TsR297z448jj42VbBeph+RzbTeWuCI/uyRSv9UrJRAKWF49jUZRMtNNv3RQOSj2haHpszrJwQzU+fyC0CJa4iyJVUSi5gKMQnli8id7dO3L9yUex59MUvvlcp994OPFHduqrQ7zAtJPd5PRXcqelQtDicHV0rXrPnszP/qOd3jrvR7FP9gkx9gTt7WA3JJSonkybK2DmueEn/ivDPX+29T6Fvn2PCMtXMjEyG2vIWemJoOmxOYmTIlvhtxfB8ogfaUGrTxWFkpLKqhpq9tdx23kjQq6m8k+zLFRzcPqvWftFHUP8q2HYlHChWjz27rBrKfqMTtBMr5ZQvCLgs91HvYcTu3pEEkZPC/dcevEntoJZ+He78V+/8bbF4Q9aPf6D9jbAzHM5wn8QdnQIF+F17E64GaAEt5XWSmlxAbOuncDTS4t4dF8fzuvxCX2PmawxCiU3qKyq4ZKHFlLvC5DntZjznQmtNy4Rh8+POIMhZfckPmBzBfzzrHDMwHlKd1c/nzndntB3umIcgYBtbYy+JFyYZ3ntyT9mTYkg21fa19u2LHyMU//QbzyxgXMTUh4CYeXhWBTejhpvaGM8s3QLdb7u/MVbyqxRgyhtoftqeqySlAff/oQ6n73Ocp0vwDNLt6Q8p00RHcwGexEih80VdtHdzshlYLGssJ//ynlw6q1w1Xy4+mU46hTiWhqh68Zbww7b4nD6O3ny7e1o5bE3WK/hxBtO+aUW1rURFm6ops4XIGCg3hdg4YbqFru3WhRKQiqrangjqkNsimTQtodTtOZWFsOmhP8OtRx3vTPiiWwDEu33L/sFbPxvbFDbuW6fqO6xzrYTL3EHm7evBFxr7619Nbx4kcYb2hROnMKx7icMaLk1Q1RRKHGprKrh509/GNGuwyO03rqJxtJvfGy/J3csI9Q3yrW4keUNxieSXDOq+yxI+JwD1djGfsD+faA68lz35H+gGpBwY/OAT+sf2ijhOMWWhkS+mgVVFEoMsxdt4pfPfhRhPVgCd10wsk3FJ9Km33i4eHaSA0zk3xFxhQTEdJ814QneyXyKF1+ITo8NVo6H1v7WeESbx45TBHh66RZmXdsyMUNVFEoE8ZQEwMi+PVp33USmiE6PBcDA0n+Fs5jiMXoaVD4WGdh22oMkqcCOSY8F29WFz/591v+LPF5rItoU8eIUqiiUFsXpDhsvDnHRsaok4hIvPRZSu4D6jY9aGEmSu5ggfnpsjyJwrAkIXyNei3FVFq2ebLUc16wnJcTCDdURa2CDHSC94aQBak0kwnn6H3elvaCQQyoX0OYKe1GkECZ+U8II4qTHBhVVACvyns2w3K2Sezgtx+d0vIcfe55k6CvfbpHmjWpRKFRW1fD00i18UFUTMRWNLyng52cNa59xiYbgPP2Pnpa4gjuaGJdVlEURj9HTYOnjtrViecOurXgLF2nfpjZLt20L8QTq7fXeW6iNhyqKdk5lVQ2XzHifOn+sJXHykMNUSWSK6B5Slje9yVw8gD/4O0i8hYu0b1ObpLKqht9XdOOflpc8fGDl4W2BhwBVFO2chRuqY5QEgMeSFs3TbvVsroCZ54Qn/g9mRfRiiqHfeDjr9zDvx9hpsGlUqGxcYDcPxNi/Uz1Jah1Fm2PhhmoW1Q/kUrmFCdYqFvuGc3Mg8xXaSRWFiIxN4xr1xpjYxQqUVsGEAYV4hIh6CZFW3kI8G2xcAP768HY6LoFtywi1Jw/4UqfUdioMNxg0gTRiGkpbY8KAQjyWsDQwmKX+wVhCi2Q+pbIo3gYWk7yzWX+gpLkEUlqW0uIC7rpgJLc+9xF+YxfV3XXBSA1eN5R02pFHE73+Rbz1MNwkK8RT2gWlxQXcOeVobnt+BQFjyG+hCu1UimKxMeaUZAeIyJvNKI+SBYb06cZF449EgKlji9SSaAyOK2nh322T7Ljvpnb7RK9/EW89DDfJCvGUdsO0445kSJ9uLNxQzYQBhdmvo0ilJNI9RsldKqtquPQfC6nzBcj3Wkxtby06movolelilkyNwyFHRa4ZcchRye+hAWolSGlxQYs+0KVVRyEiXxeRHq7tniJyQcakUlqMbHakbFMsnx21ZOrB1LUL+3cm345Hv/Ew8aeqJJQWJd2Cu9uNMaFVf40xu4DbMyKR0qIUdM7HEsESWrwjZdsiThgvVbDZ3YU23rai5AjppsfGUyiaWtvKmb1oE7c9vwJ/wOCxhNvOHaHxicYy+hKonBm57GmqYLPThXbV86lX2FOUKCqralosTpHuZL9ERP4E3Bfc/h5QmeR4Jcdx+jr5gi07AsZQs78uxVlKQravjF0bO5301XFXqoJQGkx0bDHTXWTTdT19H6gDngDmArXYykJppUT3dbJEC+yahHvVO4f1r7W8HEq7oKVjiyktChHxAPOMMZMyKkni+w8Dfgj0At4wxtyfDTnaGuu274moBb72a/3V7dQUOveK3bdzXcvLobQLWnq1u5QWhTHGDwTcWU/pIiKPiMgOEVkRtf9MEVkjIutF5OYU919ljLkB+BZwYkNlUGKZvWgTzy37LGJft055CY5W0iJexlKvQS0vh9IucFa7+8nkIS2yeFG6MYq9wEci8hqwz9lpjPlBivNmAvcCjzk7ghbKfcDpwBZgsYi8AHiAe6LOv9oYs0NEzge+CzyeprxKEl5a8XnEtiWo26mpDJsCn7hqT8UDJ/4we/IobZ6WrKUQY1I3IxORK+LtN8Y8msa5Jdiuq6OD28cDdxhjzghu/yJ4rWglEe9aLxpjzknw2nXAdQC9e/cunTt3bqrLZY29e/fStWvXrNx7fY2fJ9ccZO2u8P/97BIv3xraoUHXyeYYmpPmHMfhn73C4Z+/xsEOh7C531S+6jG0Wa6bDm3h/9EWxgCtexyTJk2qNMaMi96flkXhVggiUgD0M8Z82EhZ+gKbXdtbgOMSHSwiZcBUoAMwP4mMM4AZAOPGjTNlZWWNFC/zlJeX09LyOWtOPLF4E/5gco4IXD9xADefPazB18vGGDJB846jDMcojl4RO9O0hf9HWxgDtJ1xuElLUYhIOXB+8PhKYIeIvGuM+UkGZQPAGFMOlGf6Pm0ZJ5Wutj4qfdNobEJRWiu5WEfRwxjzlYhcCzxmjLldRBprUWwF+rm2i4L7lAzxzNItsUoCsHTNCUVplTgPfyP8qznoXU2X877J0GNPy9j90lUUXhE5HDvz6JdNvOdiYJCI9MdWEBcD05p4TSUBsxdtYk7Fppj9Atyla04oSqtk4YZqRvhX86+835KHj8D8Z6FPkoWymki6BXd3Aq8A640xi0VkAJAySVxE5gDvA0NEZIuIXGOM8QE3Ba+3CnjSGPNx48RXklFZVcOtz31EIE6+wvUnDdA1JxSllTJhQCEneFaRhw+vBLAC9Wxd9mrG7pdqhbtLgFeNMU8BTzn7jTEbgG+kurgx5pIE++eTJDCtNB633/LBtz8hziqnWGhsQlFaM6XFBSweMon6dc+C8VGPl/f9w7kwQ/dL5Xo6EnhKRPKAN4CXgAqTTk6t0uK4+79YAr7YsASW0GKrYimKkjmOnXgmV635FeeZt7EsYUzf7hm7V6qFi34H/E5EugGnAVcDD4jIKuBl4BVjzPaMSac0CHf/l2h3U5/uHfjBqYOp2V/XYqtiKYqSOUqLC7jj/BEc9dLdeE098sq70OeFjMQp0q2j2AM8G/xBRIYDZ2FXXJ/R7FIpjcLp/3KwPkC0yfeDUwdrTEJR2hjdti3EE6hHCNhL5G5ckBFFkW4wGxEZJSLni8hUYCjwqVNdreQGpcUF3HbuCLp08ETsH3hoF1USitLGqKyq4acV3ThovPiMRR1eVnccnZF7pVtw9wgwCvgYcDzfBngmI1IpjaKyqoY7XlhBXVQE++qvDciSRIqiZIqFG6qp8A3kUm5hgrWKhYFhrHyhnlmH1TS7azndOooJxpjhzXpnpdlZuKGa+iglcfrw3mpNKEobxHE1f1A/mKX+wQB4gmtTNLeiSNf19H4wLqHkMBMGFJLnCa/dnO+1uOHko7IokaIomcJpNX7JcUeS77XwZHDd+3QtisewlcU24CB2Ya8xxoxqdomURlNaXMCc647n6aVbEGDq2CLNblKUNozTavwbY4sy2vcpXUXxMHAZ8BHhGIWSg7Rkj3pFUXKDTH/v01UUXxhjXsiYFIqiKErOkq6i+EBEZgP/wXY9AWCM0awnRVGUNk66iqITtoKY7Nqn6bGKoijtgHQrs6/KtCCKoihKbpJWeqyIPCoiPV3bBcEiPCWHmL1oE5c9vIjZi2LXn1AURWks6bqeRhljdjkbxpgaERmTGZGUxjB9/ioeeGcDAAvW7QTQQjtFUZqFdAvuLBEJ5V6JyCGkr2SUDFNZVcOMBRsi9r204vMsSaMoSlsj3cn+j9gFd09hF9tdCNydMamUBvHM0i0xbcXPOvrw7AijKEqbI91g9mMisgQ4JbhrqjFmZebEUtKlsqqGp5ZsDm2LwPUTdZlTRVGaj7TbjBtjVhpj7gXqVEnkDgs3VOMLmhMCXDL+SG4+e1h2hVIUpU2RtqJwcUOzS6E0moLO+SG3kwGOPqJHVuVRFCU7VFbVcN9b66msqmn2azcmIC2pD1FaivI1O0J/W0DN/rrsCaMoSotTWVXDA29/wpurd2CMId9rMevaCc3a+6kxiuK8Zru70iSmz1/FqyvDS5Z7PJKRFsOKouQmlVU1XDLj/YjFyuoysCZFuivcdQC+AZQAXhHbqDDG3NlskigNIl5K7PDDu2vnWEVpR8RbrMyS5n9gTNeieB7YDVTiagqoZIfKqhru/M/HMSmxFx2rmU6K0p5wFitzLAqPJdw55eisLYVaZIw5s1nvrDSKeKampsQqSvukpRYrS1dRvCciI40xHzW7BEqDWLihOkJJgKbEKkp7piUWK0tXUXwNuFJEPkWXQs0qBZ3zY/ZpSqyiKJkkXUVxVkalUFJSWVXD00u38PHW3TGvaUqsoiiZJN0WHlWZFkRJTLy4hEO+19KUWEVRMkpSRSEiS40xY5t6jNI04qXAjS7qwdF9e2QseKUoiuKQyqIYJiIfJnldAHWQZ5gJAwrxWOALhPcdP6BQA9iKorQIqRTF0DSu4W8OQZTElBYXcNGxRzLLtXLdP/77KaeP6KPWhKIoGSepotDYRO4wdWwRTyzeHOoUGzCm2cv0FUVR4tGY7rFKC+J0hAS4c8rReC3BEg1iK4rScuhypjlMZVUNl/5jIXW+QKgj5BPXH8/CDdVMGFCo1oSiKC1CWhaFiAyPs6+suYVRwlRW1fCX19dS5wsQMFBbH+CBtz+htLiA700aqEpCUZQWI13X05Mi8nOx6SQi/wfck0nB2jOOJfHfdTsjGv+9tnI7s10BbUVRlJYgXUVxHNAPeA9YDHwGnJgpodyISJmILBCRB9qLFfPM0i3U1geILa+Dl1Z83uLyKIrSvklXUdQDB4BOQEfgU2NMIPkpICKPiMgOEVkRtf9MEVkjIutF5OYUlzHA3uB9t6Qpb6ulsqqGJxYnthrOOvrwFpRGURQlfUWxGFtRHAtMBC4RkafSOG8mENGeXEQ8wH3Y/aOGB681XERGisi8qJ/DgAXGmLOAnwO/TlPeVsvCDdX446jgksLO/PbrI7WVuKIoLY4YE8/BEXWQyDhjzJKofZcZYx5P49wSYJ4x5ujg9vHAHcaYM4LbvwAwxiSNeYhIPjDbGHNhgtevA64D6N27d+ncuXNTjitb7N27l65du8Z9bX2Nn+kVtfhc/xZL4PJh+ZQdmddCEqYm2RhaEzqO3KEtjAFa9zgmTZpUaYwZF70/3fTYHSIS/Sj7diNl6Qtsdm1vwY6BxEVEpgJnAD2BexMdZ4yZAcwAGDdunCkrK2ukeJmnvLycRPKVAWPG2p1i12/fQ+WmXRhjmLvOx3knj8uZbKdkY2hN6Dhyh7YwBmg743CTrqJ4ETtWINixgv7AGmBEhuQKYYx5Bngm0/fJNfr27IQAS6pqCBioz8CC6YqiKOmQbpvxke5tERkL3NjIe27FzqByKAruU4gssvNagtdj4fcHyNNKbEVRskSjKrONMUtFJKG7KAWLgUEi0h9bQVwMTGvktdoMlVU1LNxQzWe7DoSK7PwBw0Xj+9G3ZyetxFYUJWukpShE5CeuTQsYi11Lkeq8Odhu914isgW43RjzsIjcBLwCeIBHjDEfN1TwtoRjRRysD2AJeCxBAoY8r8U3dL0JRVGyTLoWRTfX3z7smMXTqU4yxlySYP98YH6a927zLNxQzcFggZ3fAAHDxeOP1EWJFEXJCdKNUbT5+oVsUtA5P6IK22/giJ6dVEkoipITpFoK9T8Qt5MEAMaY85tdonaCE5Mo6Jwfty1HQef8LEilKIoSSyqL4g8tIkU7Y32Nnz+8sTDkbpKo1wWo2V+XBckURVFiSaUoPjXGaLvSZubdrb6QkoCwySZib3gsUYtCUZScIVWvp+ecP0QkZfBaSU1lVQ0Ltvri+vMsAcsSAsZw57yPqayqaXH5FEVRokmlKNxekQGZFKS9sHBDtZ3ZhP3mDjysa+hNDgQgEDARldiKoijZJpWiMAn+VhrJhAGF5FngEeiQZ3H1if3pkGfhEcjzCHne4N9aia0oSo6QKkYxWkS+wn747RT8m+C2McZ0z6h0bZDS4gL+99iOHOxZHKq2HtKnW2gdbEDXxFYUJadIqiiMMZ6WEqQ9MbDAQ1nZwNB2aXFBhFJQBaEoSi6R7sJFiqIoSjtFFYWiKIqSFFUUWaKyqob73lqvKbCKouQ8jWozrjSN8k31/OvV9wkYg9cSvjmunzYAVBQlZ1FF0UK4ezs9vqouVEtR5zfMXrSJp5duYda1E1RZKIqSc6iiaAHcq9ZZIiEl4WDQpU4VRcldVFFkkHir1mEMHgFj7N5OliUEgosUaYGdoii5iCqKDJFs7euLB3k5tF9/LbBTFKVVoIoiQyzcUB2z9jXYJe1F7ODaSZEFd4qiKLmKKooMMWFAIflei3qfbUUcfUQP7pz3cWhd7M5HbGLacUdmW0xFUZSUqKLIEKXFBcy6dkLIrRS9LvZtz69gSJ9uak0oipLzqKLIINE9nDyW4AvYKU8BYzTLSVGUVoFWZmeIeJXXk4YehscSBMjXLCdFUVoJalFkAHfGU77X4rZzR3DnvI9DGVAnFXm56bzxak0oitIqUIuimamsquEvr68NZTzV+wK8tOLziAyowk5CaXGB9ntSFKVVoBZFM+JYEk7Q2gquVHfW0YezeOOXoQyooYd4YqwObd+hKEquooqiGXFqJwy2qXbiwF786LTBMavY7fl0eUSdhbbvUBQll1FF0YxE1044SgIiM6DKP409VgPbiqLkKqoompHo2glHMTg9n9z7Eh2rKIqSa6iiaGaiayfixSISHasoipKLaNZTMxIvi+mZpVs4WB8Zi1AURWlNqEXRTCSyHJ5ashln+QmPxwoGs7dkT1BFUZQGooqiCbhjD/GymIBQyw4BLiy1lzst/zSLQiuKojQQVRSNJF71dXQW05pte7BEAEO+1+IbY4uyLbaiKEqDUUXRSNwWxMH6AG+t2cHUsUUIMDWoEO6c9zH+gMFjCbedO0ID14qitEpUUTSSCQMK8VpCnd9ggNdWbscSu9nf1LFFEcV3xhhq9tdlW2RFUZRGoVlPjaS0uIBvjuuHuPa54xNOQZ0n2MZDC+oURWmtqEXRBKaOLeLppVuoqw8QINzbqaBzPgs3VHPbuSNY8dnuCGWiKIrS2sh5RSEiE4FLsWUdbow5IcsihXBXVxd0zqdmfx0FnfO544UV1PsNHgssy8LnD/D00i3a+E9RlFZJRhWFiDwCnAvsMMYc7dp/JvBXwAP8wxgzPdE1jDELgAUicgGwOJPyNobo6upbnv2IOr+dEusLAIEAoI3/FEVpvWTaopgJ3As85uwQEQ9wH3A6sAVYLCIvYCuNe6LOv9oYsyP49zTgmgzL22Si3UweS8AYjVMoitJqyaiiMMa8IyIlUbvHA+uNMRsARGQuMMUYcw+29RGDiBwJ7DbG7MmkvM3B1LFFPFW5JVRPccd5I6jZX6eN/xRFabWIMSb1UU25ga0o5jmuJxG5EDjTGHNtcPsy4DhjzE1JrvFr4BVjzHtJjrkOuA6gd+/epXPnzm2+QTSQ9TV+Vn/pZ+ghHgYWeGJe37t3L127ds2CZM1HLo1BROjSpQseT+x7nQpjDCKtP90gE+Pw+/3s27ePTM8RDrn0mWoKrXkckyZNqjTGjIven/PBbABjzO1pHDMDmAEwbtw4U1ZWlmmxEpLqzuXl5WRTvuYgl8bw6aef0q1bNwoLCxs8We7Zs4du3bplSLKWo7nHYYyhurqaPXv20L9//2a7bjJy6TPVFNrKONxko45iK9DPtV0U3NfmcbrLrq/xZ1uUNkVtbW2jlISSGBGhsLCQ2trabIui5ADZsCgWA4NEpD+2grgYO1DdpnH3hvIKjBlbozGLZkSVRPOj76nikOn02DnYnpheIrIFuN0Y87CI3AS8gp3p9Igx5uNMytEYnM6wTn1EqmB0vFXs3Lh7Q/kMmiqrKEqrIdNZT5ck2D8fmJ/JezcF5+n/YL3dq8np4ZSoYC7eWhTRx7nXyPYImirbxigpKaFbt254PB68Xi9LliwB4Msvv+Siiy5i48aNlJSU8OSTT1JQUIAxhh/+8IfMnz+fzp07M3PmTMaOHZvlUShKfLTXUxzcDf2AlKvTJVqLwo1Txf2TyUP432M7qjXRBnnrrbdYtmxZSEkATJ8+nVNPPZV169Zx6qmnMn26XVv60ksvsW7dOtatW8eMGTP47ne/2+j7+ny+JsuuKMlQRREH5+nfeXOcHk6JrIB0GwCWFhfwvUkD46bMKi2Lk1iwbMtXGb3P888/zxVXXAHAFVdcwXPPPRfaf/nllyMiTJgwgV27dvH555/HnH/XXXcxZMgQvva1r3HJJZfwhz/8AYCysjJ+9KMfMW7cOP76179SXl7OmDFjGDlyJFdffTUHDx4EbEtn586dACxZsiSUjXPHHXdw2WWXcfzxxzNo0CAeeuihjL4PSuumVaTHtjTxejgli1G4j9fCutwnwlXosZj1nc5N/p+JCJMnT0ZEuP7667nuuusA2L59O4cffjgAffr0Yfv27QBs3bqVfv3CyX9FRUVs3bo1dCzA4sWLefrpp1m+fDn19fWMHTuW0tLS0Ot1dXUsWbKE2tpaBg4cyJtvvsngwYO5/PLLuf/++/nRj36UVOYPP/yQhQsXsm/fPsaMGcM555zDEUcc0aT3QWmbqKJIQHQPp+Y+XskeEa5Cf/P04Prvf/9L37592bFjB6effjpDhw7lpJNOijhGRBqUSfTuu+8yZcoUOnbsSMeOHTnvvPMiXr/ooosAWLNmDcXFxQwePBiwLZf77rsvpaKYMmUKnTp1olOnTkyaNImKigouuOCCtOVT2g/qelLaHRGuQk/z9ODq27cvAIcddhhf//rXqaioAKB3794hl9Lnn3/OYYcdFjp+8+bNofO3bNkSuka6dOnSJeUxXq+XQLAxZXRNRLTS0nRYJRGqKJR2hzux4KFLRzXZmti3bx979uwJ/f3qq69y9NF2s+Tzzz+fRx99FIBHH32UKVOmhPY/9thjGGNYuHAhPXr0iHA7AZx44on85z//oba2lr179zJv3ry49x8yZAibNm1i/fr1ADz++OOcfPLJgB2jqKysBODpp5+OOO/555+ntraW6upqysvLOfbYY5v0PihtF3U9Ke0Sx1XoTPBNYfv27Xz9618H7AykadOmceaZZwJw8803861vfYuHH36Y4uJinnzySQDOPvts5s+fz8CBA+ncuTP//Oc/Y6577LHHcv755zNq1Ch69+7NyJEj6dGjR8xxHTt25O9//zvf/OY38fl8HHvssdxwww0A3H777VxzzTXceuutMW0lRo0axaRJk9i5cye33nqrxieUhKiiUJQmMmDAAJYvXx73tcLCQt54442Y/SLCfffdl/LaP/vZz7jjjjvYv38/J510UiiYXV5eHnFcWVkZH3zwQcz5EydOZO3atXGvPWrUKB577LG4rymKG1UUipLDXHfddaxcuZLa2lquuOIKLcpTsoIqCkXJYWbPnp2R695xxx0Zua7SNtFgtqIoipIUVRSKoihKUlRRKIqiKElRRaEoiqIkRRWFojQDd999NyNGjGDUqFEcc8wxLFq0iLKysohOshs3bgwV4pWXl9OjRw+OOeaY0M/rr78O2Kmz3/72t0Pn+Xw+Dj30UM4999y4977mmms44YQTGDVqFBdeeCF79+4Nvfbkk08yfPhwRowYwbRp4fXBzjzzTHr27JnwmoriRrOeFKWJvP/++8ybN4+lS5fSoUMHdu7cSV1dXcrzJk6cGLfaukuXLqxYsYIDBw7QqVMnXnvttaTtPf785z8jInTr1o2f/OQn3Hvvvdx8882sW7eOe+65h3fffZeCggJ27NgROud//ud/2L9/Pw8++GDjBq20K9SiUNonmytgwR+xPqts8qU+//xzevXqRYcOHQDo1atXk6uczz77bF588UUA5syZwyWXxF0DDIDu3bsDYIzhwIEDoZ5NDz30EN/73vcoKLBblDh9pgBOPfVUunXr1iQZlfaDKgql/bG5Ah49H968m85PXWRvN4HJkyezefNmBg8ezI033sjbb78deu3SSy8NuZbOPvvsiPMWLFgQ4Xr65JNPQq9dfPHFzJ07l9raWj788EOOO+64pDJ897vfpU+fPqxevZrvf//7AKxdu5a1a9dy4oknMmHCBF5++eUmjVNpv6iiUNofGxeAvw6MH/z19nYT6Nq1K5WVlcyYMYNDDz2Uiy66iJkzZwIwa9Ysli1bxrJly5g/P3L134kTJ4ZeW7ZsGUcddVTotVGjRrFx40bmzJkTo2Dicf/99/PZZ58xbNgwnnjiCcCObaxbt47y8nLmzJnDd77zHXbt2tWksSrtE1UUSvujZCJ48kE84Mmzt5uIx+OhrKyMX//619x7770xnVobw/nnn8/PfvazGLfTGWecwTHHHMO1114bI8PFF18cundRURHnn38+eXl59O/fn8GDB7Nu3bomy6W0PzSYHUVlVY2uVNfW6TcerngBNi5g/2GldOk3vkmXW7NmDZZlMWjQIACWLVtGcXExK1asaNJ1r776anr27MnIkSMjmgC+8sorob+NMXzyySf07t0bYwwvvPACQ4cOBeCCCy5gzpw5XHXVVezcuZO1a9cyYMCAJsmktE9UUbiIWCLTazHr2gmqLNoq/cZDv/EEmqHN+N69e/n+97/Prl278Hq9DBw4kBkzZnDhhRcmPc+JUTj86le/ijinqKiIH/zgB0mvYYzhiiuuYNeuXYgIo0eP5v777wdsy+PVV19l+PDheDwefv/731NYaC/SNHHiRFavXs3evXspKiri4Ycf5owzzmjkO6C0dVRRuIhYItPXPEtkKm2f0tJS3nvvvZj90a3AS0pKQlZGWVkZu3fvjns9dx2EQ1lZWcx6EgCWZfHuu++yZ8+emCwmEeFPf/oTf/rTn2LOW7CgaXEZpX2hMQoXEUtkeptniUxFUZTWjloULpwlMjVGoSiKEkYVRRTOEplK68IYEyo0U5oHY0y2RVByBHU9Ka2ejh07Ul1drRNbM2KMobq6mo4dO2ZbFCUHUItCafUUFRWxZcsWvvjiiwafW1tb2yYmw0yMo2PHjhQVFTXrNZXWiSoKpdXjFJQ1hvLycsaMGdPMErU8bWUcSm6iridFURQlKaooFEVRlKSoolAURVGSIm0xU0REvgCqsi1HEnoBO7MtRBNpC2MAHUcu0RbGAK17HMXGmEOjd7ZJRZHriMgSY8y4bMvRFNrCGEDHkUu0hTFA2xmHG3U9KYqiKElRRaEoiqIkRRVFdpiRbQGagbYwBtBx5BJtYQzQdsYRQmMUiqIoSlLUolAURVGSoopCURRFSYoqCkVRFCUpqihyDBHpIiJLROTcbMvSWETkAhF5SESeEJHJ2ZanIQTf/0eD8l+abXkaQ2t+/6Np7d8HEbFE5G4R+T8RuSLb8jQWVRTNhIg8IiI7RGRF1P4zRWSNiKwXkZvTuNTPgSczI2VqmmMcxpjnjDHfAW4ALsqkvOnQwDFNBf4dlP/8Fhc2AQ0ZQ669/24a8fnK6vchHg0cwxSgCKgHtrS0rM2FKormYyZwpnuHiHiA+4CzgOHAJSIyXERGisi8qJ/DROR0YCWwo6WFdzGTJo7Ddeqvgudlm5mkOSbsL/Xm4GH+FpQxFTNJfwwOufL+u5lJ+p+vXPg+xGMm6f8vhgDvGWN+Any3heVsNnQ9imbCGPOOiJRE7R4PrDfGbAAQkbnAFGPMPUCMKS0iZUAX7A/aARGZb4wJZFLuaJppHAJMB14yxizNsMgpaciYsJ/6ioBl5NCDVEPGICKryKH3300D/xddyfL3IR4NHMNmoC54TC49eDQIVRSZpS/hp1OwJ6HjEh1sjPklgIhcCezMhS9FkAaNA/g+cBrQQ0QGGmMeyKRwjSTRmP4G3Csi5wD/yYZgDSDRGFrD++8m7jiMMTdBTn4f4pHof/FX4P9EZCLwTjYEaw5UUeQgxpiZ2ZahKRhj/oY94bY6jDH7gKuyLUdTaM3vfzxa8/fBGLMfuCbbcjSVnDGt2yhbgX6u7aLgvtZGWxmHm7YwprYwBmgb42gLY0iIKorMshgYJCL9RSQfuBh4IcsyNYa2Mg43bWFMbWEM0DbG0RbGkBBVFM2EiMwB3geGiMgWEbnGGOMDbgJeAVYBTxpjPs6mnKloK+Nw0xbG1BbGAG1jHG1hDA1FmwIqiqIoSVGLQlEURUmKKgpFURQlKaooFEVRlKSoolAURVGSoopCURRFSYoqCkVRFCUpqiiUdo2I+EVkmesnnVbwGccl1xFJjrldRO6J2ndMsCkgIvKWiOwVkXGZlldp22gdhdKuEZG9xpiuzXxNb7AAqynXSCmXiAwGXjbGDHDtmw7sN8bcGdwuB35mjFnSFHmU9o1aFIoSBxHZKCK/FpGlIvKRiAwN7u8SXLimQkQ+EJEpwf1XisgLIvIm8IaIdBaRJ0VkpYg8KyKLRGSciFwtIn9x3ec7IvLnNOSZLCLvB+V5SkS6GmPWAjUi4u7k+y1gTrO+GUq7RxWF0t7pFOV6cq8It9MYMxa4H/hZcN8vgTeNMeOBScDvRaRL8LWxwIXGmJOBG4EaY8xw4FagNHjMk8B5IpIX3L4KeCSZgCLSC3sRotOC8iwBfhJ8eQ52XyFEZALwpTFmXcPfBkVJjLYZV9o7B4wxxyR47Zng70rsJVIBJgPni4ijODoCRwb/fs0Y82Xw769hr0WAMWaFiHwY/Htv0Oo4NxhLyDPGfJRCxgnYi/e8a68JRT52ryGAJ4D3ROSn2ApDrQml2VFFoSiJORj87Sf8XRHgG8aYNe4Dg+6ffWle9x/ALcBq4J9pHC/YSuiS6BeMMZtF5FPgZOAbwPFpyqAoaaOuJ0VpGK8A3w8u94qIjElw3LvY8QLEXjt5pPOCMWYR9toF00jPAlgInCgiA4PX6xIMZDvMAf4MbDDGbGnYcBQlNaoolPZOdIxieorj7wLygA9F5OPgdjz+DhwqIiuB3wAfA7tdrz8JvGuMqUkloDHmC+BKYE7QhfU+MNR1yFPACNTtpGQITY9VlAwgIh7s+EOtiBwFvA4MMcbUBV+fB/zZGPNGgvObJW1X02OV5kAtCkXJDJ2B/4rIcuBZ4EZjTJ2I9BSRtdhB9LhKIshXqQruUiEibwEDgPrGXkNRQC0KRVEUJQVqUSiKoihJUUWhKIqiJEUVhaIoipIUVRSKoihKUlRRKIqiKElRRaEoiqIk5f8DnmzVQDz3ryYAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots()\n", "ax.loglog(energies[:-1], flux500_mean, '.', color='C0', label='500 group')\n", "ax.loglog(energies_shem[:-1], flux_shem_mean, '.', color='C1', label='SHEM-361')\n", "ax.set_xlabel('Energy [eV]')\n", "ax.set_ylabel('Flux [n-cm/src]')\n", "ax.grid()\n", "ax.legend()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.1" } }, "nbformat": 4, "nbformat_minor": 4 }