uproot_example.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# High level analysis example via uproot\n",
    "\n",
    "`uproot` is a python package that allows users to read `ROOT` files without acutally compiling `ROOT`. I believe this is a step up in terms of making the experimental analysis more accessible. This is a very minimal example on how to extract info from a `TTree` object, and how to read the histograms (`TH1`, `TH2`).  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import uproot\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example on how to read a root tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[b'DataTree;1', b'HeaderTree;1']"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# There are two trees in this root files\n",
    "file = uproot.open(\"data/test_flat_v3.root\")\n",
    "file.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[b'fUniqueID',\n",
       " b'fBits',\n",
       " b'iEventNumber',\n",
       " b'iRunNumber',\n",
       " b'volumes',\n",
       " b'volumes_sName',\n",
       " b'volumes_iVolumeID',\n",
       " b'volumes_dTotalEnergyDep_keV',\n",
       " b'volumes_iTotalOptPhotNumber',\n",
       " b'volumes_iTotalThermElecNumber',\n",
       " b'primaryParticles',\n",
       " b'primaryParticles_sName',\n",
       " b'primaryParticles_dPosition_mm_x',\n",
       " b'primaryParticles_dPosition_mm_y',\n",
       " b'primaryParticles_dPosition_mm_z',\n",
       " b'primaryParticles_dDirection_x',\n",
       " b'primaryParticles_dDirection_y',\n",
       " b'primaryParticles_dDirection_z',\n",
       " b'primaryParticles_dEnergy_keV',\n",
       " b'primaryParticles_dTime_ns',\n",
       " b'primaryParticles_sVolumeName',\n",
       " b'primaryParticles_iVolumeID',\n",
       " b'tracks',\n",
       " b'tracks_sParticleName',\n",
       " b'tracks_iParticleID',\n",
       " b'tracks_iTrackID',\n",
       " b'tracks_iParentID',\n",
       " b'tracks_sCreatorProcess',\n",
       " b'tracks_dWavelength_nm',\n",
       " b'tracks_dCharge',\n",
       " b'tracks_steps',\n",
       " b'tracks_steps_iStepNumber',\n",
       " b'tracks_steps_sProcess',\n",
       " b'tracks_steps_sVolumeName',\n",
       " b'tracks_steps_iVolumeID',\n",
       " b'tracks_steps_dPosition_mm_x',\n",
       " b'tracks_steps_dPosition_mm_y',\n",
       " b'tracks_steps_dPosition_mm_z',\n",
       " b'tracks_steps_dDirection_x',\n",
       " b'tracks_steps_dDirection_y',\n",
       " b'tracks_steps_dDirection_z',\n",
       " b'tracks_steps_dTime_ns',\n",
       " b'tracks_steps_dParticleEnergy_keV',\n",
       " b'tracks_steps_dEnergyDep_keV',\n",
       " b'tracks_steps_trackIndex']"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# The DataTree has many branches. \n",
    "t = file[\"DataTree\"]\n",
    "t.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAD8CAYAAACSCdTiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAECFJREFUeJzt3X+MZWV9x/H3p6yooBYoo0VgO9AoCTWm0NGAttaCKAKBNrEJG3+sitnEplZtG11CUtP/8EesNTbiRldIS1ct4o+AVglqSRO7dhcRFhcElOIqukNI1WhTJH77xz2Lw3R+3nN25s7j+5VM5pznnHvP97nPnU/OnHvOuakqJElt+rX1LkCSdPgY8pLUMENekhpmyEtSwwx5SWqYIS9JDTPkJalhhrwkNcyQl6SGbVrLjR1//PE1PT29lpuUpA1v7969D1XV1DiPXdOQn56eZs+ePWu5SUna8JL817iP9XCNJDXMkJekhhnyktQwQ16SGmbIS1LDDHlJapghL0kNM+QlqWGGvCQ1bE2veNXjTW+/ccnl91954RpVIqlV7slLUsMMeUlqmCEvSQ0z5CWpYYa8JDXMkJekhhnyktSwZUM+yc4kB5Psm9f+piR3J7kzybsOX4mSpHGtZE/+auD8uQ1J/gi4BHhuVf0O8J7hS5Mk9bVsyFfVLcDD85rfCFxZVf/brXPwMNQmSepp3GPyzwb+IMnuJP+W5HlDFiVJGsa4967ZBBwLnAU8D/hEklOrquavmGQbsA1g8+bN49YpSRrDuHvyB4Dra+RrwC+A4xdasap2VNVMVc1MTU2NW6ckaQzjhvyngXMAkjwbOBJ4aKiiJEnDWPZwTZJdwIuB45McAN4B7AR2dqdVPgJsXehQjSRpfS0b8lW1ZZFFrxq4FknSwLziVZIaZshLUsMMeUlqmCEvSQ0z5CWpYYa8JDXMkJekho177xqt0PT2G9e7BEm/wtyTl6SGGfKS1DBDXpIaZshLUsMMeUlqmCEvSQ0z5CWpYYa8JDVs2ZBPsjPJwe5boOYv++sklWTB73eVJK2vlezJXw2cP78xycnAecADA9ckSRrIsiFfVbcADy+w6O+AtwF+t6skTaixjsknuRj4XlV9YwXrbkuyJ8me2dnZcTYnSRrTqkM+yVHAFcDfrGT9qtpRVTNVNTM1NbXazUmSehhnT/63gVOAbyS5HzgJuDXJbw5ZmCSpv1Xfariq7gCefmi+C/qZqnpowLokSQNYySmUu4CvAqclOZDkssNfliRpCMvuyVfVlmWWTw9WjSRpUF7xKkkNM+QlqWGGvCQ1zJCXpIYZ8pLUMENekhpmyEtSwwx5SWqYIS9JDTPkJalhhrwkNcyQl6SGGfKS1DBDXpIaZshLUsNW8qUhO5McTLJvTtu7k9yV5PYkn0pyzOEtU5I0jpXsyV8NnD+v7SbgOVX1XOBbwOUD1yVJGsCyIV9VtwAPz2v7YlU92s3+B6Mv85YkTZghjsm/Hvj8AM8jSRpYr5BPcgXwKHDtEutsS7InyZ7Z2dk+m5MkrdLYIZ9kK3AR8MqqqsXWq6odVTVTVTNTU1Pjbk6SNIZN4zwoyfnA24E/rKqfDVuSJGkoKzmFchfwVeC0JAeSXAZ8AHgqcFOS25JcdZjrlCSNYdk9+araskDzRw5DLZKkgXnFqyQ1zJCXpIYZ8pLUMENekhpmyEtSwwx5SWqYIS9JDTPkJalhhrwkNcyQl6SGGfKS1DBDXpIaZshLUsMMeUlqmCEvSQ1byZeG7ExyMMm+OW3HJbkpyT3d72MPb5mSpHGsZE/+auD8eW3bgZur6lnAzd28JGnCLBvyVXUL8PC85kuAa7rpa4A/HrguSdIAxj0m/4yqehCg+/304UqSJA3lsH/wmmRbkj1J9szOzh7uzUmS5hg35H+Y5ASA7vfBxVasqh1VNVNVM1NTU2NuTpI0jnFD/rPA1m56K/CZYcqRJA1pJadQ7gK+CpyW5ECSy4ArgfOS3AOc181LkibMpuVWqKotiyw6d+BaJEkD84pXSWqYIS9JDTPkJalhhrwkNcyQl6SGGfKS1DBDXpIaZshLUsMMeUlqmCEvSQ0z5CWpYYa8JDXMkJekhhnyktQwQ16SGmbIS1LDeoV8krcmuTPJviS7kjxpqMIkSf2NHfJJTgT+ApipqucARwCXDlWYJKm/vodrNgFPTrIJOAr4fv+SJElDWfY7XhdTVd9L8h7gAeB/gC9W1Rfnr5dkG7ANYPPmzeNubmJNb79xvUuQpEX1OVxzLHAJcArwTODoJK+av15V7aiqmaqamZqaGr9SSdKq9Tlc8xLgO1U1W1U/B64HXjBMWZKkIfQJ+QeAs5IclSTAucD+YcqSJA1h7JCvqt3AdcCtwB3dc+0YqC5J0gDG/uAVoKreAbxjoFokSQPzildJapghL0kNM+QlqWGGvCQ1zJCXpIYZ8pLUMENekhrW6zz5tbTcjcDuv/LCNapEkv6/Sc0o9+QlqWGGvCQ1zJCXpIYZ8pLUMENekhpmyEtSwwx5SWpYr5BPckyS65LclWR/krOHKkyS1F/fi6H+HvjXqnpFkiOBowaoSZI0kLFDPsnTgBcBrwWoqkeAR4YpS5I0hD6Ha04FZoGPJvl6kg8nOXqguiRJA+gT8puAM4EPVtUZwE+B7fNXSrItyZ4ke2ZnZ3tsTpK0Wn1C/gBwoKp2d/PXMQr9x6mqHVU1U1UzU1NTPTYnSVqtsUO+qn4AfDfJaV3TucA3B6lKkjSIvmfXvAm4tjuz5tvA6/qXJEkaSq+Qr6rbgJmBapEkDcwrXiWpYYa8JDXMkJekhhnyktQwQ16SGmbIS1LDDHlJaljfi6F+JUxvv3G9S5CksbgnL0kNM+QlqWGGvCQ1zJCXpIYZ8pLUMENekhpmyEtSw3qHfJIjui/yvmGIgiRJwxliT/7NwP4BnkeSNLBeIZ/kJOBC4MPDlCNJGlLfPfn3AW8DfjFALZKkgY1975okFwEHq2pvkhcvsd42YBvA5s2bx91cL8vde+b+Ky9co0pWZ6PWLWly9NmTfyFwcZL7gY8B5yT5p/krVdWOqpqpqpmpqakem5MkrdbYIV9Vl1fVSVU1DVwKfKmqXjVYZZKk3jxPXpIaNsj95KvqK8BXhnguSdJw3JOXpIYZ8pLUMENekhpmyEtSwwx5SWqYIS9JDTPkJalhg5wnv9Etd48YSY/nfZU2DvfkJalhhrwkNcyQl6SGGfKS1DBDXpIaZshLUsMMeUlqmCEvSQ0bO+STnJzky0n2J7kzyZuHLEyS1F+fK14fBf6qqm5N8lRgb5KbquqbA9UmSeqpzxd5P1hVt3bTPwH2AycOVZgkqb9B7l2TZBo4A9i9wLJtwDaAzZs3D7G5BXn/mcfr83p435Fhred9Xpba9nLb7fMeavHeNhs1Y3p/8JrkKcAngbdU1Y/nL6+qHVU1U1UzU1NTfTcnSVqFXiGf5AmMAv7aqrp+mJIkSUPpc3ZNgI8A+6vqvcOVJEkaSp89+RcCrwbOSXJb93PBQHVJkgYw9gevVfXvQAasRZI0MK94laSGGfKS1DBDXpIaZshLUsMMeUlqmCEvSQ0z5CWpYYPcoExtOZw3Yup7Y6r1uuHWcg7XDbdavNEX9BvHw7XdVrknL0kNM+QlqWGGvCQ1zJCXpIYZ8pLUMENekhpmyEtSw/p+/d/5Se5Ocm+S7UMVJUkaRp+v/zsC+Afg5cDpwJYkpw9VmCSpvz578s8H7q2qb1fVI8DHgEuGKUuSNIQ+IX8i8N058we6NknShEhVjffA5E+Bl1XVG7r5VwPPr6o3zVtvG7Ctmz0NuHuVmzoeeGisIidXa32yP5OvtT79qvXnt6pqapwn7nODsgPAyXPmTwK+P3+lqtoB7Bh3I0n2VNXMuI+fRK31yf5Mvtb6ZH9Wrs/hmv8EnpXklCRHApcCnx2mLEnSEMbek6+qR5P8OfAF4AhgZ1XdOVhlkqTeet1Pvqo+B3xuoFoWM/ahngnWWp/sz+RrrU/2Z4XG/uBVkjT5vK2BJDVsokN+o9w2IcnJSb6cZH+SO5O8uWs/LslNSe7pfh/btSfJ+7t+3Z7kzDnPtbVb/54kW9erT10tRyT5epIbuvlTkuzuavt494E7SZ7Yzd/bLZ+e8xyXd+13J3nZ+vTksVqOSXJdkru6sTp7I49Rkrd277d9SXYledJGGqMkO5McTLJvTttg45Hk95Lc0T3m/UmyTn16d/eeuz3Jp5IcM2fZgq/9Ytm32Pguqaom8ofRh7n3AacCRwLfAE5f77oWqfUE4Mxu+qnAtxjd6uFdwPaufTvwzm76AuDzQICzgN1d+3HAt7vfx3bTx65jv/4S+Gfghm7+E8Cl3fRVwBu76T8DruqmLwU+3k2f3o3bE4FTuvE8Yh37cw3whm76SOCYjTpGjC48/A7w5Dlj89qNNEbAi4AzgX1z2gYbD+BrwNndYz4PvHyd+vRSYFM3/c45fVrwtWeJ7FtsfJesaa3fnKt4sc4GvjBn/nLg8vWua4W1fwY4j9GFXyd0bScAd3fTHwK2zFn/7m75FuBDc9oft94a9+Ek4GbgHOCG7g/loTlv1sfGh9EZVmd305u69TJ/zOautw79eRqjUMy89g05RvzyivPjutf8BuBlG22MgOl5gTjIeHTL7prT/rj11rJP85b9CXBtN73ga88i2bfU3+BSP5N8uGZD3jah+zf4DGA38IyqehCg+/30brXF+jZJfX4f8DbgF938bwD/XVWPdvNza3us7m75j7r1J6k/pwKzwEe7Q1AfTnI0G3SMqup7wHuAB4AHGb3me9nYYwTDjceJ3fT89vX2ekb/VcDq+7TU3+CiJjnkFzp+NtGnAiV5CvBJ4C1V9eOlVl2grZZoX1NJLgIOVtXeuc0LrFrLLJuI/nQ2Mfo3+oNVdQbwU0aHAxYz0X3qjlVfwujf/GcCRzO6I+x8G2mMlrLa+ieuX0muAB4Frj3UtMBqg/dpkkN+RbdNmBRJnsAo4K+tquu75h8mOaFbfgJwsGtfrG+T0ucXAhcnuZ/R3UXPYbRnf0ySQ9dWzK3tsbq75b8OPMzk9IeulgNVtbubv45R6G/UMXoJ8J2qmq2qnwPXAy9gY48RDDceB7rp+e3rovtA+CLgldUda2H1fXqIxcd3UZMc8hvmtgndp/YfAfZX1XvnLPoscOjT/q2MjtUfan9Nd8bAWcCPun9NvwC8NMmx3Z7aS7u2NVVVl1fVSVU1zeh1/1JVvRL4MvCKbrX5/TnUz1d061fXfml3ZscpwLMYfRi25qrqB8B3k5zWNZ0LfJMNOkaMDtOcleSo7v13qD8bdow6g4xHt+wnSc7qXp/XzHmuNZXkfODtwMVV9bM5ixZ77RfMvm68Fhvfxa3VByxjfoBxAaMzVe4Drljvepao8/cZ/dt0O3Bb93MBo2NoNwP3dL+P69YPoy9cuQ+4A5iZ81yvB+7tfl43AX17Mb88u+bU7k14L/AvwBO79id18/d2y0+d8/grun7ezRqc3bBMX34X2NON06cZnY2xYccI+FvgLmAf8I+MztLYMGME7GL0ecLPGe29XjbkeAAz3WtzH/AB5n3ovoZ9upfRMfZD2XDVcq89i2TfYuO71I9XvEpSwyb5cI0kqSdDXpIaZshLUsMMeUlqmCEvSQ0z5CWpYYa8JDXMkJekhv0fYi+T8jp0zFYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "dep_ene =  t[\"volumes_dTotalEnergyDep_keV\"].array()\n",
    "h = plt.hist(dep_ene.flatten(), 40)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example on how to read histogram"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[b'hB8;1',\n",
       " b'hhep;1',\n",
       " b'hDSN;1',\n",
       " b'hatm;1',\n",
       " b'hPP;1',\n",
       " b'h2vBB;1',\n",
       " b'hKr85;1',\n",
       " b'hRn222;1',\n",
       " b'hRn220;1',\n",
       " b'hDetER;1',\n",
       " b'hDetNR;1']"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "h_file = uproot.open(\"data/All_BG.root\")\n",
    "h_file.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEWCAYAAAB8LwAVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztvXmUZNdd5/n9RURG5J5ZWYska6FkSXjBbm81JS/dXrGxuy2YdkNjwxyg8UHdPW16mzMezzBtuqH7HAM9LD24F40xGmhGNtACLxhsjrzR2JYtgWXLlo0lS1glqVR77hnrnT8iMu/3/vK9l5FZGZWRVd/POXXqxXv33Xdjyfu797daCAFCCCGEp7TXAxBCCDGcSEAIIYTIRAJCCCFEJhIQQgghMpGAEEIIkYkEhBBCiEwkIIQQQmQiASFEH5jZY2b2vfT6rWZ23sxeldP+tWb2F2a2YGbfNrPb6dqrzaxjZkv078cvxfsQYjtU9noAQuw3epP5LwP4OyGEz7lrFQAG4A8AvBPAHQCOAfiUmd0bQnig1/TJEMJ1l3DYQmwb7SCE2Aa9ncD/BeD7QgifM7OjZhbM7O1m9h0AnwQwB2AawG+HLl8C8BCA5+7dyIXYPhIQQvTPPwbw8wBeF0K4z117FYDnoCs4ngZwF4B/YGZlM3sZgO8C8N+p/REze9rMHjWzXzGziUvxBoTYDhIQQvTP6wF8AcBXM6796xDCcghhtff6LgDvBlAH8GcAfiaE8Hjv2jcAvBDANQBeC+Al6KqshBgqJCCE6J9/BOC7AbzPzMxdW5/8YWbPBvBBAD8GoArgewC808z+DgCEEE6GEL4eQuiEEB5F11bxg5fiDQixHSQghOifUwBeB+BvAfiP7hqnRX4egG+GED7eEwLfBPBHAN6U029A17AtxFAhASHENgghPImuWuiNZvYrOc3+EsAtPVdXM7ObALwZwAPAhpvrDb1r1wN4D4APXYrxC7Ed5OYqxDYJITxuZq8F8FkAV2dcf8TMfhLAf0DXOD0P4HcA/EavyYt7rw8AOAvgDwH8H5dg6EJsC1PBICGEEFlIxSSEECITCQghhBCZSEAIIYTIRAJCCCFEJvvai+nQoUPh6NGjez0MIYTYV9x///1nQgiHt2q3rwXE0aNHcd99PiWOEEKIIszsr/tpt68FhBBC7FdeX/qh5PWfdn5vV9sVte0XCQghhNhF+p3Q+528d6NdlvDoBwkIIYQYIEWT88Wu8LfzXH7W5lyT2UhACCHEENPvjmQQSEAIIcQA4QndT/b8+mJVUbJBCCGuaAYxCe42PJ433fJOvOmWd267j+3YDPq1Pezkc5KAEEJcduylIOFnV266MblWuelG/PG3fvGSjGM3kIAQQlwW7NRTp6iPi3U9bT3yaN/P2m28cNrJTkYCQghxWXIpVU95AiNLCOxkou7n2V4g7Ab7uh7EsWPHgiKphRCXgn7dVYtUTJ6iHUYWvr9+7/f3/cnDv3R/COHYlvdtY2wDxcxeDeDnAXwNwAdCCJ/e0wEJITYYtE5/N1w5t5qY90L3v9UEvj7Ofif67QoUvm8nO4yBCggzez+6tXhPhRCeR+ffCODXAJQBvC+E8B50C7cvARgFcGKQ4xJCDI6deOAUCQSvksmb6IvcSfPw7XgS3Woy3pFXUM4k7c+3D05uHJfPLiXX6jfMxWufun/j2I4/P2kXvvjV5PVOhMugdxB3Avh1AL+1fsLMygDeC+D16AqCL5nZhwH8WQjhM2Z2FYBfBvCjAx6bEKJPtjMZbidX0MWSJzz8GF5f+qG+XD79JLodgdEPRf0z5bNLsa1rU/vOuY3j5Tcf3zgePbWatLPjz8fakTEAwMRDp9MHPNzfeAcqIEIInzWzo+70cQAPhxC+DQBm9gEAPxBC+Hrv+nkAtbw+zex2ALcDwA033LDbQxZi2+xlpGs/FKmHtnPtYpPEFak4tiM8ilRJO905ZLGVSiZP4Gx13/ruwO8MGBYkfOx3CZsmfmL5OYczjwEMh4DI4VoAj9PrEwBuNbO3APg+ALPo7joyCSHcAeAOoGukHuA4hbjs2Q1h1m8fu2UD8IFo/YyFJ/P2a16ycVxZbiTteNJmNQ8A2ME4OXv1TR5ZO4b1Z/j+Ex6JhzxeVikB6Q6i9tEvJtcm6LjwWQXshYDIyhIVQgh3A7j7Ug9GiItlWHYMeSv+nY5vO9lBLzqlwzZUOf14CRXtEJJJ1q3I8/T76/gV/DpveOnPxRc0GVeQv5so2kGc/bGXxe7+/OTGcds9vzEVi4K233Jrcm387ns3jhePvzR9wL3oi70QECcAXE+vrwPw5B6MQ4grgt3wQNpKJdOvaqdfQ3K/7ViYVG668aJjDIqMvnn4dut6fwCouZ1GnkDzKqCRpc7GMa/+/c6l8+woSFhYAMAEjctf65e9EBBfAnCLmd0I4AkAbwXwI9vpwMxuA3DbzTffPIDhCbE/udhkb9thO/74F+uyyn14AZD33CKvICBOtK2JanK+nNFmHRYE3H/d9VFZbueOg3coQBRIo258vLtI7iF1EwDMfiO28zsSFjqTT7awEwbt5noXgFcDOGRmJwD8bAjhN8zsHQA+ju738f4Qwte2028I4SMAPnLs2LGf2u0xC7HOsBuf+2UQMQytRx7t29C9ft2vtO348/vW5edFC1duunFj8m/T+U39PoK+YH0/ewsBwPmbola/SrsEFghZ9yXjJZsHj5F3HUDX13+dpeujAGJ1E5AKj3PPuzq51hqN2vzSzuSDIqmF2K/0O/Fvpf7ZzRV+P8/Lou1WxnmqHj/xs8DgFbN3+WS4D/9cnsA37S5oTCuk7584sZK04/sWbkz7mHsw2+5w8uVT6Tioy1KrvzmaBQIAdOjRpdQWjwfe+y/3VyS1EPuF/ZByejsTeFGw2DpbGY4TI60jT5fPEzWvurcTb5DYIJz6hklUNgW7hPQedyLPfnDdeHpfI07o46fS3QULj7PPi9784093kMfqoWg/sLQ7VJfis1rpMJJdw+j5/P6LkIAQYpsMQhjsNIhsJ2PpN2rZG4A9626rry/90KaV/fozWHAUrf6LzrfoviLDMa/wywU2iCTQzBmHWV1UP5BOj7XzccZlD6E6uZoCQLsaV/L12dQ4XGqVkcXagbRdk0wStQvxuONm7LW5fDVSfZbHtH+M1BeNjNRCFLNV0Zp+YhL6dT3NmuhZMHziC+/eGAdTJIDaOa6iLSdkFt8W3TerR+JE3XZ2AVb1zD0Yzy89I50CeQLvVPLrNpdJKKwcSSf92oW4WvdqH6ZOG57pR1M1UmOGhUw873cJTHCyZ5R2QPUD+fcVIRuEEANkO2mg864VTdR5QWNbrf7z+st77nbu9av8In//dYoCudiA6901px9ZjuMh9Y1f/U8/cGrjmHcNPhL5/PGrNo4PfPHp5NpT33fNxjGrbPyYePU/8VSq2lm+Jrbl3YBf/Tcpyq26GI/9DiLQo73wCAXL/2/9jGwQQgwt21EN5QmFonKWlZtuzJ3AeTK2g8/fWOGv95nVDkgnet93keqnKLXEej/tG/KjhTkArLqYTrgL5FnEwsILCBYKPKHP/+1rkna1+bhgXnjBkeRanh6/Xct/ff7ZqfAY4cl+JB771X9rMo6jfigeT/512t/S0bhTGj2ddtIcu/jFvwSEEBeJN9DyhJtH1sRelGhuq3uB/o273haQF5lcPru0MSb/TG6XF2MApIKjqN3K1XG2TP15Unw6iQapmFhYbBIkN8SprrLW38Q5/8x0wmXPIt4l+NV/uZ7fZ2Mm+7wXEEypEdVNy9elY68sskosvY8N2u1JGamF6ItBxzdcbNbSIvXTVoFiWWqqnXoF2fHnbwi/lvMQYkPvpp0BCQ+OHRg/2UyaLVE6icnHox/m6f/55Um7I/fFZfcFugdwEce1OJGyQADSCbgxHdt5r6Cla+O1UefF1JiOx0bzbSl9W4nwCM42zJN4Zyx2Uqqntor2aBQElWVKp1FJBUSg181xd61Kr9377BcJCDF0DHuA2nazo25178VeY9YFCAeQZZFnF8iNTHZJ7fpN/sZqH17hF/U/91B6je/zq+Sla+PMz6qd0XPpZFmfJaFAkztP+v7a4tH0GquH2hTJ1nRbnrT/dBzVC3EczXI87qThEiivRKHQnKHZfSyd6UsL+VN4qMSBmJdUfbIvBYS8mMR22Ou4hX6FR9Y9QL5doMhtdKuCM62J7JgBVqf7RWeeUPHpI1oTcdJmfX9nJF0lc5qIc8+L4/WeP+zrv3R9eu0Zn436nNMviqNfucqvyKm/eTqfBjAnBOeAxJO4kVrJG4c7tPo3t7tozNE1Vk1ZKkiah+lih1xZR1JVUWeUXlfTa0ZvoDzrIuX6ZF8KCKXauLzZjcjenfbTDxeT6K4fT6UiWCj0m+fHZwBlAy779gPA2VfEdA2lZpy0vEcPe/uwZw4ATD8WRQsLCB8R/PStcek9+QTlL6qlMzMHio268gdP/s0oFCoUPO1VR6BH8y6kXXU6/ZX47OCuseDq1OiaM2kkwmMu/XxtjT4Po/68DYL6tHJ8EZbdlM07irYTrNNReDZXdzbV70sBIcR22Gk5y36zjW6HvIhjXzM4z2XV35OwRarqdfJcPoE0109RwRlWHa0dSHUszUnSmdOWZOyMW/1W4qzI9yxfnU50ZVr8Np1mi1f1rcTlM521S83YZ4s0Xe2xfL39yIVU8DUOZNsMvDtpe5zepxMe3D8LJ3M2ZOOdQpvG4XYJI6NRAI2MpMKo0YgDG5sssJwXIAEhLgt2a7ewmxXJPHnZQLcir623A7Bqp0LP8rr/Kc5F5Fw5n35ddPuc/VacVC7ckvpyNqZinxMn0+U6T/ZsHF44mi6TeSVfpzIxXmXDQsZ7+5TXqN0Yn3cqJs5LxBdKbgZvxquNuXQyDtS2fZDe85rT79PL0op7z+PxvlItHgevz6LXNkKurBOpqqi+Gr2/RkaSS2g347PrnSvIBiHEdihSP/W7Y9htiuIIvK6fBYHP+slwmogil0+OPuadAACM3/3VzHZ+0q7QBNwac0niciKQvbcPT/zsGurU8anu36mO2PbanIrvc2QxnRDZ26fFu4aOG/tIvFbyz5qkEyxYJrwLKe0u3A6lfCHO4mGUhEXZ9cEqJnorLBAAoDYWP9S1lWrutXZbAkKIDXZrJ7C+ei8yDjNFUcStAvtBkRciRwhznh9OCgekaaHZ/dPD7qY+Y+k5cjHlnYH3Cho9F6+de3b+NMIxApwyAkgnfu/FwzRotV4964LBZkhNQ0NsTXjdDj2L7Aclt/rnXULnQCrRAu0uWKzYUvr+jWwB7clU7cN9Vkbyv/VAK/5qjYTAcvpBlWi8s7NpVtk1UjGVvQDqk30pIOTFdOWyndQV62SpaLLa7yRFhL/XCw+fhn8jvbS3EXC9ggL1EwsFFhZAOmlzfqDKqjPE0k7Dq47Y24dVQk3nobo2F6eOoloDfN8mPTu/ZrW9UyPVKEK45VbkHdbVc4yA28RYnQ0UtMJ3KqYwRW9m0U2PZDMojVG7WjrRk+0ZFbcdGiGhUCpRuo6Gi9ugt2XUx7QTAtVy7G+tme4uqvSsRrMgEq+AfSkg5MU0XOxG3MJOJv6dkpenyK/nEjWQKziTl1rCF8Hx9QbyznMcwMJ1+RnZFm+If+iHHkgNj5yQzqeZZjhqedO4SF6wV1DFlVdokabLu9iz0bZJKSOsZa5dvpGWSVJSuECxEvdJQqBzMNXVh2r2A3wqOqPYgTDpLpIVvEO7ieqE22mQ2qrdTD+c6ngULI0meZPVnJRl9RsJEhYIQCqMJkfT30O7wx5TO0u7sS8FhLj0DFPwWtazi9JTFEUS9xvwVZRrqEjllOdN5FNEL10fx8HupUAaP8B1A3y+Idb9F6WW4HarqXNSMtknieC86YMNsc7Fvk0Tepuihc25YfLrQLr/9oSbBHn1P5leK52OD2vNxknWnG2hTPp+nrTLLq6ABcbIWDppN9fiZzpOXkF+8h2pxGctrxXozojxavohNsnvlQVEx22NRgoka5sM3dO1tdx2RUhAiItmJ8Ki3+Lyvu+8BHV+oi+qmcyTM+f2ySqJuXE/6ep9VDH35yuZsVG5yH7AQsEHlCXpGQpSULMg4JQRa4fdpLKIXFhdxEKh6SKCk92AcynlNBE8UftFbJvUOeztYy23JSEhw7sEAGhPx8mYXUM5dgBIjcC10fj9ecNupUrxGG61XiL1Fht9O+6NsUoodPK/r4mxKGSqlfRZY2TRZ9XR6Ei6W3nqXEzudGQ2/WJHSLCsNPsTVB4JCHHJyEsKNwjaByf7rne8jq8MxkVhajTeJ25LM4CyqufMC9OYAC4Wz/2d/JdpviFW7XCkL5BO2uyRVFRkhr19WqNJsyRCmI28AFBepeAtWtV7f/4O+/A7bRYLhQ4FcvEuAQCsQe6w1XwVSKKPr/kVc7yvNkn1nr13Er1utaMwqo6mu4RRMgh7dc75xfj74JiDsneVJUacIXo0uY92Bs7Ntd2i3cpIfF9+B3H0cEwY1XKurIcmor3izFpBIYkCJCDEnsBF74F8W0PRToN3BZwuGigOFOMaAEn5SRdVzDsD1u7OPpy24x1FxakU1ubiZDT/01EoeM3A9GPxhI9MBmkHVq7Kd1dkOwF7DI2d3tx2456VdMJZvYpWvzShb9Lh885g0vmvnq9mtoObcDGerR6pjrn+qIuiFXkl2SWk3xF79CRGZDemOrVbC6mtZmI8/gp4cueVOgCs0eR+7cz55NqZ5Wi1H63E9+nVVPx6tBzfy/l6qutjoTBaSd/zUit+D81Nodr9IQEhLhlbVTnbLlzzwNcWbrl2DOf5Z6HgVUcVChHgncGm/P9VWlk691IOGuNaA43JdKJbvD5/4m9RNTA2APvYqmQipflg1ReSIfWLjyTOVWm71T8omVyn7iafnGCwETdpt6l6W+LK6dQ+LDDKlXSA7ObJk+rSavol8TV2+fSr/1o1jtEbsEcoMGJuPPUmYiZo0l52qh2+r0IfdstZ+nlHcY5W/5PVem47vwth4XFwNH+8RexLASE31/3BTvMNATElRXAG5rzqat4QXVSvIC/VRCWNGcvNPsoTPQAsX00upc4gnJda2gsZntC9yob1/+zR4wPPkmuUCqI15XTp7FkzkU7agdQ+HFUcnOAoz8WJquNcKEuVbKHgV8k8GfPK/cCB9Itos8rKqVHYIMwqIe/R0yKbAe8gvFpmgozFIz5SjlhtxQ+7Vk4/Q36f/loeVfcsnuyrtbg1HK+kX3qDdgYrLefmSp/H0pVkg5Cb6/4ny7Moq9C9dwc1VucgW1gA6cTv4wWmKb6B4wXW5lJBwqv8JOBrJu2Pcwx5FRBP9qz22ZTeeS3/GqeI5vTUHae375AaKLAxeyyd3TvsCbScv/rnXUPZuWEG8kCamEkN86zv5xQS7YLJnXX/rYKoXz/x10md06T+feoK9uJhdUvVqWV4jDPO84cn4Eop33uoRHYR366Ro+opmty5P7/TYGE0XbC7qBT5ERewLwWEuDzZUDn16XoqhBgsEhBXALsRw1B0z076r9x0Y5Kemm0I7FLqbQt5cQub1FfkvrpCNY2B/OLzm9xLaSXPKqCOizPjXYOPOE7iCvivraDWQGPWhxzHwySLqLMLsCdQoBVjyXnSBHIjrRxMV8kdWr1XaddQX0nfdJncQWsuiyjvANhw7FfTvPrl1b83+jJeZTNCK+0meSeNOlUMP5t3BmvtdApkFc5CI9UD8op8jtQ+Fxqpmxgbi70Ka7Ya71tqxv4nR1L7V15gm7czFF3j9zzap6rLIwFxBXKpgt58v6w68jmAWEVU5piFoohjMjZPnEiNcCsviEKhPpv+kc49FLfi557jjAFEdSH+ka7Nsftn2o4TzfkCNJwTqMQJ7g44vfVqnNz8xA82KpMQs2Y6IXD6BxYCPJkDQKvAE2iEhMI4xQs06ulUUSX7gZ+YajRBtgsqma1SJPFkLT6raHJfdaqYAzRR82TsJ34eY9FkyX34SXu8El/zBD5bdUKWPw/39hdImBwZo1KqjXw31KIoaL42Xk7Hu0YRi16I9YsEhMhkuykusoSM91Jig3MZzuDs7l2/VnE7CF8UZx2uSAak+Yd8LqJTL4lCgVWzfi5bPZId5OUFRFp+0q2SqeA8VxrzLkhJTiEXccwupkbHpZl04ucAsNXl+B59Ggf2zfe6+irtBngncOhAmpeKPXzMyRt2AR0rkyunC6aYGc1e/Xt4wi1bvk6fV+ezblJtkMAo8vxJV92poGrlGIS9IKmSAKp4W8hIFCY1Mmx5I3WVrq1QrnI/psRzq+VSstNn6sfYLxIQVyDDUON5UzZTjkdwQiGPdNeQrsDYvZST2HVf00qeqoT5lNbVhXjMxuFNxWhY7eMMx+28kpB+EU87CHM5/zmYq0Hpnssu+naUVvXsvlmUydO7cs6MkfqFIni9MZc9fBputc6umLy69p40PBmXaAJrOGHBwWHei4f7WCG/f17tA26lXYnjO+DiOxZpkuUVOADMVuPvjSdqH7zGgqTiJn42MrdIb+mNzxV+z7QzaLii3KO0fZ2ouXQd1OdSUzuIoWOvayFfStZ3C4ldYa8GI4TYFSQgLhO2I4zyUl7kZTldh4vrMP0GvOU9CwAqOcbnBVfxjAPbiorn8EKLo5kBoD5LqiNaQDZcvQKmNUnGbJdRtMM7Cpe9M9AOwjioywV8lSkYrOUCz1gNNDu3TOfTMSbqF3r/Yy5/D6du8Dr9MVqhF/nws0qkaFXPeBfPPHWO17lP0orfr+p5hc4xDE+vpSlP5mrxc+PV9Ew1ddFld1Cv9mGmSFW06Fbn/F78zoCZb8Tf7/SIy+FF75N3Gl7FxDuKlU2fb1GVkf6QgBggl+tuoQg2InuPI57QfQjaMl2rkAHbG58vPDv2yfEIPtvowa/Hye3US1xBF5r32H6wKfCMo5b579wbkdl+MOaC0miyYwHhi8XwtZpLQcF68iZFHxdlAGXKBV5Bh8fSoDQ27vJk6dM4sIdPkdGXJ7Qj/rMhmwTr2b2xlSdZVg8B6QRZp+PDo64qX2KkjmPiSRoAJmnin3DPWqUxJvaIdqo6O0CqqJKLMMzLvtopiG+o0A/WC0jGCw+2SXj7T7/sSwGhSOq9YX2n4Cf+nWV5EUIMO/tSQCiSevfwqp68Qj39psngKmxnX3F1cu3gn5/cOD5//KrcPrjimc+wykZldi+tO+em898df9qb0lqQSojLT3rVUXmVInOpRoG51X+J7uv4mge0Mxgdo6yczj+evYdqbrXOsQSbitsTc+PRbZIjjH1EcF40L5DuBthYXHVl445wLWS3qmWVEK/+J9zOgI277O206vqbqqzRNZenicbFfYy4fCVNUtOcp/69MXuhGXcUfkXO4yhS+/BqvenK483Tqp53Q34HwTsUVmd5Fvo0Pnu3337ZlwJCbKZIndVP1lSPr43AVc682ofhFBeTjzdyr3HaagBYekb8KZ5/NqkvzuSPce0gjdepW1eu4YR03vWUEs1xniKn+25RaunadPxj5sIxQOpeyp5EQFobgCeOItfQ2dFUH51M1D4jKsGTezLRO4EwSZNikU6bJ08vIFZoojtcyy8wwRN6yalXmjnqoVmnj2+TIJmu5E+W3L8XJOxSOkUTc829LxZoPj0Fj5EndP/ZsOAryufEffjviPtcbmWrtoD0e/bqt3P1qMgt+t0UIQFxmbLT2guv/P5fAgDszClOCHE5IQExhPQb6bxbEdG+/GUW7Vq+aoNtEo2puCq6cHP685r6DqmOrnYeF5QFlYPL1g6lz+KdeP0qLjGZtuMylcEXmTlAOxvyOqq4dNSB0z9QzIFPycHF4SdrzojKqh7aGXhf/6Lkb+wlVBTIxav/PDWPv8+vSA/WUuNuVt/dMeWX3GSjOq/42wXlMjkdddOp3/L9gNK0HHzfmFdn5XgT+THxyv1wNf0szlEelQ4HTm5K1R2PvaH7fE7EtN+FsNqOP18fc3GGcrlXSqnBnWNQFEktEnjXwLEJZWdzyBMO3oW0XI8/0qdvnfLNMxk9l04cC0fjZNSYSduGMmUsnaVMlg33xzeXPbn7jKVcxGZi3HnFkCdQhewCPtAoj9nxVAXCbqOb0jGPxGfN1ci7xXmVLCW6aRcoR5MHT7gLrlTcVDJpUxCaT7Ndyvc6YjUK6/EnR+q57Xx/9RzV0ZgTaHxfnb2nClJfe9XRKjX1kydToffCz11q56fBfmot/ZHypM3R0gtu8i0S9vxe+HvhID8g/Q34awxncPV5n06vRoFmBek6ipCAGAK2sgusXx90mU4hhGAkIPYxRQKDvZPYiccnv+N6CM2JuLqZeTRdFc4/syBXDjmdcFbS5evSdmUuiXnEZRjlWALevk+6FA+UkK40mb9S49XZ5Gj+zsAXgWc4UGyyIBUExwhMOpUCGyl5J7ApwV2Zjaj5hlhe/fqdQTkxiFNJTGfBX6EvabyUfjYzVLeUE+3VnYppklzI2k4JlHgT0Wq67fNP5Xhk+Wfx+L2hm591gTyQvDGbdx7tgtX0mbX8VPN58R5FKjafH4nhWJKiPhap3ZSr+cDXfLzLQj3ubLzqs18kIPYIzmy6W7z2de/ZOFZsghDiYpGAGALs+PPxiS+8O/Maq598zIKvtsawbeE87RJ8ZlNWH3M08sINLl0y7RJ8UjtOac1Ryp2rXdnHFepz1EUSk/64MhpXheNj6QqXYwJYrzrtfP159etTRnAKhdlaak9gvM489p2ufCfJIOpXuGm7/FX3qA/jzqFW0K7Zyf5z9rr5A6VoC/HZUZnETuJsvHzfCPJdKJsdsoWU8lfJ/vNg2NbScnEFi2RP4h2U34Us0w5iAo3Me4Bid9C8tNtFq3/P+bX4R3aAnBb4PAC0SlRjpKCsKMfM8I4BSCPtF9ZkpB46imwLPs7A5zPKy3skhBCXin0pIC6HVBvhi1/dOC6yJbAg8VHF7HraGnN57WmnwMnpGtNpO3biqB+M96w+w1Uro9s64y6SmHYDfFdwUcXladoZjKe7C7YFLK1FveqIT2lNK6bZgmLunZx0yUB+Ejqv+8/T6Xe8zp30+F7fn6y0C3TpPuJxIur0AAAgAElEQVSW4fvqHXZ/TPuYq0S3zDOt6Gk2U053SYuUgMr3kfesFVcoO+99+f75mve6YvdQdkst+iyKdm/nmvHvo1QwtZ1tRO+ezTUf4pjOrKYZw4p+b0xig3ApzlfoNe+AR9zOZZF2A9xuze0guI6GT954fil+Hp2CIlFF7EsBMcypNrZbaGerfvxOQwghLhX7UkDsJzh4jdVI3p7A8LWFm+IqZvxkumpZuDGuRvzOgHcNbFtoO1Uk2w/aUxR4VnY7CNoNlJz9YGo6rqw42+hYLR0vrzrHqvn5a+bGo458U61ease+4rMjafoPXp36YCXWT7NdwMN2AV5pj1q6A+EVr1/9TpbXMq/5VTeP118boedxISD/LNbj866haKfRDOkUsEav2Y7j9epL9EM6OJKm2kjacvXNTbr6eMxBaD5orJGT4gJIdyUcL+DHW83xJjuzlu4S+L4j42mg3AVa1S/TTsB7D3FMR8t5vPEOZbWZP97lOqcoIW1AM/2+Au0MRt3f29oieVCV8j23ipCA2AHb2SWsp64ANqe4FkKIYUYCYoD4KGVOVsfCYlNmU3q9fA2VDbw29anmhZZPScGL0NYY7RJ87ADtBriIzfhUurLmVcyBiXS1zkVmOEbA62m5D+9Tzjn0OZK27CKOecV4mNJCTDs9O6/+N6d4oLQZXGuglP+eWQfvPYlqiK9Hza3iApUIJTvGVDm1dzRL8XNbc0WveRxTpfg+L7RTmxT3P0/XvBcTe8X4VT1/NkW2AF4ln8VU7jW21yy7mIA0ypqKGLlsrj7NB5NXH8F7+5xvx230FEWFn19NP0O2hfk62az/X2nE40pBSVe/M+Ad9tJK/DzaLZd6hZI+1pdot1LL97JaOJMuQUtVsg22ixKW5CMBQRTlNtot24IQQuwXJCB2wFYpL9ZdVDlwDQDqB+LH7UtpMovXUwI5WhgvX5eu9tq0M/AFo8JMXIFwXMGYq2VQo5VKUbUythlMubw8bCfg1f+cswuw91Dd5aefpje6RAUcrq7NJ+2eQcc+CpgpWgnzSrsoHTPvQvgeD9tWNsUV0PcyU46fRz24VTLr/p2XVJv64J1M3e00msYxAXEc55tpdDBHY7fd7yZv17Dkimpwmc2xTrpr4h3AheW4cvf5nNhmsGBsS0jfF+9EfbwAw0kJve6fPYHmQ3zW0lp+nqPzi+nuokIedSGp8pd+iGsrsc/qmMvNtRLfG9sPOEOAb1euxu+ytZx+Npw23lquFgnZIDYlrOyTK0JA7FbW037xgkEIIfYjV4SA6JciwcFRy7acv3JtUT3lc89JV13VJVrF3ERRpG4RE6hCWWuWcu/46GNaFdXcSmWc8g9x0fpNmT1zisx4GwF7DPlUyrWcql5e/8oeQzO1/KJDB2gVO1VKdfWJ/aAgeyd7O7He3sO7hHmn02/npIj2Uc/Lnfx8O2yT4P5X3D3s7eQ9dTiugFf4qy42oZLYVvI/G28LYNgucLoedx7Vgp2WTyf+1PL0xjH79/u4Ao5+X6Vdw4WV/F1C22UsbTTieKu0G245nXu7FV+3G+RN5lK8nzodxx5cJuHSgTje+hKtzv3KnWx59fnUbTCpRki3tZ5O33MgL0Ieh49G54/eB3p3KL9ZZX5nyXd2ZrkQQghx2XNZ7iD6NSgXtfMeSJflByWEEAVcMfPexdodWhO8RYtb+wu3pNt19shbuSbtY2GKKkOR4alUTzdy4VBUxRgFqI1Ppka+pI6xK2LDrqfTZBz0qQVYJTRNxdErLliLDcIzlfygNFaBTDpXTlZ7HK6kwVWNkJ10zcOpFljV4wO+mAnnvnq2FVUnZcpP4I3e7EbKRmpvyOXXSy4SkYPU+HPyrqysRvLqITZaswHYp65gVd88pb726iFf/5g5V4/vmdWMT5DaCABmqFrZY+fnkmu1kWxVT3Dv68x8/B7GSCXqg8Ga9fzvtkNFoxoLpHLzKjZSzRgdN844dRapdsypmJqPx/GOrJKRetqlV6GklO2p1DhcPRG/v/Zo/u98/ClyD76WxuTe18hyfO01oon6aWclqa8cASE3VSGE2B6XpYDwu4V+hQMborm2MpDWIV54CZVlTBfCWPquuGLwrmU2EVdWIxTEMjuZrshnx8jlkwqCXD2Rn9Lg8GiaFoDTCYzTLuFAZTlpl5dCYsSlk+DALm+8zXMH9SkeeAU95QLb2jm7kE0pKZC9CveG6MVOXBluCiijcZ1rRcNpUdrm+U5s591LeYx+Z8DX2FXUu41yQODJtXS1zkGEvDPwieA4keFZCgDzRuQKvf+T8+mzZifi53iaXEXbLtnbqfkYHNdxSRnXVimFBBlwzaf7phV/czz+TdlC+hlyckhruveyQOlA6Da/Yu7UKDHePN3jZkDeHHuPZd4c8t+9OZda3ohWnXG4uhCPm5Pxs/EZX7iP0dOU5sblmaQMJRg7nV6bOBk/hKKCX0XISC2EECKTy3IHIXWSEEJcPJelgOiXxbe9NHm9clX+1pOpHyQ/59Q+hzAT1S0zc6k6p0aRmAfG8jOWsgrgBQeezB0H1yEYd3tUzubJah/vw8/5h9g33+clYn9+Xw+BDczTpOrx8QFsLPWqGI53YBWTVw95g/M6rFLy93mjL+dHWiFVz4jbv7N6aJ7K6HHtYwCYpBiRFZcbiCOOOW/QQtPl1Upy/qc/vhNLM3214xiBogyga2txHObUak9zHADFDnhdg81Tn84ezCqcMEmOGZX0WSOLlMH2capr4J7VHosnvOpo9Ew8rs/G48kn0mctUXYCVt941c70Y1QPY9YZqalme22eaoU8nj6L5xHnz4ESRV2PnSFjuVNpj53JtiovX52qikZIszz7rfTNsAPN3Deya6Bsxb5UMZnZbWZ2x/z8/NaNhRBC7Ih9KSBCCB8JIdw+MzOzdWMhhBA74rJUMXkvJs6NxPEMnEobAFav4hQXbotHHknVifxUGwenolrJFxzhZGWHqvHaiCs4wvED7Jvv1T58zauO2KunqDA9xw9MVy5sHJdc9r+itBZr5OHDaiWfJuNCJ1UXMafbUbXB6if/3BVKE8EqoHPtNI1DUsTH6Qv52om1qJeYdiVHz1FpSk457dVDrBL0HkNcIGaeEs01nOcLl4RsutTPnF6h2STPIteutUyeQB32j0+/y9Jqmdohvdak+6j7ynL6PVRp8+774AJVI4vkqdNI++CPin++k0+kf3tcWre6mD5s5Ugc5Myj8b7KsktKeSF+NrXzlDTyQPrbmDjBOiGXrG+VvY7oO5lMv/Op78Rn89iB1BuS38vUXV9I2rFHZflT98ex03kAqFDaHy4uBgCHvuxcLHfAvtxBCCGEGDyXzQ6CPZc21XGeyE/pK4QQIpvLRkAUsXJ1VFmc+x4KyHJh8KUjUcUwVk2t/lxf+abZsxvHPnXFXDWqmLzKgj2NWP2y2bOoQe2iqqhq+Z4Ivk4yw0FvPsiN03UsdqLHjVdLsTJjLfhgsOxns9qo20fs03snsdqHvalWXMbSvAyrT9fTZ3FWUq/24XQSjSQYLG3HNSvYK6gw26jLIrp6gT5TyvkfXOBZUvFr1QU1UbBZkjXDBajVSAXC3j7+I+MMKBNP5ad7WJuL/XlvnMoaqVgm0nHMkSdQkTqn3KBKh1VSsTmVTamZP8YD7/98HFNhnZbobsjjSFVKwNqRsdxrnKm5UpDRuXw2VS2v14cB8l3w/aK29p1z8bnc96fuz61zf+DsUvIsZqeu/7kCwsze0sf9ayGEj+3oyUIIIYaaoh3E/wPgQ9jk5ZzwSgASEEIIcRlSJCD+OITwk0U3m9l/3eXx7BjeXp59dlpi8dStVDjjYFQB1SqpGoXz0HzX9Pnk2tWj0W3jAJXSnHE5hbispM83xCoWDvjy3j6cl4c9dcrOs2iBVDHVTdlXo9rqgsswypxtx8+Kg+F8kBurlTZlIqUgNc5TVKQe8rmIVqlqEmcp9eqhhQaVi6RcRD5TaJMylvqAMm67tBrHUXbeZGvL1D+rhFadqmSJAiyrzmOIju0ceRkhpXae1DmuvhEXlOLgqvqsK1RDw68uxHaHP3Myabf8nMPxuR/9YnKN09wf/oxL7sPtbogqG69u8SqWjWd9J33N45h+4FRm30DqxZNFP5maWcXCys1QcL9Xy5R38KzKTTcmaiDPen/+Wa2cZ72+9EM7Uhf5vsyK1v2RXC+mEML/tNXN/bQRQgixP9mxDSKEcPfuD0cIIcSwUKRiuq33/xEALwfwyd7r1wD4NIA9FxB/df+3N7ZbxR4MQgixP9nL5KO5AiKE8A8AwMw+CuC5IYSneq+vAfDeSzO8/mG95Zk3pDr9647ECOGbpqOL6uFaqiu9uhrb+epiB8ux7Vw5W8cKpO6hB0srue3YbdLbFtqkoeZEc75ewQSNcc1VVDvdzo5uLkq093jz4Maxd6ldbEc7g6+odqoZXUyfoloGFecqe6ER+1houER+VNz+3ErUEq80nB2D3EGT6GNXdSywC6izT5Qo0RwnkPOV/WoL8b4k/7+LHG7E0giYfdhFoFPkbH0mHvtoYXb59InbOGp3/O57N45X3nJrOl5yKWU3Sdb1A6ndwbtXTjyU2h3W3Sb9JMXfXuuRR5NrrD/3+vc8N8xBs5OKkjutQnkpnzVo+omkProuHHo8DeC7BzQeIYQQQ0I/gXKfNrOPA7gLXaP/WwF8aqCjEkIIsedsKSBCCO/oGaz/Vu/UHSGEPxjssIQQQuw1FkJ+GPuwc+zYsXDfffcBAF74jl/ZOH/j2x5O2h0/EHWkz6xGHevhykLSbopiE3zqilkudhPYRpCOaapAadegz7qJ/D5KpDI/2crPgLoSon5+c9xCvPZ0M6ZF9/YDjlXg+AaORQCA02tR0d7opH3M12PbZbIZLNdT+0GL0lrU59PxWivHL9udL9Hr6vmCD5vsBC4cAzUKcZl4KjZkXT+QpnyY/Ua0O3k///bBNO6G4bZsJ/O+/W2XpfOT97xr4zjPj97r/j39+vcPq/5bDA4zuz+EcGyrdlvaIMzsLWb2LTObN7MFM1s0s4Wt7hNCCLG/6ccG8YsAbgshPDTowQghhBge+vFielrCQQghrjz6iaS+z8w+COAPAWwo6RVJLYQQw0NWQN3F2pf6iaQGgBUAb6DXAUMQSc2sPCMev2jm8eTa905+feP4aCUmpKshNbbWjAKokCaka1LSuCkKWGu7ILeVEI3ZKy7nf5s2bBfIqNx0QW5cmpONz6dbU0k7Toy36BLycTK8RTI4n15LDapsYD69lG9sXV2Jlt62T1a3mB14Nnom3aBSvkN0rkk/t8oSfVZ024irmsjJ6kpNKgPrSlFyMJg3AG+V/G1jTBSdzwZhn0ityOjL18qPPLqjP9jdDi6TUVr0y5aR1EIIIa5M+vFieqaZfcTMTpvZKTP7kJkp8ZEQQlzm9GOk/v8A/C6AawA8A8DvAfjAIAclhBBi7+lHQFgI4bdDCK3ev/8KYCDRdWY2YWb3m9mbB9G/EEKI/uknDuJTZvYudHcNAcAPA/gjM5sDgBDCubwbzez9AN4M4FQI4Xl0/o0Afg3dIk3vCyG8p3fpf0N3t7JtGtfFKOhnjT6VXGPD9MFSzHJaD82kXYfk3lJIs7kud6IRdJGikRcLqqad66RG3yeaBzaOOYL5ZGM2aVdvx6+FI5qfXJlO2lWoAtqTC+m1TieOY/EcRWOvpYZ5zmDaGaHKe4vOgJ+fwDaJTB49z2uH1HA8/Uis5scF4Luv4/N8lTOGDcc+S+lnP/y/bhwnxmFX6H23KepbBmFxqRjEb60fAfHDvf//oTv/k+gKjGcW3HsngF8H8FvrJ8ysjG668NcDOAHgS2b2YXTVV18HkF8fUwghxCWjn2R9OzZIhxA+a2ZH3enjAB4OIXwbAMzsAwB+AMAkgAkAzwWwamYfCyG4zPtCCCEuFUWBci8OIfxF0c39tMngWgAcqHACwK0hhHf0+vwJAGfyhIOZ3Q7gdgC44YYbtvloIYQQ/VK0g/hNM3s1gJw0mwCA3wDwom0+M6u/DcV1COHOoptDCHcAuAPoZnPd5rOFEEL0SZGAmAFwP4oFxOmCa3mcAHA9vb4OwJM76CepSY07X1LcWFzRXKoU14NIdyDEXlEUSX10QM/8EoBbesF2T6Bboe5HLrbTkbFYv6Hkiga3qQ7D2U70pFlzWqwnyWPoidaR5Nq5dvRIOtGIef3PN9N6DafqMR3GuXp67cxK9KBaXInPql9IayMwRuk6ygupZ1EYie+rdjr1WC6Rg9bMajw++GAdeVSWo+fWwk0TybWpu74Qx/vm48m1Iq+j3LQT2NnEyX3UClJXFKXCEEL0Rz9xEDvGzO4C8HkAzzKzE2b29hBCC8A7AHwcwEMAfjeE8LVt9nubmd3RQmPrxkIIIXZEP26uOyaE8Lac8x8D8LGL6PcjAD4ybXM/tdM+hBBCFDPQHYQQQoj9Sz/J+v6umc3Q61kz+x8HOywhhBB7TT87iJ8NIWxk8g8hXADws4MbkhBCiGGgHxtElhAZqO1iK8zsNgC3jWFiy7ZCCCF2hoVQHGvWS7h3Ad38SQHATwM4EEL4iYGPbguOHTsW7rvvPgDA0d/6hY3z77z1T5J2Lxl9bOO4QYn2Hm8dTNp9ZSWGZzyylCaCe2w+JtprdmIfF06lVd5KlOSuspjK1rFTdHwmutiW6+l3MP1AbLjwguhuO373vSiCXTvf8NKfS6594gvvBtC/n778+YW4fDGz+0MIx7Zq14+K6acBNAB8EN1Mq6sA/snFDU8IIcSw00+yvmUA77oEYxFCCDFE5O4gzKxsZv/QzH7ezF7urv2fgx+aEEKIvaRIxfRfALwKwFkA/7eZ/TJde8tAR7UF65HU8/PzWzcWQgixI4oExPEQwo+EEH4VwK0AJs3sbjOroTiB38AJIXwkhHD7zMzM1o2FEELsiCIBsVETsleL+nYAXwbwSXSL+wghhLiMKRIQ9/VqR28QQvg5AL8J4OggByWEEGLv2TIOYpjhOIibfymaSF7xqjQ57PdMxnITi+2YZvvLF65L2j185tDG8doT6SapVI9atYNfjedHltKU4RMnVjaOWxNVeD55T9chrKg+waWqXSCEuDK56DgIM/sfzOxqev1jZvYhM/sPZjaXd58QQojLg628mBoAYGavBPAeAL8FYB69kp97hbyYhBBi8BQJiHII4Vzv+IcB3BFC+G8hhH8F4ObBDy0feTEJIcTgKRQQZrYeaf06dL2X1tnTZH1CCCEGT9FEfxeAz5jZGXTzL/0ZAJjZzeiqmYQQQlzG5AqIEMK/M7N7AFwD4BMhujuV0E3gJ4QQ4jKmUFUUQvhCxrm/GtxwhBBCDAuqSS2EECITCQghhBCZ7EtvpPWSozffHL1tpx+Jkc6fv/Fo0v6bU7Eq28mnZmM/i+nbHzsZ5eUt9ywm18pnlzaOW488mju2gP4in4vaKHJaCDEM7MsdhOIghBBi8OxLASGEEGLwSEAIIYTIRAJCCCFEJhIQQgghMpGAEEIIkYkEhBBCiEwkIIQQQmSyLwWECgYJIcTg2ZcCQoFyQggxePalgBBCCDF4JCCEEEJkIgEhhBAik32ZzTWLg3d8buO4Nf7y5NpSZXLj+NC5sHF8+DMnC/v842/94i6NTggh9h/aQQghhMhEAkIIIUQmEhBCCCEykYAQQgiRiQSEEEKITCQghBBCZLIvBYRyMQkhxODZlwJCuZiEEGLw7EsBIYQQYvBIQAghhMhEAkIIIUQmEhBCCCEykYAQQgiRiQSEEEKITCQghBBCZCIBIYQQIhMJCCGEEJlIQAghhMhEAkIIIUQmEhBCCCEyqez1AAbBVb/6Ofxp5/f2ehhCCLGv0Q5CCCFEJhIQQgghMtmXAkIFg4QQYvDsSwGhgkFCCDF49qWAEEIIMXgkIIQQQmQiASGEECITCQghhBCZSEAIIYTIRAJCCCFEJhIQQgghMpGAEEIIkYkEhBBCiEwkIIQQQmQiASGEECITCQghhBCZSEAIIYTIRAJCCCFEJhIQQgghMpGAEEIIkYkEhBBCiEwkIIQQQmQiASGEECITCQghhBCZSEAIIYTIRAJCCCFEJhIQQgghMhkaAWFmzzGz/2xmv29m/3ivxyOEEFc6AxUQZvZ+MztlZg+68280s2+a2cNm9i4ACCE8FEL4RwD+PoBjgxyXEEKIrRn0DuJOAG/kE2ZWBvBeAG8C8FwAbzOz5/aufT+A/w7gngGPSwghxBYMVECEED4L4Jw7fRzAwyGEb4cQGgA+AOAHeu0/HEJ4OYAfzevTzG43s/vM7L7Tp08PauhCCHHFU9mDZ14L4HF6fQLArWb2agBvAVAD8LG8m0MIdwC4AwCOHTsWBjdMIYS4stkLAWEZ50II4dMAPn1phyKEECKPvfBiOgHgenp9HYAn92AcQgghCtgLAfElALeY2Y1mVgXwVgAf3k4HZnabmd0xPz8/kAEKIYQYvJvrXQA+D+BZZnbCzN4eQmgBeAeAjwN4CMDvhhC+tp1+QwgfCSHcPjMzs/uDFkIIAWDANogQwttyzn8MBYZoIYQQe8/QRFILIYQYLiQghBBCZLIvBYSM1EIIMXj2pYCQkVoIIQbPvhQQQgghBo8EhBBCiEwkIIQQQmRiIezffHdmdhrAX+/1OHI4BODMXg+igGEfHzD8Yxz28QHDP8ZhHx8w/GPcyfi+K4RweKtG+1pADDNmdl8IYWgLHw37+IDhH+Owjw8Y/jEO+/iA4R/jIMcnFZMQQohMJCCEEEJkIgExOO7Y6wFswbCPDxj+MQ77+IDhH+Owjw8Y/jEObHyyQQghhMhEOwghhBCZSEAIIYTIRAJiFzCz95vZKTN7kM7Nmdmfmtm3ev8f2MPxXW9mnzKzh8zsa2b2z4ZpjGY2amZfNLMHeuP7N73zN5rZvb3xfbBXgXDPMLOymf2lmX10SMf3mJl91cy+bGb39c4NxXdMY5w1s983s2/0fo8vG5Yxmtmzep/d+r8FM/vnwzI+Gue/6P2dPGhmd/X+fgbyW5SA2B3uBPBGd+5dAO4JIdwC4J7e672iBeB/CSE8B8BLAfwTM3vuEI2xDuC1IYQXAHghgDea2UsB/AKAX+mN7zyAt+/R+Nb5Z+hWQVxn2MYHAK8JIbyQ/OKH5Tte59cA/EkI4dkAXoDu5zkUYwwhfLP32b0QwEsArAD4g2EZHwCY2bUA/imAYyGE5wEoo1u2eTC/xRCC/u3CPwBHATxIr78J4Jre8TUAvrnXY6SxfQjA64dxjADGAfwFgFvRjQ6t9M6/DMDH93Bc16E7ObwWwEcB2DCNrzeGxwAccueG5jsGMA3gUfScY4ZxjDSmNwD482EbH4BrATwOYA7diqAfBfB9g/otagcxOK4KITwFAL3/j+zxeAAAZnYUwIsA3IshGmNPffNlAKcA/CmARwBcCN0a5gBwAt0/jr3iVwG8E0Cn9/oghmt8ABAAfMLM7jez23vnhuY7BvBMAKcB/GZPVfc+M5sYsjGu81YAd/WOh2Z8IYQnAPx7AN8B8BSAeQD3Y0C/RQmIKwgzmwTw3wD88xDCwl6PhwkhtEN3a38dgOMAnpPV7NKOqouZvRnAqRDC/Xw6o+le+4y/IoTwYgBvQleN+Mo9Ho+nAuDFAP5TCOFFAJax9yqvTfT0998P4Pf2eiyenv3jBwDcCOAZACbQ/b49u/JblIAYHE+b2TUA0Pv/1F4OxsxG0BUOvxNCuLt3eqjGCAAhhAsAPo2urWTWzCq9S9cBeHKPhvUKAN9vZo8B+AC6aqZfxfCMDwAQQniy9/8pdHXnxzFc3/EJACdCCPf2Xv8+ugJjmMYIdCfcvwghPN17PUzj+14Aj4YQTocQmgDuBvByDOi3KAExOD4M4Md7xz+Ort5/TzAzA/AbAB4KIfwyXRqKMZrZYTOb7R2PoftH8BCATwH4wb0eXwjhfw8hXBdCOIqu6uGTIYQfHZbxAYCZTZjZ1Poxujr0BzEk3zEAhBBOAnjczJ7VO/U6AF/HEI2xx9sQ1UvAcI3vOwBeambjvb/r9c9wML/FvTYGXQ7/0P0xPQWgie4q6e3o6qjvAfCt3v9zezi+v4nulvMrAL7c+/e3h2WMAP4GgL/sje9BAO/unX8mgC8CeBjd7X5tCL7rVwP46LCNrzeWB3r/vgbgZ3rnh+I7pnG+EMB9ve/6DwEcGKYxouskcRbADJ0bmvH1xvNvAHyj97fy2wBqg/otKtWGEEKITKRiEkIIkYkEhBBCiEwkIIQQQmQiASGEECITCQghhBCZSECIKxYz+5leVsyv9LJ33to7/w4ze9jMgpkd2maf6xlVC4vIm9mdZvaDRW18OzP7HTM71899QuwGla2bCHH5YWYvA/BmAC8OIdR7gmA9RfKfo5sE7dM77P41IYQzFz/KlBDCj5rZnbvdrxB5aAchrlSuAXAmhFAHgBDCmRBTVfxlCOExf4OZTZrZb/Z2CF8xs7+31UN6O4pfsG69iy+a2c10+ZVm9jkz+zbtEszMft3Mvm5mf4ThSFwnrlAkIMSVyicAXG9mf2Vm/9HMXtXHPf8KwHwI4fkhhL8B4JN9PmshhHAcwK+jm8NpnWvQjXJ/M4D39M79XQDPAvB8AD+Fbp4dIfYECQhxRRJCWEK3KMzt6Kag/qCZ/cQWt30vgPdSH+f7fNxd9P/L6PwfhhA6IYSvA7iqd+6VAO4K3ey2T6J/ISTEriMbhLhiCSG00bUzfNrMvopukrM7C24x7CyNcsg5rru+s9oIsWdoByGuSHr1h2+hUy8E8Ndb3PYJAO+gPvqtTfzD9P/nt2j7WQBv7RVQugbAa/p8hhC7jgSEuFKZBPD/9ozBXwHwXAD/GgDM7J+a2Ql08+p/xcze17vn3wI40CsW/wD6n7xrZnYvujWt/8UWbf8A3ayhXwXwnwB8ZhvvSYhdRdlchdhFekWFjq27uWFPU7MAAABYSURBVPrXu9D/neimG//93ehPiCK0gxBidzkN4J6tAuV2gpn9DoBXAVjb7b6FyEI7CCGEEJloByGEECITCQghhBCZSEAIIYTIRAJCCCFEJhIQQgghMvn/Afxs2xl5RDKkAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "h =  h_file[\"hKr85\"]\n",
    "xmax = h.xhigh\n",
    "xmin = h.xlow\n",
    "ymax = h.yhigh\n",
    "ymin = h.ylow\n",
    "nbinx = h.xnumbins\n",
    "nbiny = h.ynumbins\n",
    "\n",
    "x = np.linspace(xmin, xmax, nbinx+1)\n",
    "y = np.linspace(ymin, ymax, nbiny+1)\n",
    "X, Y = np.meshgrid(0.5*(x[1:]+x[:-1]), 0.5*(y[1:]+y[:-1]))\n",
    "\n",
    "Z = np.array(h.values)\n",
    "eps=1e-6\n",
    "np.seterr(divide='ignore') #ignore log(0) produced warning message\n",
    "plt.pcolormesh(X,Y,np.log10(Z))\n",
    "plt.yscale('log')\n",
    "plt.xlabel('S1c [phd]')\n",
    "plt.ylabel('S2c [phd]')\n",
    "plt.title('Kr85')\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}