cms-jpsi.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Measuring the J/$\\psi$ Meson Mass\n",
    "\n",
    "This notebook will walk you through a simplified measurement of the mass of the [J/$\\psi$ meson](https://en.wikipedia.org/wiki/J/psi_meson).\n",
    "\n",
    "We will use data taken by the [CMS experiment](http://cms.web.cern.ch/) hosted on [CERN](http://cern.ch)'s [opendata portal](http://opendata.cern.ch).\n",
    "\n",
    "After importing some modules and defining helper functions we download the data and then start exploring it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [
    {
     "data": {
      "application/javascript": [
       "require(['codemirror/mode/clike/clike'], function(Clike) { console.log('ROOTaaS - C++ CodeMirror module loaded'); });"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "IPython.CodeCell.config_defaults.highlight_modes['magic_text/x-c++src'] = {'reg':[/^%%cpp/]};"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:rootpy.plotting.style:using ROOT style 'ATLAS'\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Welcome to ROOTaaS 6.05/02\n"
     ]
    }
   ],
   "source": [
    "# Required imports and setup\n",
    "import os\n",
    "import numpy as np\n",
    "\n",
    "from rootpy.plotting import Hist, Canvas, set_style, get_style\n",
    "from rootpy import asrootpy, log\n",
    "from root_numpy import root2array, fill_hist\n",
    "\n",
    "from ROOT import (RooFit, RooRealVar, RooDataHist, RooArgList,\n",
    "                  RooVoigtian, RooAddPdf, RooPolynomial, TLatex)\n",
    "\n",
    "style = get_style('ATLAS')\n",
    "style.SetCanvasDefW(900)\n",
    "style.SetCanvasDefH(600)\n",
    "set_style('ATLAS')\n",
    "log['/ROOT.TClassTable.Add'].setLevel(log.ERROR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def get_masses(filename, treename='events'):\n",
    "    # Get array of squared mass values selecting opposite charge events\n",
    "    expr = '(e1 + e2)**2 - ((px1 + px2)^2 + (py1 + py2)^2 + (pz1 + pz2)^2)'\n",
    "    m2 = root2array(filename, treename=treename,\n",
    "                    branches=expr, selection='q1 * q2 == -1')\n",
    "    # Remove bad events\n",
    "    m2 = m2[m2 > 0]\n",
    "    # Return the masses\n",
    "    return np.sqrt(m2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def plot_hist(hist, logy=True, logx=True):\n",
    "    # A function to plot the mass values\n",
    "    canvas = Canvas()\n",
    "    if logy:\n",
    "        canvas.SetLogy()\n",
    "    if logx:\n",
    "        canvas.SetLogx()\n",
    "    hist.xaxis.title = 'm_{#mu#mu} [GeV]'\n",
    "    hist.yaxis.title = 'Events'\n",
    "    hist.Draw('hist')\n",
    "    return canvas"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Fetch the data\n",
    "\n",
    "The data has been preprocessed and converted to a [ROOT](http://root.cern.ch) nTuple.\n",
    "\n",
    "These events were extracted from the [CMS Mu Primary Dataset](http://opendata.cern.ch/record/14). Thanks to http://openstack.cern.ch for providing a CernVM running SL5 where we could set up CMSSW 4.2.8 and to https://github.com/tpmccauley/dimuon-filter for generating the CSV from CMS AODs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "if not os.path.exists('events.root'):\n",
    "    !curl https://cernbox.cern.ch/index.php/s/ZE45HERahm7DeZ6/download -o events.root"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Plot the dimuon mass spectrum\n",
    "\n",
    "Calculate the dimuon invariant mass for each event and plot the full spectrum. This shows many of the well known dimuon resonances. Several of the strange features on the plot are explained by varying trigger threshold for each of the regions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4QAAAJYCAIAAAC1p7+MAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElE\nQVR4nO3dW7aqStYuUDgt6wWWTCiZYck4D5ErftILMhWIAHp/2G1uFyqiyOeIWz0MQwUAADn8v9w7\nAADAeQmjAABkI4wCAJCNMAoAQDbCKAAA2QijAABkI4wCAJCNMAoAQDbCKAAA2QijAABkI4wCAJCN\nMAoAQDbCKAAA2QijAABkI4wCAJCNMAoAQDbCKAAA2QijAABkI4wCAJCNMAoAQDbCKAAA2QijAABk\nI4wCAJCNMAoAQDbCKAAA2QijAABkI4wCAJCNMAoAQDbCKAAA2fwn1xOHEEIIVVW1bdu2ba7d2J26\nrnPvAgBwfMMwbPNE9QbP1HVd3/fX67XrunhL27b3+z1t0DRNDKZ8VNcLvGU/PsgXd//rXeZvv8gB\nOZLCD8j2u7fSMy71sL88znf3/dO9nIlfK/yAOBOXfZxDXha3/JCsXhl9ruSFEMZJtKqq+/3etq08\nCgBwNuv2GU2l0KZpUlt8unEYhmEYrtdrVVUP8RQAgDNYtwYby6IPrfDxxnGr/fMtvKOZ/vedObbC\nD4jGwQUfRzN9yQo/IM7EZR/nkJfFLT8kW4ymH0fMlErHNzZNs8FuHEY9KffeAQClm84SG8eJLcLo\neLD8y46hcQN9RmcaJuXeOwCgdNNZYuM4sUUYHafMvu+rqor9RB82MMETAMDZbNpMn1LpQ+40emlj\nP/7i+eLuf72LEu/XCj902+/eSs+41MP+8jjf3fdP9yr841Sywg+dM3HZx3FZ/NG6vVNDCJfLpfo3\nmj6WRavRMU0bVEc/0EspvFP89hwQKIEzEQqxxwFM684z2rZt0zT3f+KNt9st/jFOog8N9wAAnMHq\nzfQhhHHQvN1uz31DTeoEAHBOmRtWQgjGLf2JtrAHDgiUwJkIhdhjM72vj53xjf/AAYESOBOhEHsM\no+s207dtW9f1xwlEu66r61pLPQDA2WwxtdNMJr3nC4oxUAJnIhRijyfj8qPpu65LsTKOoP9Y8oyb\n6Tw60/QiXXv8FAIAWypq/fDlw2gI4WES+5lz2gujM4mbAMAvPmaJLdPq8mF0nCljDG2a5uO9uq4T\nRgEAzmbdoVJt297v95dzi/IdQ1YBgLUdagWmSvs7AABvKLPtjMooALC248wzCgAAE7YIo3Hq+482\n2BOgBE58AJJ1+4xWhU1kBRQitv74fgBgiwFM0fV6XfW5AADYnXV7p8ayR9M0lvpcigFMHED6GPs8\nA5Rpy+/nLcKoi82CXLw5AGEUoHBG0wMAcApbhFFt9AAAvLRuGI2DlrquW/VZzsYMWQDAL4qacHP1\nDgFxeXpjmJaijx0HoM8oQOGOszZ9CCGG0fv9Xtd10zTTG6+6MwAAlGbd2BuT6MyNFUjmUEniAFRG\nAQp3nNH040nvAQDggbLEzqgkcQAqowCFO05lFAAAJgijAABks1EYjcPq0+RVcebREIIpSAEAzmyL\nMNq27eVyeR5WH0Lo+76ua5M6AQCc0+phtK7rFEObphlPNZrG2l8uF3kUAOCE1p30PrXCX6/X+HfX\ndSmbtm07DENccupyuZxwUO3Lqa/kcgDgPNYNo33fV1XVNM1E39CUR2O/0lX3pyghhJcrApztOAAA\nZ7ZuGI0+jlJqmuZ+v58zhD0n9RMeBADgtLYIozPT1dmap+Pr7bpO+gQATmuL0fQzU+bZMtnZwjcA\nwLMtwujHZvrYdfJsYTTquu5h+lUAgPNYN4xer9eqqu73+0TMiqOXqvOF0RjB7/d7mvGq7/uzHQQA\n4OTWDaNd16WYVdd127axbTquvRTXZIpb3m63VfdkWfG1TG8zfo1t276L47fbLYQQQogzW00HdwCA\ng6k3mN2zbduXcxglaRbSXei6ru/7pmkmOn3GbZ5vnz7aIYTL5TL9yHW9xVsGq0ofY59ngDJt+f28\n0TPFMuFzJN1XDI1iNXciMsZMGf+OHRXGU4pOH/D44BPbuHhzAMIoQOG2/H7eYmqnqqpSA/1+xcb0\nl/XOBymJjt/FlFC7ruu6Lj5a27bj5v54iMYrpgIAHNsWo+kPoK7ry+UyJ4mmzP3QC7Zt29R9Nt7S\n932KreP7GsMEAJzHumE0Dt/ZXUP8L9KLfc6U6ZZYE43ZNNaMYzeGmFNPdbgAgJPbojL6MJR+j4aR\n6Wb02Df05TYpZaYpBeI6qJfLJZZdm6bRfw4AOJWNpnaqqiqmrr2n0pnmr4B6u92u1+v1eo1zPM25\nV/2bX14XAFC+H6PCxmlh3QFMKXfGITtpmvfYV3KPQ+kXMQ6dD2OY5lA9BQAm/B4VtsyjGw1gimF0\nGIY41VF0gOb7B+mFGIQEADDH1qPpu66LPS9TKk3N9xvvCQAA2WWb2imm0nGh9ADG4+XXe5ZyOnkA\nAHt0oj6j77xbkIk59BkFAH7xMUscsM9oFOfXrOv6crmkJNo0zfV6PVXA0qMUACDaojJ6pIXpZzrM\nkCwAgFWtPs/oQx20qqpYBx2G4ZBJNM6r+rIHQnq9h3zhAABfWDeMjguEcV73o2bQ5GGZpbGlyqXl\n9DgGAPboXAOYmqbpuu48vSTTK71cLuOOsKmjwu8TCJyqfy0AsLiiBjCtG0bP2XXydrvFJabquo7R\nM60+FaN53t0DAChHnqmdHsTMepjqadu2KY/2fZ9ub5rmnOkcAOCdesE23xBCTGDPjxlCiPM6PSfO\niXvtXSyCphe+SNr+WDY/3mHkeOr6v9886Q8AtjSnFX6z7+clK6MTZb/YYzJG0gWfsXArtci7eAMA\nvyiqz2i25UABAEAYBQAgG2EUAIBsihhND5zWuFuS/tAAJySM7s90n2KXc3Zk/HG1fhjAZor6yhVG\n90fcZKeK+u4DOLOiRtMLo8B2/JQC4IEBTAAAZCOMAgCQzSrN9O+WWYqLMM3cGACAw1sljMa15v/6\nTwAAnI0BTPtjaicA4BdFTW9SL5tdfmlzb9t2sf04rrpe+C2DzUx/en22Acqx5Xeyb/+dccFmv4RR\ngL3Y8jvZaHoAALIRRgEAyEYYBQAgG2EUAIBshFEAALIxz+j+mGcUAPhFUfOMCqP7I24CAL/4mCW2\nTKua6QEAyEYYBQAgG8307NLH5gOdGQBgF4RR9mp6Yckt9wQA+JpmegAAslEZZRnzi5Ea0AGARBhl\nMXNSpgZ0AGBMGIWjMboLgB0RRvfHCkwnN6e6PD26a/wIPjAAJ1RUQ6Uwuj/Sw6m8/L745TMwvu9D\nMP3xkQHYi6JWYBJGoXTrBcSHR1Y0BWB7wigUJG+7yUTRVDAFYCXCKGSzeBP8grTmA7ANYRS+93sh\ncxepTms+SXzrvenAgoRR+GwidJ7wqvyn1vyiBmwCUCBhFP7LnFlfeAimH7cBgAfCKOcica5KHRSA\nvxJGOR2JcyUOLABfEEY5JiU6ANgFYZQDUqIDgL0QRnkrV3FRURPK5NwE1iCMMmWNEuOc65nSJpRp\nGAaRFFiWMMqmpEwAYEwY3R+TEwEAvyiqiUMYPZqiPl6wEuuRAvzi4zfnlnFCGN0fl16OKn33zfmQ\nx238+sqirmtfRMBS/l/uHQCoqqoa/sm9I3zgPQKWJYwCAJCNZnpg3/7UuA9AaYRRYMdSANV5FGCn\nNNMDAJCNyihQHC3vAOchjAJl0fIOcCrCKFAuJVKAwxNGgUIpkRbFuwCsxACmUoQQ2rbNvRdQqPp/\n5d6dk1KfBtagMlqEEMLlcsm9F1AoGah88ReCdwr4gvWFizC/Y5wloeElp8baxkf4+WgLo3AwW36p\naqbPT+s8AHBawmhmXdfd7/fr9Zp7RwAAMtCwlVlq25rZyKUtEl5yaqztoZk+/hFvMQMXHI9m+tLV\ndf2xbT2E0HVd27Zx467rXj5O5esb2JthGB6+uHyPAV8zmv7PXsbK5236vk//e7/f7/d73/fj7+v4\nOLfbbfldhFNSn8tOfRr4gi+OP4sXvKZpQggvNxjP0xQ7g4YQ7vd7vCUe8LjN+EE008MinCNrePkF\nlQ718x/A3m15OquMzhVCCCGM653vpCQ6fhdTQu3+qarqfr8/TN+t4R4o0/zvJd9jwJ/oMzpLXdeX\ny2VOEk2Vzof297Ztm6apqio+iOmcYD1WaQLYEWF0YalH6XPcTLfEsU3D/4r/9DwsAPgTJ9HGXjbv\nAMwnjM4yTo2xwPlO7Bv6cpuUU991NgUWZCH7pUwcw4cW+fGva0cemEmf0VWs2gr/41e8ohFnMJ4R\nUxfG3708emmC5D/9E7CBfZ2AwmgGLyuj86+Urqkwnyrdjxw62KPfo8KW574wuqSUMo1PgtKYhXSO\nhzXhookj5mACvxNGgeN7brV//qcDew7iL0Pnw/ZnODJACQxgWtJ4vPx6z1JPWu954QBezmJxJO++\nFlJ3hfEXRfqnlxbZmd8fBFjDdJbY+ORVGd2fQ15BIZeSC6XPbeVzrhDvXsLGL80YJijZnBUft9mT\nShjNQo9SKMTz+pYltFCP92Fc0Yz/mkJeUbn5ndoCocAnwugqzCQKe/Sc/Ko3JckNRqpODCHaS7yL\nuVkeBabpM7qwON19nPr+QZr0Pv3xnXI6ecAhPfSbHJdLH258OPU+9sF62P5dZ80j9Wp92VcVyG7O\nd9RmhNGFTSyztFS59NiXLijKy5Pr5Xk3P1+e7Wx9OEoiKZSgqC8ozfQLS/1BL5fL+L3sui6WS6/X\na5YdA9ZznmT5i+Gp0px7j4AiCKPLu91ul8ulqqq6rmP0DCGkNet/bKMH2LV3JVLZFE5LGF1e27Yp\nj/Z9n25vmsbAJoDqzTwG8iick5N/RbEIGkJo//n9MT92t/KGAntUv589AFjcnN7bm52JwujOKB4A\nx6b5HkqwZd7QTA9AQcbXP6Od4AxM7QRAocwJBWegMgpA0Z7zqEIpHInKKAA7MJ6IW5UUjkRldH+m\nv4UVDIBjG/4teV/5xoNvFfWLztDsnTGaHiDSag/rMZoeAD54ObZJMIXdEUYB2LGHqaAy7gnwHWEU\ngONQJYXdMZoegIN4XvI+154A8xkNszPWpgeYyQgneMfa9HzPaHqAPzEJFHxhy7yhmR6AI7OgKBRO\nGAXg4MZLN4mkUBqj6QE4vuc8quEeCiGMAnAiz632UinkpZkegNMZhkEGhUIIowCcVKqS5t4RODXN\n9Psz/b3ptz7AX+lFytkU9RvMpJU7Y55RgMVZRBQemGcUALYz/FMVVjGCMxBGAeB/mI4UtqTPKAD8\nl+lIYXsqowDwyKJNsBlhFABeUBOFbQijAPCaUU2wAWEUAD6TR2ElBjABwJTn5ewrjfiwHJVRAPhM\n+oSVWM5nZz62E3lDAdZm1if2bk63k80+4Zrp98fXH0AJRFL26+PndstO0prpAeBvxhdyE5HCj1RG\nAeDPxrPi590T2DuVUQD4lUgKXxNGAeB744nxRVL4gjAKAL+yVhN8TZ9RAFhSyqMG2sMcKqMAsAzp\nE74gjALAYlJ7faXJHubRTA8AC3se0qRoCu+ojAIAkI0wCgArMsoepgmjALCKcf/RSh6FN2q9WPbl\n43eZNxSgQDqPUpQ5P402+6wKoztT194ygL0SSdmLLfOGZnoAALIRRgFgIw+9SIFKGAWA7RnMBIkw\nCgAZjKfEhzOzAhMAbGq8XqghTaAyCgB56EIKlTAKACXQZM9pCaMAUAR5lHMSRgEgp3FjvTzKCQmj\nAJDfy1FNcAZG02cWQgghVFXVtm3btpn3BgBgWxY6z+n5t+/tdpuOpNamBzi8eHXwbU9G1qY/hRg6\nr9dr7C10vV6rqrpcLpl3C4DcYgjQXs9JKLNl8/zDt23b+/0+XRxVGQU4CfPhk5HK6Fk0TZN7FwAo\nlAzKSSizFWROJyGVUYDzGDfT+/JnSyqj59K2beoYdLvdcu8OAKUQQDkDYfQbdV1/nIYphNB1XQya\nbdt2XTfnkeM0TwDwwGAmjkqb7591Xdf3fdM0E8ExbvN8+8cm+OrT7E6a6QHOSZM9W9oyb0g2fxa/\nDibCaAghzdAUJ2wKIdzv93jLxAGPd5yOucIowDkJo2xJn9ESxWb3Oa0kKYkOw9B1Xdd1IYTUGTS2\n18eHmtl2DwAxGZiClOMRRmep6/pyubxseX+QipoPQ5Hato0TOcUHiQ3xDxXQtC7o7zsMwPEoiHJI\n1qZfWCp2PmfKOKd9VVUhhJhN7/d7GtsUQog5VbkUgAmpOJp7R2AZKqOzDCPTM9XHuPlym5QyYwU0\nhBDz6OVySWVXUzsBAKeiMrqKmU3tIYQUTNu2nXmvH38Na+UBOIY5S6VwTvsqnAujGYy7isYA+qd+\nor53AIAJv0eFLeOsZvolpZRpEBIAq4o9xyoj69k/ldH9mf7SUTcFOI9hGOJFQZM9f1LUDxiV0SWl\nguiqq3oOk9Z7XgAKV1TCoGTTWWLjOKEyCgA7lnKDJMpOqYxmoEcpAIvTOMZOqYyuYtVmegB4xxL2\n7I7K6MLidPdx6vsHadL7H9dYqif98sgAwBlMZ4mN44QwurCHZZbGliqXltPjGICiuAowU1EDmITR\nhaX+oJfLZXx713WxXHq9XrffKwBOYpwkNJexC/qMLu92u8UkWtd1jJ4hhLRm/Y9t9AAwX13XyqUU\nThhdXtu2KY/2fZ9ub5rGwCYAtjH8mw8fCucH04piETSE0P7z+2N+/GbxhgIQC6KWZeKdOT9UNvvk\nCKM7o8EFgJnGYdTlgz/Z8gNjABMAHFOKodrrKZkwCgCHNS5uiaSUSRgFACAbo+n3Z/qnrS5BADww\nsp4HRX0ehNH9ETcBmC9eNVIeNZKJakaW2DKtaqYHgLNIeTT3jsD/EUYB4BQURCmTMAoAZ5FWrk/z\nPamSkp0wCgDnokRKUYRRADgpZVFKYEjdzlibHoBFjC8ocay9K8h5FLU2vamd9seXBQDwC1M7AQCZ\nWSmUQgijAHBeWtvIThgFgJN6SKJpsieFUrYkjALAqSmOkpcwCgD8l5oo2xNGAeDsFEfJyKRiO2Oe\nUQDW8zD5aMY9YVXmGeUnvh0AgF+YZxQAKJF6B9sTRgGA/yGSsiVhFAD4P5IoGxNGAYAXzH7PNoRR\nAGCKSMqqhFEA4NEwDLG9XhJlbcIoAADZCKMAwGsGM7EBKzDtjBWYANjYw6XHheYArMDET3wLALCl\nYRj0HD0YKzABAEBVCaMAAGQkjAIAH6SZniqTPbE0YRQA+Bt5lAUJowDALONRLzGPWjKU3wmjAMBc\nZnRhccIoAPANNVEWIYwCAJCNMAoA/Jn2epYijAIAP9Fezy8sBwoA/IGaKMtSGQUAfqU4yteEUQBg\nAfIo36kV2/fl46nuDQVgS88XJlei8s355bDZ+6jP6P44yQEoxzAMaqK78zFLbPmeaqYHACAbYRQA\nWFJatj73jrAPmukBgJ+M23xlUP5KZRQAgGyEUQBgMUbZ8lfCKACwMI31zKfPKACwpPFkT+NUqmjK\nSyqjAMAq1EeZQxgFABb2sggqm/7iwEdPGAUAKNqBk2ilz2h2IYQQQvy767qcuwIA66vrWudRxnwg\ncmrb9n6/j2+5Xq/TkdQ5DMAuvCzmjcc2uZzNFw/alkdsy7yhmT6bruvu93vTNLfbbRiG2+1WVVXf\n96lQCgAHM06odV0fu/V5cUc9XMps2Tz/yum6ru/7pmkm8qjKKAB7Ea9Z0xEqXtRc3aZtX05WGT2p\ntm2rqnpouAeAnZIvF/FQTs64JysRRrOJrfPjW2JBtGmaPDsEAJs7ZLriT4TRb9R1HauYE0IIXde1\nbRs3fh6W9PAIIYS+7ytj6gE4nGEYVEmXcrz4bmqnP5sTFmPvz/S/9/v9fr/3ff/uVEzD6q/X68eY\nCwD79bEXKWNnOFb6C/9Z/FhMDDMKIVwul/j39XqNt6SeoM9N82nj2+32MYnq4g3Afr0c0hRvcXV7\naTx0acthTFu+Iyqjc8XZ6cf1zndSuBy/iyl0dl2XaqupIDonhgLAgaWcKpWejT6js9R1fblc5iTR\nVC6N84YmbdvGkUnpQdI8o8MwSKIAnEEMmuLmd4563ITRhY2rng//lG6JgTWmUlPcA3BOKVqdoVsk\nEzTTzzL+LfK8hudY/KeX0zOlUU3jAPoysxpQD8CxpRLpcxLVUv/s2AdEGF3FnImf4h9fTHH/4y/I\nY3+gAdgjQ+yXta+DKYxmEMdCfV3+lCYBODz10V/8fui2jLP6jC4p1TsNSAIAmENlFADIZlzDS9W4\n+Ef8r/ro4Y+AyuiSHsbLr6SetN7zAsD2XNqe/X5MprPExsdcZXR/Dv8LCQB4aamRXh+zhD6jB6dH\nKQD8ifrogQmjqzCVPQAsToe0QxJGFxanu385e2iay+nHOe3L6eQBAMt62Xzs6ra4ovqMCqMLS0Hz\nuTi6VLl0mLTIUwAA2a0XCqezxMZxQhhdWOoPerlcxrd3XRfLpdfrdfu9AoDCqaecltH0y7vdbjGJ\n1nUdo2cIIa1Zb915AJgQU+lzUVBj/VEJo8tr2zbl0b7v0+1N0xjYBADvzC+O1qdZLPQML/Ms72UW\nsQgaQmj/+f0xP/4u9IYCcAzvLnlxrs0zXO9evsx6iYWp5pSZNzvCp3gvj+Qkpx8AVJOZadyaf9Qr\n47uL/gZhYMu8YQATAFCoM08Uc54+ssIoAFA0k48+O9LLF0YBAMhGGAUAdilVB49UJjwhUzvtz/Qp\nd9q+NQBwsGG+672WouK7MLo/RzrNAGCmd5PhR0Wlq/J9zBJbHk/N9ABA6RRiDkwYBQB2ZmY2VS7d\nBc30AMBuKJEej8ooAEBBZi79fZi6r8ro/hhNDwBRXKd+fEv6XxfECUUFWWF0f5xdADDTeLKnd3+/\n3DivtcOi0fQAAL8qJDjyI2EUANixd5F0XNur67qohukJe9nPBQmjAMD+KIsehjAKAJxR4TXI86Rt\nYRQA2Lfz5LZDEkYBgHOZrokWXjE9HmEUADiIL0qkRUXPonZmM+YZ3R+T3gPAfA/D6o90ofz65RSV\neoXR/TnSWQQAyxovlfm8PtODOZlsm/y6cUo26T0AwJKGYZgZ5t7FrHj7zHj6p31jmjAKABzNX6uM\nJeTLP+3DkZpJhVEA4NQ+psASouqBCaMAAK9tEEMlXWEUADigcUP2kRq1j0cYBQD4H8/VSvXL9Qij\nAACP6rpePIBKtC8JowDA8S3eUj9/KqiXd5y+5ZdH2x2T3u+PFZgAYKaVMujLG3Ndgr+YML+oCCuM\n7o+4CQBfGC/O9J2H+86ZE2riqp1xbVIrMAEAZPMcxRR6MhJGAQD+Zlw4TH+/HPP0rphqdFQijAIA\nLOZPoXC/CXJBwigAwMKeS6czc+fwz1p7Vh5hFAA4kZjzxv8lL2EUAOCz+U3qGt//RBgFAM5o7bLo\ny0FOL//3O+P933X8FUYBgFPTZJ+XMAoAQDbCKABwXl8XRNeopJ6zOiuMAgCQjTAKAEA2wigAQFWd\ntZU8u//k3gH+bHr6BicSADCtqKmghNH9ETcBYEvDMBSV3n73MUts+Xo10wMAkI0wCgDwf4Z/Hm6c\ned/524w3fn7G8xBGAQCOILat13W9r04F+owCAPzXREE0/h1z3mmrmGtQGQUAOIh91UQjYRQAYEl/\nqpsqsgqjAABHVni5VJ9RAIAvxSlIn7uWzrnjw39PS2UUAOBoZo6pL6FoKowCALzwsmD5ReHz+V7L\nlkL3XlgVRgEAyEYYBQA4pi9a4befM18YXVhd123bTm8TQui6rm3buHHXdVvsGQCwqJmLf67XjH6M\nRUSNpl/SnFjZdV3f9+l/7/f7/X7v+/4AHyYAoDRpjdCH28sJHnU5u3IA8Z1umiaE8HKDEMLlcol/\nX6/XeMv9fo+3zHkvHuaPAADyendpfrg9xcE1ruPxuf7UvJ5WN33Yn1X38yWV0QWEEEII43rnOymJ\njt/glFC7rtNkDwCcijLbr55/hbyrjKbQebvdHvqVtm0b66Mf3w6VUQAoyn4ro+NdSnu7fWXUAKbt\npKrn8windMu79n0AoEzvQtvz3KLKSS8Jo78aRpqmmdgy1j5fbpNyqjAKAJyKPqNb+zjxEwDA2sZt\n+u9G3G9DZbQsKqMAcGArtdTvugOAMLqRlDJVRgEAEmEUAIBshNGNLDhevv7Nj88OALAgA5j2Z9f9\nQgCAVf11wtHsVEbLokcpAHAqwujWjJcHAEiE0e3E6e7j1PcP0qT3c9am1yUUADgMYXQ7E8ss/alc\nOkxaZl8BgL3ZaQwQRreT+oNeLpfx7V3XxXLp9Xrdfq8AADIymn5Tt9stJtG6rmP0DCGkNevntNED\nAByJMLqptm1THu37Pt3eNM38lvrpjqE7LdEDAOdUyy5ZxCJoCKH9Z+Yd69pbBgC8FnPCIgOaN8sb\nks3OCKMAwDt7DKMGMAEAkI0wCgBANgYw7Y8BTADAYQij+yNuAgCHoZkeAIBshFEAALIRRgEAyEaf\n0f0xgAkAOAxhdH/ETQDgpRQSJqa+X2pW/KVopgcAIBthFACAbIRRAACyEUYBAMhGGAUAIBuj6ffH\n1E4AwGEIo/sjbgIAH6UpnGJyqOt6HCHKmeBJMz0AANkIowAAZCOMAgCQjTAKAEA2wigAANkIowAA\nZGNqp/0xzygA8KNC5nWqhNE9EjcBgAkpKkxkhudZSB/+3oxmegAAshFGAQDIRhgFACAbYRQAgGyE\nUQAAshFGAQDIRhgFACAbYRQAgGyEUQAAsrEC0/5YDhQA+FGKE+NckWWNUGF0f8RNAOCvxvnhZZao\n69pyoAAAnIswCgBANsIoAADZCKMAAGQjjAIAkI0wCgBANsIoAADZCKMAAGQjjAIAkI0wCgBANsIo\nAADZCKMAAGTzn9w7wJ/VdT3xr8MwbLYnAMAevcsS0xljJcLo/oibAMAvPmaJLVOpZnoAALIRRgEA\nyEYYBQAgG2EUAIBshFEAALIRRgEAyEYYBQAgG2EUAIBshFEAALIRRtm3LAuXAQ+ciVCIPZ6MwigA\nANkIowAAZCOMZlbXddu2ufcCACAPYTSnruty7wIAQE7CaE593+feBQCAnHT6ZKYAAAgPSURBVITR\nDEIIXdftcbwbAMCy/pN7B06nhAxa1/UwDFve/a93+XEPz6zwQ7f97q30jEs97C+P8919/3Svwj9O\nJSv80DkTl30cl8UfqYwCAJCNMLq1YaRpmty7AwCQkzAKAEA2wigAANkIowAAZCOMAgCQjamd9ufH\nyaEOPDcEAFCVMY/kfMLo6cQP6I8f0y/u/te7zN9+X6fcBgo/INvv3krPuNTD/vI43933T/dyJn6t\n8APiTFz2cVwWfyGM7oy6JgBwJPqMAgCQjTAKAEA2wigAANkIowAAZCOMAgCQjTAKAEA2tamCAADI\nRWUUAIBshFEAALKxAhMAwLmEEEII8e+u63Luij6jAACn0rbt/X4f33K9XjNGUs30AABn0XXd/X5v\nmuZ2uw3DcLvdqqrq+z4VSrenMgoAcBZ1XVdVNY5/Xdf1fd80Ta48qjIKAHBebdtWVfXQcL8lYRQA\n4Cxi6/z4llgQbZomzw4JoxxAXdfxVx2wgTlnXAih67q2bePG2cfqwjEscvY9PEIIoe/7KuuYelM7\nsW8ucrClOWdc7H+W/vd+v9/v977vDVGAX6xx9qVh9dfrNWNZxwAm9i12xM7Y7RpO5eMZF0K4XC7x\n7+v1Gm9JfdFcceBry559441vt1veBkaVUXYpztY7/v0HrGf+GZcub+MrX7rsdV2nNQP+ZI2zLxVE\ns8fQSBhlf+KvQ2Ab88+4VLCJMxcmbds2TRObC4VRmG+Nsy/NM1pOi6IBTAAsY1x3efindEs51z84\nkvlnXyyyFnUmCqPszzCScSoKOIn5Z1xs+Hu5TbpSFnUJhMItfvalE7B9YjQ9AAdRQhc0OKc5Ez/F\nPzJOcf9MGAVgOyqjkEusjBbYb1szPQALGDf/5dwPOJ+9n33CKAAA2QijACzAeHnIZe9nnzAKAEA2\nwigA29lpnzY4gGLPPmEUgCXttKEQDmCnZ58wCsAy4oTbL+cvTLPJFDitDBzArs8+YRSAZUwss7TT\ngg3sxa7PPmEUgGWkHmmXy2V8e9d1sWBzvV633ys4g12ffVZgAmAxt9stXgvruo4XvxBCWjW72FZC\nOID9nn3CKACLads2XRH7vk+3N01Tflsh7Np+z756GIbc+wDA0cQyTAih/SfzDsFp7O7sE0YBAMjG\nACYAALIRRgEAyEYYBQAgG2EUAIBshFEAALIRRgEAyEYYBQAgG2EUAIBshFEAALIRRgEAyEYYBQAg\nG2EUAIBshFEAALL5T+4dAADYvRBCCCH9b9u2bdtObJk2bv9ZZDe6rpt+9pdbzt/5NdTDMGz2ZAAA\nh9R1Xd/341ueI9bzNmPX6zUGxF/Udf3u2cdCCJfLJT1puteCezKfyigAwGKapqmq6rmy2Lbt/X5/\n2KyqqnRj3/d93/9YJWyaZvws74zrslVVXa/XeMuc+y5OGAUAWMy4vTsZJ9Hb7fYQVUMIXdfFDer6\np1brrutiybPruonqZirQxj1JWz6USLdhABMAwIpS0GyaZhiGl0XTEEKqlf7SRJ4efKI/QIrL1+v1\n6ydakDAKALCilAtfFk2T9K8TOXKOFGrfPd1DG312wigAsBvjgd5t29b/xOJivD2OyEm2HIvz7E9l\nyOltHl7X+CU/bPbw1A8e2ujzGwAAdiKml9vt9jLV3G63VBcci+3jq0o58uH2tD+/PPi711tV1fV6\nfd5+Iualh5q448t/Wo/KKACwM3GMzu12i2kmBcHL5XK/31P0TMEryyDx8VO/jMjzxddbVdX1eo2v\nOr20vu+fS7/pgDwXR8f14192aUHCKACwP8NoJNA4VzVNM+4TmULbdH/Ntf3SIJ7ue7vduq6L/9u2\nbSpkPvcxnWip/7FD6hqEUQBgZ54Lje9G7ZTSLfKN+o1xvE611efX8jFqv4uehYyjj4RRAGBnXs6O\nlGE/1jfdqp5e9XMYfdlSnx6knDb6yqT3AMDu7C56vqtcPlcoH2qZ6Y4Py8d/lJYe7bruYdKoHzuw\nLk4YBQBY17sRVM8VyncN69N9PV/m1Lg06PNTlxblhVEAgLXMXCz+o+leni/zZVoaNITQtm2ZbfSV\nMAoAsJ60Kv30YvHRc4GzbdvU2v7FU8c/Ykt9mW30lQFMAADrSSFyzpxKaT7RZGKIUnr8ca/QBzF6\njkuzpZVFK2EUADiP5/LkQ5L7uMEXUgt7XdcTm0135XzOqdW/UUp937+7b3ot4zVUp/Y1B2EUADiL\nGN3S/8bG63HW/LjBF7quS43jD3OIpmep6/pd19KUZR9y5JyW93SXRRaCWokwCgCwrhBCCoJ938eZ\n7du2jX+kqud4adMkZdn7/R7vFaV7TWfl8QMW2EZfCaMAABsIIaQFk6JxKfR6vcblPV82o4+z7P2f\nqqqapnl4zGfjByywjb4ymh4A2JG0IPuDibHq7+6yvbigfPgn3hI9bPN837h96sD6fMfpJ11i99ci\njAIAbGdmiHypzHb2H2mmBwAgG2EUADipj8PkfxxHzxzCKABwLm3bhhDiJJ3VqDvm/A0m1P8svtvr\nybvP+owCAOdyv98vl0vTNMMwxAk+67oej/L5uMFHZc7o+U7TNO9mOd3A344sAMB+xUwZQhgPIRr/\n78cNWJwwCgCcRV3Xt9ttIll+3IDF6TMKAJxC7Po5MSbp4wasQRgFAE7kY9VTWXRjmukBAMhGZRQA\ngGyEUQAAshFGAQDIRhgFACAbYRQAgGyEUQAAshFGAQDIRhgFACAbYRQAgGyEUQAAshFGAQDIRhgF\nACCb/w9ibZRM+Vk3IQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "Canvas('Canvas_qPpcHsJFLVvBSJ5cQhwEt8')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "masses = get_masses('events.root')\n",
    "\n",
    "mass_hist = Hist(1500, 0.5, 120)\n",
    "fill_hist(mass_hist, masses)\n",
    "plot_hist(mass_hist)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Zoom in on the J/$\\psi$ resonance\n",
    "\n",
    "The whole spectrum does not show anything weird, time to zoom in on the region of the J/$\\psi$ meson."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4QAAAJYCAIAAAC1p7+MAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElE\nQVR4nO3dT4gzyWH//+ofD/uclB5h5/TVDIMXTHIIIXswedaHbgUc7Bw3l4RkN7BOZvTHl7AG78xF\nrdM8CckhhxkNgmzADvi0zjEJBNSCxEtsSHwKDgnD4O0Q8g+piS82xvM7lJ9Ku/+UuqXuLnX3+3Wa\np9Uqlar1aD5T1VVlPT09CQAAAMCE/890BQAAANBdhFEAAAAYQxgFAACAMYRRAAAAGEMYBQAAgDGE\nUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMYRRAAAAGEMYBQAAgDGEUQAAABhDGAUAAIAxhFEAAAAY\nQxgFAACAMYRRAAAAGEMYBQAAgDGEUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMYRRAAAAGEMYBQAA\ngDGEUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMYRRAAAAGEMYBQAAgDHPTL2w7/u+7wshXNd1XddU\nNRrHsizTVQAAAO339PRUzwtZNbyS53nz+Xw2m3meJ4+4rrter9UJjuPIYIqdLKucS0Y5lEM5lEM5\nlEM5lFN1OXlU3jOa7MnzfT+aRIUQ6/XadV3yKAAAQNdUe8+o6gp1HEeNxauDT09PT09Ps9lMCBGL\npwAAAOiCavtgZbdobBReHoyO2iePIMuxdb9TDuVQDuVQDuVQDuUcoo7Z9NGIqVJp9KDjOKW8iuu6\nlmVZluW6blau9X3f0spT+M47CtT5mppEq1SocAAAgPZ4qlLyJeSgfOpBx3H2exVVZtJqtcp/claD\npJ6WVdvVapWzJlJqENc0RVmXjHIoh3Ioh3Ioh3Iop+py8qhjaSfZ8yd/ns/nQohYIlRrPO1XuCxT\nFit7FtUcqeFw+PTTaVL1O6amwGQd1BF126t8ufV67XlerNfT9/3hcKgqI3+Q5ydrIiKrCqjCZc2Z\n0QUAALqi0qgrX0L186lew1g3oTw4m832fonkG1GvFetlLPRaKlDGzldBNk9lsmqSdVwVntqZWvUl\ng0Q714N2rgftXA/auR60cz3qbOdqXymat6K9ockT9nvP6ump4TI1MhYKo5q6JcvRVCb15oSsRPuU\nCPGxh/LUHAeinetBO9eDdq4H7VwP2rkedbZztROYXNeVkWu9XqvBdBXaUge1C0mdDqWog8nx7kK3\nBGjqpt6UvjLq5aI1UQP0yWJVo+WvJAAAQBNVPpve9/1omFutVskgWNuiTioL5gmj+pOTCVVzi2dy\njdXkQzsPAgAAtE9NSzupnthozHJdV94WuXcSVSVnPapeKPZQdDUlKVmHQslVaHs6k/SFp/akAgAA\ntE8dYVSjui5ANaE+mg5VthsOh7FB8Pl8nrXIaKqieTEWUvXPomcUAAB0ROX3jFqWlWeJ+NS+yb1F\n70aNvnr059lstlqtVqtVdMA9mgL365XU58jkbaD68+kZBQAA7VbHOqM5lRK8fN/3PE9lvqwl6KMj\n+3KTJLnk53q9jq6KepwK9eAmZd3VgChaqR60cz1o53rQzvWgnfM4MCrUrPww6nmeipUyFO7s8pSn\nHR4B1RryQgjHcZLpVpN35TahIlL/aGlHhf+HAABA4/CoUGecLT+Mqps1lZyR7pAwGh2XFxlz9nMq\nGkBjL6Tv31V3jrquK5eFSu2ILTp3CgAAoKHKD6PR/JR/grkcKN/vFT3Pi+4Iuve9p47jRJOoyoup\nit5UcJydrAAAAGZVMkyvfpYj3YcEzZ3U3vSp4/JlvYSm81KRWVYfOlU50cn49IwCAIDOqmMHpkoT\nlRydz5NE5V2hmtn9sX5c/eJNyU5fzdvUhMvUwplEDyA/Bl4ANFq1YTQ6mami8uUPeV5Fsw2SyMiL\nMmsmB+tTT1bFJstPXX5fLimV+luk0Pr5ALopCILr6+t+v++6br/fv7q6CoLAdKUAoDCr0VOz1YR3\n/db20Syozo9Gxuj8p2iDqOOx89UUs1jrqfLzFKLKiXXrphYSfUqjLxmAUgRBcHl5eXZ2dn9/L4+M\nx+PHx8flcjkYDMzWDUAL1Jk36nilnGsk7VGTnOsORCfXR58iux6jdUtOw1eVV/cbqI7S1Dn7qnyZ\nj6NrCyTfYGzqVbTwrJlYhFEAQojr6+vNZqOSqDQej23bvrm5MVUrAK3RqjCaf52qesKo5llZC0Kl\nhumsk2OLTClZ7y6aRxXNmgCEUQBCiH6/v91uk8dt25bH1+s1t/oA2Ft7wmg0xulH0kWOtfFL5L8i\nhHBf0T8luhj+zqqq8ksvnDAKIAzDk5OTrEd///d//8/+7M+22+3JycloNJpOpwzcAyiqPWE09Z5I\nHIIwCkBk94w+e/bsd3/3d7mRFMCB2hZGCU8lIowCEEJcXV2FYbhYLKIH33jjjR//+Mff+c53oge5\nkRTAHgijyEQYBSCECILg4uLi/Pxc5dHJZLJcLn/0ox8lT1Y3kgJATnXmjWrXGZUYoweAcg0Gg+Vy\nadu2bdtCCNu2X3vttdQkKoQIw5AwCuBoVRt75VRx7hktET2jAGLUrsI7p9gDQE7tGaYXrybUk0fL\nQhgFkCX1RtLJZNLr9V6+fGmqVgCaqD1hVC5vpNbR1C96R1rNgzAKIEvqjaQPDw/L5fL09NRs3QA0\nS3vCaM69lyQyVh471/mnGYEuC4Lg9vZ2sViEYWjbtlxnlCQKICbPtkEtCaOp2wtlIUXlQc8ogDzU\njaQAsIf29IyidIRRAABQtbYt7QQAAACkIowCAADAmJrCqLx7yXrF8zx5UP4AAACAbqojjLquOxwO\nk9Pq5apPlmWxqBMAAEA3VR5GLctSMdRxnOhSo2qm53A4JI8CwE75F8sDgKaoNoyqUfjZbPb09BRb\nasR1XTVRazgcVloTAGiuIAiur6/7/b7ruv1+/+rqKggC05UCgHJUG0blIqOO42juDVV5lM5RAEgK\nguDy8nKz2cj95bfbbRiGFxcX5FEA7VDHPaM7ZynJsXvCKAAk3d3dnZ2d3d/fqyOLxeL8/Pz29tZg\nrQCgLM9qeI2cu4AQRgEgabFYyD7R2EHbtm9uboxUCQBKVEfPaM6Uyc51ABAThmEyie58CAAa5CiG\n6eX8UMIoAMTYtn1yclL0IQBokGrD6Gw2E0Ks12tNHrUsS/5AGAWApNFoNB6PYwcnk8loNCpaFCtD\nAThClS/tJCcnycXtXdeVQ/Zy7yW5J5M8c7VaVVoTAGio6XT6+PgYzaOTyeTh4WE6neYsgZWhAByz\nyofpfd9XC92v12v5d/l6vZ7P5+pv9NlsRrcoAKQaDAbL5dK2bdu2hRC2bfd6veVyeXp6mufprAwF\n4MhZapnPSsmu0OQI0Ww2Y3v6QiyrpksG4AjFtg7J4/r6erPZRFeGEkKMx2Mm4wPQqDNvkGwaRt3Y\nkIULCiCq3++nTrq3bZvJ+EBn7YwTosZEQRhtGHpGAeQXhqFmxv1ms2E+PoBUdeaNau8ZlVOUGIgH\nACNYGQrA8atjndHYVHoAQG1KXBkKAKpQ09JOQoj1ej0cDkmlAFCnw1eGAoBKVT5M7/v+09PTbDZL\nplKG7wGgageuDAUAVat7NoznefP5PHrEcRy5AH6d1WguJjAB2NseK0MB6KZOLO2UTKVkrDwIowAA\noGrtmU2v4XmeHL43VQEAAAAY98zIq2ZtyAQAAIBOqbVnVN6uZFnWcDhUSdRxnNlsxtAzAABAB9XR\nM8rG9AAAAEhVbRhNzlISZFAAAAC8Um0YjS5uP5vNXNdlVREAAAAolQ/Ts4woAAAAsrBoZcOwzigA\nAKhaJ9YZjfJ9n93qAQAAOqjMMOr7vmVZlmWlPuR5Xmri9H1/OBwOh8MSawIAAIBGKPOeUU3vplza\niR5QAAAARJnZgQmHSO17VrijFAAA6OmzRM0Io81D3AQAAIfYmSXqTKtHMYEJAAAA3UQYBQAAgDGE\nUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMYRRADg66/XadBUAoCaVrDOatc2S3IQp58kA0DVBENzd\n3S0Wi+12e3JyMhqNptPpYDAwXS8AqFAlYVSz0Tx70ANAqiAILi8vz87OttutEGK73YZheHFxsVwu\nyaMAWoxhegA4Cnd3d2dnZ/f39+rIYrE4Pz+/vb01WCsAqJpV7t6Sh4y5u65bWj3ay7JKvmQAjkS/\n35d9ojG2baceB4Dq1Jk3SDYNQxgFWikMw5OTk6xHN5uN5lEAKF2deYNhegAwz7btrLipeQgAWqCS\nCUwAgKJGo1EYhovFInpwMpn0ej1TVQKAGtAzCgBHYTqdPj4+jsdjdWQymTw8PEyn00pflzVNAZhF\nGAWAozAYDJbLpW3btm0LIWzb7vV6y+Xy9PS0ipcLguD6+rrf77uu2+/3r66ugiCo4oUAQI/ZMA1j\nWZb+BC4o0AK+71e6wIha01StJDUejx8fH1nTFOiInXFC1JgoCKMNw2x6AIe7vr7ebDbRNU2FEOPx\n2Lbtm5sbU7UCcDxY2gmZCKMADseapgD0WNoJAFCVMAyzEqfmIQCoCGEUALqFNU0BHBXWGQWAzmFN\nUwDHg55RAOgcU2uaAkASYRQAOqfmNU0BQIOp2Q3DbHoA5ap6TVMATcRs+n14nue6rmVZlmW5rut5\nnuZk3/fV+TtPThbu+37O8/MULn8T5C8cAEpEEgVgVhu62TzPm8/nqQ+tVqvk92zW+VlNkbpLgeM4\nqanR9/3hcJizJkII13WTG0NnFS7oGQUAANVj0fsCouFvNpvJnkXf91XCi73B2PnySNbJIhIWHceR\naVIF2dlsFuv1TBYePT9P4aoyWXmUMAoAAKpWa954arisN7JareRxx3F2nq9Ons1m0eMqUMaOO46T\n+qL6wmM1yTquCl+tVqnvN3kQAACgRHXmjWYnm6wQKSUjozo/mfNS86UmsidfV1MZFWp3vmK08FhI\nVQ8lDwIAAJSozrzR7AlMaiA7dZKQOpg8LXn7pjqSHBxXUTIpeu+ppjKphavh+GSx8mDyXlIAAICW\naXYYLUqT/5LJVf2QOvEomVA1s+BVCZqcuvMgAABA+zQ7jHqepzqTUx+VP8SyXc6opw+jSZqkW7Rw\nTTctAABAmzQ7jGpEp6UXfWKe04rmxVg19M+iZxQAAHREO8NodImlnMPuqYXs8dL6wpO3gerPp2cU\nAAC0W9vCqNzNSCVRNcO9TazDmK4+AACo1oFRoea00KowKmOoGp1/enqK9jsWHVg/2rHyAxdQMF19\nAABQrQOjQs1poSVh1Pd9y7LUIPhqtaptgDuWWfWvq+4c1SfjoncUAAAANFQbwqjnedFNOGMdokWp\n55Z7NyeLhgIAACQ1Poz6vi9Xnpfj8qmr3yefUvQl8hzMs1J9MunSMwoAALqs8WFU9ok6jpMnYmry\nokqxydVJU0tOrhulCY6acJkz6QIAALRSs8Noctukvc9PLUFmzeien7GTo+EyGWeTD0XPl3s4pSbj\n/VZIBQAAaByr0dOrXdeVuU2zfbz46XSoViuIvnHP82TinM1m0ZPVeqWx46mFROsTPZ5ViCon1q2b\nWkj0KY2+ZAAA4PjVmTeanWxyroO1Wq2iY+7R2U7ip/dqSvaPqmjoOI4sRHWURotNVilWuEgLlyoE\nq/Oj/0y9/5UwCgAAqkYYzWuPMCp+Oo8qmrtOVR7VlKkvXGR0c4qfzqNKVhIVhFEAAFA9wmgdZNqT\nOzZJO8+XadV13Z1z9v1XSi+cMAoAAKpGGEUmwigAAKhanXmj2bPpAQAA0GiEUQAAABhDGAUAAIAx\nhFEAAAAYQxgFAJNSt2EDgO4gjAKAAUEQXF9f9/t913X7/f7V1VUQBKYrBQAGEEabx9IyXTsAuwVB\ncHl5udlsttutEGK73YZheHFxQR4FUA99lqg5TrBoZcOwzijQAtfX15vN5v7+PnpwPB7btn1zc2Oq\nVgCgsOg9MhFGgRbo9/uyTzTGtu3U4wBQMxa9B4DWCsMwK3FqHgKAtiKMAkCtbNs+OTkp+hAAtNUz\n0xUAgM4ZjUZhGC4Wi+jByWTS6/VMVQkATKFnFADqNp1OHx8fx+OxOjKZTB4eHqbTqcFaAYARhFEA\nqNtgMFgul7Zt27YthLBtu9frLZfL09NT01UDgLoxNbthmE0PtIzv+67rmq4FAPwUZtMDQFccbRJl\nn1IA9SCMAgD+D/uUAqgZYRQA8BPsUwqgfoRRAMBP3N3dnZ2dRfcpXSwW5+fnt7e3BmsFoN2YDdMw\nTGACUB32KQUgMYEJAFA39ikFYARhFAAgBPuUAjCE7UABAD/BPqUA6kfPKADgJ9inFED9CKPNY2mZ\nrh2ABmOfUqAj9Fmi5jjB1OyGYTY9gHqwTynQZXXmDZJNwxBGAQBA1VjaCQAAAJ1AGAUAAIAxhFEA\nAAAYQxgFAACAMYRRAAAAGEMYBQAAgDGEUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMYRRAAAAGEMY\nBQAAgDGEUQAAABhDGAUAAIAxhFEAAAAYQxgFAACAMc9MVwCFWZalefTp6am2mgAAgCbSZ4maEUab\nh7gJAAAOsTNL1JlWGaYHAACAMYRRAAAAGEMYBQAAgDGEUQAAABhDGAUAAIAxhFEAAAAYQxgFAACA\nMYRRAAAAGEMYBQAAgDGEUQCoyXq9Nl0FADg6hFEAqFYQBNfX1/1+33Xdfr9/dXUVBIHpSgHAsSCM\nAkCFgiC4vLzcbDbb7VYIsd1uwzC8uLggjwKARBgFgArd3d2dnZ3d39+rI4vF4vz8/Pb21mCtAOB4\nWE9PT6brgAIsi0sGNEm/35d9ojG2baceB4BjUGfeoGcUAKoShmFW4tQ8BACdQhgFgKrYtn1yclL0\nIQDolGemK4DCLMvSPMogPnBURqNRGIaLxSJ6cDKZ9Ho9U1UCAH2WqBk9o83zpGW6dgB+ynQ6fXx8\nHI/H6shkMnl4eJhOpwZrBaDj9Fmi5jhBGAWACg0Gg+Vyadu2bdtCCNu2e73ecrk8PT01XTUAOApM\nzW4YZtMDzeX7vuu6pmsBALvVmTdINg1DGAUAAFVjaScAAAB0AmEUAAAAxhBGAQAAYEwLw6hlWZop\nAr7vW1qpz/I8z3VdeYLrur7v6+ugzndd1/M8/clyTkP+wgHAlPV6bboKANqmbbNhPM+bz+eO42RF\nOnmCpoRkg6Qm1KyX8H1/OBwmj69Wq9SI7Lpu8stdU38mMAGoXxAEd3d3i8Viu92enJyMRqPpdDoY\nDEzXC0BV6swbbduBSR80hRAq5DmOk3w0mRfVEcdx5M/yJdbrted5sV7PaBKdzWbRKg2Hw+RFVUlU\nFe77/nq9Xq/XdJECOBJBEFxeXp6dnW23WyHEdrsNw/Di4mK5XJJHAZRg5xL8jbBarVT4E0I4jpN1\npjxhNpvlKVaVGTtfBdnUwmPHV6tVaq2yjqvCV6tVav3z1BwAynJ1dTUajWK/O8bj8fvvv2+6agCq\nUmfeaMOYb3IYXT/MLYSYzWY7b+WMlpxspWQ5qls0Wbi6NyBajuoWzSo89V0wTA+gZv1+X/aJxti2\nnXocQAuwzmjlCm2CEu1zjYneFaCCYzLmqpeLhks1QJ8sVh5kogAA48IwzEqcmocAIL82hNFoT29q\ntlNUFswTRvUnJxOq5hZPVYImp+48CAD1s2375OSk6EMAkF8bwugeoqspScmYWCi5Cm1PZ2oFNIWn\n9qQCgBGj0Wg8HscOTiaT5I2kALCHboVRle2Gw2FsEHw+n2ctMpqqaF6MhVT9s+gZBXA8ptPp4+Nj\nNI9OJpOHh4fpdGqwVgBao6NhVAgxm81Wq1VsGn40Be7XK6nPkcnbQPXn0zMKwLjBYLBcLm3btm1b\nCGHbdq/XWy6Xp6enpqsGoA3ats5oTrFZ7XLDJLnApxzBN1e13Qr14CYxGR9AUYPB4Obm5ubm5vi/\nIQGIg6NCzTrXM6pWz0o+JH9QN48e7Rfugat5ma4+gAY72i9GAFEHRoWa00K3wmgeRRdUin016wfW\n1Z2j+ltOi86dAgAAaCjC6P+JzTEq925OFg0FAABIIozupu+8VPKsVK8CLj2jAAAAolNh1Pd9uaRo\nVqdmbKFQfV5MriqqCY6acJkz6QIAALRSh8KoZhskkZEXZdaM7vmpOVkVmyw/dVKUXFIqtSe10Pr5\nAAAAzdWhMCoiI+mxvOj7/nA4lD9HH8rKl6kni4zw6vu+DJexHUSzpu2z/RIAAOgOq2Vr/cjlQh3H\nyUpy0ZW3knd5rlarZDRU/ZTyIZU1kydHy5fRUyVRkbbAp+d5qjR5fvSfqT24ltW2SwYAAI5NnXmj\nbclmZxgV2SvBpoZLEcmjeU6OdrJGZbVzNI8qWUlUEEYBAED1CKOV818RQriv6J/ieZ46PyspJssv\nvXDCKAAAqBphFJkIowAAoGp15o1uTWACAADAUSGMAgAAwBjCKAAAAIwhjAIAAMAYwigAAACMIYwC\nAADAGMIoAJQvuVMGACAVYbR5LC3TtQM6LQiC6+vrfr/vum6/37+6ugqCwHSlACBOnyVqjhOE0eZ5\n0jJdO6C7giC4vLzcbDbb7VYIsd1uwzC8uLggjwI4NvosUXOcIIwCQDnu7u7Ozs7u7+/VkcVicX5+\nfnt7a7BWAHDk2FuyYdgOFDha/X5f9onG2LadehwAjhbbgQJAw4RhmJU4NQ8BAAijAFAC27ZPTk6K\nPgQAeGa6AgDQEqPRKAzDxWIRPTiZTHq9nqkqAcDxo2cUAMoxnU4fHx/H47E6MplMHh4eptOpwVoB\nwJEjjAJAOQaDwXK5tG3btm0hhG3bvV5vuVyenp6arhoAHC+mZjcMs+mBRvB933Vd07UAgD3VmTdI\nNg1DGAUAAFVjaScAAAB0AmEUAAAAxhBGAQAAYAxhFAAAAMYQRgEAAGAMYRQAcJD1em26CgAajDAK\nANhHEATX19f9ft913X6/f3V1FQSB6UoBaB7CaPNYWqZrB6ATgiC4vLzcbDbb7VYIsd1uwzC8uLgg\njwKNoM8SNccJVlBvGBa9B3AMrq+vN5vN/f199OB4PLZt++bmxlStAJSFHZiQiTAK4Bj0+33ZJxpj\n23bqcQDNwg5MAIDjFYZhVuLUPAQAqQijAIBibNs+OTkp+hAApHpmugIAgOYZjUZhGC4Wi+jByWTS\n6/VMVQlAQ9EzCgAobDqdPj4+jsdjdWQymTw8PEynU4O1AtBEhFEAQGGDwWC5XNq2bdu2EMK27V6v\nt1wuT09PTVcNQMMwNbthmE0P4Nj4vu+6rulaACgTSzshE2EUAABUjaWdAAAA0AmEUQAAABhDGAUA\nAIAxhFEAAAAYQxgFAACAMYRRAAAAGEMYBQAAgDGEUQAAABjzzHQFUJhlWZpHWRIfAADo6bNEzQij\nzUPcBAAAh9iZJepMqwzTAwAAwBjCKAAcZL1em64CADQYYRQA9hEEwfX1db/fd1233+9fXV0FQWC6\nUgDQPIRRACgsCILLy8vNZrPdboUQ2+02DMOLiwvyKAAURRgFgMLu7u7Ozs7u7+/VkcVicX5+fnt7\na7BWANBEFlOzm8WyuGSAef1+X/aJxti2nXocAJqlzrxBzygAFBOGYVbi1DwEAEhFGAWAYmzbPjk5\nKfoQACAVi94DQGGj0SgMw8ViET04mUx6vZ6pKgFAQ9EzCgCFTafTx8fH8Xisjkwmk4eHh+l0arBW\nANBEhFEAKGwwGCyXS9u2bdsWQti23ev1lsvl6emp6aoBQMMwNbthmE0PHBvf913XNV0LAChTnXmD\nZNMwhFEAAFA1lnYCAABAJxBGAQAAYAxhFAAAAMYQRgEAAGAMi943j2VZmkeZ3gQAAPT0WaJm9Iw2\nz5OW6doBgBBCrNdr01UAkEmfJWqOE4RRAEBpgiC4vr7u9/uu6/b7/aurqyAITFcKwFEjjAIAyhEE\nweXl5Waz2W63QojtdhuG4cXFBXkUgAZhFABQjru7u7Ozs/v7e3VksVicn5/f3t4arBWAI8d2Pg3D\nDkwAjla/35d9ojG2baceB3C02IEJANAwYRhmJU7NQwBAGAUAlMC27ZOTk6IPAUALw6hlWa7r6s/x\nfd/zPNd15cme5+nPVyfL833fz3l+nsJ93y9UOAAcp9FoNB6PYwcnk8loNDJSHwCN0LYbED3Pm8/n\njuNoIp08J3k8qylSF4bNegnf94fDYfL4arVKjciu6yZX49PUn3tGARytIAguLi7Oz88Xi4U8MplM\nHh4elsvl6emp2boBKKTOvNG2ZCODoybMRcPibDaTR1QcTLaGCouO48g0qYLsbDaL9XomC4+en6dw\nVZmst0AYBXDMgiC4vb1dLBZhGNq2PRqNptMpSRRonFrzxs4l+BthtVqp8CeEcBwn68zUN75areTB\n2WwWPa7KjB13HCe19fSFx2qVdVwVvlqtUuuf9dYA4HikfoMBaIo680YbutmSw+g7x9CTg+aqkzLa\nIKrkZCvJh6Kdo6rwZI+pujEgWk7qK0YLT30X9IwCAICqsbRTVVRGTN6+qY4k81+0zzUmeu+pemJy\nxlJq4Wo4PlmsPMjOzgAAoPXaEEajPb2p2U7R5L9oB2fsh9SJR8mEqpkypUrQ5NSdBwEAANqnDWG0\nqJxRTx9GkzRJt2jhmm5aAACANuliGNXLmf+K5sVYSNU/i55RAADQER0Ko0V7OvfrldQXnrwNVH8+\nPaPA8eA2bgCoQofCaGtYhzFdfaBhgiC4vr7u9/uu6/b7/aurqyAITFcKAHQOjAo1p4UOhdGiA+tH\nO1Z+4GpepqsPNEkQBJeXl5vNZrvdCiG2220YhhcXF+RRAMfswKhQc1roUBitSCyz6pOuunNUn4yL\n3lEAoCJ3d3dnZ2f39/fqyGKxOD8/v729NVgrAGgTwmicyn/l3s3J3WZAEy0Wi2gSVQfV3usAgAN1\nMYwWzZH6zkslz0r1yaRLzyhwtMIwlKPzhR4CABTSrTCqyYtqOfrkLk2peTG5qqgmOGrCZc6kC6B+\ntm2fnJwUfQgAUEi3wmhymyUlNf/JrBnd8zN2cjRcJuNs8qHo+XIPp9RkXHJ/5DQAACAASURBVGj9\nfADVGY1G4/E4dnAymYxGIyP1AYD26VYYVVlwOBxGj3ueJ/NfbJPPrHypnh47nhpefd/XFx7rLmX7\nJeB4TKfTx8fHaB6dTCYPDw/T6dRgrQCgTZ6ZrkDdVquVjJKWZcl0qMKi4zixcOm6ruM46/V6Pp/7\nvi9josqaq9UqVrjv+3JprmThIq3HdDabzefz9XqtzleFx5IrACMGg8Fyuby9vbVtOwxD27Z7vd5y\nuTw9PTVdNQBoCatlC0+6rrterx3H0XQr+r4f6xkVQmieIsuMHVytVqk3iaYWLoTIamfP85K3Acxm\ns2RylSyrbZcMaAr1FykAtF6deaO7yUamPfnbRdp5vkyrrutmJUXFf6X0wgmjAACgaoRRZCKMAgCA\nqtWZN7o1gQkAAABHhTAKAAAAYwijAAAAMIYwCgAAAGMIowAAADCGMAoAAABjOrcDUwvITZ6ysPAT\nAADQ02eJmhFGm4e4CQAADrEzS9SZVhmmBwDUJLm1MgAQRgEA1QqC4Pr6ut/vu67b7/evrq6CIDBd\nKQDHgjAKAKhQEASXl5ebzWa73QohttttGIYXFxfkUQASYRQAUKG7u7uzs7P7+3t1ZLFYnJ+f397e\nGqwVgONhMRumWSyLSwagSfr9vuwTjbFtO/U4gGNQZ96gZxQAUJUwDLMSp+YhAJ1CGAUAVMW27ZOT\nk6IPAegU1hkFgLj1eu04julatMRoNArDcLFYRA9OJpNer2eqSgCOCj2jAPATrEBUhel0+vj4OB6P\n1ZHJZPLw8DCdTg3WCsDxIIwCgBCsQFSZwWCwXC5t27ZtWwhh23av11sul6enp6arBuAoMDW7YZhN\nD1Tk+vp6s9lEVyASQozHY9u2b25uTNWqZXzfd13XdC0A7FZn3iDZNAxhFKgIKxABgMLSTgBQK1Yg\nAgBTCKMAwApEAGAMSzsBgBCsQAQAhtAzCgBCsAIRABhCGAUAIViBCAAMYWp2w1iWpT+BCwocjhWI\nALTbzjghakwUhNGGYWknAABQNZZ2AgAAQCcQRgEAAGAMYRQAAADGEEYBAABgDGEUAAAAxhBGAQAA\nYAxhFABg2Hq9Nl0FAMYQRgEAZgRBcH193e/3Xdft9/tXV1dBEJiuFIC6EUYBAAYEQXB5ebnZbLbb\nrRBiu92GYXhxcUEeBbqGMAoAMODu7u7s7Oz+/l4dWSwW5+fnt7e3BmsFoH7sLdkwbAcKoB36/b7s\nE42xbTv1OIA6sR0oAKDNwjDMSpyahwC0EmEUQKcxj9sI27ZPTk6KPgSglQijALqIedzGjUaj8Xgc\nOziZTEajkZH6ADCFMAqgc5jHfQym0+nj42M0j04mk4eHh+l0arBWAOpHGAXQOczjPgaDwWC5XNq2\nbdu2EMK27V6vt1wuT09PTVcNQK2Ymt0wlmXpT+CCAjsxj/vY+L7vuq7pWgAdsjNOiBoTBWG0YVja\nCThQGIaa+TGbzYbZMwDA0k4AUBXmcQPAUXlmugIAULfRaBSG4WKxiB6cTCa9Xs9UlQCgs+gZBdA5\nzOMGgONBGAXQOczjBoDjwWyYhmECE1Au5nEDQFKdeYNk0zCEUQAAUDVm0wMAAKATCKMAAAAwhjAK\nAAAAYwijAAAAMIYwCgAAAGMIowAAADCGMAoAAABjCKMAumK9XpuuAgrgegEdQRgF0HJBEFxfX/f7\nfdd1+/3+1dVVEASmK4VMXC+gawijANosCILLy8vNZrPdboUQ2+02DMOLiwvyzXHiegEdxN6SDWNZ\nlv4ELigQdX19vdls7u/vowfH47Ft2zc3N6ZqhSxcL6AeO+OEqDFREEYbhr3pgUL6/b7sY4uxbTv1\nOMziegFHgr3pAaAEYRhmJRjNQzCF6wV0E2EUQGvZtn1yclL0IZjC9QK66ZnpCgBAhUajURiGi8Ui\nenAymfR6PVNVggbXC+ggekYBtNl0On18fByPx+rIZDJ5eHiYTqcGa4UsXC+ggwijANpsMBgsl0vb\ntm3bFkLYtt3r9ZbL5enpqemqIQXXC+ggpmY3DLPpgb35vu+6rulaIC+uF2BQnXmDZNMwhFEAAFA1\nlnaqlqXl+37yKb7ve57nuq5lWa7rep6nfwl1sjw/tczU8/MUDgAA0Bpd7GbT7zqwWq1iA0Oe583n\n8+SZWU2XWr7jOFkxdzgc5qmGKryDlwwAANSpzrzRuaWdVCJ0HCfn+SqJzmYzeWS9XouM66QSpOM4\n8mf59PV67XlerNczmkRl4er84XBI6AQAAK3XuW421c2Z842rbs7o+SpEzmazaL5UhceOu64r82vs\nRfWFp3am0jMKAACqxj2jx0JlwdVqFT3uuq7sWI0N36t/JntAk8fVQdUnqgqXR2R+BbAH/vu0FVcW\naJ/OhVEZAXOO0avsmLx9Ux1Jdl7GwmVUNLymJtSdhQPQCILg+vq63++7rtvv96+uroIgMF0plIAr\nC7RY58Ko/Ks659p18uTU5KoSpMqL6ofUwpMJVRM0VQnMrAfyC4Lg8vJys9lst1shxHa7DcPw4uKC\n1NJ0XFmg3ToXRhXP86IrOmniac7kqg+jSZqkC2APd3d3Z2dn9/f36shisTg/P7+9vTVYKxyOKwu0\nW0fD6Hw+j93uuV6vsxYZ1cv5lKIj74RUoKjFYhHNK+rgYrEwUh+UhSsLtFu3wmg0CDqOs1qtnp6e\nVquVSn7RJT+L9nTud3+nvnBu1QdyCsNQjuEWegjHjysLtF63wqgym83UrsdyhyQ1X/74t0LW7yC1\nk+nqA5Wwbfvk5KToQzh+XFlgDwdGhZrTQrcWvXddN2vRLLla03q9Vp2R0YH1PAlVLSZaNdYZBVKN\nRqMwDGNDt5PJpNfrmaoSSsGVBYo6PCrUmUc72jOaqp4FlWK5Vv9a3DkK5DedTh8fH8fjsToymUwe\nHh6m06nBWuFwXFmg3Qij/2e/AXr1LP3TiwZc7hYFihoMBsvl0rZt27aFELZt93q95XJ5enpqumo4\nCFcWaDfC6G5Fc2Tq+cmDstdTHzqP/wZW4KgMBoObm5vtdrtarbbb7cuXL8kr7cCVBVqsW2FU3pOb\ntZJ8cvq8Ji+qQpK7NKWG0eSqopqgWXQiP4AY/u+0FVcWaJ9uhdHUDeWV5PHkNktKauLMKj81XCbj\nbPIhvnYBPW5oAYCm61YYzbPXvFrjKXowuv6oLEf+Coxt8pmVL9XTY8dTw6vv+6mFA1DYqRwAWsPq\n2jpB0QWY1KB59Eisy9P3fRUlZTpUYTF5crR8x3FkllVZc7VaJUOwWjohVrjIWJfBsjp3yYAYuVN5\ndH/I8Xj8+Pi4XC4Hg4HZusGI9XrN2iNAuerMG11MNp7npY7Uz2az1NtJo3lUSU2iUuqCo6lJNKtw\nkb1CGGEUuL6+3mw2sf0hx+Oxbds3NzemaoX6BUFwd3e3WCy22+3JycloNJpOp/xBApSCMFoHmTvV\ngvZZs5pSz5d2ni/Tquu6Owv3X9lZOGEU6Pf7qZtA2rbN5pDdQQc5UCnCKDIRRtFxYRhqdoDcbDbs\nD9kRdJADlSKMIhNhFKBnFIKPAVCxOvNGt/amB9AC7FSOMAyzEqd8iA5yoEG6tbQTgBZgp3LYtp0V\nNzUPAThOhFEADcNO5RBCjEaj6B8k0mQyGY1GRuoDYG/cgNgw3DOKzkpdS1ItiIGuCYLg4uLi/Pxc\n3bAhO8j5swQoRZ15g55RAEdNv9kSSbSz6CAHWoNutoZROzZl4YKiTVhLEnnQQQ4UtTNOiBoTBWG0\nYRimR6ewliQAGME6o8hEGEWnsJYkABjBPaMAsHstyZrrg2ZZr9emqwAgF8IogCPFWpLYg37GG4Aj\nRBgFcLxYSxKFyBlvm81Gdpxvt9swDC8uLsijwDEjjAI4OmqAlc2WUMjd3V107QUhxGKxOD8/v729\nNVgrAHrMhmkYJjChxYIguLu7WywWcm/x0WgkQ+ft7e1isQjD0LZteZC1JJGKGW9AWZjABKBzsgZY\nhRA3Nzfb7Xa1Wm2325cvX5JEkYoZb0BD0c3WMPSMoq1YUhSHo2cUKAvrjCITYRRtRYzA4a6ursIw\nVLvVS5PJpNfrvXz50lStgCZimB5AtzDAilIw4w1oIsIoAPNYUhSlGAwGy+XStm3btoUQtm33er3l\ncsl9xsAxY8y3YRimR8us12vHcQQDrCib7/uu6yaPq48cAA2G6QG0XHKbnLfeeosBVpQolkTZmQk4\nWoRRAHVLXcVpNpvN53MGWFEFdmYCjhljvg3DMD1aYOcqTlkDrMB+WDgMKIqlnZCJMIoWYBUn1IyP\nHFAUYRSZLMvSn8AFxZELw1AzO36z2TB3HuXiIwck7YwTosZEwT2jzfOkZbp2QKb1ei1YxQm14yMH\nJOmzRM1x4lmdLwagg4IguLu7WywW2+325ORkNBr91m/91o9//OPUVZxMVRLtNhqNshYOkz+z3hNg\nED2jACqUOov5n/7pn7773e+yihNqk7Uz01tvvcV6T4BxhFEA+5Mj75qDd3d3Z2dn0VnMi8Xi537u\n537+53+eVZxQm9Sdmebz+Xw+Z70nwDgmMDUMs+lxDJIj77JTM3nwF37hF/SzmFnFCTVTHznWewI0\nmE2PTIRRGCdH3qP9nePx+Lvf/a5lWZ/+9KejB//lX/7lb/7mb7LKYRYzzGK9J0CDMIpMhFEYl9qf\n9MYbbzw9Pf3jP/5j9OB4PP7ggw9+8IMfJAvh9z3MyrPeE7Oa0GWEUWQijMK4rP6kVM+fP3/33XdT\nZzG/fPmygtoBeWV9knu93pe+9KXYDSeDwaD+GgIG1Zk3mMAEIBc5LSkMw0I9mj/4wQ/+9V//lYnz\nOEKj0Sj6yZTeeeedn/3Zn2VWE1AnwigAnSAIomvfvHz5Us5Hzsm27Q8++ICJ8zhCqes9/e3f/u2b\nb74ZW//h/Pz89vZW/jN1BQkAhyCMAp2m/82aukroJz7xibfffjt25htvvPFLv/RLsYOTyWQ0Gg0G\ng5ubm+12u1qtttvty5cvSaI4BqnrPf3P//zPn//5n8fOXCwWt7e3rEgKVIQwCnRCLHTG+jtjv1nV\nyamrhH72s5/9u7/7u1h/0s/8zM+cnJzoh+NZwgnHJvaX0tXVVRiGqWf+7//+73/8x38wdg9UgTAK\ntFlq6Ezt77y4uPj2t78dO/nu7i42a14I8bWvfe2//uu/Yv1JX/3qV7/61a8yHI+Gkn8pabaqf+21\n1z744AP1z9jYPYBDMDW7YZhNj6jUpWfUwdQFQR8fH19//fUf/ehHsZT59ttvf/Ob3/zc5z6njr/7\n7rvR374xcu2b1CXrWccezSU7R2PrPzx//pwVytA1LO2ETITRLtBHTJF7A6Tvf//7ydCpWfvztdde\n++EPf6g/IvE7GG0VBMHFxcX5+bnKo1/84hf/9E//NOt8ViRFW7G0E44Ls0cL2bu5sobUYwe//e1v\nJwfZ33nnnd/5nd+JHby/v08Osi8Wi9QkKoRI5s4f/vCHr732WuygnJa033sEjlxyVtMnP/nJrBUk\ner3eH/7hHzKrCTgQYbSpUhNPoRi0s4Scc1wOqcMetaq6hHKjZP5apd7HmRoxf+M3fuMTn/hEbFJR\nGIaxXZEWi8WPfvSj/d5L1A9/+MN3331X/ZNVQtF6yfUfxuMxK5ICFXpCowghrq6u5C32Jycn77//\n/scff/zxxx8nD2aVkHpy8uC3vvWtX/u1X4t2gI3H4y984Qvf+ta3Ul+rUB32rtXh7ytagu/7+qdH\ny4+enDz48ccf52+u1Ne6urpKdjemrpeU7Kos6vnz58nfrM+fP089udfrvf/++6qX6Ctf+cr3vve9\n/FcWaIGPP/74C1/4QmyxiE996lO//du/Hfv/Mh6P33//ffXE1O8NoBFEjRGRntEmkX9wyz/ZhRBq\nycbT09PoweRf56oH7sDut09+8pO/+Zu/mToLO3V2tqqDvl8wf60Of1+pc8Y19c85dP4Hf/AHySWQ\nUpsrq1ap89b/4R/+Ibbbu0gbTC9kMpm8++67ybW+/9//+3/J1UMnk8lkMmGVUHRcoRVJF4uFfpwE\nOAbWLrXWprbYi8NdXV0lr2Bq55n86zzZA/elL33pwO635PHxePzixYtksePxeDKZ5OkXzF+rw9/X\n22+//frrr8e6MD/1qU8lc5isf7K/czgc/sqv/Ers4LNnz0SaZHOl1io6Dl6uZ8+exULn5z//+e99\n73sff/xxrL/z7//+75N9P/Jk0x984IisVqunpyf9BL5f/dVfTY6TqMEWukvRCKLGiEgYbZKsBfBS\n9Xq9ZJDKykwVefbs2SFJLtXh7ysZELOS9/Pnz3Nm3JrlzLiTyUTGd80gu/zNKiUTKkkUyJL1hZz6\nvZH1x7kqjYSKYyMIo0jaYyWd6vrb9lZKkjvC91WK5PtKba533nnnU5/6VKwL03Xd4XCo6deMhs6d\nCp0MdNP777+fvPd6Mplk3X6d/ONccxe+REKFQYJ7RpGk2Roki2bF8jyS3W9ZX7L5pd4EWdSB76s2\nhZqr1+v9+7//e549Nv/zP//z61//etENkAqtQs+S9cBO0+k0ee/1P//zP2ctnRZb9Fd/F/4eK5kA\nDVZb7MXh3n///eQVTO08++IXv5j/M5C/+y313srJZPLLv/zLyR6CA1X0vrICYupqmoXSZPLuzNTm\nyhpP/8pXvpI6Sq4fOk/twqRfE6hH6n/PQr0GyW+e1Pvad/ahZqFvFXsT9IwilVzcMU/n2b/9279l\nrdKczEz5u9++/vWv//d//3fszIeHhz/5kz9J9hAUuo8zZ60Of1+pc8bfeeedwWCQfF9f/OIXkyE7\nK03+3u/9Xp7mSq2VXLkzubrh6elp6kH13NQuTPo1gXqk/vccjUb5/zhPLo7xta997eOPP87Zh5rV\nY7rfKtGAKYTRJhkMBkKInOOzqas0p2amrBI+85nPxL5nP/OZzyTXN5FnJo9fXFwckuQqel+pATE1\neS+Xy6urq2TIzkqTV1dXeZpr53i6IGICTRP975k6fF/oj/OcCfX8/Pz29jbnJm2mRv8Jvsiltj5Y\nlEJdsp3js6mrNO+c0ZJ/hDfrTHU8tQL7zbMp930VGvjeY+hc3yw7DwJouuRXhFzaIvb79/C78FNX\nF0m9QWi/0f+sUX79JiBPB2+DgmMgmE2PLIU+HMZX6ikryeUpdo/qHT7BnDQJQEP/x3lqaiyaUJOr\ncORfJbpoQj1wAz/9SquFbm/NX8Lhr5W/AnWWUMPdwIIwiiz7fTiMZ6aKkpzx9wUAOeXcaaKUhJpf\noYSacxMQzTYie2z7vLPD9cCDWa+VvwL5n76zhDyfotq6nAVhFFnq/HAAAEqnv2sof0IttLpIIfl3\n2iu0gV/ObVBkH2rODtfUEvIfzHqtQj2+qR3JWWUe0me839N3FptF1Jg3LPl6aArL4pIBQNv4vq9m\nQQVBcHt7u1gswjC0bXs0Gr311lue552fny8WC3mOnDf50UcfhWGYLO21116LzYJ6/vx51gKoBr3x\nxhtPT0+xxafH47FlWY+Pj2dnZ2rO1ng8/uu//uvPfvazX/va13aWkP9g1mt997vftSzr05/+9M4K\nvP3229/85jc/97nP7Xz64+Pj66+/HltuVtXBtu3FYrHdbk9OTkajkVw85+7uLnrw+9//fv6nDwaD\nIAhiJaQWK0+Wpa3Xa8dx5M+15o3aYi9KwSUDgI7YOfMydReo/KtEVzf6f6DULVWzOlyreK1CPb45\n92cej8f5d+cqtHV21uZeOTuSszp3BT2jyELPKAB0WawP9eLiItljOp/P/+Iv/mJn3+pf/dVfJbv6\nSulDTXbNopDUftz8xuPxd77znV/8xV+MdaOmFpvaufv4+PiXf/mX9eWN2mIvSsElAwAo+deqO+T+\n1Kyd9vJv4FdopVXULPUuYcEOTAAAYKf8m7Qlz0zdmKPQTnv5N/ArtA1K/n2bU0vIf7Dots/JChR6\n+vPnz0vfOrsUyW5s1YNeD8IoAACNl3+Ttj0SatZOe/k38Cu0oV3qVsyp+zanlpD/YKFtn1MrkLrF\ndFbwfffdd0vfOrvSLme5iVcNCKPNY2mZrh0AoHn0CVX2tqYez9k1OxgM8mfZ1OSa2uGaWkL+g4W2\nfU6tQGpHsmbX6NK3zk59eqHbKjTklKY61HZDAErBJQMANFqebVAK7dt84ME9NgvcudDBzp0Cy9o6\nO2uL7NQbglOLzbpLuM68QbJpGMIoAKA76txp78DNAvd++oFbZ2edmbPY1Nj6+c9/vs68wTpBDcPS\nTgAAtFJ03S79wfxPz1Nscp+F6XR6dnZWW94g2TQMYbQetHM9aOd60M71oJ3rQTtXJ5pQ62xnrmjD\n8J+wHrRzPWjnetDO9aCd60E716POdmY2PQAAAIwhjAIAAMAYwigAAACMIYya9+Uvf/nNN9+0LOvN\nN9988803v/GNb5iuEQAAQE24C9iw1D2T3nvvvT/6oz/KOp9LVgPauR60cz1o53rQzvWgnevBbPqu\nUEn0xYsXb775phDij//4j+WRDz/88K233kp9SimXjHIoh3Ioh3Ioh3Iop+pycr0WYdSUL3/5yzJ6\nxnKnSqipl+bYPmSUQzmUQzmUQzmUQzmH4J5RY2QSffHiRawH9MMPP5Q/cPMoAABoPcKoGSpoytH5\nKJVNs24bBQAAaA3CqBnf/OY35Q+pifPFixdCiI8++qjWOgEAANSOMGqGCqOpkt2lAAAArUQYBQAA\ngDGEUQAAAJjzBBNk47948SL1UTWh/sMPP0x9IgAAQKWqT0M/8cz0O+2oFy9efPTRR1lTlNQdpcl1\n78mjAACgTRimBwAAgDGEUTOYLw8AACAIo8dJv/ATAABAaxBGzVA9o1/+8peTj8p7SeXS9wAAAC1G\nGDVDzUxKdoKqnUJTcyoAAECbEEaNee+994QQH3300Te+8Q3P81zXtSzLsqxf//Vflyckp9Lr+b6v\nCnFd1/O80uvcdNF23ruJYoX4vl9yLdvO933Vhq7r0oaHKL0xZVEl1a4Nymphvp/1ymrnUr7kW6yi\nr1/Z4AcVUdsiUkhSA/EvXrx47733ZDyVPvzww9Vqlb+orOtbqJAWm81mhzfRarXi/9GBHMdJbUDH\ncUxXrXmqaEw+z1GltLDme2M2m1VV9UYppZ1L+ZJvt4q+flU5BxVyyJNxoNVqlXpjqFr0Pmc50Y/U\nbDabzWbRz1ylb6ERor8MZrPZarXar4lo5wNFm0tdiLK+ELumisZUZZZd2UYqq4X13xvk0VLaOfol\nz/dzqoq+fqOFHFI9rpBJ8vp9+OGH7733nuwclVsuqf9XeT4f6qMQ++OvUCHtlvVfpVATqf/JsXZW\n7c8vlZ2yWjurbaFRbmOuVit+c8eU0sKa7wdaWyqlnbNOVl/yfD9X8fUb6/U/qHqHPBmH0P8nyd9F\noTmTfo6n8tq5lEK6TP1WTj7EL4yiSmzM1JG7kqvbQGW1sKZJ+dg/ld3OqWdmhbBOqejrt8TvDSYw\nGaPuGk69w1od3Hlz8Xq9FkKk/kZxXXff2rVHKe2csxDshw9qiWjMqu3Rwqm3M6pymMCXKn87qwZM\nfYr85Sh/USJp728M+UR5U8Th1WBv+pZI/Z/Gd1xZ9C0Z/aVCFNBwHIf2KUtZjRn9bHueN5/PDy+z\nHQ5vYX1IglTKJ1kmTtpZo9yvX8/zZOqQ0/NLKPGQblVUJ//Ib9Y9SYwB5ZGznTVjHE809cGy7nvG\nHg5pTP3nHFJZH1e+N/RKaWcaeac92lm1qnxKKd8bfOkco6Jzj2Kz5KKzCDt+o4xe/naO/d+LYQ7T\nHlarVWxZA1pvb2U1JmE0SxUfV7qEkspqZ1kOc/KyHNjOsacQRtspOj0t/7NS7xkliWoUbWfNyTT4\nHmKfVZLoIcpqTMJoltI/rnzyU5XSzrFbGPlaTjqknZPDiaV8bzCB6YjIOw6Hw6H8p2ap5BjXddU9\no47jqM/Ker3mHpqk/dpZ/X+zLEveCibvlWG7mv3E/nyaz+fMA9sbjVm1Els4+qUxm824UlFVfJL5\nPZi0dzv7vi/DRjLRHuqQJIsSRT8chf6Sy+pmV+GJvwuj9m7np+ztKximP0R0NI17Rg90YGPSM7rT\nIS0c+7uXT7tGWV8L/B7UK9rO+rBxSE340jHvkG8o/d3ZzAuJKuU3QbIbVZZDGD2Qak/TFWmDvRuT\nMJrTHi3Mrkt7KOVrgd+DO+VsZ/kZTsZ6wmgbRO9u2eMbaueHgD8KpQPbWU8/vQl5EINKtHdjchVy\nKtpQ6suHr+JCyvpA8jeAXs52FrntUQfuGTXJ9325qp/8htrj5hhWEs3j8HbeWb78gTuTNFzXVbfb\n4kA0ZtVKbGF1h+hqteKSxZTSzvI2XL5+NY7/G4NF702Sc2gcx9n7IxKduoQsh7ezeLXNkuu6ya+8\nY/4ffjzkB9XzPJrrcDRm1cpqYfV1UahjqTtK/CTv/FXY5bRaSjtnTfZVC+Dnn3WdYo/eVJSilAEI\n/T2jrPf7VPZAT+ooG42ch35/Ab6RCqmoMRmmV8pqYb4c9Epp55zrQB9e2+aq9Ou3lBamZ9QY/Xbn\nSvRR1Q+qPj2u6zqOs16v5/N5rNPO9321elGXVw8ppZ2FELPZbD6fr9drz/OiJ6sxuC43ch6e58kP\npOu6sb/O1ee2lD2Ou2CPxkz9VCNLKS0c22o167VSx1s6opR2VmcOh8PVahX7PcgOt6IR3xiHJFkc\nIucFiv6pl3XVoudH1xlNltBB9bRzxxs5p9i6WrE2ZG5HIUUbM+tTHUUfUtThLZzzj6uO95uW8kmO\nDhAnCxF8RVfzjSExm77Z8nxJiXwh6Sl7Byb+B9bQzjRyfpq1Wk1XrXkKNabmU60QRmMObOGsp/Ph\njynlk5wV/fkrVyn9G0Mq5XvDemLIpi1831c7A7mvGK5TG/mv0Mh7k0OWtGEpaMyq0cL1KKWdo9/P\notu3QGQ5zs8zYRQAAADGsM4oAAAAjCGMAgAAwBjCKAAAAIwhjAIAHscfAgAABR9JREFUAMAYwigA\nAACMIYwCAADAGMIoAAAAjCGMAgAAwBjCKAAAAIwhjAIAAMAYwigAAACMIYwCAADAGMIoAAAAjCGM\nAgAAwBjCKAAAAIwhjAIAAMAYwigAAACMIYwCAADAGMIoAAAAjCGMAgAAwBjCKAAAAIwhjAIAAMAY\nwigAAACMeWa6AgAAAI3n+77v++qfruu6rqs5U53svlJKNTzP07966pn5K18F6+npqbYXAwAAaCXP\n8+bzefRIMmIlz4mazWYyIB7CsqysV4/yfX84HKoXVc8qsSb50TMKAABQGsdxhBDJnkXXddfrdew0\nIYQ6OJ/P5/P5gb2EjuNEXyVLtF9WCDGbzeSRPM8tHWEUAACgNNHxbiWaRFerVSyq+r7veZ48wbIO\nGrX2PE92eXqep+ndVB20sibqzFgXaT2YwAQAAFAhFTQdx3l6ekrtNPV9X/WVHjJErgrX3A+g4vJs\nNtv7hUpEGAUAAKiQyoWpnaaKelSTI/NQoTbr5WJj9MYRRgEAQGNEJ3q7rmu9IjsX5XE5I0epcy5O\nUqFuSP05sfcVfcux02IvHRMbozfvCQAAoCFkelmtVqmpZrVaqX7BKDk+XimVI2PHVX0OKTzr/Qoh\nZrNZ8nxNzFNFaZ6Y+lB16BkFAAANI+forFYrmWZUEBwOh+v1WkVPFbyMTBKPvnRqRM5Pvl8hxGw2\nk+9avbX5fJ7s+lUNkuwcjfYfH1KlEhFGAQBA8zxFZgJFc5XjONF7IlVo09+vWbVDBsTVc1erled5\n8p+u66qOzOQ9ppqR+gNvSK0CYRQAADRMsqMxa9bOsdwWmcHKEI3Xqm81+V52Ru2s6Hkk8+glwigA\nAGiY1NWRDNSjevpRdfWuk2E0daReFXI8Y/SCRe8BAEDjNC56ZvVcJnsoY32Z6omx7eN3UluPep4X\nWzTqwBtYS0cYBQAAqFbWDKpkD2XWwLr+Xs/UnCq3Bk2+9LFFecIoAABAVXJuFr+T/i7P1Hyptgb1\nfd913eMcoxeEUQAAgOqoXen1m8VLyQ5O13XVaPseLy1/kCP1xzlGL5jABAAAUB0VIvOsqaTWE1U0\nU5RU+dG7QmNk9Ix2zR5bt6ggjAIAgO5Idk/GktzOE/agRtgty9Kcpr+VM5lTxatZSvP5POu56r1E\n91DV1dUEwigAAOgKGd3UP+XgdTRr7jxhD57nqcHx2Bqi6lUsy8q6tVRl2ViOzDPyrp5SykZQFSGM\nAgAAVMv3fRUE5/O5XNnedV35g+r1jG5tqqgsu16v5bMk9Sx9Vo4WeIRj9IIwCgAAUAPf99WGSVK0\nK3Q2m8ntPVOH0aNZdv2KEMJxnFiZSdECj3CMXjCbHgAANIjakD1GM1c96yn1kxvK+6/II1LsnORz\n5fnqBtbkE/UvWkb1q0IYBQAAqE/OEJnqOMfZD8QwPQAAAIwhjAIAgI7aOU3+wHn0yIMwCgAAusV1\nXd/35SKdInI7Zv4TNKxXSq92dczWmXtGAQBAt6zX6+Fw6DjO09OTXODTsqzoLJ+dJ+x0nCt6ZnEc\nJ2uV0xoUa1kAAIDmkpnS9/3oFKLoP3eegNIRRgEAQFdYlrVarTTJcucJKB33jAIAgE6Qt35q5iTt\nPAFVIIwCAIAO2dnrSbdozRimBwAAgDH0jAIAAMAYwigAAACMIYwCAADAGMIoAAAAjCGMAgAAwBjC\nKAAAAIwhjAIAAMAYwigAAACMIYwCAADAGMIoAAAAjCGMAgAAwBjCKAAAAIz5/wET7n2pWgf5/wAA\nAABJRU5ErkJggg==\n",
      "text/plain": [
       "Canvas('Canvas_AughySWzwtsLXcz4a6C6TN')"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mass_hist_zoomed = Hist(100, 2.8, 3.4, drawstyle='EP')\n",
    "fill_hist(mass_hist_zoomed, masses)\n",
    "plot_hist(mass_hist_zoomed, logx=False, logy=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Measure the mass\n",
    "\n",
    "The final step of the analysis is to build a fit model using [RooFit](https://root.cern.ch/roofit) and use it to extract the mass of the J/$\\psi$ meson."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:rootpy.stl:generating dictionary for deque<RooAbsArg*> ...\n",
      "INFO:rootpy.stl:generating dictionary for stack<RooAbsArg*, deque<RooAbsArg*> > ...\n"
     ]
    }
   ],
   "source": [
    "def fit(hist):\n",
    "    hmin = hist.GetXaxis().GetXmin()\n",
    "    hmax = hist.GetXaxis().GetXmax()\n",
    "\n",
    "    # Declare observable x\n",
    "    x = RooRealVar(\"x\",\"x\",hmin,hmax)\n",
    "    dh = RooDataHist(\"dh\",\"dh\",RooArgList(x),RooFit.Import(hist))\n",
    "\n",
    "    frame = x.frame(RooFit.Title(\"Z mass\"))\n",
    "    # this will show histogram data points on canvas \n",
    "    dh.plotOn(frame,RooFit.MarkerColor(2),RooFit.MarkerSize(0.9),RooFit.MarkerStyle(21))\n",
    "\n",
    "    # Signal PDF\n",
    "    mean = RooRealVar(\"mean\",\"mean\",3.1, 0, 5)\n",
    "    width = RooRealVar(\"width\",\"width\",1, 0, 100)\n",
    "    sigma = RooRealVar(\"sigma\",\"sigma\",5, 0, 100)\n",
    "    voigt = RooVoigtian(\"voigt\",\"voigt\",x,mean,width,sigma)\n",
    "    \n",
    "    # Background PDF\n",
    "    pol0 = RooPolynomial(\"pol0\",\"pol0\",x,RooArgList())\n",
    "    \n",
    "    # Combined model\n",
    "    jpsi_yield = RooRealVar(\"jpsi_yield\",\"yield of j/psi\",0.9,0,1)\n",
    "    model = RooAddPdf(\"model\",\"pol0+gauss\",RooArgList(voigt,pol0),RooArgList(jpsi_yield))\n",
    "\n",
    "    result = asrootpy(model.fitTo(dh, RooFit.Save(True)))\n",
    "\n",
    "    mass = result.final_params['mean'].value\n",
    "    bin1 = hist.FindFirstBinAbove(hist.GetMaximum()/2)\n",
    "    bin2 = hist.FindLastBinAbove(hist.GetMaximum()/2)\n",
    "    hwhm = (hist.GetBinCenter(bin2) - hist.GetBinCenter(bin1)) / 2\n",
    "    # this will show fit overlay on canvas\n",
    "    model.plotOn(frame,RooFit.LineColor(4))\n",
    "\n",
    "    # Draw all frames on a canvas\n",
    "    canvas = Canvas()\n",
    "    frame.GetXaxis().SetTitle(\"m_{#mu#mu} [GeV]\")\n",
    "    frame.GetXaxis().SetTitleOffset(1.2)\n",
    "    frame.Draw()\n",
    "    \n",
    "    # Draw the mass and error label\n",
    "    label = TLatex(0.6, 0.8, \"m_{{J/#psi}} = {0:.2f} #pm {1:.2f} GeV\".format(mass, hwhm))\n",
    "    label.SetNDC()\n",
    "    label.Draw()\n",
    "    \n",
    "    return canvas\n",
    "\n",
    "import ROOT\n",
    "_=asrootpy(ROOT.RooFitResult)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4QAAAJYCAIAAAC1p7+MAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElE\nQVR4nO3du7XrSGI2UOBfE4BCmDR6HIApSPKVgLwZn6AveQpApoxpWfIJOpo0JgRlgN+ovtW4eBEk\n3sDeq1evc3nAIlgkwe/UM62qKgEAgC38v61PAACA6xJGAQDYjDAKAMBmhFEAADYjjAIAsBlhFACA\nzQijAABsRhgFAGAzwigAAJsRRgEA2IwwCgDAZoRRAAA2I4wCALAZYRQAgM0IowAAbEYYBQBgM8Io\nAACbEUYBANiMMAoAwGaEUQAANiOMAgCwGWEUAIDNCKMAAGxGGAUAYDPCKAAAmxFGAQDYjDAKAMBm\nhFEAADYjjAIAsJk/bH0CfCZN061PAQA4v6qq1nmgZcNoWZZlWT4ej/avsizLf1j0HM5nljdHmqbK\nUY5ylKMc5ShHOX3lTC9k7GMtEXvLsiyK4vV6jTz+fr8XRTH7aZzSDt+sylGOcpSjHOUoRznfP9a8\nj1QURaMdNLaA1m8syzL8vx5YsywLtzNgb28y5ShHOcpRjnKUo5xJjzXjI+V5HsJllmVFUYzsf290\n5WslHba3N5lylKMc5ShHOcpRzhTzzKYvyzJN09frdb/fq6oqy3L8SNA8z4uiqKrq+XxmWfZ4PL4b\nphDib5qmaZqGMgdOdcCYwt+24MbjB86kfkofFQ4AcB7VHO73e4ih04VI+umj9z275/M5/uC+Cuk8\nrO8kn8/nyDMJsiwbX3g4mVGV8o5ylKMc5ShHOcpRztLljLFeG+xCyrK83W7h5/v9HloW66NRG0+w\nPpagXVq7IbN+fGjuHRhR0DiZ8EM8vl3V7cLjmfeNoN1b87tylKMc5ShHOcpRziRLp92l9T2R2ELZ\naGUMN45sx42BsnF8DLJjTqbvTPpuj4V3Nqae4CU7BPW8DvW8DvW8DvW8DvW8jjXr+dg7MMW2w3bn\ne57nIdWNX2GqLTZqtltA27f3nUye5+GWxpnE+zZaQDsLBwA4pZOE0c7c1pf2kiT5aKX9gWGm9XWs\nBk4mPlz9TAZGC0yP0QAAh3DsMPqpmAXHhNHhg9sJdWAWfCxhIKe+vREA4HyOHUbDmlBVzwDbmPza\n2a6+mlLQjokfJddksKWzbbjwzpZUAIDzOXYYHVCfll6/Mfxwu90aneCfrm/6aV5shNThe2kZBQAu\n4pxhtL7EUj321X++3+/P5/P5fNY73Osp8LtWyeEc2R4GOny8llEA4Nz+sPUJzKwsy6IoYubrW4K+\n3rMf1hYNS36+Xq+Pto9ayD//8z8P/Pa7HaqivlEN1KmldajndajndajndajnMSZGhZWdqmU0z/PY\n/x4W72zEyrIs+8aYtufCbxhJ//rXvw78duJqXqs9CwBgExOjwsppYZEw2tm5XN/evXPC0MRHTNO0\n3iD6dQf3pwsqtfPuwMFx5OjwkNNP504BABzUzGG0KIo0TW+3W2OkZpqmj8ejnvPChKFZxkQWRVHf\nhLNqNYiO1JhjNO9oTouGAgC0zRlGi6KoLwIfxaSYJEmWZfXM14itXyjLMjxo6JdfYtei4cbLaMxK\n9THgahkFAEhmDKMxFCZJ8nw+26kr7LRelmUYuBmnFtWj6hfC3bMsextqQwPtQHNsYymo4bzYXjdq\nIDgOhMuRSRcA4JymD3ENYizrHAAbes8bYh4NOfULcVWmkceHg0Mbat/J1E+170l1Hlz1P9nOcgZO\nfuA8Z3zJgANLko7/AGayZt5Iq5kmTIVFBO73e72jPK732fconfcaL6zHlAxuH5/8PEE+Ht95no1T\njbc3jo8rJjSeVyx/TCGxnEazbmch9bvM9ZIBB9a5bouLAzCTNfPGxmE0BK8xnewDD/pWfdhA/S7t\nUZ71I+tnGA4Ov+ocjdA+pZCP40ZQSVcl1EfZhuPr/+wM6MIokCTCKLCsNfPGbIveZ1n29YTxNafp\nVFUV82LjhDvDZVgD//XD8MHh9pC/G3O5Ol/REDfDkfXjv24qBgA4ltlib2zkaxQ40BEf2037gt1y\nyh+SJMl/GL5LURTx+LdJMZY/e+FaRoEk6WkZ7eSKAXzukN30yc891DFR9SXOvmGaDBNGgSQRRoFl\nHTWM1vNlUhtkGcdNtm9JtmgWPTRhFEgSYRRY1lHDaBBn/IwhiX5KGAV6mdUEzGTNvDH/3vRhTfv7\n/d7YXbMuy7Ip+3YCMFLY7GPrswDotVLsjbNzVnisc9MyCvTSMgrM5Njd9CxKGAV6CaPATI7dTd/J\nZusAi6uqjv8A9m2RMBrW10xrbrdb/Z95nounX0sHbX12AMDeDWeJlePEzG2wjdWd3jKb/lO66QGA\npR11zGh9p/WktqpokiT1ptD6IqOJPPohYRQAWNohw2i9TXRMvrQD03eEUaCu0Znm8gDM4pBhNK51\n/1GBYVCCxtHxhFEgMnseWMghZ9OHJHq/3z+6VzjeZCaAT/VNMDCPETiWlZZ2GiaMAnxkuHdeHgUO\nZLYwGjb//DRW2pkJ4FOdSVQeBQ5qtjAaAuXr9RqfR+O0emEUoFea1v8baBM1WhQ4ojlHp8ZrZFjU\nqSiKviOLoohJNMsy3fTjmcAEl/Nz/EyT368AnReD+uGuFsB3Djmb/rfienqGsiyrry1av10S/Ygw\nCpdTu66+TaKJMArM4ZCz6YOqqjon1Hcm0fv9LokCzKv+9WHkKLB/S8Xe8ofwz9frFWY4heGheZ4b\nJ/odLaNwOT8SZb1ZdPjqrXEUmOjA3fQsTRiFy+kKo28vA/IoMMWBu+kBWFqV6H0HzuMPyxVd76Yv\nyzJ20CfWcgIYLzZOpK1bBu8UG0fTtJVfNZYCuzF/G2xZlkVRdM5YqsuyrCgKqfRTuunhmr7odv/p\nLsIo8IkDjxnN8/xtDK2ztNOnhFG4pu/GgPbmUZcRYNBRx4ymaRqT6P1+fz6fVZfn8xlm1idJ8nq9\nBtbGp1M6aOuzA5YlRgLTDWeJlePEbLG3LMvb7ZYkyf1+H5kvi6J4PB5JkmjqG0/LKFzQlKnx8b71\nltH0t6JcTIBuh2wZDQE0jAQdf5fQRKpxFGBNoZ9q67MASJIZw2jooP90QlI43rBRgDEESOB8Zguj\ncRgoADOaa+xWmlRJ9eM/gN2YedH7MAb00+Mt8AQAcE0zjxlNPkmW8UhjRgHe+q5BUzMosHOz7cCU\n53mWZa/X6/V6hRUB7vd7jJt5ntd3Y6o3oN7v97nOAeBkLNcGnN7M8/bjak0jPZ9PffQfsbQTXMqU\nRZ3mLQS4lEMu7RQURVFV1ZjGzvv9XlWVJAowhhAJnNWysbe9ZpP0OZGWUbiOGVs0f1/63vUDGOHA\ne9OzNGEUrkMYBbZy4G56AGY34zeCGVHA3mwcRsuytP0SAMBlzba00xfSH3+h63cGaJi3CbOqtIkC\nO6WbHgCAzWzZMmo7e4C3Zu86SlPTmIAdMTX7YMymh4uYff67pe+B8dbMG1u2jALQaYnxnYaNwheK\noqjPtM5/mF5ynMNdFMVcRx5VxaF4QeEKkuT3/5YolqU9n8/7/R5Go2VZdr/fn8/nd+XEIW1Zlg0X\n0jg47HS4rUYlTCxkTCWEg6c/aFVVA4MJ357GGHG7yrdFxTOZ+Ih1u4qIazzS8/mc/poRiJtwBQul\nxoUyLg19ISbLsunl9BXS96BT0tjXGbqqqufz2Xk+nxb40ZPqC1WfVkLj5LOaueq2fsJv3xjfvX8m\nSo4eRut/xHSa5U+Ka1rzzQFsRRg9rvrXXwhzsQHso/jSaONstA42Dm48RBBv+foL9+u8VQ9z7fP5\n6AQalTBQk+2Dv6uE+sl35r9Z6rZxzmPOZ+Wm7jXzxsyPVH+FxhBJP7XmmwPYxKKRURhdVMwN7RDz\nUU9rXzmxkPq358CDfpH/Gnf/LgB1PtmB8+wUE0XjHDoL7zu4+rwSxkSUemAdWWynMT31S/TRj3HU\nMNpIovFPk+cPjb/t5NEvCKNwesLoccXvwfavYnwZ8603kD/aeW7g4I8etPOxvgujfblwoH76Chmo\nhHr5A5Xw0YMOhNqvjxz2NqC/PWAhhwyj9b8Sxrzp5/qr4mpUF5zeonlxnZ76MLou/lxvpIhfEI32\niz3MtpluuBFr/DMdyB/thxgOK1OizHevy5hEPr4SOo9sP6mBp/lRIv8omQycYX0y2fBTHn7PTPxz\nYopDhtHvmpHPdA1ahzAKp3eCMBq/PpMu7e/pKYFpV0LvX186GfmVN5w/2lFvStYZ9t0X9JhE/va1\n/rQSQs0PHzw+jE6MJQNDFtsHD78xtuqjrw4aRr97/cILJoyOJ4zCua0QFlcLo40E0B7KFW6/SEfZ\n+FQxnMPavx2Od+uH0THn8zaMDvetf9TzPlfNj9Se3tSYxNa+y0CNDVfmotb8PO5ib/r6crIAnEZV\nVXF58Pp63VmWxSt/nucxBJzs6yCsVV4URZ7nr9crSZIx03zrNTP+sUL5429vyPM8bUmS5PF4dN4+\n0cizGjjb8EPfGybWfJqm4bH62ukb92qU/4XH45H8iI+hnDzPi6KIJ9BeuD7E5XadxPM57Vr3P8y2\nA1OWZa/X69PrSDh+lp0MAI7tt+/439pCqiRN0mU37lx6k/p2X3z4pkhaGeKs3wK3263+z+fzOf2Z\n1nNY+Pl+v4cAFG9pHzy+2LrwYg2v1biQiX+WfFfzww9aFEWo57bqxwcppsZ2UXme9yWloijCCRdF\nUc+d8eezfkCi2VpGQ019lEeLoghv9NPXMsAFta/tu73ad7YLtk08/9vttkTTb4wst9st7l0Zgun4\nBshwx4YkSe73e+ftb+3qtV6o5tvCo/TF95iUOm9PfrSqRhv+MbCy2VpG418Mt9sty7LQKN15ZHgr\nxxoPB891GgAMW22T+gNd2yf2GveJDWZlWYamr9vtNkv7aMPz+QzlPx6PRqCJrdErmxj+PkrSbZ01\nX73rCBh+0ParVpZl4/jwz9fr1RmBBuokNm+3Dz59H32SzD069dNHP8HcyZXN/pIBe/FbD/2P2UVL\nTjJafo5UknTNAhmYd9J3lzMZ+c0bBxcO/7ZdV40mtPANu88JTG/PZ3iK0nAV9R08ftbUyPdh38oG\nb7XLb8+p33AefbDmQ8/WMhqrbGBQRfvFONDfzQBLS5Pml1mYLFItOrSTtXzaSNkeA/r2+E9PaQnD\nT3OWruePnmmsw7eVHwNMHO3w3Wm87fIdGMHyeDxCU+h1+uiTGbvpozD8tnNkSZxWJoMCjLF0DF16\nDtMsQudpo/0i5Ia9xfTwx0NfU8vIruf2FKW6T2d8b5Vpps+X75uVlbRSYOyL73s/fPpnwMiD24fF\nB/qib73eU3+tPvolwmggcQJ8rUpWGdR5ZGVZzhiwRi5XVF+RakBfi+Z3C868vTE2/bSDyyaZJmbu\ndj2MP5/hRN5I2MMHt48fEBNhnudvJ9f3/aovyw68UklPu+xVotRqAwKYhZcMTmzNjeMX3ucpSeYb\nMxoyRONXfQd/Z2SufTvocJbNh4bLCbfXh3KOGUX69hE73e/3rze1TybvTd938p1bFg3U7Ud7kFa1\nttWB4xurlrZv76y3t2+AkTuIruPrt803j7XaIzGLNd8cwMrWC6NJ0pwpNeuj9n0fzxVGw1f+jCc8\nl/gE24Ej5ox2qu6sk85y+vJZZyHxZNbPNJ1PdmAiUef594XIzoPjI7afbF8N96kHzc4tRts7F3Q+\nXOOOfW+AgZK3nc+35kcsrXY24IZhb/uSvKBwXPHzvfjnOE3T+ur6cz9q39DJ2BHZvlINjLZsjxkN\nP8Qez9A5G7s+w45HnQevoH6Vjvkj9tu2O/rj8Y06qc8GDhml/s9GP2/94MZ2PiOHFnzUjz/m4Pi8\nwmye+ipIfa9++1exxz9OCYpPs/1WGa75dqUNiINQGwXW+99Dge0pho37Nu749jTqz2LRL/QxQ1PW\nSxSrxV5m4SWDs1phV/r6gx26ZbTxz1BsbBJr/Lz+ZbOv07+zhXLg67hz79C+Zs7Og8evn/jRANyR\nZY6/78BvO0+sr8mw74S/a2IcqJNYYOdp9209OqaJeqCJd2XjX+jptIweTJp6yeCc6u0Ua7aMJrFx\ndK/XlkbLaGh5qrraosKNaZqG9qdw+8imwXnFqSqh1XbKpN44nWVgN5nGweGx9jAROy6tM6US6uvz\nvK2EGWs+CA8XH318gfXnnuzj5fjImnlDsjkYYRTOar0++qQ7jKa/PfrurjCNMNrZ7V7vLY1hNHRe\nbxJG4ejWzBtLLe0EwJ5VSdpYY3+HMbRTCKZbnwUwm/+39QkAwO/izimdv429pQO/HXkwsBPCKMD1\nVNVPowEa/9xUWZaNPaXrU5j71rofmGsijMLOzRNGwzoasxRVluXInTAATsNlr6HdMhpX9un8uqnP\nMolCihVGYefmCaNhdHmapm+3zxoQYmiY+TjLWQEwxq6icIiVj8ejPhM5fC8MbKWY53mWZfUm1XCY\n0aVwADMuE1UPkePXx3o+n/U7fr3fwP1+j+VkWTZ8As/nMx7/9uB24W9P8qPC6zXwtvB5XzJgJ9bc\nCLT9oHu7rvTtcBNuH3nHMddqoM+aeWPmeftlWRZFUR/fE3dNiB0lcQGw+mHJh7sj1NV3nmjo3Myj\n7/i+qugcNtC3Vkh724aBM0lq20uMKTyxtBOc1KrrOu3gcceor08Zb0lGdLvHpZ0WPT04t8OvM9qO\npMOmXDXq4e9+v4dxAgPbjjWOT36Oxe3aGNiLrH3a7cLrx48pPJ5MXx4VRuGUhNEZCaMw3ap5Y9F2\n13pveF2WZXN1oPQ9kThOqLEZWufx8eBGl3oMlI3b4zMaczJ9Z9J3eyy8s36WfsmA9W3YXb7PbvqJ\n2hdt4FNr5o1jX4H6QmTQjozx+HbO68yXA5G9/bgDJ9O5HXNfoo2Fd+4pLIzC+Qij8xJGYbo188ax\n1xmNHdmd3THxxvZh7SFHjTFJde2h9FF97OnAyXQWHrvj28WGG8ePcwDOYcO+8l1NqJ/o+Xzqo4cD\nOXYY/dRA/msn1/hD52D5dkIdWNMqljCQU9/eCMAYLqFwLMcOo0VRxMbkzt+GHxoXppHXqeEw2jaQ\ndD8tfKCZFmBGZ5q3BBzUscPogPq09E/vOOawT/Ni4zSG7+XPeriOM/WPA3zhnGG0vsTSyG73zkK+\neOjhwtvDQIeP1zIKAJzb2cJoWCE5JtGDbgT3z//8zwO/TadZ7VkAH9m8x9zlAU5jYlRYOS2cKoyG\nGBp756uqqrc7ftqxvmFf+V//+teB305cQGG1ZwEAbGJiVFg5LZwkjJZlmaZp7AR/Pp+rdXA3Muvw\n48aRo8PJ+NMRBQBf8ycqsK0zhNGiKOqbcDYaRD8V7zvvaE6LhgJtesYB5gmjRVGkaZrn+frrDJdl\nGVaeD/3yY07g0xw53HgZjVmpvp10tYwCAFc2Wxi93++v1+vxeIRxr6ul0tAmmmXZmIg5kBfjCbdX\nJ+0sub1u1EBwHAiXI5MucG76yoHLmq2bPq4/f7/fsyyrp9LlolV726Svj+8sIWTN+p6fjYPr4bId\nZ9u/qh8f9nDqTMbfrZAKMJFhA8D60uUmTIUYGsNWlmVFUczb75zneSh/YPv45Od0GFcrqD/xoihC\n4rzf7/WD43qljds7C6mfT/32vkJiOY1m3c5C6ncxIx5OI4a/bT/W9QzqAgMk6+aNNR6pLMuiKJZI\npSPXwXo+n/U+9/psp+TnvZra7aMxGmZZFgqJDaX1Ytun1Cg86QqXMQTH4+v/7BzqIIzCmewkjCZ7\nOhNgD1bNG9NXohrv+XzWu56zLHs+n1MKHPkcG4/SuRJ+mP/UqbO7vO/M+5bZ7yu8s003rAnQ95Tf\nVApwEEny+3+b28+ZAHuwZt7YppmtLMs4Cz7pbwVcVHjEsGNT8Pb40G46ZtGA8ofZC9cyCqexq85x\nLaNA3dm66QeExJZ0TfqhkzAKp7Gr/LerZAxs7kJhlE8Jo3AawiiwW2vmjTPswATARAIosBVhFACA\nzQijAABsRhgF2MCe9zra87kB5yOMAqyuEffSVAAELusPW58AHxved8pceziWKtlLDK0qkRiuYuQe\nlusQRo9H3AQApnibJdZMqxuE0bjQ/ZitiQAAOLFlx4wWRZGmaX13pTzPb7fb4/F4PB63200YBQjS\nNN1VxxnAOhZsGW1fVcuyfL1e9Vter1ee56GhFODKdjUCJ00tgw+sZKmW0dgammVZbP6MN1ZVVVXV\n/X5PkqQRTwEAuI6lNh4NzaJZltVbPcON9/s9ptL2LQyzNz2cwJ43go/ntrcTA9Z0nr3p6xEzptL6\njVmWLXoCAADs2bJhtD4/qXNgaDjAmFHgmrQ+AiwbRusp8/F4JEkSxok2DjCnHmBvzOwH1rFSN31M\npY3cafYSwI7YpxRY3VKjU8uyvN1uyY/Z9KFZNKmtXRIPSHa2oMnOmcAEJ7DfSUJpmiRJmvx2Wr9t\nVbq7swQWt2beWGqd0TzPsyx7/RBufD6f4Yd6Em103AOcm6ZGgLoFu+nLsqwHzefz2R4balEnAIAr\n26zPtyxL85a+oJsejm6/ffRJdzd9ajAVXM+aeUOyORhhFI7uQGE0CXl0jycKLOsMi97neZ6m6dsF\nRIuiSNNUTz0AwDUtNYHpIxa9B9iFqkqSpKrPstIsCixszjBaFEWMlWEG/dsmz3CYwaMfSQfn4urE\nhz0zlR7Yg+EssbI5w2hZlo1F7EeuaS+MfkTcBACmeJsl1kyrc4bReqYMMTTLsrf3KopCGAWuxh+V\nAMFSU6XyPH+9Xp1rizKF2fRwaLueSl9zlPMEFnKSHZgS/e8AB5da2QlYmGa2g9EyCod2lBbH+mix\nnZ8qsIQzrDMKQMOeZq++IYACq1k2jIal799a9BwAANitBRe9lzIBOml3BIiWncAU3O/3hR4FAIBD\nW2p0amgWzbLMVp/zMoEJjusos5eCY50tMK8188ayYVRsmp0wCsd1rHh3rLMF5mU2PcDZHHcU/XHP\nHDiEZcOoPnoAAAYsFUbDpKWiKBYqH+CgjtLrfZTzBI5uwQEBYXt6c5jmZcwoHNQRh2Ae8ZyBWZxh\nb/qyLEMYfb1eaZpmWTZ88EKnAQDAni0Ve0MSHXmwpr7xtIzCER10q3cto3BZZ5hNX1/0HgAAOmlm\nO5i3m6x6QWGHDtrEeNDTBt4as2f74ceMshxxE1hZmsqjcCpvs8SYtDoXi94DALCZxcNomFaf/hBW\nHi3L0hKkADunNRRYwbJhNM/z2+3WnlZfluXj8UjT1KJOwOnZThNgwIJhNE3TGEOzLKsvNRrn2t9u\nN3kUAOCylgqjsRf+fr9XVRU66+Nv8zyPI2dvt9tC5wCwK3q9AdqWCqOPxyNJkizLBsaGxjyqcRQA\n4JqWHTP6dpZS6LsXRgF2zshXYCGLT2Aac5gwCgBwTcuG0ZEp096hwFkdvUHROFdgaRt304fp9sIo\nAMA1LRVG7/d7kiSv12sgj8adpoRR4PSqJE3SH/8B8EO63EbneZ7X1xlNkuT1emVZlud5WZbxV8/n\nUxgdL00XfMmA2cXkWSW1DHqoT/HvT+FIZw1MsmbeWPaR6nm00/1+ty/oR4RROJYzhdHkYCcOfG/N\nvLH4BKbn81nfeykKi+FLogAAV6aZ7WC0jMKB/NSmWGsZDT8d6LOspx6uZs288Yd1HgaA6EAxFGBp\nS3XT53mepqle+CWkg7Y+O6DDTwNGAbY2nCVWjhPLjhl9PB5pmobp84s+0KVUg7Y+O6BLVf30H8Cm\nhrPEynFiqTBaFEWct/R6vW63m1QKAEDD4qNTi6KoryoaWNHpayYwwYGcZt6P1Z3gas6zzmhdURSP\nx6N+S5ZlRVFY8f4jwigcxZkC3JmeCzDGOcNo1E6l0tV4wigcxckC3GlaeYExzrPofaeiKKqqCpvX\nA5ye9AYwYO11RsuyLIpieI9QAAAuYqWW0bIsw8qjt9stJtEsy8KmoPM+Vpi2P3AmXyysFca2hgPG\nrAkQj8/z/O1UrVg5FhwAAC7n7UJTUwxvTL+E0PufZdnwAR9VSOdhfQ/xfD47j38+n53Hd9bPwPkv\n/ZIBc0mS3/47vCSpkuT3p3OSZwUMWTNvLNVN356llKyyolP7QRtiu2NnCmw3qcZbsiwLP4eHeL1e\nRVE0nk5ZlrfbLfwcU284/na7Va1cm+d5aCeOhYdlsF6vlyZSYJ/SpLKhFDCnhUJuPerd7/e+dsG5\nPJ/PepPncMtiMrppNpbZOD4+u87CG7fHttLGWfXdHgvvrLTlXjJgRrEd8Qwf2Z9bRn9rHAVObc28\nseCY0SzLQpxaejHRMBT1bZvoF2KZ7RbQ9u3xxsZIgDzPwy2NaVvxvo0W0M7CAbalNRRYyFJhtCzL\nMC9nofIn+ujEBoaZ1hPwQIiMD1fPnbGDvl1suNGaA3ACZ13XqW+6J8CnNlhntC5k1omF1Ft6O7Nd\n/eHCD2PC6PDB7YQ68ERiCQM59e2NAPtRnTVlA6ubJ4zG9ZI6fxW2p+/81e12izN+1lRfTSlox8SP\nkmsy2NLZeQIDhXe2pAIAnM9sYbTvV2Fa/U6GP8bzrC93Gjwej496nT7Ni42QOnwvLaPAjlTVb//9\nkCaaRYHZbNxNv7J6BAxz/BvT8Osp8LtWyeEc2R4GOny8llE4IsMpAcZbezvQnaiPdgqbJIUlP1+v\n167mXf39739v3zhx3oCRXsB3qkrOhmM41hTDy7WMxtWz2r8KP8QRBXuIpH/84x/bN05czWv1JwHX\n5QMHbGJiVFg5LVwrjI7x6YJKjcw63LEeR44ODzn9dO4UAMBBCaO/a8wxmnc0p0VDAQDahNH3hhsv\nozEr1ceAq2UUACC5VBiNi6H2NWo2FgodzovtVUUHguNAuByZdAF25VCzI4Bdu1AYHdgGKenJiyFr\ntne97zw4Ftsuv3NSVOeG9cFH6+cDuyKlAXzkQmE0qfWkN/Ji2Asq/Fz/VV++7Dw46QmvZVmGcNnY\nQbRv2r7tl4Ddsj4AMLuZ1xkd6AFv/2r9sBV66pMkeTwej8ejPcrz+XzWjx+qCYAAACAASURBVM/z\nPMuy1+v1eDzi+qMxazYOrpefpmmInjGJJl0tpvf7/fF4vF6veHwsvJFcgcOR2wDGSGdZSirs+fn1\n3WdcziqsXZ9l2UDS7VsJ9vl89u0U3+5M7zu43sha1/ccO6vufr/3baCapvO8ZMBC4gXmrJ/U0z9B\nIFk3b5wtjI5U/pAkSf7D8F2KoojH9yXFdvmzFy6Mws6dPqud/gkCyRHDaDKtz90CRuMJo7Bzp89q\n9b6lsz5H4JBhlHUIo7BzwihwAmvmjWvNpgdY1BXWdRJAgXkJowAAbEYYBZif5kOAkYRRAAA2I4wC\n8KUrjJEFliaMAgCwmZm3A2UFfTtIBRZ+gq389NE89QpPVaVNFI5tOEusTBg9HnETAJjibZZYM63q\npgeYWZXsqMkBYOeEUQAANiOMAixoVwOzlnD25wcsThgFWJBB3gDDhFEAPiZjA3Mxmx5gBj/1Vktq\nAKNpGQUAYDPCKMCctIoCfEQYBWASE+qBKYRRAAA2I4wC8A0DEoBZCKMAAGxGGAWYyqBJgK8JowAA\nbEYYBZjNZYdRahsGviaMAvC5NG0mUIEU+IrtQI8nHbziV5dtmQFWVyVpmrjmwPEMZ4mVCaPHI24C\nAFO8zRJrplXd9AAAbEbLKMAke+rs2lIaOu313gAfEkYBmIcYCnxBNz3APK6ZxKpEyzAwiZZRAD5X\nj97iKDCBllEAADYjjAIAsBlhFIB5WFgA+IIwCgDAZoRRgO9pC0yuuowAMBdhFACAzQijADPQOgjw\nHWEUgNkYtwB8ShgFAGAzwigAUxmlAHxNGAX4ki5pgOnsTX886eAXYKWBAgAYNJwlViaMHo+4CXvj\nQwkcy9sssWZa1U0PwJz21OACHIAwCgDAZoRRgM+lP3dhpan2QGMVgO8IowAAbEYYBZikSq7eJgow\nhTAKAMBmhFGAeexq3b5tqQlgPGEUYB7WAAb4gjAK8LE0kTs7SOPAF+zABDCNCAYwgZZRgO8JogAT\nCaMAAGxGGAUAYDPCKMBnrFs0hloCRhJGAQDYjDAK8CWzl9rUCfApYRQAgM1YZ/R4hrcctAcMADBs\nV9sXC6PHI27ChvZ0AQf40tsssWZa1U0PwEzStJ7W01R4B94TRgG+oYsCYBbCKABzqhKtocAHhFEA\nADYjjAKMZQDkF9I03dW8XWBvzKYHYEEWAAGGnbBlNE3TPM+HjynLsiiKPM/DwUVRDB8fDw7Hl2U5\n8vgxhZdl+VHhwObkq5HSRE0Bb6Qn+5u1KIrH45Fl2UCkC8e0b++ris4Opr6HKMvydru1b38+n50R\nOc/z1+s1svBwMid7yeBA4sXAp3BY/aqpruCI1swbZ0s2ITgOhLl6WLzf7+GWGAfbtRHDYpZlIU3G\nIHu/3xutnu3C68ePKTyeTN9TEEZhKwLWRwR3OLRV80Z1Cs/nM4a/JEmyLOs7svOJP5/PcOP9fq/f\nHsts3J5lWWftDRfeOKu+22Phz+ez8/z7nhqwqCT5/T/eUldwaGvmjTM0s7W70d/2obc7zWMjZb1C\nYsntWgq/qjeOxsLbLaZxYEC9nM5HrBfe+Sy0jMJWNPV9RHXBoa2ZN044gWlAzIjt4Zvxlnb+q7e5\nNtTHnsY7tmcsdRYeu+PbxYYb22NJAQ7Hsk7AsDOE0XpLb2e2iwbyX72Bs/FD58SjdkIdmDIVSxjI\nqW9vBAA4nzOE0U+NjHrDYbRtIOl+WvhAMy2wCc17n9I7D4x0xTA6bGT++zQvNkLq8L20jAIAF3Gh\nMPppS+d3rZLDhbeHgQ4fr2UU9kaD3xe0KwMDbAe6a3//+9/bN07c5dlkfPjA7x+36vdbfIiAfZsY\nFVZ2oZbRTzvW99BX/sc//rF948TVvFZ/EsBFud7AViZGhZXTwoXC6EIamXU46caRo8PJ+NMRBcCi\nbLAOsBxhtCnmv3lHc1o0FACg7Yph9NMcOdx4GY1Zqb6ddLWMwoFUSZokSZqmxxqPtQcqDOhzrTA6\nkBfjcvTtXZo682J7VdGB4DgQLkcmXWA/DL8GmNG1wmh7m6WoM/+FrFnf87NxcD1ctuNs+1f148Me\nTp3J+KP18wF2S2gH3rpWGI1Z8Ha71W8viiLkv8Ymn335Mt69cXtneC3LcrjwRnOp7ZdgR6rqpzzV\n+CcAk6Un62zK8/z1emVZ1pfkyrKMUTKkwxgWO+8VCgy/DTExZs3n89nudo8jyRqFJ10LfBZFEUsL\nx9f/2W5hDeWf7CWDnasPdvTh+0KsQLUHB7Jm3jhbsnkbRpOf82g0cJeYR+s6k2hf4Un/UvP1PBr1\nJdFEGIXVCaMTqUA4ImF0DSHtlWWZ//D2+JBW8zzvS4pR+cPshQujsDINexMJo3BEwii9hFFYmTA6\nnTqEw1kzb1xrAhMAALsijAIAsBlhFKCXfYMAliaMArAS4R5oE0YB3jPzBmAhwigAyxLlgQF/2PoE\n+Fg62NFl4ScAYNhwlliZMHo84iasY0/XaoA5vc0Sa6ZV3fQAAGxGGAV4Q2/EJGlab2ROU23OwE+E\nUQAANiOMAnTQeDevKlGhQDdhFACAzZhNDzDEgNHZpUkVWp6tDQIkwigA6xNDgUg3PUCTAaNLMGwU\n6KRlFCBJkkYC1W43n3ojqDgKtGgZBeilM3kh2p6BSBgFAGAzwijAT1J99IvR0gy0CaMAvdI0TfUo\nAyzJBCaAblWSaspbTqp2gSRJtIwCALAhLaMASZLUxjOmrVuYT1WZSg/8RMsowO/kJICVaRk9nuHp\nFDbZA47CsFHYyq6mZgqjxyNuwgp8zoATe5sl1kyruukBWJWgD9QJowC/2VO31VWoc0AYBQBgM8Io\nQJN+5KWpYSASRgEA2IwwCsCWDBuFixNGAZJEJALYiDAKwAaqpBb/0/S3/4DrEUYBfmJuzfrSRKXD\ndQmjAABsRhgF0D+8jZ966oGrEkYB2Is1t8MGdkIYBfidAaNbCcNGKy8AXM8ftj4BgNX93Pxm9gzA\nhoRRALZQVUmSVPU/DTSLwiXppgf4jSwEsD4to8czPMDfiCvgoNLU3wOwkl1NFhRGj0fchBnVB4ym\naerzBVzB22vdmmlVNz3AbyTRTah1uDhhFCBJRCKAjeimB67nR/Dc06ApksSwUbgkLaMAAGxGGAXQ\nGrcx9Q9XJowCF6WPHmAPhFEAdsQfCXA1wihwdfqIATYkjAJXpPltb/xJAJcljAIAsBlhFLg0DXI7\npN0aLkUYBS5H1gHYD2EUgF2ot1L7gwGuQxgFrksfPcDm7E0PnFqrhS1NJND9qpL09xcovnb+aIBT\n0zIKwB75swEuQsvo8aSDY6kqTQgwjs8KcFnDWWJlwujxiJvwNY1t+1fvqU+Tqkp29JUJp/E2S6yZ\nVnXTA7Bru2rCAWYnjAJXpINhzxqtobqD4NyEUeAq9NEfkVcNTs+YUeDUflpIfbvTYKT4enmx4DK0\njAKXo9d3/7xGcB3CKHAJ5sAA7JMwCsCu+UMCzu2KYTQdVJZl+y5lWRZFked5mqZ5nhdFMfwQ8eBw\nfGeZncePKRyYQv8vwK6kF1wyY3jJuufzmed5/ZaiKB6PR/vIvqrrLD/Lsr6Ye7vdxpxGLPyCLxlM\nVP9Q+gAdiN3pYStr5o3LzaaPiTDLspHHxyR6v9/DLa/XK+l5nWKCzLIs/Bzu/nq9iqJotHrWk2go\nPB5/u92EToAgTeVROK3LNbPFZs6RTzw2c9aPjyHyfr/X82UsvHF7nuchvzYedLjwzsZULaPwKc2i\nx+W1g62smTeuOGZ0vJgFn89n/fY8z0PDaqP7Pv6z3QLavj3eGNtEY+HhlpBfAS5LAIUruFwYDRFw\nZB99zI7t4ZvxlnbjZSNc1tXDa2dCfVs40CtN2/9pWjuNNP35JQbO4nJhNDQ3ds4N6ju4M7nGBBnz\nYvyhs/B2Qh0ImrEEM+sBgHO7XBiNiqKor+g0EE9HJtfhMNo2kHSBWdS3NdcselA/7+fqVYQTumgY\nfTwejeGer9erb5HRYSPv8mnPu5AKAFzBtcJoPQhmWfZ8Pquqej6fMfnVl/z8tKXzu/Gdw4WbwwSz\n0Cx6aFXy+wjR2DgaOrU2OiNgTtcKo9H9fi/LMgTBsENSnC8/Mnou6l//9V//7//+r++3wztIvbXm\nE4EN6dI9t6qqrHMHfSZGhZXTwrXCaJ7n4frVOYE9tI/GxshPO9ZnTLH/8R//8Q//8A99v62mmesk\nAVbT2TgK9JkYFVZOC9cKo8PWWVCpkVmHH8vIUfhAVcX/TF06lfCy1sijcCbC6O++a9qM9xq++6cB\n12hRgDp/VMBZCaPvfZojO49v39gYFdBpDwNY4dAkmDP5aZknA+DhLK4VRsOY3L6V5NvT5wfyYiyk\nvUtTZxhtryo6EDQ/ncgP1IkpAAdyrTDauaF81L69vc1S1Jk4+8rvDJftONv+lTAKUKdxFM7nWmF0\nzF7zcY2n+o319UdDOaGls7HJZ1++jHdv3N4ZXsuy7Cwc+E3XHvT1X0b66AH2L73aWj95nsdu99hp\nXr+l0eRZlmWMkiEdxrDYPrhefpZlIcvGrPl8PtshOC7l1Sg8SZLOlyZNL/eSQVNng9iPz4UwegVe\nZVjamnnjismmKIrOnvr7/d45nLSeR6POJBrU827UmUT7Ck96kmgijEIyFEZllIv46YVOvOowP2F0\nDSF3xn2Y+mY1dR4fvD0+pNU8z98WXv7wtnBhFPrCaONmH5Rz686jXnWYiTBKL2EUOsNomiSJhe6v\nRBiFRa2ZN641gQk4L0n0WmwQCqchjAKHJ4sQhMWktz4L4DN/2PoEAD7Ubvk0g+WSqiSNf4ekSVUZ\neAbHpGUUODYz6Ak0kMNBaRkFDkyX7HVVVZIklfcAHJ8wCuzY4Pr2427m5Krq97dJmnobwPEIo8cz\nPDzfkCmuQwc9bfIojLGrqX7C6PGImwB19cZRYIy3WWLNtGoCE3BImkWpq78HBFM4FmEUOBhLSfKW\ndwgciDAKHEyjd0mzKEHjnZCmSZKmYinsnzAKHIwOevrU9whNrDwKB2ECE7BjrbCpnYth9W2ZgEMQ\nRoF9GLGkaOMQzaJ0am0TCuyabnrgGCRRvpOm+1pSEWgQRoEDkET5SGPwaKLjHnZMGAX2ThLlC83J\nTNpGYa+MGQX2K03TRpuWJMqQxrJfP2dQO4XCPmkZBfZMEmWSjsVHgZ0RRoE9SpMqlUSZgzwKO6eb\nHljdu1Wc2r+XRJmisfhomiaVRZ9gN7SMAvsiibIEmzPBbgmjwF60NxKvKkmU2ZhfD/skjAK7oEGU\nFcijsEPGjB7P8FYilS9wDqjdZ+qNzELa40cT7zeuZ1fbkgmjxyNuciZiKOtr5NHEEqRcz9sssWZa\nFUaBJQ1eziRRtiKPwn4Io8DqqioxSJQ1db23qjRNfv6LSJc9bEIYBdYmhrIfnU2kifckrEgYBT40\nfiBR1/e5rnn2pp1HE5EUViSMAivpXGbclz17EJZ8EklhE8IosDgxlEOofoTPxu0iKSxKGAUW1Lfp\nou919qlKkkQrKaxLGAUWkTY3u/mN73IOQcc9rEYYBfp9teix1lCOofMd+fN7XiSFFQijwGzEUE5p\nOJIm3uEwjTAK19PZ3jn+67R15ED7qS9pTqMvkiY/fwS85+FTwiic10eh88Me+beH+0rmlAYiaaC5\nFD4ljAJjpb9PNR7iO5jTq5K0b1fbSHMpjCSMwtF8Nano+0frbwFq8HXLGYyY1fTjtjRJkqp2/Jhg\n6mMCbcIo7MO6EfOtkRnUNyuX9dt7v/bJrZJEcyl8QRilx8Q5LldznOpK07TvtARQmEX9MyKYwlvC\n6KGES9eYJrSZJqmMuntVTY1i489q+vMa3Qc3w5fD9Oe1gKr2eumCh1HGXzf67zp8uGDKlQmjxzN8\n/XMRe+OjOL6zrvMpfs+daTL+beJLET7W88ftyObS9m99DFlCuqcvOGH0eMZelzZ/ny10Aps/r90b\n397ZyTcffG1gGEzd+GDad4DPKRNV795Da6ZVYRSOamLojHyrwVyqvmFL/Ynyo2A6ojw4HmH0nEb+\nde6xjvJYciecUCtRNsLppy1TEioHJYyeU+9f5x5rt4/1w1yP6RsIlrXQZ6yxVlTrsWZJqH1cN9iE\nMEqPCbNHr+hddb1p2pxQr7484PzmTqgjHud3LjIsTRg9lNBWN+XCMPGi0nf3icl1+qVu4jJSX53A\nT/my+7nOfAn3lQB0+uni8PnabcN2uUgdpyKMcjGdbQm9l9rNrqyu6XBUa3YrDRZbve1zmbttdXoh\nLn3XJIye1+af6YVO4EexP131DjJ8YPPXBOB3a/X+jyfOXpMwypwuOKbUhQ94b4Hm0qVW/OhLqINn\nO9eQgOnW/Bpy/Z+LMHo8Fwx8i3I1AbYxLaGuueLHW++HBLyznzg73m6q/72df9MJo0dyoPf9anb+\nAQP4zFVXMvk9zn4bso8YZ1ez83eQMHoke/ozeIiACDAnV9URvmmdFXwHTVy/Zzxh9HhclADodaDx\nqYe1ZvD9Ru2xpuTm1fKGMAoAxzQ+LBxnfOqutmI+gS+H81bVmj2xwigAXNJC41OPE3zPGrIPRxgF\nAH74aKc9RhN8BwijAMDC1tz2eSHTG5KX2Dq76+7Tg2+67nRpYfR4ht8i1eYfVwCYYs0vsomPddjv\n3F2tzSOMHo+4CQAnsVAaflfsm11fq2rNxlFhFADgqnbQwvX/tj4BAACuSxiFDiuP3b4s9bwO9bwO\n9bwO9Xw+wigAAJsRRgEA2IwwCgDAZoRRAAA2I4xu7y9/+cuf/vSnNE3/9Kc//elPf/rv//7vrc8I\nAGAlqRXUt9U5K/DPf/7zv/3bv/Ud7yVbgXpeh3peh3peh3peh3pex5r17BXdUkyiv/zyy5/+9Kck\nSf793/893PLrr7/+4z/+Y+ddZnnJlKMc5ShHOcpRjnKWLmfUYwmjW/nLX/4Somcjd8aE2vnS7O1N\nphzlKEc5ylGOcpQzhTGjmwlJ9Jdffmm0gP7666/hB4NHAYDTE0a3EYNm6J2vi9m0b9goAMBpCKPb\n+N///d/wQ2fi/OWXX5Ik+dvf/rbqOQEArE4Y3UYMo53azaUAAKckjAIAsBlhFACA7VRsIVT+L7/8\n0vnbOKH+P//zPzvvCACwqOXT0G/+sPUzvahffvnlb3/7W98UpTii9F/+5V8av5JHAYAz0U1/GGVZ\n/s///M/WZ3F+6nk1ZVlufQrnV5blf/3Xf219Fuenntfh+ryO9d/PdmDaRtx+qbP+h38LAHAaWkb3\naHjhJwCA0xBGtxFXEv3LX/7S/m0YSxqWvgcAODFhdBtxz892I2jcKbQzpwIAnIkwupk///nPSZL8\n7W9/++///u+iKPI8T9M0TdN/+qd/CgfEwDpSWZaxkDzPi6KY/ZyPrl7PX1dRoxCzcD5VlmWswzzP\n1eEUs1dmKGqmszuDuWrY9XnYXPU8y0X+xBa6/IYKn1TEaotI0RY74n/55Zc///nPIZ4Gv/766/P5\nHF9U3+v7USEndr/fp1fR8/n0OZooy7LOCsyybOtTO54lKtP7uW6WGh64btzv96VO/VBmqedZLvLn\nttDlN5YzqZApd2ai5/PZOTA0Lno/spz6W+p+v9/v9/p7btGncAj1L4P7/f58Pr+rIvU8Ub264gsx\n1wXxapaozFjm3Cd7SHPV8PB1Qx6dpZ7rF3nX504LXX7rhUw5Pa/QlsLr9+uvv/75z38OjaO//vpr\nVftcjXl/xLdC44+/jwo5t76PykdVFD/JjXqO9e9L5a2+2u6rWwbMW5nP59M3d8MsNTxwfVDbwSz1\n3HdwvMi7Pi9x+W20+k86vSl3ZorhD8n4JoqBI7VzVPPV8yyFXFn8Vm7/yhfGp2aszM6eu5lP94Dm\nquGBKvW2r+au584j+0LYpSx0+Z3xumEC02biqOHOEdbxxreDi1+vV5Iknd8oeZ5/e3bnMUs9jyyE\n73ijzkhlLu2LGu4czhjLMYGv0/h6jhXYeZfw5Ri+KGn7+ooR7hgGRUw/DXvTn0TnJ801bi7DNVn/\nUhEFBmRZpn7mMldl1t/bRVE8Ho/pZZ7D9BoeDkkEs7yTQ+JUzwPmvfwWRRFSR5ieP0OJU5pVWc74\nnt++MUn6gMYYWc8DfRyVqp6sb9wzX5hSmcPvc4K53q6uG8NmqWeV/NYX9RxrNdxlluuGi84efTr3\nqDFLrj6L8OIDZYaNr+fGZ6/BHKYvPJ/PxrIGau9rc1WmMNpniberJqG2ueo5lGNOXp+J9dy4izB6\nTvXpaePv1TlmVBId8Gk9Dxyswr/QeK9KolPMVZnCaJ/Z367e+Z1mqefGEEaX5bYp9dzuTpzlumEC\n046EEYe32y38c2Cp5IY8z+OY0SzL4nvl9XoZQ9P2XT3Hz1uapmEoWBgrY7ua7zT+fHo8HuaBfU1l\nLm3GGq5fNO73u1eqbol3su/Btq/ruSzLEDbaiXaqKUmWGdXfHB/9JdfXzB7Dk78L676u56p/+wrd\n9FPUe9OMGZ1oYmVqGX1rSg03/u71bh8w12XB9+CwT+t5OGxMORMXne1NuUINj842L6Rulm+CdjNq\nKEcYnSjW59YncgZfV6YwOtIXNWzXpS/MclnwPfjWyHoO7+F2rBdGz6A+uuWLK9TbN4E/CoOJ9Txs\neHoTY4hBM/q6Mr0KI31aUfHi41L8kbnekP4GGDaynpPRvjgHY0a3VJZlWNUvXKG+GBxjJdExptfz\n2/LDD0YmDcjzPA63ZSKVubQZaziOEH0+n16yhlnqOQzDdfkdsP8rhkXvtxTm0GRZ9vVbpD51iT7T\n6zn5sc1SnuftS96eP+H7Ed6oRVGorulU5tLmquF4ufioYek6Znwnv/0qvHJanaWe+yb7xgXwx8+6\n7vBFayqzmKUDYnjMqPV+q7k7ejp72VTyGMP7C7gifWShytRNH81Vwy4Ow2ap55HrQE8/2+Na9PI7\nSw1rGd3M8HbnUf23sR00vnvyPM+y7PV6PR6PRqNdWZZx9aIrrx4ySz0nSXK/3x+Px+v1KoqifnDs\ng7tyJY9RFEV4Q+Z53vjrPL5vZ9nj+Aq+qMzOdzV9ZqnhxlarfY/V2d9yEbPUczzydrs9n8/G96Ad\nbpNDXDGmJFmmGPkC1f/U63vV6sfX1xltl3BB69TzxSt5pMa6Wo06NLfjI59WZt+7uk4bUt30Gh75\nx9XF201neSfXO4jbhSQu0ctcMQKz6Y9tzEUqGReSqv4dmHwCV6hnlTzewFqtW5/a8XxUmQPv6kgY\nbZhYw3139+ZvmOWd3Bf9/ZUbzX7FCGa5bqSVLpuzKMsy7gyU/7DxOZ1R+YNK/lroslSHs1CZS1PD\n65ilnuvX5+TaQyD67PP9LIwCALAZ64wCALAZYRQAgM0IowAAbEYYBQBgM8IoAACbEUYBANiMMAoA\nwGaEUQAANiOMAgCwGWEUAIDNCKMAAGxGGAUAYDPCKAAAm/nD1icAAHB4ZVmWZRn/med5nucDR8aD\n8x9mOY2iKIYfvfPI8Se/hLSqqtUeDADglIqieDwe9VvaEat9TN39fg8BcYo0Tfseva4sy9vtFh80\n3mvGMxlPyygAwGyyLEuSpN2ymOf56/VqHJYkSbzx8Xg8Ho+JrYRZltUfpU+9XTZJkvv9Hm4Zc9/Z\nCaMAALOp93dH9ST6fD4bUbUsy6IowgFpOqnXuiiK0ORZFMVA62ZsoA1nEo9sNJGuwwQmAIAFxaCZ\nZVlVVZ2NpmVZxrbSKV3ksfCB8QAxLt/v968faEbCKADAgmIu7Gw0jeJvB3LkGDHU9j1co49+c8Io\nAHAY9YneeZ6nP4TGxXB7mJETrTkXp+2jZsjhYxrPq/6UG4c1Hrqh0Ue/vQoA4CBCenk+n52p5vl8\nxnbButA/vqiYIxu3x/OZUnjf802S5H6/t48fiHmxqIE7dv5qOVpGAYCDCXN0ns9nSDMxCN5ut9fr\nFaNnDF6bTBKvP3RnRB4vPN8kSe73e3jW8ak9Ho9202+skHbjaL39eMopzUgYBQCOp6rNBKrnqizL\n6mMiY2gbHq+5tCkd4vG+z+ezKIrwzzzPY0Nme4zpQE/9xAGpSxBGAYCDaTc09s3a2cuwyB5pj3q8\njm2r7efyNmr3Rc+dzKMPhFEA4GA6V0fa4DyWN9yrHp91O4x29tTHQvbTR59Y9B4AOJzDRc++lst2\nC2WjLTPesbF9/Ftx69GiKBqLRk0cwDo7YRQAYFl9M6jaLZR9HevDYz07c2rYGrT90HuL8sIoAMBS\nRm4W/9bwKM/OfBm3Bi3LMs/zffbRJ8IoAMBy4q70w5vFB+0GzjzPY2/7Fw8dfgg99fvso09MYAIA\nWE4MkWPWVIrriUYDU5Ri+fVRoQ0hetabZvfWLJoIowDAdbSbJxtJ7u0BX4g97GmaDhw2PJSznVOT\nH7OUHo9H333jc6nvoTp0rlsQRgGAqwjRLf4zdF7Xs+bbA75QFEXsHG+sIRofJU3TvqGlMcs2cuSY\nnvd4l1k2glqIMAoAsKyyLGMQfDweYWX7PM/DD7HVs761aRSz7Ov1CvcK4r2Gs3K9wB320SfCKADA\nCsqyjBsmBfWm0Pv9Hrb37OxGr2fZ1w9JkmRZ1iizrV7gDvvoE7PpAYADiRuyNwzMVe+7y/rChvLl\nD+GWoHFM+77h+DiAtX3H4Qed4/SXIowCAKxnZIjstM9+9ol00wMAsBlhFAC4qLfT5CfOo2cMYRQA\nuJY8z8uyDIt0JrXhmOMPGJD+MPtpL2fbczZmFAC4ltfrdbvdsiyrqios8JmmaX2Wz9sD3trnip59\nsizrW+V0BZ/VLADAcYVMWZZlfQpR/Z9vD2B2wigAcBVpmj6fz4Fk+fYAZmfMKABwCWHo58CcpLcH\nsARhFAC4kLetnppFV6abHgCAzWgZBQBgM8IoAACbEUYBANiMMAoAwGaEUQAANiOMAgCwGWEUAIDN\nCKMAAGxGGAUAYDPCKAAAmxFGAQDYjDAKAMBmhFEAADYjjAIAsBlh/aUxsAAAADBJREFUFACAzQij\nAABsRhgFAGAzwigAAJsRRgEA2IwwCgDAZoRRAAA2I4wCALCZ/w+iq7kh97FDcgAAAABJRU5ErkJg\ngg==\n",
      "text/plain": [
       "Canvas('Canvas_eaGb6FuDUeecZuDP5pDNf4')"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[#1] INFO:DataHandling -- RooDataHist::adjustBinning(dh): fit range of variable x expanded to nearest bin boundaries: [2.8,3.4] --> [2.8,3.4]\n",
      "[#1] INFO:NumericIntegration -- RooRealIntegral::init(voigt_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n",
      "[#1] INFO:Minization -- RooMinuit::optimizeConst: activating const optimization\n",
      "[#1] INFO:NumericIntegration -- RooRealIntegral::init(voigt_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n",
      "[#1] INFO:Minization --  The following expressions have been identified as constant and will be precalculated and cached: (pol0)\n",
      "[#1] INFO:Minization --  The following expressions will be evaluated in cache-and-track mode: (voigt)\n",
      " **********\n",
      " **   13 **MIGRAD        2000           1\n",
      " **********\n",
      " FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.\n",
      " START MIGRAD MINIMIZATION.  STRATEGY  1.  CONVERGENCE WHEN EDM .LT. 1.00e-03\n",
      " FCN=-236203 FROM MIGRAD    STATUS=INITIATE       79 CALLS          80 TOTAL\n",
      "                     EDM= unknown      STRATEGY= 1      NO ERROR MATRIX       \n",
      "  EXT PARAMETER               CURRENT GUESS       STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  jpsi_yield   9.00000e-01   5.00000e-02   0.00000e+00  -5.87205e+02\n",
      "   2  mean         8.59915e-01   5.00000e-01   8.52460e-01   1.10649e+02\n",
      "   3  sigma        2.05619e+00   2.50000e+00   0.00000e+00   1.10265e+04\n",
      "   4  width        1.00000e+00   5.00000e-01   0.00000e+00   2.71010e+03\n",
      "                               ERR DEF= 0.5\n",
      " MIGRAD MINIMIZATION HAS CONVERGED.\n",
      " FCN=-499021 FROM MIGRAD    STATUS=CONVERGED     654 CALLS         655 TOTAL\n",
      "                     EDM=4.52326e-07    STRATEGY= 1  ERROR MATRIX UNCERTAINTY   2.1 per cent\n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  jpsi_yield   7.30934e-01   1.61929e-03  -1.58744e-06   6.55353e-02\n",
      "   2  mean         3.09195e+00   7.46992e-05  -2.60707e-07  -4.16278e-01\n",
      "   3  sigma        2.22841e-02   2.32619e-04  -1.48699e-06  -9.50308e+00\n",
      "   4  width        3.28755e-02   4.94770e-04   1.89099e-06  -3.20096e+00\n",
      "                               ERR DEF= 0.5\n",
      " EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  4    ERR DEF=0.5\n",
      "  2.622e-06  4.547e-09 -2.466e-07  6.185e-07 \n",
      "  4.547e-09  5.580e-09  1.701e-10  3.120e-10 \n",
      " -2.466e-07  1.701e-10  5.411e-08 -1.058e-07 \n",
      "  6.185e-07  3.120e-10 -1.058e-07  2.448e-07 \n",
      " PARAMETER  CORRELATION COEFFICIENTS  \n",
      "       NO.  GLOBAL      1      2      3      4\n",
      "        1  0.78488   1.000  0.038 -0.655  0.772\n",
      "        2  0.06053   0.038  1.000  0.010  0.008\n",
      "        3  0.92333  -0.655  0.010  1.000 -0.919\n",
      "        4  0.94637   0.772  0.008 -0.919  1.000\n",
      " **********\n",
      " **   18 **HESSE        2000\n",
      " **********\n",
      " COVARIANCE MATRIX CALCULATED SUCCESSFULLY\n",
      " FCN=-499021 FROM HESSE     STATUS=OK             29 CALLS         684 TOTAL\n",
      "                     EDM=2.56515e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                INTERNAL      INTERNAL  \n",
      "  NO.   NAME      VALUE            ERROR       STEP SIZE       VALUE   \n",
      "   1  jpsi_yield   7.30934e-01   1.67135e-03   1.09005e-03   2.66149e+00\n",
      "   2  mean         3.09195e+00   7.41372e-05   1.48553e-05   2.39051e-01\n",
      "   3  sigma        2.22841e-02   2.36863e-04   2.89734e-05  -1.60065e+00\n",
      "   4  width        3.28755e-02   5.18970e-04   4.30546e-05  -1.60706e+00\n",
      "                               ERR DEF= 0.5\n",
      " EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  4    ERR DEF=0.5\n",
      "  2.793e-06  1.782e-09 -2.637e-07  6.848e-07 \n",
      "  1.782e-09  5.496e-09  5.696e-10 -3.126e-10 \n",
      " -2.637e-07  5.696e-10  5.610e-08 -1.133e-07 \n",
      "  6.848e-07 -3.126e-10 -1.133e-07  2.693e-07 \n",
      " PARAMETER  CORRELATION COEFFICIENTS  \n",
      "       NO.  GLOBAL      1      2      3      4\n",
      "        1  0.80526   1.000  0.014 -0.666  0.789\n",
      "        2  0.06724   0.014  1.000  0.032 -0.008\n",
      "        3  0.92733  -0.666  0.032  1.000 -0.922\n",
      "        4  0.95128   0.789 -0.008 -0.922  1.000\n",
      "[#1] INFO:Minization -- RooMinuit::optimizeConst: deactivating const optimization\n",
      "[#1] INFO:NumericIntegration -- RooRealIntegral::init(voigt_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n"
     ]
    }
   ],
   "source": [
    "fit(mass_hist_zoomed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}