{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Coffea Processors\n",
    "This is a rendered copy of [processor.ipynb](https://github.com/CoffeaTeam/coffea/blob/master/binder/processor.ipynb). You can optionally run it interactively on [binder at this link](https://mybinder.org/v2/gh/coffeateam/coffea/master?filepath=binder%2Fprocessor.ipynb)\n",
    "\n",
    "Coffea relies mainly on [uproot](https://github.com/scikit-hep/uproot) to provide access to ROOT files for analysis.\n",
    "As a usual analysis will involve processing tens to thousands of files, totalling gigabytes to terabytes of data, there is a certain amount of work to be done to build a parallelized framework to process the data in a reasonable amount of time. Of course, one can work directly within uproot to achieve this, as we'll show in the beginning, but coffea provides the `coffea.processor` module, which allows users to worry just about the actual analysis code and not about how to implement efficient parallelization, assuming that the parallization is a trivial map-reduce operation (e.g. filling histograms and adding them together). The module provides the following key features:\n",
    "\n",
    " * A `ProcessorABC` abstract base class that can be derived from to implement the analysis code;\n",
    " * An interface to the arrays being read from the TTree, either [DataFrame](https://coffeateam.github.io/coffea/api/coffea.processor.LazyDataFrame.html#coffea.processor.LazyDataFrame) or [NanoEvents](https://coffeateam.github.io/coffea/api/coffea.nanoaod.NanoEvents.html#coffea.nanoaod.NanoEvents), to be used as inputs;\n",
    " * A set of accumulator types (value_accumulator, list_accumulator, set_accumulator, dict_accumulator, defaultdict_accumulator, column_accumulator) as described further [here](https://coffeateam.github.io/coffea/modules/coffea.processor.html#classes) to be used as output; and\n",
    " * A set of parallel executors to access multicore processing or distributed computing systems such as [Dask](https://distributed.dask.org/en/latest/), [Parsl](http://parsl-project.org/), [Spark](https://spark.apache.org/), and others.\n",
    "\n",
    "Let's start by writing a simple processor class that reads some CMS open data and plots a dimuon mass spectrum.\n",
    "We'll start by copying the [ProcessorABC](https://coffeateam.github.io/coffea/api/coffea.processor.ProcessorABC.html#coffea.processor.ProcessorABC) skeleton and filling in some details:\n",
    "\n",
    " * Remove `flag`, as we won't use it\n",
    " * Adding a new histogram for $m_{\\mu \\mu}$\n",
    " * Building a [JaggedCandidateArray](https://coffeateam.github.io/coffea/api/coffea.analysis_objects.JaggedCandidateMethods.html#coffea.analysis_objects.JaggedCandidateMethods.candidatesfromcounts) for muons,\n",
    " * Calculating the dimuon invariant mass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from coffea import hist, processor\n",
    "from coffea.analysis_objects import JaggedCandidateArray\n",
    "\n",
    "class MyProcessor(processor.ProcessorABC):\n",
    "    def __init__(self):\n",
    "        self._accumulator = processor.dict_accumulator({\n",
    "            \"sumw\": processor.defaultdict_accumulator(float),\n",
    "            \"mass\": hist.Hist(\n",
    "                \"Events\",\n",
    "                hist.Cat(\"dataset\", \"Dataset\"),\n",
    "                hist.Bin(\"mass\", \"$m_{\\mu\\mu}$ [GeV]\", 60, 60, 120),\n",
    "            ),\n",
    "        })\n",
    "\n",
    "    @property\n",
    "    def accumulator(self):\n",
    "        return self._accumulator\n",
    "\n",
    "    def process(self, df):\n",
    "        output = self.accumulator.identity()\n",
    "\n",
    "        dataset = df['dataset']\n",
    "        muons = JaggedCandidateArray.candidatesfromcounts(\n",
    "            df['nMuon'],\n",
    "            pt=df['Muon_pt'].content,\n",
    "            eta=df['Muon_eta'].content,\n",
    "            phi=df['Muon_phi'].content,\n",
    "            mass=df['Muon_mass'].content,\n",
    "            charge=df['Muon_charge'].content,\n",
    "        )\n",
    "\n",
    "        cut = (muons.counts == 2) & (muons['charge'].prod() == -1)        \n",
    "        dimuons = muons[cut].choose(2)\n",
    "        \n",
    "        output[\"sumw\"][dataset] += df.size\n",
    "        output[\"mass\"].fill(\n",
    "            dataset=dataset,\n",
    "            mass=dimuons.mass.flatten(),\n",
    "        )\n",
    "\n",
    "        return output\n",
    "\n",
    "    def postprocess(self, accumulator):\n",
    "        return accumulator"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If we were to just use bare uproot to execute this processor, we could do that with the following example, which:\n",
    "\n",
    " * Opens a CMS open data file\n",
    " * Creates a lazy data frame (roughly equivalent to [uproot.lazyarrays](https://uproot.readthedocs.io/en/latest/opening-files.html#uproot.tree.lazyarrays) but with some specializations needed for othere execution environments)\n",
    " * Creates a `MyProcessor` instance\n",
    " * Runs the `process()` function, which returns our accumulators\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'sumw': defaultdict_accumulator(float, {'DoubleMuon': 10000.0}),\n",
       " 'mass': <Hist (dataset,mass) instance at 0x1074336d0>}"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import uproot\n",
    "\n",
    "filename = \"root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012B_DoubleMuParked.root\"\n",
    "file = uproot.open(filename)\n",
    "df = processor.LazyDataFrame(\n",
    "    tree=file[\"Events\"],\n",
    "    entrystart=0,\n",
    "    entrystop=10000,\n",
    ")\n",
    "df[\"dataset\"] = \"DoubleMuon\"\n",
    "p = MyProcessor()\n",
    "out = p.process(df)\n",
    "out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One could expand on this code to run over several chunks of the file, setting `entrystart` and `entrystop` as appropriate. Then, several datasets could be processed by iterating over several files. However, [run_uproot_job](https://coffeateam.github.io/coffea/api/coffea.processor.run_uproot_job.html#coffea.processor.run_uproot_job) can help with this. One lists the datasets and corresponding files, the processor they want to run, and which executor they want to use. Available executors are listed [here](https://coffeateam.github.io/coffea/modules/coffea.processor.html#functions). Since these files are very large, we limit to just reading the first few chunks of events from each dataset with `maxchunks`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "e8cfeaac1dd14b4a8ba2a9c46c64faa1",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(IntProgress(value=0, description='Processing', max=8, style=ProgressStyle(description_width='in…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'sumw': defaultdict_accumulator(float,\n",
       "                         {'ZZ to 4mu': 399752.0, 'DoubleMuon': 400224.0}),\n",
       " 'mass': <Hist (dataset,mass) instance at 0x12337e110>}"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fileset = {\n",
    "    'DoubleMuon': [\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012B_DoubleMuParked.root',\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012C_DoubleMuParked.root',\n",
    "    ],\n",
    "    'ZZ to 4mu': [\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/ZZTo4mu.root'\n",
    "    ]\n",
    "}\n",
    "\n",
    "out = processor.run_uproot_job(\n",
    "    fileset,\n",
    "    treename=\"Events\",\n",
    "    processor_instance=MyProcessor(),\n",
    "    executor=processor.iterative_executor,\n",
    "    maxchunks=4,\n",
    ")\n",
    "out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, if we want to use more than a single core on our machine, we simply change [iterative_executor](https://coffeateam.github.io/coffea/api/coffea.processor.iterative_executor.html) for [futures_executor](https://coffeateam.github.io/coffea/api/coffea.processor.futures_executor.html), which uses the python [concurrent.futures](https://docs.python.org/3/library/concurrent.futures.html) standard library. Optional arguments to these executors can be provided via `executor_args` parameter of `run_uproot_job`, such as the number of cores to use (2):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "7167c2d0bec14069857a18d004f616d6",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(IntProgress(value=0, description='Processing', max=8, style=ProgressStyle(description_width='in…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'sumw': defaultdict_accumulator(float,\n",
       "                         {'ZZ to 4mu': 399752.0, 'DoubleMuon': 400224.0}),\n",
       " 'mass': <Hist (dataset,mass) instance at 0x106bb6890>}"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "out = processor.run_uproot_job(\n",
    "    fileset,\n",
    "    treename=\"Events\",\n",
    "    processor_instance=MyProcessor(),\n",
    "    executor=processor.futures_executor,\n",
    "    executor_args={'workers': 2},\n",
    "    maxchunks=4,\n",
    ")\n",
    "out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Hopefully this ran faster than the previous cell, but that may depend on how many cores are available on the machine you are running this notebook."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Getting fancy\n",
    "Let's flesh out this analysis into a 4-muon analysis, searching for diboson events:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import awkward as ak\n",
    "\n",
    "class FancyDimuonProcessor(processor.ProcessorABC):\n",
    "    def __init__(self):\n",
    "        dataset_axis = hist.Cat(\"dataset\", \"Primary dataset\")\n",
    "        mass_axis = hist.Bin(\"mass\", r\"$m_{\\mu\\mu}$ [GeV]\", 600, 0.25, 300)\n",
    "        pt_axis = hist.Bin(\"pt\", r\"$p_{T,\\mu}$ [GeV]\", 3000, 0.25, 300)\n",
    "        \n",
    "        self._accumulator = processor.dict_accumulator({\n",
    "            'mass': hist.Hist(\"Counts\", dataset_axis, mass_axis),\n",
    "            'mass_near': hist.Hist(\"Counts\", dataset_axis, mass_axis),\n",
    "            'mass_far': hist.Hist(\"Counts\", dataset_axis, mass_axis),\n",
    "            'pt_lead': hist.Hist(\"Counts\", dataset_axis, pt_axis),\n",
    "            'pt_trail': hist.Hist(\"Counts\", dataset_axis, pt_axis),\n",
    "            'cutflow': processor.defaultdict_accumulator(int),\n",
    "        })\n",
    "    \n",
    "    @property\n",
    "    def accumulator(self):\n",
    "        return self._accumulator\n",
    "    \n",
    "    def process(self, df):\n",
    "        output = self.accumulator.identity()\n",
    "        \n",
    "        dataset = df['dataset']\n",
    "        muons = JaggedCandidateArray.candidatesfromcounts(\n",
    "            df['nMuon'],\n",
    "            pt=df['Muon_pt'],\n",
    "            eta=df['Muon_eta'],\n",
    "            phi=df['Muon_phi'],\n",
    "            mass=df['Muon_mass'],\n",
    "            charge=df['Muon_charge'],\n",
    "            softId=df['Muon_softId'],\n",
    "            tightId=df['Muon_tightId']\n",
    "            )        \n",
    "        \n",
    "        output['cutflow']['all events'] += muons.size\n",
    "        \n",
    "        soft_id = (muons.softId > 0)\n",
    "        muons = muons[soft_id]\n",
    "        output['cutflow']['soft id'] += soft_id.any().sum()\n",
    "        \n",
    "        twomuons = (muons.counts >= 2)\n",
    "        output['cutflow']['two muons'] += twomuons.sum()\n",
    "        \n",
    "        dimuons = muons[twomuons].distincts()\n",
    "        \n",
    "        twodimuons = (dimuons.counts >= 2)\n",
    "        output['cutflow']['>= two dimuons'] += twodimuons.sum()\n",
    "        dimuons = dimuons[twodimuons]\n",
    "        \n",
    "        opposite_charge = (dimuons.i0['charge'] * dimuons.i1['charge'] == -1)\n",
    "        \n",
    "        dimuons = dimuons[opposite_charge]\n",
    "        output['cutflow']['opposite charge'] += opposite_charge.any().sum()\n",
    "        \n",
    "        mass_20GeV = (dimuons.mass > 35)\n",
    "        dimuons = dimuons[mass_20GeV]\n",
    "        \n",
    "        exactlytwodimuons = (dimuons.counts == 2)\n",
    "        output['cutflow']['== two dimuons'] += exactlytwodimuons.sum()\n",
    "        dimuons = dimuons[exactlytwodimuons].compact()\n",
    "        \n",
    "        leading_mu = (dimuons.i0.pt.content > dimuons.i1.pt.content)\n",
    "        pt_lead = ak.JaggedArray.fromoffsets(\n",
    "            dimuons.offsets,\n",
    "            np.where(leading_mu, dimuons.i0.pt.content, dimuons.i1.pt.content)\n",
    "        )\n",
    "        pt_trail = ak.JaggedArray.fromoffsets(\n",
    "            dimuons.offsets,\n",
    "            np.where(~leading_mu, dimuons.i0.pt.content, dimuons.i1.pt.content)\n",
    "        )\n",
    "        \n",
    "        near_z = np.abs(dimuons.mass - 91.118).argmin()\n",
    "        far_z = np.abs(dimuons.mass - 91.118).argmax()\n",
    "        \n",
    "        output['mass'].fill(dataset=dataset,\n",
    "                            mass=dimuons.p4.sum().mass)\n",
    "        output['mass_near'].fill(dataset=dataset, \n",
    "                                 mass=dimuons.mass[near_z].flatten())\n",
    "        output['mass_far'].fill(dataset=dataset, \n",
    "                                mass=dimuons.mass[far_z].flatten())\n",
    "        output['pt_lead'].fill(dataset=dataset,\n",
    "                               pt=pt_lead.flatten())\n",
    "        output['pt_trail'].fill(dataset=dataset,\n",
    "                                pt=pt_trail.flatten())\n",
    "        return output\n",
    "\n",
    "    def postprocess(self, accumulator):\n",
    "        return accumulator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "9b7d0eb45c0e456aa7a9a52889536aeb",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(IntProgress(value=0, description='Processing', max=20, style=ProgressStyle(description_width='i…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "{'mass': <Hist (dataset,mass) instance at 0x106bb6dd0>, 'mass_near': <Hist (dataset,mass) instance at 0x106bb6710>, 'mass_far': <Hist (dataset,mass) instance at 0x106bb6b50>, 'pt_lead': <Hist (dataset,pt) instance at 0x106bb6f10>, 'pt_trail': <Hist (dataset,pt) instance at 0x106bb6910>, 'cutflow': defaultdict_accumulator(<class 'int'>, {'all events': 1999940, 'soft id': 1803995, 'two muons': 1394358, '>= two dimuons': 584764, 'opposite charge': 576563, '== two dimuons': 158843})}\n"
     ]
    }
   ],
   "source": [
    "import time\n",
    "\n",
    "tstart = time.time()    \n",
    "\n",
    "fileset = {\n",
    "    'DoubleMuon': [\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012B_DoubleMuParked.root',\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/Run2012C_DoubleMuParked.root',\n",
    "    ],\n",
    "    'ZZ to 4mu': [\n",
    "        'root://eospublic.cern.ch//eos/root-eos/cms_opendata_2012_nanoaod/ZZTo4mu.root'\n",
    "    ]\n",
    "}\n",
    "\n",
    "output = processor.run_uproot_job(\n",
    "    fileset,\n",
    "    treename='Events',\n",
    "    processor_instance=FancyDimuonProcessor(),\n",
    "    executor=processor.futures_executor,\n",
    "    executor_args={'workers': 6, 'flatten': True},\n",
    "    chunksize=100000,\n",
    "    maxchunks=10,\n",
    ")\n",
    "\n",
    "elapsed = time.time() - tstart\n",
    "print(output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What follows is just us looking at the output, you can execute it if you wish"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.0, 3000.0)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAEPCAYAAABlZDIgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3dd3xUZfb48c9JL4TQexekSZNiQcWKCIriWtBVgdVFUXRdt1l2bbv+1l0LtkXFlWJBRVe/otiQsqIizQUERKVKKAESSEJCSJnz++O5SSaQhBBmMpPkvF+vec29zy1zMknmzL1PE1XFGGOMOV4RoQ7AGGNM7WAJxRhjTEBYQjHGGBMQllCMMcYEhCUUY4wxAWEJxRhjTEAELaGISJyILBWRVSKyVkQe8so7isgSEdkgIm+JSIxXHuutb/C2d/A71z1e+Q8icmGwYjbGGFN1wbxCOQScq6p9gL7AMBE5FfgHMElVOwP7gBu9/W8E9nnlk7z9EJEewGigJzAMmCwikUGM2xhjTBUELaGoc8BbjfYeCpwLvOOVzwAu85Yv9dbxtp8nIuKVv6mqh1R1M7ABGBSsuI0xxlRNUOtQRCRSRFYCu4G5wEZgv6oWeLukAK295dbANgBvewbQ2L+8jGOMMcaEiahgnlxVC4G+ItIAeA/oFqzXEpHxwHiAxMTE/t26Be2ljDGmVlqxYsVeVW1a1eODmlCKqOp+EVkAnAY0EJEo7yqkDbDd22070BZIEZEoIBlI8ysv4n+M/2tMAaYADBgwQJcvXx6sH8cYY2olEdl6PMcHs5VXU+/KBBGJBy4AvgcWAFd4u40B3veWZ3vreNvnqxu5cjYw2msF1hHoAiwNVtzGGGOqJphXKC2BGV6LrAhglqp+KCLrgDdF5G/A/4CXvf1fBl4VkQ1AOq5lF6q6VkRmAeuAAuA271aaMcaYMCK1cfh6u+VljDHHTkRWqOqAqh5fLXUoxhhTJD8/n5SUFHJzc0MdSp0VFxdHmzZtiI6ODuh5LaEYY6pVSkoKSUlJdOjQAdfVzFQnVSUtLY2UlBQ6duwY0HPbWF7GmGqVm5tL48aNLZmEiIjQuHHjoFwhWkIxxlQ7SyahFaz33xKKMSYsRUZG0rdvX0466SSuvPJKcnJyytzv9NNPr+bIynf22WdztAZBTz31VLk/S1UtXLiQr7/+OqDnrApLKMaYsBQfH8/KlStZs2YNMTExvPDCC6W2FxS4EZyC8UFadO5gsIRijDEhdOaZZ7JhwwYWLlzImWeeyciRI+nRowcA9erVA9yH6pAhQ7j00kvp1KkTd999N6+//jqDBg2iV69ebNy4EYAPPviAU045hX79+nH++eeTmpoKwIMPPsj111/P4MGDuf766znrrLNYuXJlcQxnnHEGq1atKhXXwYMHGT16NN27d2fUqFEcPHiweNuECRMYMGAAPXv25IEHHgDgmWeeYceOHZxzzjmcc8455e4HcPfdd9OjRw969+7N73//ewD27NnDL37xCwYOHMjAgQP56quv2LJlCy+88AKTJk2ib9++LFq0KKDv/TFR1Vr36N+/vxpjwtO6desqtV9iYqKqqubn5+vIkSN18uTJumDBAk1ISNBNmzYdsd+CBQs0OTlZd+zYobm5udqqVSu9//77VVX1qaee0t/85jeqqpqenq4+n09VVV966SW96667VFX1gQce0JNPPllzcnJUVXX69OnFx/zwww9a1ufKE088oePGjVNV1VWrVmlkZKQuW7ZMVVXT0tJUVbWgoECHDBmiq1atUlXV9u3b6549e4rPUdZ+e/fu1RNPPLE4zn379qmq6jXXXKOLFi1SVdWtW7dqt27dimN/7LHHKvW+Finr9wAs1+P47LUrFGNMWDp48CB9+/ZlwIABtGvXjhtvdFMnDRo0qNzmrgMHDqRly5bExsZywgknMHToUAB69erFli1bANds+cILL6RXr1489thjrF27tvj4kSNHEh8fD8CVV17Jhx9+SH5+PlOnTmXs2LFHvN4XX3zBddddB0Dv3r3p3bt38bZZs2Zx8skn069fP9auXcu6devKjLms/ZKTk4mLi+PGG2/k3XffJSEhAYDPP/+ciRMn0rdvX0aOHElmZiYHDhwo87yhYP1QjDFhqagO5XCJiYnlHhMbG1u8HBERUbweERFRXC9y++23c9dddzFy5EgWLlzIgw8+WOa5ExISuOCCC3j//feZNWsWK1asqHTsmzdv5vHHH2fZsmU0bNiQsWPHltlMt7z9oqKiWLp0KfPmzeOdd97hueeeY/78+fh8Pr755hvi4uIqHUt1sisUY0ydkpGRQevWbkqlGTNmVLjvTTfdxB133MHAgQNp2LDhEdvPOussZs6cCcCaNWtYvXo1AJmZmSQmJpKcnExqaioff/xx8TFJSUlkZWVVuN+BAwfIyMhg+PDhTJo0qbjuZujQoTz77LPF5ypKuP7nDCVLKMaYOuXBBx/kyiuvpH///jRp0qTCffv370/9+vUZN25cmdsnTJjAgQMH6N69O/fffz/9+/cHoE+fPvTr149u3bpx7bXXMnjw4OJjxo8fz7BhwzjnnHPK3S8rK4uLL76Y3r17c8YZZ/Dkk08CrlJ/+fLl9O7dmx49ehS3fLvkkkt47733Ql4pb4NDGmOq1ffff0/37t1DHUal7Nixg7PPPpv169cTEVG7vn+X9Xs43sEha9c7ZIwxAfLKK69wyimn8Mgjj9S6ZBIsVilvjDFluOGGG7jhhhtCHUaNYmnXGGNMQFhCMcYYExCWUIwxxgSEJRRjTNi7+sXFXP3i4lCHYY7CEooxps4pGhq/Z8+e9OnThyeeeAKfz1fl8xUNUHm4sWPH8s4771R47NixY0lISCjVMfHOO+9ERNi7d2+VYwoFSyjGmDqnaFiXtWvXMnfuXD7++GMeeuihkMXTuXNn3n//fQB8Ph/z588v7s1fk1hCMcbUac2aNWPKlCk899xzqCq5ubmMGzeOXr160a9fPxYsWADA9OnTmThxYvFxF198MQsXLixe/+1vf0vPnj0577zz2LNnzxGvs2LFCoYMGUL//v258MIL2blzZ/G20aNH89ZbbwFuGP7BgwcTFeV6dWzZsoWTTjqpeN/HH3+8ePyxlStXcuqpp9K7d29GjRrFvn37ADfR15/+9CcGDRrEiSeeWG29560fijEmZB76YC3rdmQeUb5uZ+mynENuYMdeD35aqrxHy/pHHNujVX0euKTnMcXRqVMnCgsL2b17N6+99hoiwnfffcf69esZOnQoP/74Y4XHZ2dnM2DAACZNmsTDDz/MQw89xHPPPVe8PT8/n9tvv53333+fpk2b8tZbb3HfffcxdepUAE488URmz57Nvn37eOONN7juuutKjf9VnhtuuIFnn32WIUOGcP/99/PQQw/x1FNPAW6SsKVLl/LRRx/x0EMP8fnnnx/Te1IVllCMMcbPl19+ye233w5At27daN++/VETSkREBFdffTUA1113HZdffnmp7T/88ANr1qzhggsuAKCwsJCWLVuW2ufyyy/nzTffZMmSJbz44otHjTMjI4P9+/czZMgQAMaMGcOVV15Z6nzgxiMrGro/2CyhGGNCprJXEkUtvN66+bSgxLFp0yYiIyNp1qxZuftERUWVqrgvazj6IiJSal1V6dmzJ4sXl99S7eqrr6Z///6MGTOm1FAvx/K6/oqG7o+MjAzqlMb+rA7FGFOn7dmzh1tuuYWJEyciIpx55pm8/vrrAPz444/8/PPPdO3alQ4dOrBy5Up8Ph/btm1j6dKlxefw+XzFrblmzpzJGWecUeo1unbtyp49e4oTSn5+fqmJvQDat2/PI488wq233lqqvHnz5uzevZu0tDQOHTrEhx9+CEBycjINGzYsrh959dVXi69WQsWuUIwxdU7RbJD5+flERUVx/fXXc9dddwFw6623MmHCBHr16kVUVBTTp08nNjaWwYMH07FjR3r06EH37t05+eSTi8+XmJjI0qVL+dvf/kazZs2KK9iLxMTE8M4773DHHXeQkZFBQUEBd955Jz17lr5Cu/nmm4+INTo6mvvvv59BgwbRunVrunXrVrxtxowZ3HLLLeTk5NCpUyemTZsWyLfpmAVt+HoRaQu8AjQHFJiiqk+LyIPAr4GiZhD3qupH3jH3ADcChcAdqvqpVz4MeBqIBP6tqo9W9No2fL0x4asqw9cH+5ZXXRSM4euDeYVSAPxOVb8VkSRghYjM9bZNUtXH/XcWkR7AaKAn0Ar4XERO9Db/C7gASAGWichsVS17gmZjTK1jiaRmCFpCUdWdwE5vOUtEvgcq6qlzKfCmqh4CNovIBmCQt22Dqm4CEJE3vX0toRhjTBiplkp5EekA9AOWeEUTRWS1iEwVkaKJmlsD2/wOS/HKyis3xhgTRoKeUESkHvAf4E5VzQSeB04A+uKuYJ4I0OuMF5HlIrK8rF6qxhhjgiuoCUVEonHJ5HVVfRdAVVNVtVBVfcBLlNzW2g609Tu8jVdWXnkpqjpFVQeo6oCmTZsG/ocxxhhToaAlFHE9e14GvlfVJ/3K/buHjgLWeMuzgdEiEisiHYEuwFJgGdBFRDqKSAyu4n52sOI2xoShaSPcw4S1YLbyGgxcD3wnIiu9snuBa0SkL64p8RbgZgBVXSsis3CV7QXAbapaCCAiE4FPcc2Gp6pq6R5BxhhTCe+9994RowqvXr2aOXPmsGnTJl566aXi8oKCAtauXcu6detKNa/dsmULX3/9Nddee+0xv35mZiY9evTgsssuKzXWV20RzFZeXwJSxqaPKjjmEeCRMso/qug4Y4ypjFGjRjFq1Kji9SlTpvD6669z4YUXEhERwW233Va87d5776Vv375H9NXYsmULM2fOrFJC+ctf/sJZZ51V9R8gzNnQK8aYOunHH3/k4Ycf5tVXXy01dhbAF198waxZs5g8efIRx919990sWrSIvn37MmnSpHKHuz/cihUrSE1NZejQoaXK69Wrxx/+8Ad69uzJ+eefz9KlSzn77LPp1KkTs2e7u/tHGzo/XNjQK8aY0Pn4btj13ZHlu1aXXs/Lds9/b1u6vEXvI49t0QsuqnAwDfLz87n22mt54oknaNeuXalt+/fvZ+zYsbz66qvUr3/k8PiPPvoojz/+ePGYWk888USZw93HxcUVH+Pz+fjd737Ha6+9dsQw8tnZ2Zx77rk89thjjBo1ij//+c/MnTuXdevWMWbMGEaOHFnhzxJO7ArFGFPn/OUvf6Fnz57FQ877u+WWW7j++usZPHhwpc715Zdfct111wHlD3c/efJkhg8fTps2bY44PiYmhmHDhgHQq1cvhgwZQnR0NL169aq2YecDxa5QjDGhc5QriWJFLbzGzTnul1y4cCH/+c9/+Pbbb4/YNmPGDLZu3cprr7123K/jb/HixSxatIjJkydz4MAB8vLyqFevHo8++ijR0dHFw91HREQUDzsfERFRPOx8VYewr26WUIwxdca+ffsYN24cM2fOJCkpqdS2TZs2ce+997Jo0aLi6XfLkpSURFZWVvF60XD35557bqnh7v0VDYcPrj5k+fLlPPpoJZMp0KFDByZPnozP52P79u2lhs4PJ5ZQjDF1xgsvvMDu3buZMGFCqfJ77rmH+fPnk5OTc8Rsi88++yxnnnlm8Xrv3r2JjIykT58+jB07ttzh7gOpoqHzw0nQhq8PJRu+3pjwVZXh6wN5y8s4NW34emOMCQxLJDWCtfIyxhgTEJZQjDHGBIQlFGNMtauNdbc1SbDef0soxphqFRcXR1pamiWVEFFV0tLSSvXkDxSrlDfGVKs2bdqQkpKCTYQXOnFxcWX22j9ellCMMdUqOjqajh07hjoMEwR2y8sYY0xAWEIxxhgTEJZQjDHGBIQlFGOMMQFhCcUYY0xAWEIxxhgTEJZQjDHGBIQlFGOMMQFhCcUYY0xAWEIxxhgTEJZQjDHGBIQlFGOMMQFhCcUYY0xAWEIxxhgTEEFLKCLSVkQWiMg6EVkrIr/xyhuJyFwR+cl7buiVi4g8IyIbRGS1iJzsd64x3v4/iciYYMVsjDGm6oJ5hVIA/E5VewCnAreJSA/gbmCeqnYB5nnrABcBXbzHeOB5cAkIeAA4BRgEPFCUhIwxxoSPoCUUVd2pqt96y1nA90Br4FJghrfbDOAyb/lS4BV1vgEaiEhL4EJgrqqmq+o+YC4wLFhxG2OMqZpqqUMRkQ5AP2AJ0FxVd3qbdgHNveXWwDa/w1K8svLKD3+N8SKyXESW29SixhhT/YKeUESkHvAf4E5VzfTfpqoKaCBeR1WnqOoAVR3QtGnTQJzSGGPMMQhqQhGRaFwyeV1V3/WKU71bWXjPu73y7UBbv8PbeGXllRtjjAkjwWzlJcDLwPeq+qTfptlAUUutMcD7fuU3eK29TgUyvFtjnwJDRaShVxk/1CszxhgTRqKCeO7BwPXAdyKy0iu7F3gUmCUiNwJbgau8bR8Bw4ENQA4wDkBV00Xkr8Ayb7+HVTU9iHEbY4ypAnHVGLXLgAEDdPny5aEOwxhjahQRWaGqA6p6vPWUN8YYExCWUIwxlTdthHsYUwZLKMbUdZYkTIAEs1LeGFNbFObD/16FjBSQCDh0AGLrhToqE2bsCsWYuuyHT2DPevdY8qIrK+uKZfu38OFvYf8W2LcJfvqs2kM14c8SijF11bQR8H8TICcNDu6DLyeVv6+vwD037OitFwY/PlPjWEIxpq6LSYDEJi6xPDsQtn4FW7+ESSdBwaHS+4p9ZJjy2V+HMXVR6lrI2gUFuaXLCw4CCnENIGObqyuZNgI+/pPbftpt1R6qqTksoRhTV/jXjcz5HaRvgPwcaD8Yul4EhXmu0h2B+EZuvzeuCVm4puaxhGJMXaBa8gCXPGKTofVAGD0TBo2H+q2hfito3NnvOB9k73V1LP6++Kc1NTZHsGbDxtQFL54Fu1a75UXeWK0SAVGxEBkNzXuWVLiPmwNLpsDHf4Dty0qfJ6Fx6fWipDJuTvBiNzWGJRRjapuyPuT3/gixSZCXA3t/cmVt+sN1/yn7HO1OcfsXXdE07gLXvgWHstz6vs3u9lijzhBX3xKLASyhGFM7VOYDPbY+FOTBhs9d3Um7U0tv9z+2ZR9o0af09qTm7gpl4K9hzbtwMA0yt8GhZEhs6q52TJ1mCcWY2kDV9RV5bwLs+d7dznrp/JIP+cI8/52P7dz+iSYyCkY87lqJbV/u6lYO7nMdHhMaw6b/Qqchx/vTmBrKEooxNdm0Ee421K7VgELKkpJt8Q0hqlnJeq8rYN37sG+LW5fIis9d0dWORLgKffVB2o/uiicnHVa/ZQmlDjvmVl7ezIm9gxGMMaYKCg4BCkkt4cK/Q8t+rjwvG9I2QPpG98EPMOIJSG7nHufcU/XXHDcHfvURRERC0+7Qqj9Exhz3j2JqtkolFBFZKCL1RaQR8C3wkog8ebTjjDFBlLnTXW1kp7r1ei3gtFshItrbvh0O7IKsnRARBS16Q+fzoUE792jV7/hjGDen9JXMT3OtOXEdVtlbXsmqmikiNwGvqOoDIrI6mIEZY45i/YeQ6Y3+GxlbcoXw67muJdc7vwIEYuqBiLvlFSzj5pQ9VIupUyqbUKJEpCVu/vf7ghiPMaayipr0th4INx02+m+TLq5VV1mC2bQ374C7atq2DNoODN7rmLBU2YTyEPAp8KWqLhORTsBPwQvLGHPcqrtPSPOT4MdPID/F9aT/5dvV+/om5CpbKb9TVXur6q0AqroJsDoUY8LBNW+EOgLn2jfduGAx9UoaAZg6pbIJ5dlKlhljjKmjKrzlJSKnAacDTUXkLr9N9YGjNGI3xgTNtBGQuSPUURxp3Bx46dxQR2FC5Gh1KDFAPW+/JL/yTCCITUaMMcbUNBUmFFX9L/BfEZmuqlurKSZjTEXWvgfbFoPPq6cQCW08ZUlZ4a6ibLDIOqWyrbxiRWQK0MH/GFW1a1tjqlvqWjene1IrN4NiQqNQR3QkLYT8XNcvxQaNrDMqm1DeBl4A/g0UBi8cY0yFpo2A/d7Ngkad4PSJoY2nLJGxcCgTdix3Mz5e/26oIzLVpLKtvApU9XlVXaqqK4oeFR0gIlNFZLeIrPEre1BEtovISu8x3G/bPSKyQUR+EJEL/cqHeWUbROTuY/4JjTHV65Kn3Pwp0YmQvSfU0ZhqVNkrlA9E5FbgPaB4bAVVTa/gmOnAc8Arh5VPUtXH/QtEpAcwGugJtAI+F5ETvc3/Ai4AUoBlIjJbVddVMm5jaoeiUYUPpEJelhtuJVzrJ5p2hXrNISct1JGYalbZhDLGe/6DX5kCnco7QFW/EJEOlTz/pcCbqnoI2CwiG4BB3rYNXkdKRORNb19LKKbuOZDqBnuMiIYOZ4Y6moqNmwMzR7sBKk2dUamEoqodA/iaE0XkBmA58DtV3Qe0Br7x2yfFKwPYdlj5KQGMxZjw9/mDsPVLtxwZA20GwZjZIQ3JmLJUKqF4CeAIqnr47ayjeR74K+7q5q/AE8CvjvEcZRKR8cB4gHbt2gXilMaEhz0/uKuSpBZuWBNjwlRlb3n5DxsaB5yHmxflmBKKqqYWLYvIS8CH3up2oK3frm28MiooP/zcU4ApAAMGDDjGOU6NCXORMdCgffjWm5QnfaP1R6lDKnvL63b/dRFpALx5rC8mIi1Vdae3OgooagE2G5jpTdrVCugCLAUE6CIiHXGJZDRw7bG+rjE11rQRsLsWVBkWTbpliaVWq+qc8tlAhfUqIvIGcDbQRERSgAeAs0WkL+6W1xbgZgBVXSsis3CV7QXAbapa6J1nIm7o/EhgqqqurWLMxoSvsj5wM3fAwXQozIfGnWrmh3Fhvvs5YhIhLjnU0Zggq2wdyge4JADug707MKuiY1T1mjKKX65g/0eAR8oo/wj4qDJxGlNrTBsBu76DQxluPbYGfhgnt4HCPNi3ya0nNoUvJ8EZvw1tXCZoKnuF4t9vpADYqqopQYjHmLqlMB/+NQjSN7v1lTOh60Vu5kNfPsQkuR7xV78a2jirYvhjsHM1ZO92Vyk56TDvr5ZQarFK9ZT3BolcjxtxuCGQF8ygjKkz8g9C+iaIqw8ozP8bPNkTdq6E/BxXGR+bFJ7jdR2NCERGQ/3W0GagG3vM1GqVSigichWukvxK3LzyS0TEhq83JlDiG7ne74cOQGGuu8XVtBtM+Kpm1p0UGTenZsdvjkllx/K6DxioqmNU9QZcL/a/BC8sY+qI1690zwNvgsadXZ2Jr9BdsSQ0gfotQxtfoIybA32uDnUUJsgqW4cSoaq7/dbTqHwyMsZUxoTFMN1r7XXjp6GNxZgqqGxC+UREPgXe8NavxlpeGVN1L18Iu1a7Snnw6huiIKKqLflrCC20jo612NHmlO8MNFfVP4jI5cAZ3qbFwOvBDs6YWungfsjLdpXucckQnQDdLnbb7IPW1GBHu231FG7+eFT1XVW9S1Xvwg1j/1SwgzOm1tn7E/yzE6R+59brtYRGJ0DD9qGNqzqlbYAvHge1EZJqm6MllOaq+t3hhV5Zh6BEZExtlr3H3fZJauUmoYpvGOqIqk+rvq4ZcfYemP9XOLgv1BGZADvaDdsGFWyLD2QgxtQp8Y0gvkHdusXV/RJoc4rr5FjUe97UKkdLKMtF5Neq+pJ/oYjcBFQ4BbAx5jDTRkCuN5TKsEeg09mhjMaYgDvaLa87gXEislBEnvAe/wVuBH4T/PCMMbXKuDlw6i1ueebo0MZiAq7CKxRv/pLTReQc4CSveI6qzg96ZMbUFtNGuCFW0n4saSZsHFXXZNrUCpWdD2UBsCDIsRhTOxXkugroQ1kQW981EW7RO9RRhY54N0ZSlsA/2sPt/4PExqGNyQSE9XY3JpjWvAvbl5dUQjfuAle8XDMHewyUnqPc7JMJTVydUvbuox9jaoRa3i3XmBDLSXPPjTpBZBxExYU2nnCQ2ASS20L2XsjZG+poTABZQjGmOiQ0hZs+C3UU4WPcHFj7Hrw9NtSRmACyW17GVIdrZoY6AmOCzq5QjAm0ovnhAbJ2hC6OmuK9CXDzwlBHYQLAEooxgVLUPDhnrxv80Vfg5lQ3po6wW17GBIKqmxgrMwX2b3VJJXe/ayLbdQTEVTSKUV3l9T9J/Q6e7gv7fw5tOOa42RWKMVVVdGtr3Bz49D7YttitR8ZC6wFuWcTqT8rT8Sw3SGZBLuzbDOmboUG7UEdljoMlFGOqIn0T5KS75axU94EYGeM+IGMSrfd3ZSQ0cs2pczPgYHqoozEBYAnFmKp44xrYs94tf/Q79xwRDclt6tYIwsdr3BzY8hVMHx7qSEwAWEIxpiryc9xcJgWHYPMiV1fSuJMlE1OnWaW8MVUVEQ0S6a3Y7IPGWEIxpqo6nwfNe7rWXAf3QYRd8Ju6LWgJRUSmishuEVnjV9ZIROaKyE/ec0OvXETkGRHZICKrReRkv2PGePv/JCJjghWvMVUy9G/QoIN7DHs01NHUbJ/cW7pTqKlxgnmFMh0YdljZ3cA8Ve0CzPPWAS4CuniP8cDz4BIQ8ABwCjAIeKAoCRkTFtoOdBXxyW2g3amhjqZmS98Iu1ZD2sZQR2KqKGgJRVW/AA5vC3gpMMNbngFc5lf+ijrfAA1EpCVwITBXVdNVdR8wlyOTlDHVa9oI11S4yLg5Vhl/PJr3gB6XQWQ0HMqEHf8LdUSmiqq7DqW5qu70lncBzb3l1sA2v/1SvLLyyo0JjfUfQcY20IJQR1J7xDeEq2ZA486hjsQcp5BVyquqEsCmMSIyXkSWi8jyPXv2BOq0xpT27q/d0Cq+QtcpzwTO5VNCHYE5TtWdUFK9W1l4z0VTtW0H2vrt18YrK6/8CKo6RVUHqOqApk2bBjxwY9zgjzmuN3y702HIH0MdkTFhpboTymygqKXWGOB9v/IbvNZepwIZ3q2xT4GhItLQq4wf6pUZU31UYd5fYe8PoD43rIpYi3tjDhe0hvMi8gZwNtBERFJwrbUeBWaJyI3AVuAqb/ePgOHABiAHGAegquki8ldgmbffw6pqg/6Y6lHUhPWaN2DR466fSVQcDH8MulnzVmMOF7SEoqrXlLPpvDL2VeC2cs4zFZgawAqYsSUAABulSURBVNCMqVhuJnz0B9i9zl2JZHp3WZPbQv3WlkyMKYd17TXmcKlrYfWbbvTgwjx485eufOBNcPrE0MZmTBizG8HGlCfZaw+Slx3aOIypISyhGFOeIX8CBA7scutx9UMaTp3x33/aECw1lN3yMsbftBFuwieAes3gjm9h1ljXsqvvdSENrc7ISXPTAmTvhcQmoY7GHAO7QjGmyM5V8PNiSPXGM5UI13kxNgli6kGE/bsEVUJjiEuGnL2wfwus+U+oIzLHyK5QjCmSvhm0EOo1h0G/hraDXLmN01U9EhrBHzfD1IsgZYkbjcDUKPaVy5jDJbWCs/4A0fGhjqTuiYiEX84KdRSmiiyhmLpt2ogjK4AvmxyaWIyp4eyWl6m7svdC1i5A4blT3LAqJnwsnQLrbWqAmsQSiqm7lk+D9A1HlttUvuEhL9u1+MrNtCbbNYTd8jJ1V2Gee251csk0viOfhSZdQhmViYp1Y6Zl74Y938M3z4c6IlNJ9lXM1C3+9SX7twIC4xeUlJ98Q0jCMn6i4+E3q2DmaNi50k0ZYGoESyim7sjLgUNZJY98vyFV7D59eElqATcvhL82C3Uk5hhYQjF1w7QRsPdHdxulSFQcnDgsdDEZU8tYQjG1X8pyyEiBvCx3f77hCe62SnQ8XPtmqKMzptawhGJqt2kjYNdqOJTp1hMaux7ZdovLmICzhGJqr7kPlFTqxjWApt3hVx+7gR6NMQFnCcXUXmvfdU2D45Jh+OPQ64pQR2RMrWYJxdRucQ2gyYmWTIypBtax0dRunc+z+pKa7n+vw/NnwKEDoY7EHIUlFFNz+Xwwdbh7mNopJhFy9kDqdzDv4VBHY47CEoqpmfZtgb+3gZ+/gp+/hufPtGlja6Nfz4fmJ0FENOTZFUq4s4RiaqbMna6ne3xDQCHjZ9i32fU3MbVHo44w4SvXc96EPUsopmar19z1eD+UBZnbYe3/ufJpIyArNbSxGVPHWEIxNduwR+HPqdBmoFv/+ln4W3PYtgQKDwHW56TW+Gmu3dYMc5ZQTO1ww/sQGQMHdkFBLvjy3dXLgHGhjswESu4+Vzm/6q1QR2LKYQnF1A7RcXDLV9DsJFeJ2/Y0aNwF2g4KdWQmEPpcA1Hx7tbmd2+HOhpTjpB0bBSRLUAWUAgUqOoAEWkEvAV0ALYAV6nqPhER4GlgOJADjFXVb0MRtwkDhfnwwpllt/hpeiLc+lX1x2SC79z7YOvXbigdE7ZCeYVyjqr2VdUB3vrdwDxV7QLM89YBLgK6eI/xgE3fVpdtnO9m8cvY5tYTGoU2HlN9xs1xox6YEtNGhFW9Ujjd8roUmOEtzwAu8yt/RZ1vgAYi0jIUAZoQmzYCPrvfLTftDr9dCy37hDYmU/1SlofVh2i18k8gC/7uRtLetRpWvlH2/lMvcg9fYeXOfZxCNZaXAp+JiAIvquoUoLmq7vS27wKae8utgW1+x6Z4ZTsxdceOle52R2G+W7/iZUhuE9qYTGjk50DaBvj+A+h+SaijqR6q8OWT7udG4cmergEKAupzA6G2Pw1e88as+8W/IXWN6/QL8GQPuHO1mw+oyNThbsK5wnzXjyvi+NNBqBLKGaq6XUSaAXNFZL3/RlVVL9lUmoiMx90So127doGL1IRWYQG8dJ4bfiPvgBvs8aRfQOPOoY7MhEK702DXd3AgFb54vPYnFFXXCGH+I7B/C0gkqN/VRv02kLsftnwFT/tdrb881Ov0i3s+sAumXwKR0VBwyLWETF2D+24fOCFJKKq63XveLSLvAYOAVBFpqao7vVtaRXO1bgfa+h3exis7/JxTgCkAAwYMCOy7ZEJnyfOwy68itmk3GGXVaHXWsP8HO1fB7rUE+sMwpHw+eOkcd7WhCihc/JSrI3z31yX7Ne4C8Q3cFy2Am+bCzKtgw1y3Hp3gxj/L3uOSSEw9GHI3fPIn2L68dDKKioXYZEhu65rc3/gJ3H98/baqPaGISCIQoapZ3vJQ4GFgNjAGeNR7ft87ZDYwUUTeBE4BMvxujZmabuMC+OA3bvmqV6BV35Jt00bA/q1uuVkPuOz50ttN3TRuDsy8GrJq0cfA6reObMH28vkw6kW33KgzJDaBGz898tirX4Np3gCp4z6C3AyYfrFbv+xfLqnEN3SJqjAP4hvBKTe7h/8tsAAIxRVKc+A91xqYKGCmqn4iIsuAWSJyI7AVuMrb/yNck+ENuGbD1lOtNplzV0nSWPAI/NLrY7D2PdjxrVdnInDr4pCFaMJU2gb3paMmTk9QWAD/Pt8tX/+uSwIATbq5uoz8bDc2XVHz+HPuKX9On+g4GD+/9Prty0rv06ynew7ye1XtCUVVNwFHNM1R1TTgvDLKFbitGkIz1Wn39/D2ONi/DRKbQv5Bdx/4uVMAhaxdriy+EQz8VaijNeEoLwdSlsH/3QqXTQ51NMdm3kOw839u+V+D4Mzfu+Vffexuc235EqaPgIWPBub1qinp2oyNJjSmDYeD6W75pF+4b5sbPoe9fu0zmp/kRpo15nCDxsP2b+FQpuubVNMc3OeuRGKTXX3H/MPmeklo7J6z97jnxCbVG18VWUIx1WfaCK/S0edapiQ2haTWcMHDruXJjEsAARH3PPaDUEdswlXn8+APP8Hs292gkVDSjyIcb4HlZsD0kW557Gz3nNgMLnkaZt3g+ol0HAKx9d22Zt3h9xtg5miQCOh0diiiPmaWUEzgqbqmvloIV05z/zgvXwB7fwRfQcl+/a6H8x9wy1GxcPN/QxOvqflWz4Kt3tXsU71h4rKAVzgfl1k3lLRWfLSdm3IhvhGcOBRa93flY2aXPqZeUxg/r3rjPE6WUEzgbf4Cdqxwy8/0g1+8DLvXucv7qFiIiHT/UIN+XfF5jKmMrF3w/kRAXT+l/VshLzu8EkpuhmttFRHlOmaqwgnnuG3heEVVRZZQzPEp6zZD7n73XK+Fawv/iTcs2/DHoM/V1Rufqd26jiiZVO3kG6BJF/j4j/DKpfDLdyCpecXHB8vqt+HzhwB1X54yt0PHM10rxqL/mZrWkKASLKGYqpk2Ag7shoytgLh+AdceNk/FRf+At8fWuIpFU4N0HQb3+I3MtGK6e961Gp7sBr9dB/VDMPTf7ImuN7q/mET3XIuuSA5nCcVU3revwNz74eB+V1FY1Os2OgF+/ASe7usu5fOzXXmTLvD7n7yKRXEVqcYEU6+rYPHzcCjDdXzMzajehLJ4sps1tCDXDYuS1MI1fwcY/kT1xREillDMkb57B+Y9jGtp9SHEJsErl8HeH1xrrJgE1yolvhmc94D75/nwTtcRq0hsfVcZXwMrFk0NFpMAE5fAmnfhnXHw9hi45k1o1DF4r5mX7QZl3L7Cm3YaV094wUPld0aspSyhmNLysuE/N5asL3/ZtfPfucqtdx0B18w88riVb5Q09xUvEUVEVkvIxhwhOt4971kPz/SFm78I3lQHGdvdqL4x9dwXqQbt3OvXsWQCllDM4dTnnht0cENa+wrcLa7Y+lC/FYx4vOzjbvqs2kI05qg6n+/Gfzt0wP0dv3GN+8C/ZdHxt/7KzXCzhqoP+o+B5dNc+SVP18kk4i+cJtgy4WTQTa51SpGuw+G2JS6pGBPuIqPd+G+/fBsQ18pq7w/wYxmDKx6uMB+mnOseedlHbt+4wDVNztgG8//mzg1Qv3VAf4SayK5QTMWWT4X83KPvZ0w4atbNNQx55TLYvcaNtpubAa9e7oY/KZqw7br/QFNveuEP7izpRzX1Qrjly5LzvXQ+HExzyxMWu3l5Zng94NufVj0/UxizhGLK12aAG7ARbKpdU3PVawpXTYfnBsBHfygZQ87fvwa6K4wrpkLOXoiMdbfGUtfBM/3hkkmQvRe2+43iGx0HUd48IgawhGIqMsbG0jK1REJjVw9YlEziG7mycR/B7Dtcs/fM7fDaL9z2Zt2gy1D44jFI3+DGmYv2+pE07AQXPgINg9hyrIayhFJXHcqC5we7S/4mnWHXGlcWGR3qyIwJvIRG8MfNrpGJRLgriyLXvumeX7kMNi1wy3HJcM59MPAmNwr2nN8D6ir7r3kLIu2jsyz2rtQ12WluqtHsvSUdELN2uLmqi5r5xjd2/zjG1CaRURUngqte8Ua8Bq5+3TV/T2oB/a5zD3NUllDqiuy98PKFrkIyxxsKJS4ZJq5w95iNqevi6tuI18fJEkpdsXyauxcM7l7y735wvYqNMSZALKHUNpv+66ZEBRj5NLQ73fVyT9/kyn7/kxukzpKJMSbALKGEO9WSZZEjyzfMg4+8+aivfxcW/D/ITHHrRS1WisTUg3rNgherMaZOs4QSznb8D1461w3xEJ0A4//rOl8VFriJqzJ+Lr3/M/0AgU7nwOA74N2bXXm3EdDzMjf6qTHGBIkllHCUmwlTzoG8LJdM4hu6Xr3TR7hB56571yWTuAZuJODRr7tZEr962h0/aDyccK6bc9sYY6qJJZRwcijL1YFsWlBSgR4R7dq9vzoKsne7sue8OaiH/AlO8+pLWvaB02+v/piNMcZjCSWcLPs3fP6gW46MgduWulkOY5PgnhQ3odXqtyAn3c1N3dum0zXGhA9LKNXp8Ar2xZNLEoivoGQOhwlfu2EhklqU7B8RAURYBytjTNiyhBJMB/e7uo/P/gybFkLegZJt0YluKHgtcL3Uwc182H4wNO8ZknCNMeZ4WEI5Xmkbvcl2CiG5DVzyDKx6A374pKRHepG4BtC6v6tAz8+GtJ+g5+Vw5bTQxG6MMQFkCeVY7dsCK6bDnh9cYtBCyD8IiU0hbQNMH16yb0S0q+foeBZ0ucANUAelb30ZY0wtUWMSiogMA54GIoF/q+qj1fLCP3wCq2bC1sXu9pUvv/T2tqe4zoKXveBGJT2U5WY67H5xSZ3I4fw7KJbhUEEhX29MI7/AR3J8NKd0alzl8FWVtOw8VCEpLoq4aJvn3ZhAyz5UQE5eIVERQsPEmKPuX1DoY/PebPIKfRQUKgU+pUFCNO0blYxgISJERhz5WbF+VyY/p+UAMKhjIxokHP31qkuNSCgiEgn8C7gASAGWichsVV13TCdSJXfjl+TnZBAZISTERLnK8AbtIX2ju9L47M+Q483IVjS/OnCwQWeiC/LZ3+4Ccht05uvm1/LOihRiicCXpfw+NZ+3fzyBSC9Z/LJZPsnxUFCovLMihfTsPHyqxEdHEhUZQUyk269Vg3gSYqNYsimNCBGiIoX6cdE8Pa+kD8m/rj2ZLWnZZSaD2at2UD8uikU/7UUEBLhlyAnceEZHVqXs57n5G/j25/3F+088pzMvfrERwb3+lQPakBATyf6cfN5ekYIAESL8e+wAXvzvRpZsdvNH1IuJ4rzuzYiIECK9P/SICKFLs3oM7dmC9TszSc/O46fdB4iLcjNLZ+YWsHRzOo3rxZCTV8jIPq3o0CSRUzo2Yu66VH5KzQKR4tcsil/E/TO5dfcc4bd84FABXZol8fn3qcRGRfDtz/tonBjL4k1pLlcrnNy+IRed1IKmSbG09fsn7dA4kYSYSNbuyCTf+2dOz8lj+76DJMREEhUpZOUW8OrirTRIiCYuOpJ/XtGb1MxcdmXk8vXGNBJj3O8wKsL9vhJjo+jZKplIESIiKH5/oiIjaJgQTYS3HhkhRIgQHxNJvdgoFvywm33ZeRT4lB4t65MQ436/q1L2k57tvric0bkJXVskAZBxMJ9vt+5DUVo3SKBNw3jmrd9NQaGPBgnRnNO1GeL9/W3ac4D9B/OJihB6tKzPqpQMtqZl4/O+WDRNiqVLs3os3ZzOe/9zU9h2aZbE6EFtyS/0kZ6dxzeb0oiMiKBBfDSXn9y6+NxF9mXnsS8njwgR2jdOKN6+M+Mg3+/MBKBri/o0Tozh+52ZKNAyOY6WyeV80QL2ZB0iKzef6MgI2jSMJ6/QR+bBAgp87ndV6FNaNYgnJqr07OXb0nN4Z0UKqkqL5HiuGdQW4IiYAdZsz2B3Vi4Fhe5cCTGRtG+cSGpmLut3lcSdduAQby93o060SI7jiv5tKPApO/cf5OUvN7vfcYQw5vQOjJr8dfH5n7u2H4mxUWzak42qEhsdybodmWTlut9pVm4BX2/cS37h0e9UDO/VgraNEvD5lEIf+FSZ/vWWUvv0bFWfXq2TuWqg+5mXbk7nu5QMRKDQp1zatzWFPmXD7gO8v3I78d7/+4nN67F2RyZNk2LJL/Txy1PaHzWeoxGtAbdfROQ04EFVvdBbvwdAVf9e1v7dO7XR6f/vNmLzM4gszOUgMRzMK6RR5np6HfiyrEOOsN7XllyJw6fK54X9mVx4aaB+nEqJELh3eHf+Nuf7Su3ft20DVJWdGblk5uaTm+8rtf2Epols3FMyP3anJols2uvWE2IiESDf+2BbuW1/qWOb1IulXmwkPnV/oD51/9i7sw5VKrb46EgO5hdWat+qiI2KIDk+mmb1Y2mYEMOin/YG5LxN6sWy98AherSszzrvA7JIvdgocvIK8FXx36dv2wZHvM/lSY53c9RkHMyvcL9uLZJYvyuragEdQxzg/g4OHCooc/vR4uzTtgG5eYVERgiFPpfQ91Tyb6lIl2b1vKQtZBzMZ1v6wTL3u6BHc+auS6W+d3XuU2Xvgbxjeq3KGtihIcu27Kv0/kmxUTx+VR+iIoTM3Hw27D5AXJT7UvFzeg5vr3DJTATioiK9LyRQ4FPGDe5A2gH3vs1bv/uY4myZHMeuzFwSoiPJ9ylNEmPYkeGm+d76j4tXqOqAYzqhn5qSUK4AhqnqTd769cApqjrRb5/xwHhv9SRgTbUHeuyaAIH59AsuizOwLM7Aqglx1oQYAbqqalJVD64Rt7wqQ1WnAFMARGT58WTZ6mJxBpbFGVgWZ+DUhBjBxXk8x0ccfZewsB1o67fexiszxhgTJmpKQlkGdBGRjiISA4wGZoc4JmOMMX5qxC0vVS0QkYnAp7hmw1NVdW0Fh0ypnsiOm8UZWBZnYFmcgVMTYoTjjLNGVMobY4wJfzXllpcxxpgwZwnFGGNMQNT4hCIiXUVkpd8jU0TuFJFGIjJXRH7ynhuGQay/FZG1IrJGRN4QkTivocESEdkgIm95jQ5CGeNvvPjWisidXllYvJciMlVEdovIGr+yMmMT5xnvfV0tIieHMMYrvffTJyIDDtv/Hi/GH0TkwuqIsYI4HxOR9d779Z6INAjTOP/qxbhSRD4TkVZeeUh+5+XF6bftdyKiItIkHOMUkQdFZLvfZ+hwv23H9ntX1VrzwFXY7wLaA/8E7vbK7wb+EeLYWgObgXhvfRYw1nse7ZW9AEwIYYxFHUITcA02Pgc6h8t7CZwFnAys8SsrMzZgOPAxbjSXU4ElIYyxO9AVWAgM8CvvAawCYoGOwEYgMoRxDgWivOV/+L2X4RZnfb/lO4AXQvk7Ly9Or7wtrjHRVqBJOMYJPAj8vox9j/n3XuOvUA5zHrBRVbcClwIzvPIZwGUhi6pEFBAvIlG4D+2dwLnAO972UMfZHffHnaOqBcB/gcsJk/dSVb8A0g8rLi+2S4FX1PkGaCAiLUMRo6p+r6o/lLH7pcCbqnpIVTcDG4BBwY7Ri6msOD/zfu8A3+D6e4VjnP7j4CQCRS2LQvI7Ly9OzyTgj34xQnjGWZZj/r3XtoQyGnjDW26uqju95V1A89CE5KjqduBx4GdcIskAVgD7/f6JU3BXMqGyBjhTRBqLSALum1Rbwuy9PEx5sbUGtvntF+r3tizhHOOvcN+iIQzjFJFHRGQb8Evgfq84rOIUkUuB7aq66rBNYRWnZ6J3+22q3y3tY46z1iQUr+5hJPD24dvUXb+FtH2090u6FHfp2Ar3zWpYKGM6nKp+j7vV8RnwCbASKDxsn5C/l+UJ59hqEhG5DygAXg91LOVR1ftUtS0uxolH27+6eV/I7qUk2YWz54ETgL64L7tPVPVEtSahABcB36pqqreeWnQZ6T0f25CcgXc+sFlV96hqPvAuMBh3uVvUwTTkQ8qo6suq2l9VzwL2AT8Sfu+lv/JiqwnD9YRdjCIyFrgY+KWXoCEM4/TzOvALbzmc4jwB9+VxlYhs8WL5VkRaEF5xoqqpqlqoqj7gJUpuax1znLUpoVxDye0ucEOzjPGWxwDvV3tEpf0MnCoiCSIiuPqedcAC4Apvn5DHKSLNvOd2uPqTmYTfe+mvvNhmAzd4LWpOBTL8bo2Fi9nAaBGJFZGOQBdgaaiCETeJ3R+Bkaqa47cp3OLs4rd6KbDeWw6b37mqfqeqzVS1g6p2wN0uOllVd4VTnFD8RazIKEpGaj/233t1tS4I5gN3+ygNSPYrawzMA37CtVZqFAZxPoT7418DvIprPdHJ+yVtwN2uiw1xjItwiW4VcF44vZe4Lww7gXzcP+iN5cWGa0HzL1zLlO/wa10VghhHecuHgFTgU7/97/Ni/AG4KMTv5QbcPfOV3uOFMI3zP97/0GrgA6B1KH/n5cV52PYtlLTyCqs4vc+i77z3czbQsqq/dxt6xRhjTEDUpltexhhjQsgSijHGmICwhGKMMSYgLKEYY4wJCEsoxhhjAsISijHGmICwhGJMNRCRDiJyUERW+pU1F5GZIrJJRFaIyGIRGVXBORYcPoS4uKkanheReG/o8byiYdKNqW6WUIypPhtVtS+4OTGA/wO+UNVOqtofN7hpmwqOf8Pbx99o4A1VPeide0cQ4jamUiyhGHMYEXlbRJ4TkS9FZKuInCEir4rIjyLycoBe5lwgT1VfKCpQ1a2q+qwXw3UistS76nhRRCJx0xyM8AZCRUQ64AYaXRSgmIw5LpZQjDlSL2CTqp4BvAi8jBvjqgfuAz02AK/RE/i2rA0i0h24GhjsXXUU4gZrTMcN03ORt+toYJbacBcmTEQdfRdj6g4RiQMaAE95RQq8rN7gfSJSCOQF4XX/BZzhnXsG0B9Y5u6MEU/JKMpFt73e955vDHQsxlSVXaEYU1pP3DQIPm+9D7AEQETa4OoouovI3V7ZsyKSJCI9Di87yuusxU3FCoCq3oYbgbopbvDAGara13t0VdUHvV3fB87z5iFPUNUVAfiZjQkISyjGlNYLN9Jykd64UVjBJZfVwEDcaLzgRrjOKqesIvOBOBGZ4FeW4D3PA67wm0qgkYi0B1DVA7gpD6ZSeroGY0LOEooxpfXCSwze7a94Vd3nbStKLgOBdSKS6HdcWWXl8uo9LgOGiMhmEVmKu9X1J1VdB/wZ+ExEVgNzAf85K97AJTdLKCas2PD1xhwjEZmDm0siE+ilqsPKKjvsmA7Ah6p6UpBj24KbX2NvMF/HmLJYpbwxx0BEooE0Vb25orIyFALJIrKyqC9KgOOKBxYD0YDvKLsbExR2hWKMMSYgrA7FGGNMQFhCMcYYExCWUIwxxgSEJRRjjDEBYQnFGGNMQFhCMcYYExCWUIwxxgSEJRRjjDEBYQnFGGNMQPx/+9rZOJ59wCsAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "ax = hist.plot1d(output['mass'], overlay='dataset')\n",
    "ax.set_xlim(70,150)\n",
    "ax.set_ylim(0, 3000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.1, 7500.0)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAELCAYAAAD+9XA2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3dd3xV9f3H8dfnJoGEFTYikeFiDyGAynLiLKitFesAan/UXWtbq22dra2tA+vWVgW3lmpBq1VUqKMogqIiyBBBNmFvyPj8/viehAQSCHBvbhLez8fjPO4537M+hxvyyfd7zvl+zd0RERHZX7FkByAiItWDEoqIiMSFEoqIiMSFEoqIiMSFEoqIiMSFEoqIiMRFarIDSITGjRt769atkx2GiEiVMnXq1JXu3mRf96+WCaV169ZMmTIl2WGIiFQpZrZgf/ZXk5eIiMSFEoqIiMSFEoqIiMRFtbyHIiKVU25uLosWLWLr1q3JDuWAlp6eTlZWFmlpaXE9rhKKiFSYRYsWUbduXVq3bo2ZJTucA5K7s2rVKhYtWkSbNm3iemw1eYlIhdm6dSuNGjVSMkkiM6NRo0YJqSUqoYhIhVIySb5EfQdKKCJSKaWkpNCtWzc6derEueeey+bNm0vd7thjj63gyMp23HHH7fEduHvvvbfMa9lXEydO5H//+19cj7kvlFBEpFLKyMhg2rRpTJ8+nRo1avDII4+UWJ+XlweQkF+khcdOBCWUfWBmbc1sWrFpvZldY2YNzWy8mc2JPhtE25uZ3Wdmc83sCzPrXuxYQ6Pt55jZ0ETFLCKVU79+/Zg7dy4TJ06kX79+DBo0iA4dOgBQp04dIPxSHTBgAIMHD+bQQw/l+uuv59lnn6VXr1507tyZb775BoBXX32V3r17c9RRR3HSSSexfPlyAG655RYuuugi+vTpw0UXXUT//v2ZNm1aUQx9+/bl888/LxHXli1bGDJkCO3bt+fss89my5YtResuu+wysrOz6dixIzfffDMA9913H0uWLOH444/n+OOPL3M7gOuvv54OHTrQpUsXfvnLXwKQk5PD97//fXr27EnPnj358MMPmT9/Po888ggjR46kW7duvP/++3H9t98r7p7wCUgBlgGtgL8A10fl1wN/juZPB94ADDga+DgqbwjMiz4bRPMNdne+Hj16uIhUPjNmzCj3trVr13Z399zcXB80aJA/9NBDPmHCBK9Vq5bPmzdvl+0mTJjgmZmZvmTJEt+6dasffPDBftNNN7m7+7333us/+9nP3N199erVXlBQ4O7uf/vb3/zaa691d/ebb77Zu3fv7ps3b3Z391GjRhXtM2vWLC/t98rdd9/tw4cPd3f3zz//3FNSUvyTTz5xd/dVq1a5u3teXp4PGDDAP//8c3d3b9Wqlefk5BQdo7TtVq5c6UceeWRRnGvWrHF39/PPP9/ff/99d3dfsGCBt2vXrij2O++8s9z/tu6lfxfAFN+P3/UV1eR1IvCNuy8ABgOjo/LRwFnR/GDgqei6PgLqm1lz4BRgvLuvdvc1wHjg1AqKW0SSZMuWLXTr1o3s7GxatmzJJZdcAkCvXr3KfNy1Z8+eNG/enJo1a3LYYYcxcOBAADp37sz8+fOB8OjyKaecQufOnbnzzjv56quvivYfNGgQGRkZAJx77rm89tpr5Obm8sQTTzBs2LBdzvfee+9x4YUXAtClSxe6dOlStO6ll16ie/fuHHXUUXz11VfMmDGj1JhL2y4zM5P09HQuueQSXn75ZWrVqgXA22+/zZVXXkm3bt0YNGgQ69evZ+PGjXvxr5pYFfUeyhDg+Wi+mbsvjeaXAc2i+RbAwmL7LIrKyiovwcxGACMAWrZsGbfARSQ5Cu+h7Kx27dpl7lOzZs2i+VgsVrQci8WK7otcddVVXHvttQwaNIiJEydyyy23lHrsWrVqcfLJJzN27Fheeuklpk6dWu7Yv/32W+666y4++eQTGjRowLBhw0p9TLes7VJTU5k8eTLvvPMOY8aM4YEHHuDdd9+loKCAjz76iPT09HLHUpESXkMxsxrAIOAfO6+Lqlgej/O4+2Punu3u2U2a7HPvyyJSza1bt44WLcLfpKNHj97ttj/5yU+4+uqr6dmzJw0aNNhlff/+/XnuuecAmD59Ol988QUA69evp3bt2mRmZrJ8+XLeeOONon3q1q3Lhg0bdrvdxo0bWbduHaeffjojR44sunczcOBA7r///qJjFSbc4sdMpopo8joN+NTdl0fLy6OmLKLPFVH5YuCQYvtlRWVllYuI7LVbbrmFc889lx49etC4cePdbtujRw/q1avH8OHDS11/2WWXsXHjRtq3b89NN91Ejx49AOjatStHHXUU7dq140c/+hF9+vQp2mfEiBGceuqpHH/88WVut2HDBs4880y6dOlC3759ueeee4BwU3/KlCl06dKFDh06FD359r3vfY9XXnkl6TflLVQSEngCsxeAN939yWj5TmCVu99hZtcDDd39OjM7A7iScHO+N3Cfu/cys4bAVKDwqa9PgR7uvrqsc2ZnZ7vGQxGpfGbOnEn79u2THUa5LVmyhOOOO46vv/6aWKx6vWVR2ndhZlPdPXtfj5nQfyEzqw2cDLxcrPgO4GQzmwOcFC0DvE54gmsu8DfgcoAocfwe+CSabttdMhERiYennnqK3r17c/vtt1e7ZJIoCa+hJINqKCKVU1WroVRnVa6GIiIiBw4lFBERiQslFBERiQslFBGp1M57dBLnPTop2WFIOSihiMgBpbBb/I4dO9K1a1fuvvtuCgoK9vl4hZ1T7mzYsGGMGTNmt/sOGzaMWrVqlXgp8ZprrsHMWLly5T7HlCxKKCJyQCns0uWrr75i/PjxvPHGG9x6661Ji+fwww9n7NixABQUFPDuu+8Wvclf1SihiMgBq2nTpjz22GM88MADuDtbt25l+PDhdO7cmaOOOooJEyYAMGrUKK688sqi/c4880wmTpxYtPzzn/+cjh07cuKJJ5KTk7PLeaZOncqAAQPo0aMHp5xyCkuXLi1aN2TIEF588UUgdMHfp08fUlNDN4vz58+nU6dORdveddddRX2PTZs2jaOPPpouXbpw9tlns2bNGiAM8vXrX/+aXr16ceSRR1bom/MV1TmkiEgJt776FTOWrN+lfMbSkmWbt4VOHTvf8maJ8g7N6+2yb4eD63Hz9zruVRyHHnoo+fn5rFixgmeeeQYz48svv+Trr79m4MCBzJ49e7f7b9q0iezsbEaOHMltt93GrbfeygMPPFC0Pjc3l6uuuoqxY8fSpEkTXnzxRX7729/yxBNPAHDkkUcybtw41qxZw/PPP8+FF15You+vslx88cXcf//9DBgwgJtuuolbb72Ve++9FwgDhE2ePJnXX3+dW2+9lbfffnuv/k32lRKKiEjkgw8+4KqrrgKgXbt2tGrVao8JJRaLcd555wFw4YUXcs4555RYP2vWLKZPn87JJ58MQH5+Ps2bNy+xzTnnnMMLL7zAxx9/zKOPPrrHONetW8fatWsZMGAAAEOHDuXcc88tcTwIfZEVdttfEZRQRCQpyluTKHzC68WfHpOQOObNm0dKSgpNmzYtc5vU1NQSN+5L64q+kJmVWHZ3OnbsyKRJZT+pdt5559GjRw+GDh1aopuXvTlvcYXd9qekpCR0OOOd6R6KiBywcnJyuPTSS7nyyisxM/r168ezzz4LwOzZs/nuu+9o27YtrVu3Ztq0aRQUFLBw4UImT55cdIyCgoKip7mee+45+vbtW+Icbdu2JScnpyih5ObmlhjUC6BVq1bcfvvtXH755SXKmzVrxooVK1i1ahXbtm3jtddeAyAzM5MGDRoU3R95+umni2oryaQaiogcUApHgszNzSU1NZWLLrqIa6+9FoDLL7+cyy67jM6dO5OamsqoUaOoWbMmffr0oU2bNnTo0IH27dvTvXv3ouPVrl2byZMn84c//IGmTZsW3WAvVKNGDcaMGcPVV1/NunXryMvL45prrqFjx5I1tJ/+9Ke7xJqWlsZNN91Er169aNGiBe3atStaN3r0aC699FI2b97MoYceypNPPhnPf6Z9os4hRaTC7EvnkIlu8jpQJaJzSNVQRKRSUyKpOnQPRURE4kIJRURE4kIJRURE4kIJRURE4kIJRUQqtyfPCJNUekooInLAeOWVV+jWrVuJKRaL8cYbb/Dggw+WKO/UqRNmxsyZM0scY/78+Tz33HP7dP7169eTlZVVoqPJ6iShCcXM6pvZGDP72sxmmtkxZtbQzMab2Zzos0G0rZnZfWY218y+MLPuxY4zNNp+jpkNTWTMIlJ9nX322UybNq1ouvzyy+nXrx+nnHIKV1xxRYl1gwYN4oILLtjlXY39SSg33ngj/fv3j8elVEqJrqH8FfiPu7cDugIzgeuBd9z9COCdaBngNOCIaBoBPAxgZg2Bm4HeQC/g5sIkJCKyr2bPns1tt93G008/XaL/LID33nuPl156iYceemiX/a6//nref/99unXrxsiRI8vs8n5nU6dOZfny5QwcOLBEeZ06dfjVr35Fx44dOemkk5g8eTLHHXcchx56KOPGjQP23H1+ZZGwFxvNLBPoDwwDcPftwHYzGwwcF202GpgI/BoYDDzl4dX9j6LaTfNo2/Huvjo67njgVOD5RMUuIhXgjeth2Ze7li/7ouTy9k3h80+HlCw/qMuu+x7UGU67Y4+nzs3N5Uc/+hF33303LVu2LLFu7dq1DBs2jKeffpp69XbtIv+OO+7grrvuKupX6+677y61y/v09PSifQoKCvjFL37BM888s0tX8ps2beKEE07gzjvv5Oyzz+Z3v/sd48ePZ8aMGQwdOpRBgwbt8Xoqi0TWUNoAOcCTZvaZmf3dzGoDzdy9cHSZZUCzaL4FsLDY/ouisrLKSzCzEWY2xcymlDbAjYhIoRtvvJGOHTsWdTtf3KWXXspFF11Enz59ynWsDz74gAsvvBAou8v7hx56iNNPP52srKxd9q9RowannnoqAJ07d2bAgAGkpaXRuXPnCu16Ph4S2fVKKtAduMrdPzazv7KjeQsAd3czi0tnYu7+GPAYhL684nFMEUmgctQkgB1PeA3/d1xOO3HiRP75z3/y6aef7rJu9OjRLFiwgGeeeSYu5yo0adIk3n//fR566CE2btzI9u3bqVOnDnfccQdpaWlFXd7HYrGirudjsVhR1/P72o19RUtkDWURsMjdP46WxxASzPKoKYvoc0W0fjFQvE6bFZWVVS4islfWrFnD8OHDeeqpp6hbt26JdfPmzeM3v/kNzz77bNEQvKWpW7cuGzZsKFouq8v74p599lm+++475s+fz1133cXFF1/MHXeUM6HCbrvPr0wSVkNx92VmttDM2rr7LOBEYEY0DQXuiD7HRruMA640sxcIN+DXuftSM3sT+GOxG/EDgRsSFbeIVF+PPPIIK1as4LLLLitRfsMNN/Duu++yefPmXUZcvP/+++nXr1/RcpcuXUhJSaFr164MGzaszC7v42l33edXJgntvt7MugF/B2oA84DhhFrRS0BLYAHwQ3dfbaHO9wDhhvtmYLi7T4mO82PgN9Fhb3f33Xb8r+7rRSqnfem+Pt5NXhJUue7r3X0aUFpwJ5ayrQNXlHGcJ4An4hudiFQJSiRVht6UFxGRuFBCERGRuFBCEZEKVR2HHa9qEvUdKKGISIVJT09n1apVSipJ5O6sWrWqxJv88aIx5UWkwmRlZbFo0SLUm0Vypaenl/rW/v5SQhGRCpOWlkabNm2SHYYkiJq8REQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLhKaUMxsvpl9aWbTzGxKVNbQzMab2Zzos0FUbmZ2n5nNNbMvzKx7seMMjbafY2ZDExmziIjsm4qooRzv7t3cPTtavh54x92PAN6JlgFOA46IphHAwxASEHAz0BvoBdxcmIRERKTySEaT12BgdDQ/GjirWPlTHnwE1Dez5sApwHh3X+3ua4DxwKkVHbSIiOxeohOKA2+Z2VQzGxGVNXP3pdH8MqBZNN8CWFhs30VRWVnlJZjZCDObYmZTNLyoiEjFS/QQwH3dfbGZNQXGm9nXxVe6u5uZx+NE7v4Y8BhAdnZ2XI4pIiLll9Aairsvjj5XAK8Q7oEsj5qyiD5XRJsvBg4ptntWVFZWuYiIVCIJSyhmVtvM6hbOAwOB6cA4oPBJraHA2Gh+HHBx9LTX0cC6qGnsTWCgmTWIbsYPjMpERKQSSWSTVzPgFTMrPM9z7v4fM/sEeMnMLgEWAD+Mtn8dOB2YC2wGhgO4+2oz+z3wSbTdbe6+OoFxi4jIPjD36ne7ITs726dMmZLsMEREqhQzm1rsFY+9pjflRUQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQREQkLpRQRCqTJ88Ik0gVlOjOIUVkT548A7ZvhMwsWDEjlH3yOPS8JLlxiewl1VBEkmXbRnjzt7Bqbpi+fg3ytsHWdTDliWRHJ7LXVEMRSZbFU2HSAxBLBQxaHgPEIGdmsiMT2SdKKCJJE/Wj16Q9pGfC8H+H5RcugDXzkxaVyL5SQhFJttP/Aq2OTXYUIvtNCUUkGZ48A7asTXYUInGlhCJS0RZMgk05sH1TsiMRiau9TijRqImHuPsXCYhHpHrL3QqjzgDP31GW0SB58YjEUbkSiplNBAZF208FVpjZh+5+bQJjE6l+CvJCMqnXAuo0gyHPQr2Dkx2VSFyU9z2UTHdfD5wDPOXuvYGTEheWSDXX+1IYMUHJRKqV8iaUVDNrThj//bUExiMiIlVUeRPKrcCbwFx3/8TMDgXmJC4sERGpasqbUJa6exd3vxzA3ecB95RnRzNLMbPPzOy1aLmNmX1sZnPN7EUzqxGV14yW50brWxc7xg1R+SwzO2VvLlBERCpGeRPK/eUsK83PgOJ9SfwZGOnuhwNrgMIe8C4B1kTlI6PtMLMOwBCgI3Aq8JCZpZTz3CIiUkF2m1DM7Bgz+wXQxMyuLTbdAuzxl7qZZQFnAH+Plg04ARgTbTIaOCuaHxwtE60/Mdp+MPCCu29z92+BuUCvvbhGERGpAHuqodQA6hAeF65bbFoP/KAcx78XuA4oiJYbAWvdPS9aXgS0iOZbAAsBovXrou2LykvZp4iZjTCzKWY2JScnpxyhiYhIPO32PRR3/y/wXzMb5e4L9ubAZnYmsMLdp5rZcfsRY7m4+2PAYwDZ2dme6POJiEhJ5X1TvqaZPQa0Lr6Pu5+wm336AIPM7HQgHagH/BWob2apUS0kC1gcbb8YOARYZGapQCawqlh5oeL7iIhIJVHem/L/AD4Dfgf8qthUJne/wd2z3L014ab6u+5+ATCBHc1lQ4Gx0fy4aJlo/bvu7lH5kOgpsDbAEcDkcsYtIiIVpLw1lDx3fzhO5/w18IKZ/YGQpB6Pyh8HnjazucBqQhLC3b8ys5eAGUAecIV78Y6QRKqhVXPh0QEw+EFo0hZS0pIdkcgelTehvGpmlwOvANsKC919dXl2dveJwMRofh6lPKXl7luBc8vY/3bg9nLGKlK1paZD3lZYOg0e6QN9fgYn35bsqET2qLwJpbApqngzlwOHxjccEeGUP8Ly6WF+w1KNmyJVRrkSiru3SXQgIhKp2wyu+DjM39U2ubGI7IXydl9/cWnl7v5UfMMRqaZWzIQXLwYv2PO2IlVUeZu8ehabTwdOBD4FlFBEyuPlEbBq9o7l2o2TF4tIgpS3yeuq4stmVh94ISERiVRHedsgvQE0PBTOfRLqt0x2RCJxt69jym8CdF9FZG8cOgB+OHrP24lUUeW9h/Iq4akuCJ1CtgdeSlRQIiJS9ZS3hnJXsfk8YIG7L0pAPCIiUkWVq+uVqJPIrwk9DTcAticyKBERqXrKlVDM7IeE/rPOJYwr/7GZlaf7ehEROUCUt8nrt0BPd18BYGZNgLfZMVCWiIgc4Mrb23CsMJlEVu3FviIicgAobw3lP2b2JvB8tHwe8HpiQhIRkapotwnFzA4Hmrn7r8zsHKBvtGoS8GyigxMRkapjT81W9xLGj8fdX3b3a939WkI39vcmOjgRAWa9AY+fGt62F6nE9pRQmrn7lzsXRmWtExKRiOwQS4VNK2DhJLi9OXzzbrIjEinTnhJK/d2sy4hnICJSisH3Q/3WUC8LPB/WLkx2RCJl2lNCmWJm/7dzoZn9BJiamJBEpMhhJ8A1n8MlbyU7EpE92tNTXtcAr5jZBexIINlADeDsRAYmIiJVy24TirsvB441s+OBTlHxv91dDbki5fXkGbDuO2jaPtmRiCRUecdDmQBM2JsDm1k68B5QMzrPGHe/2czaEMZSaUSo9Vzk7tvNrCZhwK4ehBcnz3P3+dGxbgAuAfKBq939zb2JRSQp5n8AU0dBzizIU/d3Uv0l8m33bcAJ7t4V6AacamZHA38GRrr74cAaQqIg+lwTlY+MtsPMOgBDgI7AqcBDZpaSwLhF4uOzZ2H6y7B9A6SkQZv+yY5IJKESllA82BgtpkWTAyewow+w0cBZ0fzgaJlo/YlmZlH5C+6+zd2/BeYCvRIVt0hc1WsBLbLD1POSPW8vUoXt64iN5RLVJKYChwMPAt8Aa909L9pkEdAimm8BLARw9zwzW0doFmsBfFTssMX3KX6uEcAIgJYtNbyqVCLD/53sCEQqREI7eHT3fHfvBmQRahXtEniux9w9292zmzRpkqjTiIhIGSqkx2B3X0u4qX8MUN/MCmtGWcDiaH4xcAhAtD6TcHO+qLyUfUREpJJIWEIxsyZmVj+azwBOBmYSEkvh4FxDgbHR/LhomWj9u+7uUfkQM6sZPSF2BGGwLxERqUQSeQ+lOTA6uo8SA15y99fMbAbwgpn9AfgMeDza/nHgaTObC6wmPNmFu39lZi8BMwjj2V/h7vkJjFtERPZBwhKKu38BHFVK+TxKeUrL3bcShhgu7Vi3A7fHO0YREYkfjbooIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiUpX853q4ozX88yeQuyXZ0YiUoIQikghPngFz347f8eo0hbang8Ugbwt8+Q9YPS9+xxeJAyUUkXha+gX8/SRY+hlsWR2/46akwfnPQ/NuYUhgkUpICUUknhZNhkWfQCwNataF3iPie/zh/4YTfhPfY4rESUJ7GxY5YF35SWimEjmAqIYiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoYQiIiJxoa5XRPaXe+hZ+M3fwpY1FXfef10BNWqH/r1EKoGE1VDM7BAzm2BmM8zsKzP7WVTe0MzGm9mc6LNBVG5mdp+ZzTWzL8yse7FjDY22n2NmQxMVs8g+yfkanv0BrJwFm1ZAajqk1UrgCS18LPscvpsEo85M4LlEyi+RTV55wC/cvQNwNHCFmXUArgfecfcjgHeiZYDTgCOiaQTwMIQEBNwM9AZ6ATcXJiGRSiFva/hseBgc3B1+/hXUrJO487XpD32vhToHQVpG6N1YpBJIWEJx96Xu/mk0vwGYCbQABgOjo81GA2dF84OBpzz4CKhvZs2BU4Dx7r7a3dcA44FTExW3yD475XYYMQFqN07seTLqw0k3wy9mQs9LEnsukb1QITflzaw1cBTwMdDM3ZdGq5YBzaL5FsDCYrstisrKKt/5HCPMbIqZTcnJyYlr/CIismcJTyhmVgf4J3CNu68vvs7dHfB4nMfdH3P3bHfPbtKkSTwOKSIieyGhCcXM0gjJ5Fl3fzkqXh41ZRF9rojKFwOHFNs9Kyorqzx+1i2GR4+DRwfA0s/hkQFh2rQyrqcREanOEvmUlwGPAzPd/Z5iq8YBhU9qDQXGFiu/OHra62hgXdQ09iYw0MwaRDfjB0Zl8bF+KYzsGMYAXzoNHu0Py6aF6fGTy97PHea+Aw/0DtN3H4fyZdNhZKcwTfgjLP8qlC/5DB48OkzzPwj7i4hUI4l8D6UPcBHwpZlNi8p+A9wBvGRmlwALgB9G614HTgfmApuB4QDuvtrMfg8UPspym7uvjluU29YDDnUPhvR6cMKN8O4fYMOyUHO568iwXY06kF4fOv8APnoYcjfB5lU7jvPEQMhoCFuKhfbfP8Mnf4f6rULNx/ND+agz4NIPoUlbWDINXrs2lA8dC7Uaxu3SREQqUsISirt/QNED87s4sZTtHbiijGM9ATwRv+hKccofoNP3w3z7M+HzF+DtW2HzSsjfDiwP65ZMLblf43ZQkBceHd2UAyk1ofvFcMZd8GDv8I5CYeKplwUpNWDNPMjdApP/Bm/esONYdx4OTTvA6X+B5t2gRiLfZRARia8D6035ac/Df6LXXlr0gH7XQu0ybuB3HRImgJVz4PkfhfkB14VEEEuFYa9BbDethmc/CmOixzr7/zIkqu8+hme/D/+6DLZtCOuadoRVsyE/F5Z/CU+eBn2ugZNv3f9rFhGpIAdWQvn8+fBLPK0WzJsATdpBj3K8eN/4CLiq2MtjXc4t3/kO7gZX71SjScsIn6vmhM/UDLj0A9i6FkZ9L5St+w4+exq+mQCDH4BmnXafuEREKoHq/1tqwf9C1xR/OgQWfABZPeE3iyCtNkx5HP6+mxvvidDyGLh4HDRpH6afjA/JolZDuPzDMKXXD81kyz6HR/vBrNfLPl5+LvztpDCtX6Kb/SKSNNWzhuIOX46BCX+CLatCh30164UkUli76HRO2KYgD7J6wUFdKya2WAwOHQBXfFT2Nhf/C168MCSLVXPgrd/Bu7fD+kWwfSNg4RgXvQJjhsPiqPZ0T3s47S/Q+6cVcikHvA/vg//dDwW5yY3DC8Ij7wDDXw8dRookgXk1/Is2u0WaT/m/Yje067WAa2ckL6B9tX5JSBLFZTSEvG3hibH6rWDtAjCDui3Czf6+P4cTb0pOvAeaWxsAFu6nNWgF578AjQ6r2Bgm/An+e8eO5bZnwPnPVWwMUm2Y2VR3z97X/atnk5c71GoUnpRq0ROurKKd59U7GM57JnQ62PLI6nYAABVGSURBVPAwOHc0XDcP+v08PFW2clb4PPK0cK/GqufXWan1uxZuXBF+xio6mQAccwUMeS7cD0yrFe7FiSRJ9WzyOrgbXDcl2VHER/vvham4vr+ADmfDy/8Xls9+ZNf98raF5rxYKsTSdFO/ukqvB+3OCNOoM0Pzl0iSVM+EUt3FYtD48NCzbWnWL4H7ukPelp1WGLQ9Dc5/PuEhisiBR3+2VicWg/fvCV3J5G2B2k3DfZbjfgOZh4RHlpd+kewoRaSaUg2lOjnjHpj4pzDf7QI49qrQJAJw3K/DkLHzJiYtPBGp3pRQqpPuF4VpdzYugz8eHOZP+wscdWHi4xKRA4KavA4knb8Ph58c+hvL3Ro6pbzzCBh3Naz6JtnRiUgVpxrKgeSwE8IE8Pp1MHVU6B3509FhqpkJrY6Fsx5Sr8cistdUQzlQnf6X8P7EjSuh24WhJ4GCXJj9Bqyogi+BSrDwY/hza3huyI4xekQqiBLKgc4MznoQblgIP3oxlP3rCnjomB2Dg8kO056DP7cJU2V756PdGZCaHrrsmf0f+Pq1ZEckBxglFNmhZvRE2Nr5oZby2PGQtz0MOHZLZjTVh/89kNQwk2rKk7B1XXhxtHlXaNUn2RHtcPRl8JvFYUpNT3Y0cgDSPRTZoXlXuPxjeHkEbFwBG5fCPe1Cl/+xNKh7UOhoM+frZEeaXG36hw48RaQE1VBkBzNo2g4ufQ9+8PcweuT2zWApcOyV8PPpoRYz/Z9wd3sYfxNsWJ7sqEWkklANRUrXui9cPmnX8qwe8PXrsHE5fPjX0GlleQYpq8pyZsFnz4REuimncjVziVQiSiiyd857Jnyuj5rD1i2CZdMhPRPqH5Lc2BJl6mj46MEdvTk3r6Cxc0SqmIQ1eZnZE2a2wsymFytraGbjzWxO9NkgKjczu8/M5prZF2bWvdg+Q6Pt55hZNf9TuApJrRk+3/sLPNIH7u0Ea79LbkyJ4gXhHZ2b14Tp5FuTHZFIpZTIeyijgFN3KrseeMfdjwDeiZYBTgOOiKYRwMMQEhBwM9Ab6AXcXJiEJMlqNYQfvxXG4aiXFcruOwpuawRPn5Pc2EQkKRKWUNz9PWD1TsWDgdHR/GjgrGLlT3nwEVDfzJoDpwDj3X21u68BxrNrkpJkadkbrvgYLn0fjr4cajUOY7B88w78vkmY3v1DsqMUkQpS0U95NXP3pdH8MqBZNN8CWFhsu0VRWVnluzCzEWY2xcym5OTkxDdq2b1aDeHUP8EvZ8Flk6D7UMhoAHh4Z+WB3jDzVdhShUYT3L453Ii/r3t4mbEqmvQg/L4pPNAzPFggkmBJe2zYw2D2cRvQ3t0fc/dsd89u0qRJvA4re6tZBxh0H/xydhiaOG8LrPwaXrwQJlWhFyLnvAVjr4DV38C2dVCvebIj2jsn3hiGwa5RG1bOhhUzkx2RHAAq+imv5WbW3N2XRk1aK6LyxUDxR4SyorLFwHE7lU+sgDglHs5+FI6LbpP97QTI3XkEyUosPzd8/vgtaNgmPMVWlRxzRZiWz4CHj0l2NHKAqOgayjig8EmtocDYYuUXR097HQ2si5rG3gQGmlmD6Gb8wKhMqoIataBZxzBZCkz+G9zRGl64IDxuXBXUagR1mu54qk1EypSwGoqZPU+oXTQ2s0WEp7XuAF4ys0uABcAPo81fB04H5gKbgeEA7r7azH4PfBJtd5u773yjX6qCrkPgixchf1votLBJu/DyZO3GcFDnZEcXbF0PDx8L6xeDR62xZsmNSaQKMfe43caoNLKzs33KlCnJDkNKs25RGPO+iMGv5obEkiwF+aE5bs388E5NRgNIqwW9RsCxV0OsCvdQVNjkde5o6HjWnreXA5qZTXX37H3dX2/KS8XKzAodUG5ZA3PehA9GwqP9IZYKZ9wDR5xU8TH97QRYOm3H8ul3QecfVHwcIlVcFf7TS6qspu2g1TFh9EhLCU1MaxfAuCuTE8/aBVCzLtRvDafdCUeekpw4EqGwyW7Mj8NLpw8dEx6JFkkA1VAkedr0h98uC12bvPErmPVGKN+0ElbOCfOxFDj4KEhJS2wsXc+H0+9M7DmSodER0P9X8OkzkLs5jHOzeVV4YEIkzpRQJLlSa4TPWCps2wj3dAg1luJO/TMcfWl8z/vt+/Dy/4X56vwXe0oqnPC7MH32THi35r5uYV27M+CHTyU3PqlW1OQllUOzjqF5ZkPUkUKrvnDhP8P8+JvgT4fAM9/fv/FX1iwIXcHc0xFGnxnOtSkn/LV+IPQgfPjJcPQVUKsJpNSA7z5KdkRSzegpL6m83OHf18KX/4TtG8HzQ3larfDI8fkvhOaxBR+G8vfugtVzw3ztpjDgOsjqFe6PZLaAiX+GiX8EDPDQRcyg+5JxZcn36s9CE+MvZyc7EqlE9JSXVF9mcObIMG3bCG/8Gma/ER7xnfNWGJr4ndtg2jMl98toCOsWwrirSjlmLHRBL+Hf7+2oK/52Z4bB00T2gxKKVA0168BZD4b5SQ/BmzfAvZ1DcslsCRf8I6xrdFi4gb9gUmjSyt0My74MtRQIN6kFmnYMPUN/MBLw0N/XkGeTHZVUcUooUvW07hOeEFsyDdIyoMfF4VHk4lqp/6rd6j0iTACP9A1P2onsJyUUqXqad4WhryY7iurDYjDrdbilPqSmwzmPQYdByY5KqiAlFJED3Yk3w/z3Qxc0/7sPvv536AyzRm1oeWzV7npGKpSe8hKRoCA/PJ6du2lH2cXjoEWPkGAS/XKpJJ2e8hKR+IilwJWfwMZlsOJrGHs5PBU1fdVpBtdM3/EiqkgplFBEZIfMFmE6qCvkbw+PFi/4EGb/B/K2KqHIbimhiMiuUlIhe3iYj6WGhHJHNKhq637wgyfCfO0mGjNGiiihiMjudTon1FQKcuG9O8MN/LuKvc/T7szwPlCrY6DhYWG45MNOUKI5AOmmvIiU3+pv4Zt3wvwnj4dBybwgNIcV1++X0LQ9NGgNWft8j1cq2P7elFdCEZH9t3EFbF4NiyaX7PImNQMueSskni2rAQt9sqXXhzpNw3svBx8VHgiQpNNTXiKSfHWahqlpu9DclbsFpo6CSQ/Ao/12v2+n70Pb08NDAE2iHg8y6kOtRmG+Zj01n1URqqGISGJs27jjhUkIjx7XOzg8lrx9Uxg6YMyPy3estqeH3qe9APBw7B5Dw0MBNepA4+ieTs164YEC2SeqoYhUIms3b2fM1EXkFzgpMSNmFn1CLGakmPHc5O/IzAgvCa7ZvJ3ebRrx1ZJ1zFy6gRqpMdydbbkFnN29Be5wWueD6NQik9SYUatGFfovW7MOtD1t1/LMFjvmW/eDresgPxdWfwOxNNiwJCSMWAp8/jwUFMDahaGWYgbrFsPmlfDd/0o/b6s+sOgTqNU4JKCCPNi6NkpsDs06hdrUhuVweDQMdSwlfEJ4qi0zK+xXkBf2S68HXc4LSW3t/PC5cUV4WAHCdnWahU5ILRZ6F7CUaD76zGgYHrvesibciyo8V7OO1abJTzWUBMrLL2DjtjwAUmJG3XS9aZwM+QVO4c95zIxYrOzmkxUbtrJ+Sy4A9dLTaFovna25+UxbuJb8Ag+Th+PNWLKe9LQUUmNGSkqM1JgxdcEaxkxdVK64mmems3RdyZvZP8zO4q0Zy8nPd1JTjDWbc0usr5kao+PB9cjZuI22zeqRUSOF9NQY5/duCUCbRrWpm57K6k3bydm4DYD0tBQObVwbq07NRitmwtb1sG19mE+tCfMmhuRksfALf/sGaN4t/NKOpYZ7OIs/hYwGYSjk3C1hPy8IzW17EksLT7rtjzb94dv3SpZltoTuF4U43GHlrDBMcyxtR42ssHbWoHVoUpz56o4etAvX4+DRMoQkfUgvmDM+JHcIibrxEVD3oHC8w04ICbMgHwrysFZH66b8zjJbtvNGF9xNeloKKWYUuFPgUOCOA+5Obn647s4tMnGcgoJofbTdnBUbqV0jhbTUGIX/Dc2s2Hzh2Qx3Z9Wm7XTJyuTLxeuoXSOV1BRj7U6/DGIGXQ+pH2KJfjmF2Jy8Amf9llz6HdGEd2Yup0ndmqTGYphR9Jfuyo3bOKRhLWqkhL6VPPrBKfwKiz53Kg9lMHPpeurWTCUlxSgoCL9o+x3RmPEzl9M8M4NY9AdgLLrOZeu3MuDIJkz+djUZNVJJiUGKGRb91b09r4C66alkpKWwatN2atVICddTAPnRL9yGtWsQi0FevlM3PZVOLTIZP2M5zeqlkxb9Ek6JpjWbt9P38MZ8PG81ddNTKXAn38P3lV/grN2cS5vGtVmydguL126J/pqHjdvy6Nm6Qfg/F32/hde8aM0WcjZsK/E9HN60DgUF4d88v8BZvHZLmT9LrRvVYv6qvRsiOD0txvvXnUDNtBgFBeFnr/C7zi9wYmYclJkOhJ+DQuEP8JK/9P/3zUpmLFlPfoHzymeLaVynJsvXb2Xtllxq10jZ69gADmmYQVosxurN2zmlw0G8+sUSatVIwSx8py3qZ9CrTcOQfKOaVeHPRazw5yOKc3teAVkNMvhqyTrqZaSRn7/j33XRms3USU8jNUrgm7fn0e2QBqXGtHz9Vr5duYnaNVOibfM5sV3T6Pss+Z3m5TvzVm6kYe3wkuXqTdvpcHBm0Q984b/ojv8XZZRTcj1AWt5G0retpsBSKbAYBZZCzdx1dFvwJBbtkZK/jXlNTmDVpu2kp2ewvkYz6mxdQsyMBhkpfLdyI3VrGIazZdt2WjZM5+CVk8jYllP0i2NjrSwWNjuJAZ/+rMzvaVVmR9xSAMPNaLJm2i7bbKnZJLoWKzq2Y9TeuqzUbTO25ZR5PgC7db0Sys7MbAMwK9lxJFBjYGWyg0ggXV/VVp2vrzpfG0Bbd6+7581KV4UaZPfKrP3JspWdmU3R9VVdur6qqzpfG4Tr25/91S+1iIjEhRKKiIjERXVNKI8lO4AE0/VVbbq+qqs6Xxvs5/VVy5vyIiJS8aprDUVERCqYEoqIiMRFtUgoZlbfzMaY2ddmNtPMjjGzhmY23szmRJ+lv1FVyZlZWzObVmxab2bXVKPr+7mZfWVm083seTNLN7M2Zvaxmc01sxfNrMoOE2hmP4uu7SszuyYqq7LfnZk9YWYrzGx6sbJSr8eC+6Lv8Qsz6568yMunjOs7N/r+Cswse6ftb4iub5aZnVLxEe+dMq7vzuh35xdm9oqZ1S+2bq+ur1okFOCvwH/cvR3QFZgJXA+84+5HAO9Ey1WOu89y927u3g3oAWwGXqEaXJ+ZtQCuBrLdvROQAgwB/gyMdPfDgTXAJcmLct+ZWSfg/4BehJ/LM83scKr2dzcKOHWnsrKu5zTgiGgaATxcQTHuj1Hsen3TgXOAEn2mmFkHws9rx2ifh8yssnfKNYpdr2880MnduwCzgRtg366vyicUM8sE+gOPA7j7dndfCwwGRkebjQbOSk6EcXUi8I27L6D6XF8qkGFmqUAtYClwAjAmWl+Vr6098LG7b3b3POC/hF9MVfa7c/f3gNU7FZd1PYOBpzz4CKhvZs0rJtJ9U9r1uftMdy+t543BwAvuvs3dvwXmEv54qLTKuL63op9PgI+ArGh+r6+vyicUoA2QAzxpZp+Z2d/NrDbQzN2XRtssA5olLcL4GQI8H81X+etz98XAXcB3hESyDpgKrC32A74IaFH6ESq96UA/M2tkZrWA04FDqAbf3U7Kup4WwMJi21Xl77I01fH6fgy8Ec3v9fVVh4SSCnQHHnb3o4BN7NSE4F7YFWfVFd1HGAT8Y+d1VfX6orb2wYQ/Cg4GarNrdbzKcveZhOa7t4D/ANOA/J22qZLfXVmq2/UcSMzst0Ae8Oy+HqM6JJRFwCJ3/zhaHkNIMMsLq9fR54okxRcvpwGfuvvyaLk6XN9JwLfunuPuucDLQB9C00hhP3NZwOJkBbi/3P1xd+/h7v0J94NmUz2+u+LKup7FhBpZoSr9XZai2lyfmQ0DzgQu8B0vJ+719VX5hOLuy4CFZtY2KjoRmAGMA4ZGZUOBsUkIL57OZ0dzF1SP6/sOONrMalnoD73wu5sA/CDapqpeGwBm1jT6bEm4f/Ic1eO7K66s6xkHXBw97XU0sK5Y01h1MA4YYmY1zawN4eGDyUmOaa+Z2anAdcAgdy8+HsLeX59HgwVV5QnoBkwBvgD+BTQAGhGeOJkDvA00THac+3F9tYFVQGaxsmpxfcCtwNeE+w1PAzWBQ6Mf3LmEJr6ayY5zP67vfUKS/Bw4sap/d4Q/apYCuYTWgUvKuh7AgAeBb4AvCU/zJf0a9uH6zo7mtwHLgTeLbf/b6PpmAaclO/59vL65hHsl06LpkX29PnW9IiIicVHlm7xERKRyUEIREZG4UEIREZG4UEIREZG4UEIREZG4UEIREZG4UEIRqQBm1trMtpjZtGJlzczsOTObZ2ZTzWySmZ29m2NM2LkL8Wgog4fNLCMa3mC7mTVO5LWIlEUJRaTifONhGAKingH+Bbzn7oe6ew9C559Zu9n/+Wib4oYAz7v7lujYSxIQt0i5KKGI7MTM/mFmD5jZB2a2wMz6mtnTZjbbzB6P02lOALa7+yOFBe6+wN3vj2K40MwmR7WOR6NxKMYAZxQOOGZmrQmdar4fp5hE9osSisiuOgPz3L0v8ChhrJ3rgA6EX+g143COjsCnpa0ws/bAeUCfqNaRT+i0bzWhS5rTok2HAC+5uruQSiJ1z5uIHDjMLB2oD9wbFTnwuEedGppZPrA9Aed9EOgbHXs0YXTOT0LLGBns6MG3sNlrbPRZJUezlOpJNRSRkjoShgkoiJa7Ah8DmFkW4R5FezO7Piq738zqmlmHncv2cJ6vCMMsAODuVxB6W25C6FRxtEdDP7t7W3e/Jdp0LHBiND57LXefGodrFokLJRSRkjoTegYu1IXQizWE5PIF0JPQKyuEHqA3lFG2O+8C6WZ2WbGyWtHnO8APinV939DMWgG4+0ZC9/5PUHI4A5GkU0IRKakzUWKImr8y3H1NtK4wufQEZkRDTRcqraxM0X2Ps4ABZvatmU0mNHX92t1nAL8D3jKzL4DxQPGx2J8nJDclFKlU1H29yF4ys38TxpJYD3R291NLK9tpn9bAa+7eKcGxzSeMO7IykecRKY1uyovsBTNLA1a5+093V1aKfCDTzKYVvosS57gygElAGlCwh81FEkI1FBERiQvdQxERkbhQQhERkbhQQhERkbhQQhERkbhQQhERkbhQQhERkbhQQhERkbhQQhERkbhQQhERkbj4f637S7APS4z2AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = hist.plot1d(output['mass_near'], overlay='dataset')\n",
    "ax.set_xlim(60,120)\n",
    "ax.set_ylim(0.1, 7500)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.1, 8000)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAELCAYAAAD+9XA2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2deZhU1Zm436+q96a72XcQUFFAEKTF4IYat2g0YxJH47ihiWiiiTGZ32gWRcdMzKjRxCVq3DUqxsRxj9EoEQ2RxYBsioggzdK0DfS+1HJ+f5y6tXVVdXV3VVd19fc+Tz1V99xzz/1u36771TnfJsYYFEVRFKWnuDItgKIoipIbqEJRFEVRUoIqFEVRFCUlqEJRFEVRUoIqFEVRFCUlqEJRFEVRUkJepgVIJSJyBnBGWVnZdyZPnpxpcRRFUfoMK1eu/MIYM6wnY0guxqFUVlaaFStWZFoMRVGUPoOIrDTGVPZkDF3yUhRFUVJCTikUETlDRB6oq6vLtCiKoij9jpxSKMaYl4wxl1VUVGRaFEVRlH5HThnlFUXJbjweD1VVVbS2tmZalH5LUVERY8eOJT8/P+Vj55RCcby8DjjggEyLoihKDKqqqigrK2PChAmISKbF6XcYY6itraWqqoqJEyemfHxd8lIUpddobW1lyJAhqkwyhIgwZMiQtM0Qc0qhKIqS/agyySzp/PvnlEJRLy9FyQ3cbjczZ87kkEMO4eyzz6a5uTlmvyOPPLKXJYvPcccdR2fxb3feeWfca+kuixcv5h//+EdKx+wuOaVQdMlLUXKD4uJiVq1axdq1aykoKOC+++6L2O/1egHS8iB1xk4HqlAURVEyyDHHHMOmTZtYvHgxxxxzDGeeeSZTp04FYMCAAYB9qM6bN4+vfe1rTJo0iWuvvZY//OEPzJkzh+nTp/Ppp58C8NJLL3HEEUcwa9YsTjzxRKqrqwFYuHAhF1xwAUcddRQXXHABxx57LKtWrQrKcPTRR7N69eoIuVpaWjj33HOZMmUKZ511Fi0tLcF9V1xxBZWVlUybNo0bbrgBgN/+9rfs2LGD448/nuOPPz5uP4Brr72WqVOnMmPGDH784x8DUFNTwze+8Q0OP/xwDj/8cN577z22bNnCfffdxx133MHMmTNZsmRJSv/2XcYYk3Ov2bNnG0VRso/169cn1a+0tNQYY4zH4zFnnnmmuffee83bb79tSkpKzObNmzv0e/vtt01FRYXZsWOHaW1tNaNHjzbXX3+9McaYO++80/zgBz8wxhizZ88e4/f7jTHG/P73vzfXXHONMcaYG264wRx22GGmubnZGGPMo48+Gjzm448/NrGeKbfffruZP3++McaY1atXG7fbbZYvX26MMaa2ttYYY4zX6zXz5s0zq1evNsYYs99++5mamprgGLH6ffHFF2by5MlBOffu3WuMMeZb3/qWWbJkiTHGmK1bt5qDDz44KPutt96a1N/VIdZ9AFaYHj57c2qGkvU2lAdPgt8dBb70TakVJRdoaWlh5syZVFZWMn78eC699FIA5syZE9fd9fDDD2fUqFEUFhay//77c/LJJwMwffp0tmzZAli35VNOOYXp06dz6623sm7duuDxZ555JsXFxQCcffbZvPzyy3g8Hh5++GEuvvjiDud75513OP/88wGYMWMGM2bMCO579tlnOeyww5g1axbr1q1j/fr1MWWO1a+iooKioiIuvfRS/vznP1NSUgLAm2++yZVXXsnMmTM588wzqa+vp7GxsQt/1fSTU3EoxpiXgJcqKyu/k2lZYlK1zL6/fDUc/UMwfhhyAKjXi6JE4NhQoiktLY17TGFhYfCzy+UKbrtcrqBd5KqrruKaa67hzDPPZPHixSxcuDDm2CUlJZx00km88MILPPvss6xcuTJp2T/77DNuu+02li9fzqBBg7j44otjuunG65eXl8eyZcv429/+xnPPPcfdd9/NW2+9hd/v55///CdFRUVJy9Lb5NQMpc/wryfgrsPg7kpY9+dMS6Mo/Ya6ujrGjBkDwGOPPZaw77e//W2+//3vc/jhhzNo0KAO+4899lieeuopANauXcuHH34IQH19PaWlpVRUVFBdXc1rr70WPKasrIyGhoaE/RobG6mrq+O0007jjjvuCNpuTj75ZO66667gWI7CDR8z06hC6W1mnR+5/dwlcPccqF4Xu7+iKClj4cKFnH322cyePZuhQ4cm7Dt79mzKy8uZP39+zP1XXHEFjY2NTJkyheuvv57Zs2cDcOihhzJr1iwOPvhgzjvvPI466qjgMZdddhmnnnoqxx9/fNx+DQ0NfPWrX2XGjBkcffTR/PrXvwasUX/FihXMmDGDqVOnBj3fzjjjDJ5//vmsMMprPZTeZGEFzLsWjrsW/F5YfAssuS20/ztvw5jDMiefoqSZDRs2MGXKlEyLkRQ7duzguOOO46OPPsLlyq3f3rHug9ZD6auIgDsfvvxzuPSNUPvvj4f/PUCN9oqSYR5//HGOOOIIfvGLX+ScMkknOfWXynovr1iMmwPX7wEJ3IrmGvjvIeD3ZVYuRenHXHjhhWzbto2zzz4706L0KXJKoZi+GinvcsO122DAyFDbTYPhjung92dOLkVRlC6QUwqlT1M4AH78sVUsDnWfw/9OypxMiqIoXUAVSrZRVA7/uRnOut9ut+6FP307szIpiqIkgSqUbKR0CBx6Lly22G6v+aN1L1aUfsg59y/lnPuXZloMJQlUoWQzo2fB/EBQ1No/wZ2HZlYeRckBnNT406ZN49BDD+X222/H3wNbpZOgMpqLL76Y5557LuGxF198MSUlJRGBiVdffTUiwhdffNFtmTKFKpRsZ78j4fuBFBT7tsCvJkIOxg4pSm/hpHVZt24db7zxBq+99ho33nhjxuQ54IADeOGFFwDw+/289dZbwWj+voYqlL7A4Ilwrk3xQMse+OW4zMqjKDnC8OHDeeCBB7j77rsxxtDa2sr8+fOZPn06s2bN4u233wbg0Ucf5corrwwe99WvfpXFixcHt3/4wx8ybdo0vvzlL1NTU9PhPCtXrmTevHnMnj2bU045hZ07dwb3nXvuuSxatAiwafiPOuoo8vJsmsUtW7ZwyCGHBPvedtttwfxjq1at4ktf+hIzZszgrLPOYu/evYAt9PVf//VfzJkzh8mTJ/dq9HxOJYcUkTOAMw444IBMi5J6Dj4drl4Ldx4C7Q3w0g/gjN9kWipF6TY3vrSO9TvqO7Sv3xnZ1txmA32nL3w9on3qqPIOx04dXc4NZ0zrkhyTJk3C5/Oxe/dunnzySUSENWvW8NFHH3HyySezcePGhMc3NTVRWVnJHXfcwU033cSNN97I3XffHdzv8Xi46qqreOGFFxg2bBiLFi3ipz/9KQ8//DAAkydP5sUXX2Tv3r08/fTTnH/++RH5v+Jx4YUXctdddzFv3jyuv/56brzxRu68807AFglbtmwZr776KjfeeCNvvvlml/4m3SWnZih9Ng4lWQaOg+++bz+vfBRuHqkBkIqSQt59991gSvqDDz6Y/fbbr1OF4nK5OOeccwA4//zzeffddyP2f/zxx6xdu5aTTjqJmTNncvPNN1NVVRXR5+tf/zrPPPMM77//Psccc0ynctbV1bFv3z7mzZsHwEUXXcQ777wTMR7YfGRO6v7eIKdmKP2C4QfDFf+A++eBtwXuPhy+/0GmpVKULpPsTMLx8Fq0YG5a5Ni8eTNut5vhw4fH7ZOXlxdhuI+Vjt5BospRGGOYNm0aS5fG91Q755xzmD17NhdddFFEqpeunDccJ3W/2+1Oa0njaHJqhpLV/GamfXenQIePmAZXBZJf7vkU9m1L3F9RlJjU1NRw+eWXc+WVVyIiHHPMMfzhD38AYOPGjXz++eccdNBBTJgwgVWrVuH3+9m2bRvLli0LjuH3+4PeXE899RRHH310xDkOOuggampqggrF4/FEFPYC2G+//fjFL37Bd7/73Yj2ESNGsHv3bmpra2lra+Pll18GoKKigkGDBgXtI0888URwtpJJdIbSWxx0GvzzHjjs4tSMN2gC5BWCtw3uOwau3ZKacRUlx3GqQXo8HvLy8rjgggu45pprAPjud7/LFVdcwfTp08nLy+PRRx+lsLCQo446iokTJzJ16lSmTJnCYYeFsoKXlpaybNkybr75ZoYPHx40sDsUFBTw3HPP8f3vf5+6ujq8Xi9XX30106ZFztAWLFjQQdb8/Hyuv/565syZw5gxYzj44IOD+x577DEuv/xympubmTRpEo888kgq/0zdQtPX9xZ/+Ql88Dj8pKrzvl3htslQPAj+7Xea+l7JerqTvj7dS179kXSlr9cZSl/H74Oaj2zq+5/ttrMWRckhVJH0HdSG0tf58SdQPNh+/uXYzMqiKEq/RhVKX8flgms22M++dmjek1l5FEXpt2S9QhGRKSJyn4g8JyJXZFqerCS/CAbuZz//+TuZlUVRlH5LRhSKiDwsIrtFZG1U+6ki8rGIbBKRawGMMRuMMZcD/w4clQl5+wQ/WG3fN70J7c2ZlUVRlH5JpmYojwKnhjeIiBu4B/gKMBX4lohMDew7E3gFeLV3xexDiMBx19nPt07SSo9K7vDI6falZD0ZUSjGmHeA6MX+OcAmY8xmY0w78AzwtUD/F40xXwH+o3cl7WPMDPx5PC2w8S+ZlUVRspDnn3+emTNnRrxcLhevvfYa99xzT0T7IYccgoiwYcOGiDG2bNnCU0891a3z19fXM3bs2IhEk7lENrkNjwHCQ76rgCNE5Djg60AhCWYoInIZcBnA+PHj0ydlNjNwHFz8Kjx6Gjy/AK7TCHpFCeess87irLPOCm4/8MAD/OEPf+CUU07B5XLxve99L7jvJz/5CTNnzuwQr+EolPPOO6/L5//5z3/Oscce2/0LyHKy3ihvjFlsjPm+MWaBMeaeBP0eMMZUGmMqhw0b1psiZhcTjoIhB9r0LIqixGXjxo3cdNNNPPHEExH5swDeeecdnn32We69994Ox1177bUsWbKEmTNncscdd8RNeR/NypUrqa6u5uSTT45oHzBgAP/5n//JtGnTOPHEE1m2bBnHHXcckyZN4sUXXwQ6T5+fLWTTDGU7EF7oY2ygLWlyOn19V2jYCbWfwK41MHJ6pqVRlNi8dq39H41m14eR2+1N9j26DtDIGR2PHTkdvnJLp6f2eDycd9553H777R1WNPbt28fFF1/ME088QXl5xxT5t9xyC7fddlswr9btt98eM+V9UVFR8Bi/38+PfvQjnnzyyQ6p5JuamjjhhBO49dZbOeuss/jZz37GG2+8wfr167nooos488wzO72ebCGbZijLgQNFZKKIFADnAi92ZYCcT1+fLCcutO+bF2dQCEXJXn7+858zbdq0YNr5cC6//HIuuOACjjoqOafSZFLe33vvvZx22mmMHdsx+LigoIBTT7U+StOnT2fevHnk5+czffr0Xk09nwoyMkMRkaeB44ChIlIF3GCMeUhErgReB9zAw8aYdQmGiTWuzlAA9j/Bvv/1Z1B5CRSUZlYeRYlFEjMJIOThNf+VlJx28eLF/OlPf+KDDzqWfXjsscfYunUrTz75ZErO5bB06VKWLFnCvffeS2NjI+3t7QwYMIBbbrmF/Pz8YMp7l8sVTD3vcrmCqee7m8a+t8mIQjHGfCtO+6v0wDXYGPMS8FJlZWX/ju4bsj8UDYTWfbDqKZjTv/8ciuKwd+9e5s+fz1NPPUVZWVnEvs2bN/OTn/yEJUuWBEvwxqKsrIyGhobgtpPy/oQTTohIeR+OkxIfrD1kxYoV3HJLkgoVmDBhAvfeey9+v5/t27dHpM/PJrLJhqKkkqtWwq37w6s/VoWiKAHuu+8+du/ezRVXRCbduO6663jrrbdobm4OVjt0uOuuuyKqKM6YMQO3282hhx7KxRdfHDflfSpJlD4/m8ip9PVhS17f+eSTTzItTiTpSl+fiPuOtkbPH66HijG9d15FiUN30teneslLSV/6+mwyyvcYNcpH0d5i3/92U2blUJSeMP8VVSZ9hJxSKEoU319p3z98JuR6qSiKkiZySqGIyBki8kBdXV2mRckeDv+2fb/niMzKoShKzpNTCkWXvGIwM5Aeom4b7NmcWVkUBcglu21fJJ1//5xSKEoMxsyGU39lP+/7PLOyKP2eoqIiamtrValkCGMMtbW1EVH8qSSn3IY1sDEOww+278+c37teZooSxdixY6mqqqKmpibTovRbioqKYkbsp4KcUiga2BiHCceCuKG9AWo/tYGPipIB8vPzmThxYqbFUNKELnn1B1wuKB9tP9/7Jdi6FLztmZVJUZScQxVKf2H+a/bd1w6PnAo3Dwddx1YUJYWoQukvDBwH318V1mCgcXfGxFEUJffIKYWicSidMGgCVIyH/GK7fftkuPconakoipISckqhaBxKJ4jAD9fAKf8Tatu9FjZ0qeyMoihKTHJKoShJUjYqcvvZC+HmETpTURSlR+SU27CSJAd9Bf7fZ1A8CH4xEryt9nXjQEBg8P6hPGCKoihJojOU/krJYLsE9rNqOO66sB0G9nwKD54Mfl/GxFMUpe+RUzMUjZTvJuXRtVIMVL1vXYtHzwZX4N/E0wRNNVA6DBb8vdfFVBQlu8kphaKR8t3ksAtg6plQVAG/OwaqP7Ttfq9VLNHU77D2lkAdbEVRFNAlL8WhKOAZd8USGDCik84GbhoM9x6ZdrEURek75NQMRUkR5WOhudbOUFz59h2AMC8w44fd6+B/xsDACVYhXfJqJqRVFCVLUIWidOSyt2wU/UMnQ/EQuOxvtv3OQ6H5C2hvDPVtb7SxLAC/GGWTUI6cocpFUfohqlCU2AwYDj9YFdl29WrwtMB9x9iCXX4f+D2h/Z5m+75tKTxyutYBV5R+hioUpWvkF8NVK0LbvxgVUiQOxq95whSlH6JGeaVn/HSndS0Wd2R77UabJ0xRlH5DTikUTQ6ZIS57C27YY43z4oL8Etu+ey08cELy4/g8djntoVPTIqaiKOklpxSKJofMMBVjYfyRdtbizFj2bOrY75HT7Su67Z4jYNeHULXM5hZ74HjNL6YofQi1oSipI9wIP+5LsO2fkQrhkdNh52pbitiVBw+fZoMjz3sGtr5H0C3Z+MDrgx0fwLo/w5JfQ8EAuPT1Xr0cRVG6Rk7NUJQs4pJXoXBAWAwL8MUnVpmAbf/8Pdj2Pjx4EhExLuEsvReq18Lez9IusqIoPUMVipI+jLEeYJ8H0rc4yiQcvwdqNsQfo73Jvvs88fsoipIVqEJR0odTd+WZ/7Dvnpauj+Eom7b6jnYXgF+Os6+HT4u9X1GUXkMVipI+LvmLfW/dBx+9QnBZS9zgLoSSIeAusG3h+cPEHcpw7OD3QtVy+M1MaNlr23weG6nfVm+Xzxp2pvVyFEVJjBrle5lz7l8KwKIFczMsSS9QMhhKh0Nrnc0NBtalePRhHfte+AL89xD7eewcq0C2L4/s42uztpRfTbAKqWKsDaJ02POp9Qy77O1QmzNr0ah9RUk7qlB6Eb8xfFLdwPDyokyL0nsceDJsXhzavnK5VQTRhBfzuvD/rAK6Y6rdHn9UwGMsrI+vzSqQaPZ8Bg9/BRCoXgNtDVBQCg27oGxkqJ8qGkVJOVmvUETk34DTgXLgIWPMXzMsUrdp8fjY0+qJ58/Uv3G5YcQMwEB+kbW/lI22S2KXvAq/PwG2xylLXDwotAzWVgef/yOwQ+x47Y3w5Ddhz2bbbPy25HGnafoVRekKGVEoIvIw8FVgtzHmkLD2U4HfAG7gQWPMLcaY/wP+T0QGAbcBfVahqCLphCuWhD67XDB4Umh7x7/su7hsjEv9Dti3xbYNGAHeNutRFr4EFv4X370+coYD0Liroww6c1GUbpOpGcqjwN3A406DiLiBe4CTgCpguYi8aIxZH+jys8D+PosT43fs5GGZFSQTvPubrh8T/lDPL7YuxMOnwiWv2bYbB1kFcuov4Z3bbUyLP457cbQycXCyIjuKpLXOzpYcHj7NLq+5C20QpioaRYlLRhSKMeYdEZkQ1TwH2GSM2QwgIs8AXxORDcAtwGvGmA96VdA08cqanby5vpqDR5bxp+/2kwSKddt6dvzwaeBpjaxlP6bSuiJPPA72PwF+d7S1m3SFre9at2O/D8pGhJbFwLa3NWBnOgKjZvbsGhQlx8kmG8oYIPypUwUcAVwFnAhUiMgBxpj7Yh0sIpcBlwGMHz8+zaL2DK/P4PX5WLO9PtOi9B7OzGHAyMT94vHtNzpvW/B3W0HSGyfepWSoLRAGkFdk7SgQUhrhyuTBUwL7nWUzA027dUlMURKQTQolJsaY3wK/TaLfA8ADAJWVlX3CXNHu83feKWeQgJE9jf9yLje480OKoGCAndnsXmeXxoZOhs8DCmX41IBdxhDTulWzHnztkW1+n01cGZ2qHzoGVarCUfoh2RTYuB0YF7Y9NtCWNH0hff2Q0oKI7U9rGuP0zBE++Ss01dj4k2lnpf98122zMxGw79/+q12qGn1YyPYC8PXfw7gjQqn2o2mLmj2Ky7b52u0M6Pcn2vZHTrdLY1v/EXrtXJ2crLGyLitKHyabFMpy4EARmSgiBcC5wItdGaAvpK8vyndRnB/6s1/66PIEvXMEX1vsPF7pwpkFHf0D+z7/ldCMYchkqBhvPchceVbRDD0oELHvCkXuR5NXFJk6pnqNVQY+j/Uwwx96eZoij42nOIwBn7fzforSR8iU2/DTwHHAUBGpAm4wxjwkIlcCr2Pdhh82xqzLhHzp5Juzx/HPzbWs3VFHU5uPLbXNTF9o07KvWXhKhqVLA/FmAOlk4AQoHQazLui4b8Bw++5yhZTMI6fbtC1jZ8DY2fBePI+0sKUxb6uNdzExli2N385gvvjYLq217IWiGD9yatbbfVUr7Xnbm+y4xliPMkXpY2TKy+tbcdpfBV7t7rgicgZwxgEHHNDdIXqFRQvmYoxh4nX2UhtavbgEvnHve+S5XbmVlqVsZChepLdwuaGw3NpTooll2whv+/XU0Gd3oZ1dgY1xiaaDMnFhZylYu42nxUb4g42VeeR0W0Bs5Ax7ztbA0uwL34Pv/RN2rbZjVq+DkYegKH2NbFry6jF9YcnLQUQoClv68hto8+agkf4/nu39c4YvcXWX8rEwfErHdle+tadEM+pQ2O9IyCu2255mImY0TV/A50tDnmWAjeTHZlS+b15IQXnbkpNRl8iULCOnFEpfMMqHc+jYgRHba3fU8/GuXrQ1KB0pHgLFg+HfH7PZkkuHwaAJ1mMMsbEv44+0MyDClqV2rbVK7Mw4DonGZ1++dhv7cvecyP27kjTkh9NcC639yPVcyXpySqH0pRlKPPa1aCGpjFJUbu0eYyttdP7Qg+xsZdRMKCwLRdGPnGG3HaXiROKPOjRyPHEToXgcvvg4Kno/gad7rJlI/Q47s6n+sAsXpyjppcsKRUQGiciMdAjT31i0YC5FeR1vwTn3Lw2muVd6mejlsvBtx/bhvK7bZuNcwhl2EIyYHtoeP5cuZ3F7+Yf23SkeFovoGBlFyQKSUigislhEykVkMPAB8HsR+XV6Res6fW3JC2DG2ArmTBgU0bZmex1+0ydiM/sP8ewy5z5p3Y/DMxcXDLC2liGT7TFDD4o/boQ9JjCTqd1kFUl7k31VLYet78H980Jd/3hpYnmjZzVqb1F6gWRnKBXGmHrg68DjxpgjsOlQsoq+uOQlIkiUi2hzu4/VVXU6S+kLDJ5osx8POTDUdulfbNCk46JcOszaXKKrUEIoYyjYWBewcSzGH3r52gETSBETwBdmuFdloWQJySqUPBEZBfw78HIa5el3LFowl0UL5jK8rDCi3eP1s36nGlxzgvmv2OWyvGIbl5NfGjYzCVMoP9oQ+tzeSCgtTKCPU1Csel3iJJjbP7DR+ttXhpTNjg9CecwUJU0kG4dyIzbg8F1jzHIRmQR8kj6x+h8Th5bS0Oqh3Wfw+Q0G8Pr8waDHqaPKgX5SOriv0Vlsi7PtxKGAVS6eJjtzcdK8FA+yiiZWsKTDwkG2AFk4296HYQeHxi8ZYjMTuPLstt9n3ZhrPoIHT7bpaBQlDSSrUHYaY4KGeGPM5my1odAHAhvjked2kee2gY4ALR4/pS5BQGcruUC4knnkdDuLGD41FPwIIUUTlxi2Nb/HenxFZyXwe+0rIijThM4fLZOi9JBkl7zuSrIto/RFG4rDogVzg7OQcFrbfTS2+Wj3+mlq9apdJVeY/4p1MXbyiY2utO2jZ2G/li4YHPbDyLGvYGIHPvp90NZoX/VhOVWjI/z3bgVvmIeY2l+UFJJwhiIic4EjgWEick3YrnJsvi0lSV5es4N5nfRxlrMmXBv61egL/KB0ouh1ppJDhOcSC6dwgLW5HHouvHilXQZz5YO0Bwz1TvyKEFGvxfmcyKW4qRpuHmETaJYOt0GbndGV2YzOfPo1nS15FQADAv3KwtrrgW+mSyglPk1t3oicX86MpU/YVtY9D6f+T6alyD6i414cPngi9NkJmNz6HkHFIRLpJdYVfO1QXwVlo2D3evifsbETUo4MrHTv+jBULhlUcSgxSahQjDF/B/4uIo8aY7b2kkzdJqttKF343s8YW8Gepnaq9nasPOg31q24vNjFOfcvZcWWPZQUZn2dNKU7jP8SFA20MS7hde8dykbbLMkREffJEGb03x4on5BfGrtw2K4PA7OiKNuLQ3TCS5/H5ivz+21G50SKR5VSzpGsDaVQRB4Qkb+KyFvOK62SdYNstqE0tllD+4AkHv7F+W5KC+KvKG6pbaKl3YcxBmPAaBBkbjL0QBhxiI1jAfvgDY93KRkCBaWJx4hX3yUaT7NVBu1N9t3vs0rE77MuzJ4mO6MB60yw41+xvdFqNti8ZJ/9PbnzKjlFsj9t/wjcBzwIdPXnkELooX/hkft16/g8l+D3G/xYe8qH2+sYWJKPH2hs89HY5k1KWWWUWFl6la7x5Z/BsxcCYiPyR86AvZ+FGeLFKpn2QCVQV761ldRXWeXi8xB7umxsJUoIvUezbxs01VqvMm+rPe/g/a3Sqd0E9x8XCr50zt9Z9cpELtJKnyPZJ5DXGKLKlTsAACAASURBVPO7tErSTyjM69yXYdGCuXxe28xxt72NAQ4ZXU5JQR51LR4+2tWAP/A82NccSiS5ets+jjpgaJqkThEVYzMtQd8jejlo8lesa7HLDZcESgc9crrNOmx8gEQq7qGTbQnmIOGG/G4QXnnT22aViqfJvmJlAvA0W7n2boVBYT+mHjk9ZA9a8xxMV5NsLpDsT8aXROS7IjJKRAY7r7RK1s8ZP6SE0sI8BhTmUVJgv6h3nzcrbv8+seoV64GjdI28Als4LHq2N+rQgPvxLDtryS+xbsffeQuO/ZHtUzHOZlDuLtG2GnFDw67Qtj9s9vPWzZHHNFbbd8dNedeHob57t3RfJiWrSPYbflHg/T/D2gwwKbXi5C4Thw1AdnXeLxEHDC/D7RJ8/o7a4+ZX1lNRnN83vL2UnjGyk2Tf4QZ8lwsqL4HVz9pZzb8/DvcdFXWAiwi342gKyiJnJg7RXmGO0ghu7w59fvxrVuEpOU1SCsUYMzHdgqSCrPby6gHJKImttU3MiCrYlXWoN09qiJXWpTOcOi7lo+0Mp6AMWvbYtoISO3OJWBoLI1F9+3geZnu3wsa/hLY9zSFvMIiqXBkg3OsrngdYsp5hPT1e6RZJKRQRuTBWuzHm8dSK0zOMMS8BL1VWVn4n07KkgjULT+nQlhc1QyktcNPUHvtL3adiVJTUEu+BWTIYxh5hP9dvt0b6gftB7Ubbll8SiK4Ps7WMnGGrTEbTsjf++Y0P3v1NVJvfFgbLK9R6LjlKskteh4d9LgK+jK2LklUKpT9w8MgyVlfZei9F+S7GDCxm4+5GfH7DxuoG1m6v479fXh/z2IwomLxiu94fK8ZB6T1iKZjyMaF9f74MPlwU2idi793Qg+z+hTFc8RPFv/hiVB5tbwx5f9mTEJFbzEmcGeuzE+eiZDXJLnldFb4tIgOBZ9IikZKQArcLt0twCeS7XeS57XJEu8/Q3uzhzQ12HdtJ0RIrP1ivkldg64XEKoOrZA+fBsLKCgZY761BE6FsZGj/wP2grsoul7nyrGJI6PJroLWzQncBhfLB43a209Zgt3f8CzwtoRib9ib4fKktOmZMbBtS9FKWMaFjrtsW6uf3JV7CU3pEd91umoA+YVfJNVwuoSRB0GM4fr+hzZsFYUMaf5JdxMshBjD3e/DJGx3bK8aFlrjiLYFF0xywyXSWkj9IYLZifIDfOgKEL435fdY9+fP3ImWPNXOpWmbH8HtDbS37YNtScBd27K+khGRtKC8RcgFxA1OAZ9MllBIfJ39XoiSRrR4fHp+fVo+fVdtClR89Pj8NrV5aPT6K8nUJqt/TW0tIeUUdsx6DTd2/e22gKmXYEpk/7EdQ1Qrso0cilZKzDObkGHOWyByc5bhwf3pnxuRUu1QDfcpJdoZyW9hnL7DVGFOVBnmUJJk6qpxFC+ayfMsezr4vMqX9x9UNtHo6/iLcWttMbVM7L6zazjmHj++wP202Fv3C9k2i71uspaa8Ilvca+eqQIMEDPthNV3izVD/7R54YJ51DghXFv5w+4vTbiJnG8niHPPI6ZGeZeG2mWhiKZruKp9+prSStaH8XURGEDLOZ2W1xlx1G44m0QP/jyuqYiqT9TvraQ54g8Xar/Rz3IF6K2Ujk6tAefNIm6LFlWftLg6uvI4KZNgUWy0yPJZF3HaGAsmnXwlXNI5dBewsxOextpHgUpiEjnnwZGv7qfmo83OEK5poZRAr47IzM1KnASDJSHkR+XdgGXA2tq78+yKSdbkSsjk5pKJkNQPHw9g5tgZLMjiFwKKVR14hDJ/WycGujtUlu4qTAdkYO/PwtoSW1YyJ9EBr2h17DCXlJLvk9VPgcGPMbgARGQa8CTyXLsGU1LH8MxvAJi77q+3h9z7j1TU7O8x0tHhXP0Yk+czEYBWJUwgsckfHvhf82S6N3Tws8nydjR8xc3ERkXY/PLJ/d5ib/NZ/QH4RERi/nUm0h9lxdn0YqHBZH0q177S3JyrBHIfo2Uu6ydKltGTdb1yOMglQ24VjlTRSWtD5bwI/yaUDtOnw+0JSMCUrcJZ55r9iXcMLymxOsQufD/UZewQUlYfcxwvKknTbjeqT6BhfeEnk8P/fwCOqYZdVHuEKyRi7vetDuP/YyPGcJbRw/N5IA38sRwCnPdwDLdZYvUGGSjsnqxT+IiKvi8jFInIx8ArwavrEUpJlyqgypo8pp7TQem3VNMaoN479moWnt1+/s75DffrmQP36zTWNKEqXcOWFHvoFJTB6trWtuPMj+4w6FMbPtcoolpJwltDKR9tZjUNE3Zeo4zqzwfi9RM5uiHQa6HC8gar3ofkLu1m/E7b90x7jzESS4ZHTrfty1fuw8fXk+oePnSGl0BMSKhQROUBEjjLG/CdwPzAj8FoKPNAL8imdICKUFOQFv2Lt3vhfrrqWyF9KK7bsYfrC0D+6k9Hlyqc+CLadc//SDopHyUGcmUZP+ocbpvOL43t3Oce68+PbW2adH5ldIWJpLdEMx4RsKc7SV6yI/nhKyNtGcJZT8xE014Zynhl/aGyn8FhbvX21N9naL9ER/s55/vKT+MohliL55bjYM6B0kSLF1dkM5U5s/XiMMX82xlxjjLkGeD6wT+mDVO1twR+VsXhvk+ZWUjJA8aDI7YpxtiDYlDMj2zsoO6HT7AuD97cFxpKlrSFq+QyrLMJx3JB9bZFKyfisc0BvLRk3fWGLpvmzy2OzM4UywhizJrox0DYhLRIpPUYAd6IlZ7+hqd2HL+x//28fqSeMkiHCZyIlQ21RsBFTO/YLVz4Sw8MsFuWjO++zZ7P9hf7Zko77mmrgwZM6tsfSG/5ACeVYs5/6GGF7vxxnX4lwltiiZzFffGTryDx4Yuzj/D776uVls84suonyofegUo+SaqaOKmfZllD2VxEBYyjKc9Hu9UevIAdpbvNyzv1L+cbsUDXFz2qbmb7w9WAeMGdpLFb2Y0UBYi+XRcdmxFtSc+wjceu8BH4dhS+POctB4TOIwvJQPjCA0/4Xnr+8M8mty/H2lfC3mzru87V3zIy868PYkf8A+GH7Cmv/aW8KLbd52+xxvxyXOGbl86X2GHGH/i7hsS4diKHZ2hptihlXXiCPHr3mFdaZil8hIh1SwYvIt4GV6RFJ6SqLFszt4AJcmOdi0tBSDhlTQUF+/NusTl1KWuiqTSYW11XZbMfRD9L5r4S1BZSNuDsvPOb0i0kXvgj+TvLjOQqoq1+u8IwBxlgl1B0PMcftuTuZBXpIZzOUq4HnReQ/CCmQSqAAOCudgjmIyCRsHEyFMSbrgimzFZdLGFbWeRI8P3YGUtPQ0TssPC6lodWrsxQl9TjFtGIZoF0uKB3WsT2CJB/arvxAMGQCm0OypYg9zUmeN+pc8eJbnJiY6P6OPcedb4+NZ6SPXtLKYK2ZhDMUY0y1MeZI4EZgS+B1ozFmrjGm2wVtReRhEdktImuj2k8VkY9FZJOIXBuQYbMx5tLunqs/0tXk3NE/pFrafcEiXq2eSFuLomQleVmWQdjTEtUQiN6PNWvxttFBmcRdpE4DO1bZZboUkGwur7eBt1NyRsujwN2EFegSETdwD3ASUAUsF5EXjTGxq0UpEYgI4wYVU93QxqxxkaavPJeQ7xY8CTTDjn2RXwC/MbjpeIzjZpxNMxWtTJlDJLtM5vSLLvyVVxiZBHL2RSHbSEGp3RfrF7y3jdg/xcKKgHUFf4KlKsfQbkxg1tKJ8nBmL8Z0nI3UfGzLOzvLfeGuyhBQYCa2A4Mz1szzYpdk7gbdrYfSI4wx74jIhKjmOcAmY8xmABF5BvgaoAolSUYPLGZfVKyJ85A95/6lrN1RR1NbFtRHUZRwuqpEYjFkst3/ryfhhe913F88GIYeGKiTEotAivyY7WmgYVdkssyEBBRErIe+pwnySzu2O1R/aB0VxhzeK7O4bEqfMgYIK61GFTBGRIaIyH3ALBG5Lt7BInKZiKwQkRU1NTXplrXPE/7VifeV0TQsStYzYrr1qCpI8FCNR/Hg2O0lQ5Mfo3wslI0O5EFLcrHZSWjZVeKWXE7wPW0LKC1fe6944GRkhtIVjDG1QKe+f8aYBwhE71dWVvb7J2Fnyz8uIWgbMYA3KtCx1ePvkBmjqdWLCJQUxv+30eUnJaV0NnspqrAPcwnz9ko01u0HQ8NOu10+JhQFDwRnAgNGhNKuxMT5rggMmmA/tjda7y/TWea8QCVKb7SNJQmcJJfxvNmcxJaxFMeu1dat+pHTrWtyeHnlxuquyxKHbJqhbAfCo3zGBtqSRkTOEJEH6uo6q2WtROOP8U9oTGS6pX6vpZXs56CvJNcvv8Q+YDstB9xJmpekiBHVn2wNmFQSHfWfBrJphrIcOFBEJmIVybnAeV0ZwBjzEvBSZWVlh9iZ/oJTybEzSgrzaGr1Bs2B/jjfjXA9YwLbzW3eiMBHiCxNHN7eGZ3NaLo744l1nM6ecpTwQMHiRLHYYZQMtb+W8goDRvo4XwBxJVhqSgYJGcSTGcdd2DH9SzhtMRK3epph67tdE6s7KfqTICMzFBF5Gptg8iARqRKRS40xXuBK4HVgA/CsMWZdF8fVGUoChJDHV3fRWYrSpykebJfFEime/FhJQFxQlKSyCifb3JnTTKa8vL4Vp/1VepAWv7/PUBL98g739nICFg+fOJj1O+tpbPUmVBQuiZzBZJutvt3b0d6jKDEpHmTtB4VxZtGuvNh2GBEYcUjimUBBqY0/cVyGndovPo/d5/OEbCeJvkTitnLEm6nETfsSPU50kbIo0vBFzqYlL6UP4NTJa2r18v5ne3ALVE6I9JZxaq0ku7S0fmd9zCW0ZGj1+PjXtn3B7fBU+12VI1Uks7QWvjyoS3CpJhA7klfUac+Yxyb966Qnv2LiLbEFxowb7e6HHR/E2deV04ctv9Vti9+vi+SUQhGRM4AzDjjggEyLktWEP8SmL3y9u6FbWUGbJ7vSdyu9TCwvsIJS+8t81EzYnCAe20n7snO1zXvlbbVLVJe+Dr87BvZ8aoMGYyK2iFj0+W87CBp3RRUcK7V2nu2J0h+67LldeR3zd7nyA3m5sv9bmk1eXj3GGPOSMeayioqKzjsrQVK5XOQ3JqbHmKL0GuKyD2ZXko83EZsvKzwlflF5AmUCIB1r10fsTnRsIhniVbHsxpc01vewrbFjewq/rzmlUJSuM3VUecK4EoCJQ0sZUNj5F8TnNzS1+TQaX8ksI2eEYjXyS+y7K6r6ozOzcDIXj5zRMThy5Ay4bpt9Rcd+FJTEPrdznvLRoXGv22bPM3hSHIHDlMV12+CGPYFAyQAVYxMrr7jEikdZRdwklClAl7yUmAwdUMAXjaF1XBHBRciGEr3QdM79S9mwK+TnHs9+0Vk54fAMx+Ftjo0lfKkO4N3/d0Kwn9PWXbflaNtHIltIOm0giWTq7jjd2Z8TzL4YPngCCst653wV46BlLwwY2XFf0UCr4DytBL9Bwe04lI6wgZZOMGaWk1MzFF3ySh1F+R1nJOpJpfQ58outZ5erl347u/ICrzieYtFyhNtaYlE8MLnKlB3OlZlHe07NUJTOif416mxPuDZxios1C0/hnPuXsuyzPRHtfgNbvmhi3OASWrthIN9d30pLu4+iOEXAPt/TTKvHR2Fe8l8QY0zwuPC2rbW2joUxxla0VPofiVK5JCrQFZ3deNjBXT+X4wSw9R+R7Y7RPhGDJ0H12sR9wskrSt69OIXk1AxFAxu7T3eDHQ1Q3dBGY5s3WEOlK3xW24zXb+Ieu7OuNWHa/Vj4Deyqb8PjM8EElx6fYVd9K7vq22jzqmeY0k0GjISiQZE2jt4gv4SEhnnHEcHpkx/HxtPxwB4KFklOKRRd8uo+RfluihOUClYUBRhyAIyY1vtLSu6CxBmVCwZYJSIu62HmROj3spy65KUA1pBd1+Lho10dazRsDyu+tf+wUjbVdMwDtLW289xA4VH60bR5/KzfWU9zmzfYN9EYzW1efAaO+tXfgvua27wdPNYc54BwV+YLHnwfl0uCssQy4oefJ7z08Tn3L2XFFrvsV1KYF+F8EH59zpjRBvZ4uc6iAzKj+4SPE+18EMvxwZExOug0+m+TDNH3Il1G/JQ5CQSWmoLjdTaZiLUM1lmW486WzhLtLyix2Ymdvp2N5+z7nzGJZYKQ0mlKsoSHu0C9vJTsp6uxKH6gzevHb7o/CTcGPD4/ze1htpNujtVdjDF4/SbCfhOLxjYv7V4/Le3qYq3kDjm1xqE2lJ7h2FHcAgXunv1reLphp2j3+oNuyd2l1eNnY3VoltXey/YSq0z8bI4xiwvn05pG2rx+tu5JT9ZXJZtJonaL4+bcY1tN7zqf5JRCURtKzygpyGP2+IEUF7gZVlbI+MGdG/ay0VcqulhYJvD4EysyZwKXBaIqvY2ITRw5dk78PkMPsnaRojjPsoIyyIuVFTkKlxv7mI/1TU394z+nFIrSc/LcLqaNrmDRgrkdPL8WLZjL1SdNzpBkitJN5r+SfN36dJ9r/is2Gn7UofZhH6+/iM1FFr5v6EGhz6MOhSEHRh7jROaPnAGlw2xb6XDY78iO8S8FA6BwABzzo65dXyeoQlEURVFSghrllZQQnfW33WeCHlAOjudRMoR7fAE0tvl4PxBUGT5xag7LGxZrkcnjM7z/2Z6IX07Lt+4NfnZLyBNs/c56Glq97H+d/VUYHv4yfeHrNLR6I45z5HNSw0Dk3yHcq8s5Pjrcp7ndF/wbrdiyJ3jOFVv2MH3h6zE91+x1e4NeWtFpacKJTlvTmadZeDoZR65Y3nA9rYgZy/vNua5wr7pkxnCId33dTWWTKL1OtNedc5/i3YfO0vhE71u309qBp4XvD/vfWbezjmmDQg1N7T5iORXXNLaxu7WOaZ1ebWrQGYqSEmI9zJvavHh8kXv8xj6onVcsbKnh+MYFk+DYruIz9qHu8flpClNgnTmphV9DIlm7iz/MxtLU5u2y15zzN4oXxLmvuZ3GGPdHUXpCTikU9fJKP0dMHJJ0VL0/xgMt2cdiImN1qh/fPr+hud3XbQN5OrL1m7B3v+m+8T6el1tjm9cqHfUK6DvkFYG7gDYCnl8FpeDKD20nwp0fO5/Z1DNtvZUUudfklEJRL6/us2jB3KSWAYaVFXLY+EGUFeVZ9+Iu5NhSlD5Pbxr4o7nkVRg7B6/k222XG/KL2VRwEJeM/FNItjD53i0OZON25UWmYxk2xRrvh0+xCTRVoSiKoijZhCoURVEUJSWoQlEURVFSgioURVEUJSVoHIrSa6TDGypdhIuabEGuWO7DxpgIl18Tx83Y35m3lTHBv58JfHbGiZbNGJNSTzjnPCZMhnSQ6vH9xiB0/Pso6UNnKEoE4d5ePzjRplnJc3X8N5k6qpxPf3k6Pzt9SlLj1rV4YsaqxKM9kHk4FuHNscZ0xXh+dOU5ZUxk/8Y2H94E8RptXj/GGJrafMHjRlcUB9tWbt0bfCj7gTXb6xlebutVlBS42V3fyvKtexO6Brd4bBblf26upbHNhz8gV3OMbMWrq+qCcnT2KPX4DI1h8TexaG730djm49OaJpZv3cuST5JMjd4F/H4TvK5U4PX5Wb5lL2u2J5eiv9fpxFvspiG3ctOQW2Pu25I/ye6b/0qCSo/2zvvEbftG9/vWU6Hzj5wBw6d29QpiklMKReNQUssJBw9nvyElTBwaP0nkGYeOjtguL4o96U1UJdEdpgGcUsCdKYB0/uaMde5kwjVidXHiSMJp8/qitv1xj4+WYWddS4e2aML/1rGUazT1LZ5Ozwuwr9n227GvJUHv7pHq8EonQWhLJ2UEcpayUTBoIm+WnB7ZnlcMQybDgOFpOW1OKRSNQ0ktAwrzGFleRElB/JXRgSWRQVUjK4roygqDAEVhsSz5bldSD8HoPskc058I/jl0uad/klcI5WPY6x4S2e7OT5sygRxTKIqiKErmUIWiKIqipARVKIqiKEpKUIWiKIqipARVKIqiKEpK0MBGpcskk5W4u0wdVc7aHXW0ehI7klZOGMya7XXBOIzCPHevuoi6JPvrwbsFSgvcKU/3r2SQ+a9wU3hhsejtBMcB8Mjpsfdt/jtwXI/F0xmKoiiKkhKyfoYiIqXAvUA7sNgY84cMi6QoiqLEICMzFBF5WER2i8jaqPZTReRjEdkkItcGmr8OPGeM+Q5wZq8LqyiKoiRFppa8HgVODW8QETdwD/AVYCrwLRGZCowFtgW69dM8CoqiKNlPRpa8jDHviMiEqOY5wCZjzGYAEXkG+BpQhVUqq1CbT5/j05pGhpcVUtPQnvQxPr/pUa3zWId2NtqnNY0J+ybKReb1mw5OBF6/4bMvmoLbLXGcDBpavTTFSc4Yq/2jnQ0d2prbfdTFycdV3+rF7RK27Wnmi8Y2PtndSEu7j4bWyLGNMWytbaKh1cvmmka8Pj/1LR42h12DwyfVjWysbmB4WVGw7dMap60wphzh1DV72LS7kfKiPDbsrGdjdQOtYQ4VLe0+1m6v45Ax8VMordq2j43VDQwbUMgXTe0MLM5n7fY61myvo6Xdx56m0P/b9n0ttHl8tHp8bK5pojCJstX/+nwvG6sbaGn3sWl3I83t3ogURLWNbexpaudvG6ppaffhM1bujdUNbNvTzLjBkfnvPt/TjN9vOmSurmmw92RgcX5E/9rGNvY0e9hZ18L2vS20eW2C0jy3i621zUBorIZWDzvrWllTVcf0saG/2aftA5O+J6kimx7QYwjNRMAqkjHAn4FviMjvgJfiHSwil4nIChFZUVOT+myoSnKU5LsZWhr6B271+Pl8TwutUQkRXWJfNr24TRCZ75bgMeGUFrojcnVJ2BgOAwrzguMFZSlwJy23x9c1BVbgdlFaGBrfG6XFRKCmMfRQS6Qg4+2K1b74447/2z6/obYxvsL2+Q0761r5sGofe5raO8gK9m++q74NsHLv2NdKbVN7TEW65JMv2NvsYXdDa7DtvU22zRkjEet21FHb1M5ntc0s/riGvc2eiLT1Xr/hzQ3VCcd4c301e5s9bNvbwp6mdrbUNvHmhmqa223G5x11Idmq9rZQ09jO8i17qGlsoyqJ5JZvBMb3+g21Te1srG6M2F+1r4U9zR7+9y8fB/+eXr9hb7OH9z/b02G8nXWtVDe0Ud8Sqcide7J1T3NEu3NdK7bspWpfS/CeeHx+dtW3squ+Lfg92dPUzt5mT4e/2XtNY5O+J6kimxRKTIwxTcaY+caYKxIZ5I0xDxhjKo0xlcOGDetNEZUAJQVuCvPdTBpWSnF+5L9WdIrC0sI8SgvzELH1KkoK3BTld1QAhXkuDhldQWlhHgF9E1Qk4b/0SgrclBbmURymRKaPqQge41AWlQ05+teqC+tu67zifUHyA3JFX6M7oCiL8mIrs2RSNbqlo5x9h8w4KUfryFg1Z1JJqof3Rw2Y2uF7755kk0LZDowL2x4baEsaTV+vKIqSObJJoSwHDhSRiSJSAJwLvNiVATR9vaIoSubIlNvw08BS4CARqRKRS40xXuBK4HVgA/CsMWZdF8fVGYqiKEqGyJSX17fitL8KvNqDcV8CXqqsrPxOd8dQFEVRukc2LXkpiqIofZicUii65JV6Fi2Ym3QyyElDSxP2DXfznTqqnKmjym2Sx4WnBLejjx89sChi2y02MeSiBXNZs/AUBpZY//2LjpzA1FHlHc5ZUpgX8U8+dVR50IOrrCiPEeVFHY5xjispzGPm+IEdzn/ERCtzNG6XUDlhMMUxvNXAftlGVhRRGMcDzOlTUtizhQMRuOqEA5i93+AejRNNV9ywlQwz/xWYe0XsdidRZBrIKYWiRnlFUZTMkVMKRVEURckcOaVQdMlLURQlc+SUQtElL0VRlMyRUwpFURRFyRw5pVB0yUtRFCVz5JRC0SUvRVGUzJFTCkVRFEXJHKpQFEVRlJSgCkVRFEVJCTmlUNQoryiKkjlySqGoUV5RFCVz5JRCURRFUTKHKhRFURQlJahCURRFUVJCTikUNcoriqJkjpxSKGqUVxRFyRw5pVAURVGUzKEKRVEURUkJqlAURVGUlKAKRVEURUkJqlAURVGUlKAKRVEURUkJOaVQNA5FURQlc+SUQtE4FEVRlMyRUwpFURRFyRx5mRZAyR1+9c0ZACxaMJe3PqrmkkdXBPeVFOQxZmAxH1c3BPuEE73tsODY/Tl3zvi456zcbzBvbqgOjrFpdwMn/vqdiDG31jYx79bFMc/z4JLN3PzKhuD2QxcfzvEHDw9u721qZ9Z/vwHAlltO7yDzhGtfCW5XFOezaMFczrjrXdq9/g6y/uWHx/LcyiqeWLo17vX88ydfZnh5EUDE2J3x5SnDE+7Pd3f9t+OXp4zg6WWfx91/+ITBHdpmjB3Y5fP0hHGDi9lY3YhLIttLCvJobPOm7bwjyov4fE9zcHvyiAFsrG7s9ngFUfdn2IBCtu9rSerY02eM5qF3P4u7vzfvic5QFEVRlJSgCkVRFEVJCapQFEVRlJSgCkVRFEVJCapQFEVRlJSgCkVRFEVJCVmvUERkkog8JCLPZVoWRVEUJT5pVSgi8rCI7BaRtVHtp4rIxyKySUSuTTSGMWazMebSdMqpKIqi9Jx0BzY+CtwNPO40iIgbuAc4CagClovIi4Ab+GXU8ZcYY3anWUZFURQlBYgxJr0nEJkAvGyMOSSwPRdYaIw5JbB9HYAxJlqZRI/znDHmmwn2XwZcFtg8CPi4x8JnJ0OBLzItRBrR6+vb6PX1XQ4yxpT1ZIBMpF4ZA2wL264CjojXWUSGAL8AZonIdfEUjzHmAeCBVAqajYjICmNMZablSBd6fX0bvb6+i4is6LxXYrI+l5cxpha4PNNyKIqiKInJhJfXdmBc2PbYQJuiKIrSh8mEQlkOHCgiE0WkADgXeDEDcvRVcn1ZT6+vb6PX13fp8bWl1SgvIk8DsU4YgQAABoVJREFUx2ENWdXADcaYh0TkNOBOrGfXw8aYX6RNCEVRFKVXSLuXl6IoitI/yPpIeUVRFKVvoAolyxGRLSKyRkRWOW59IjJYRN4QkU8C74MyLWeyxMqeEO96xPLbQEaFD0XksMxJnhxxrm+hiGwP3MNVgSVfZ991gev7WEROyYzUySEi40TkbRFZLyLrROQHgfacuH8Jri9X7l+RiCwTkdWB67sx0D5RRN4PXMeigG0bESkMbG8K7J/Q6UmMMfrK4hewBRga1fa/wLWBz9cCv8q0nF24nmOBw4C1nV0PcBrwGiDAl4D3My1/N69vIfDjGH2nAquBQmAi8CngzvQ1JLi2UcBhgc9lwMbANeTE/Utwfbly/wQYEPicD7wfuC/PAucG2u8Drgh8/i5wX+DzucCizs6hM5S+ydeAxwKfHwP+LYOydAljzDvAnqjmeNfzNeBxY/knMFBERvWOpN0jzvXF42vAM8aYNmPMZ8AmYE7ahOshxpidxpgPAp8bgA3YQOWcuH8Jri8efe3+GWOMU/g+P/AywAmAk3w3+v459/U54MsiIonOoQol+zHAX0VkZSC9DMAIY8zOwOddwIjMiJYy4l1PrKwKib7g2cyVgWWfh8OWKPvs9QWWP2Zhf+Xm3P2Luj7IkfsnIm4RWQXsBt7Azqr2GWO8gS7h1xC8vsD+OmBIovFVoWQ/RxtjDgO+AnxPRI4N32nsfDRnXPVy7XoC/A7YH5gJ7ARuz6w4PUNEBgB/Aq42xtSH78uF+xfj+nLm/hljfMaYmdiA8jnAwakcXxVKlmOM2R543w08j/0nqHaWDgLvfT0jc7zryYmsCsaY6sAX2Q/8ntCySJ+7PhHJxz5s/2CM+XOgOWfuX6zry6X752CM2Qe8DczFLkU6abjCryF4fYH9FUBtonFVoWQxIlIqImXOZ+BkYC02s8BFgW4XAS9kRsKUEe96XgQuDHgLfQmoC1ta6TNE2Q3Owt5DsNd3bsCbZiJwILCst+VLlsD6+UPABmPMr8N25cT9i3d9OXT/honIwMDnYmwJkQ1YxeJkco++f859/SbwVmAGGp9Mex7oK6FXxiSsF8lqYB3w00D7EOBvwCfAm8DgTMvahWt6Grts4MGu114a73qwXin3YNd51wCVmZa/m9f3RED+DwNf0lFh/X8auL6Pga9kWv5Oru1o7HLWh8CqwOu0XLl/Ca4vV+7fDOBfgetYC1wfaJ+EVYSbgD8ChYH2osD2psD+SZ2dQyPlFUVRlJSgS16KoihKSlCFoiiKoqQEVSiKoihKSlCFoiiKoqQEVSiKoihKSlCFoiiKoqQEVSiKkmZEZIKItARyKDltI0TkKRHZHMjTtlREzkowxtvR6dFF5GoR+Z2IFAfSqreLyNB0XouiJEIViqL0Dp8am0PJicj+P+AdY8wkY8xsbHrwsQmOfzrQJ5xzgaeNMS2BsXekQW5FSRpVKIoShoj8UUTuFpF3RWSriBwtIk+IyEYReShFpzkBaDfG3Oc0GGO2GmPuCshwfqAQ0ioRuV9E3Nj04aeHFT+aAIwGlqRIJkXpMapQFCWS6cBmY8zRwP3Y3E7/D1tM6XQRKUzBOaYBH8TaISJTgHOAowKzDh/wH8aYPdj0F18JdD0XeNZoqgsli8jrvIui9A9EpAgYCNwZaDLAQyaQ0FBEfEB7Gs57DzaPVDu2oNFsYHmgllExoey9zrLXC4H3S1Mti6L0BJ2hKEqIacAHxqYpBziUQIElERmLtVFMEZFrA213iUiZiEyNbuvkPOuwZYIBMMZ8D/gyMAybUPExY8zMwOsgY8zCQNcXsFXzDgNKjDErU3DNipIyVKEoSojp2MzODjOwmVnBKpcPgcOxWWgBKowtFRurLRFvAUUickVYW0ng/W/AN0VkOICIDBaR/QCMLd/6NvAwdraiKFmFKhRFCTGdgGIILH8VG2P2BvY5yuVwYH2gPo1DrLa4BOwe/wbME5HPRGQZdqnrv4wx64GfYcs+f4gt0xpej+NprHJThaJkHZq+XlG6gIi8gq1zUg9MN8acGqst6pgJwMvGmEPSLNsWbM2RL9J5HkWJhxrlFSVJAuVha40xCxK1xcAHVIjIKicWJcVyFQNLgXzA30l3RUkbOkNRFEVRUoLaUBRFUZSUoApFURRFSQmqUBRFUZSUoApFURRFSQmqUBRFUZSUoApFURRFSQmqUBRFUZSUoApFURRFSQmqUBRFUZSU8P8BDCe9jsIg/RoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = hist.plot1d(output['mass_far'], overlay='dataset')\n",
    "ax.set_yscale('log')\n",
    "ax.set_ylim(0.1, 8000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.1, 5000.0)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAELCAYAAAD+9XA2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deZxT9bn48c8zmY1lREQ2RTYtsgiCjLiLtlatC15t/aGtC9SKYrX12voqthWh99p6r1ptXaG3uFWtXqvXBb3eulBFqQgWUFAsIpvs2wwwzJLk+f1xTjJJJskkmZNlMs/79eI1yTkn53zPZMiT7/Z8RVUxxhhj2qok3wUwxhhTHCygGGOM8YQFFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT5TmuwBeEpHzgfOrqqquHjJkSL6LY4wx7cbixYu3q2rPtpxDinEeSnV1tS5atCjfxTDGmHZDRBaranVbzmFNXsYYYzxhAcUYY4wniiqgiMj5IjK7pqYm30UxxpgOp6g65VX1ZeDl6urqq/NdFmNMS01NTWzYsIH6+vp8F6XDqqyspF+/fpSVlXl+7qIKKMaYwrZhwwaqqqoYOHAgIpLv4nQ4qsqOHTvYsGEDgwYN8vz8RdXkZYwpbPX19fTo0cOCSZ6ICD169MhaDdECijEmpyyY5Fc2f/8WUIwxBcfn8zF69GiOOuooLr74Yurq6uIed+KJJ+a4ZImddtpptDb/7d577014L5maN28e77//vqfnzFRRBRQb5WVMcejUqRNLlizhk08+oby8nIcffjhqv9/vB8jKB2no3NlgAaUdUdWXVXVKt27d8l0UY4xHTjnlFFatWsW8efM45ZRTmDBhAsOHDwega9eugPOhOn78eC644AIGDx7MtGnTePLJJxk3bhwjR47kiy++AODll1/muOOOY8yYMZxxxhls2bIFgBkzZnD55Zdz0kkncfnll3PqqaeyZMmScBlOPvlkli5dGlWu/fv3c8kllzBs2DAuvPBC9u/fH943depUqqurGTFiBLfddhsAv//979m4cSOnn346p59+esLjAKZNm8bw4cMZNWoUP/3pTwHYtm0b3/72tzn22GM59thjee+991izZg0PP/ww99xzD6NHj+bdd9/19HefNlUtun9jx45VY0zhWbFiRUrHdenSRVVVm5qadMKECfrggw/q22+/rZ07d9bVq1e3OO7tt9/Wbt266caNG7W+vl4POeQQnT59uqqq3nvvvfrjH/9YVVV37typwWBQVVX/8Ic/6E033aSqqrfddpsec8wxWldXp6qqjz76aPg1K1eu1HifKXfffbdOnjxZVVWXLl2qPp9PP/zwQ1VV3bFjh6qq+v1+HT9+vC5dulRVVQcMGKDbtm0LnyPecdu3b9chQ4aEy7lr1y5VVb300kv13XffVVXVtWvX6tChQ8Nlv/POO1P6vYbEex+ARdrGz96iqqEYY4rD/v37GT16NNXV1fTv35+rrroKgHHjxiUc7nrsscfSt29fKioqOPzwwznzzDMBGDlyJGvWrAGcYctnnXUWI0eO5M4772T58uXh10+YMIFOnToBcPHFF/PKK6/Q1NTEnDlzmDRpUovrvfPOO1x22WUAjBo1ilGjRoX3PfvssxxzzDGMGTOG5cuXs2LFirhljndct27dqKys5KqrruL555+nc+fOALzxxhtcf/31jB49mgkTJlBbW8vevXvT+K1mn81DMcYUnFAfSqwuXbokfE1FRUX4cUlJSfh5SUlJuF/khhtu4KabbmLChAnMmzePGTNmxD13586d+eY3v8mLL77Is88+y+LFi1Mu+5dffsldd93Fhx9+SPfu3Zk0aVLcYbqJjistLWXhwoW8+eabPPfcc9x///289dZbBINB/v73v1NZWZlyWXLNaijGmA6jpqaGQw89FIDHHnss6bE/+MEP+NGPfsSxxx5L9+7dW+w/9dRTeeqppwD45JNPWLZsGQC1tbV06dKFbt26sWXLFl577bXwa6qqqtizZ0/S4/bu3UtNTQ3nnHMO99xzT7jv5swzz+S+++4LnysUcCPPmW8WUIwxHcaMGTO4+OKLGTt2LAcffHDSY8eOHcsBBxzA5MmT4+6fOnUqe/fuZdiwYUyfPp2xY8cCcPTRRzNmzBiGDh3Kd7/7XU466aTwa6ZMmcLZZ5/N6aefnvC4PXv2cN555zFq1ChOPvlkfvvb3wJOp/6iRYsYNWoUw4cPD498O//883nhhRcKolPe1kMxxuTMp59+yrBhw/JdjJRs3LiR0047jc8++4ySkuL67h3vfbD1UGLYPBRjjBcef/xxjjvuOG6//faiCybZVFS/KbV5KMYYD1xxxRWsX7+eiy++ON9FaVeKKqAYY4zJHwsoxhhjPGEBxRhjjCcsoBhjCtrEWQuYOGtBvothUmABxRjToYRS448YMYKjjz6au+++m2AwmPH5QgkqY02aNInnnnsu6WsnTZpE586doyYm3njjjYgI27dvz7hM+WIBxRjToYTSuixfvpy//vWvvPbaa8ycOTNv5TniiCN48cUXAQgGg7z11lvh2fztjQUUY0yH1atXL2bPns3999+PqlJfX8/kyZMZOXIkY8aM4e233wbg0Ucf5frrrw+/7rzzzmPevHnh5//6r//KiBEj+MY3vsG2bdtaXGfx4sWMHz+esWPHctZZZ7Fp06bwvksuuYRnnnkGcNLwn3TSSZSWOmkW16xZw1FHHRU+9q677grnH1uyZAnHH388o0aN4sILL2TXrl2As9DXz372M8aNG8eQIUNyOnvekkMaY/Ji5svLWbGxtsX2FZuit9U1OIkdR854PWr78L4HtHjt8EMO4LbzR6RVjsGDBxMIBNi6dSt/+tOfEBE+/vhjPvvsM84880w+//zzpK/ft28f1dXV3HPPPfzqV79i5syZ3H///eH9TU1N3HDDDbz44ov07NmTZ555hl/84hfMmTMHgCFDhvDSSy+xa9cunn76aS677LKo/F+JXHHFFdx3332MHz+e6dOnM3PmTO69917AWSRs4cKFvPrqq8ycOZM33ngjrd9Jpgo+oIjIMODHwMHAm6r6UJ6LZIwpUvPnz+eGG24AYOjQoQwYMKDVgFJSUsLEiRMBuOyyy7joooui9q9cuZJPPvmEb37zmwAEAgH69u0bdcxFF13En//8Zz744ANmzZrVajlramrYvXs348ePB+DKK6+MmoQZKsPYsWPDqftzIS8BRUTmAOcBW1X1qIjtZwO/A3zAf6nqHar6KXCtiJQAjwMWUIwpAqnWJEIjvJ655oSslGP16tX4fD569eqV8JjS0tKojvt46ehDRCTquaoyYsQIFixIPFJt4sSJjB07liuvvDIq1Us6140USt3v8/myuqRxrHz1oTwKnB25QUR8wAPAt4DhwKUiMtzdNwGYC7ya22IaY4rZtm3buPbaa7n++usREU455RSefPJJAD7//HPWrVvHkUceycCBA1myZAnBYJD169ezcOHC8DmCwWB4NNdTTz3FySefHHWNI488km3btoUDSlNTU9TCXgADBgzg9ttv57rrrova3rt3b7Zu3cqOHTtoaGjglVdeAaBbt25079493D/yxBNPhGsr+ZSXGoqqviMiA2M2jwNWqepqABH5M3ABsEJVXwJeEpG5wFO5LGtW+BthzyZY8SIsfx6GXwAn3AC+gm+BNKbdC60G2dTURGlpKZdffjk33XQTANdddx1Tp05l5MiRlJaW8uijj1JRUcFJJ53EoEGDGD58OMOGDeOYY44Jn69Lly4sXLiQf//3f6dXr17hDvaQ8vJynnvuOX70ox9RU1OD3+/nxhtvZMSI6BraNddc06KsZWVlTJ8+nXHjxnHooYcydOjQ8L7HHnuMa6+9lrq6OgYPHswjjzzi5a8pI3lLX+8GlFdCTV4i8h3gbFX9gfv8cuA44DngIqACWKaqDyQ43xRgCkD//v3Hrl27Ntu3kJ7fj4Vdq0ETjHc/6zcwbooFFVPUMklfn+0mr44oW+nrC/7TS1XnAfNSOG42MBuc9VCyW6o0Ne6DnauSH/P6Lc6/g4+EQAP8eGluymZMgbNA0n4UUkD5Cjgs4nk/d1vKROR84PwjjjjCy3Jlbvn/wH9fmd5rtq90fu7bDl2SryhnjDGFpJAmNn4IfE1EBolIOXAJ8FI6JyiY9VAeOA7+Y2DyYNJzaOJ9AHe3st8YYwpMXgKKiDwNLACOFJENInKVqvqB64HXgU+BZ1V1ebLzFKSGPbDtM9i/K/7+6TthRg388AOo6gudE9RCgk3w8KnZK6cxxngsX6O8Lk2w/VXaMDS4IJq8Hjwx8b5r3oUSX/Pzn3zm/Ny9HnxlcPeR0cdvXgoPnQxT53tfTmOM8VghNXm1WUE0edWsa7lt6vtOraTvqPivOfAwqOrjHFN1SPS+LR/D6nmeF9OYduORc51/puAVVUARkfNFZHZNTU3uL77xH/DrOBlCb/gIeqeRW+gnn0JJTMXx8Qtgz5a2lc8YwwsvvMDo0aOj/pWUlPDaa6/xwAMPRG0/6qijEBE+/fTTqHOsWbOGp57KbDpcbW0t/fr1i0o0WUyKKqDkrYaiCrNPg8a9zdu+9Z9w+QvQ4/D0zzf1/Zbb5pzdcpsxJi0XXnghS5YsCf+77rrrOOWUUzjrrLP44Q9/GLVvwoQJfO9732sxX6MtAeXWW2/l1FOLt2+0qAJK3sw8sOW26u/D4V/P7Hw9j4SfrQVfefO2XaszO5cxJq7PP/+cX/3qVzzxxBNR+bMA3nnnHZ599lkefPDBFq+bNm0a7777LqNHj+aee+5JmPI+1uLFi9myZQtnnnlm1PauXbty8803M2LECM444wwWLlzIaaedxuDBg3npJWega2vp8wtFIc1DKR7Xvud0srdFpwPh+kXwu4h+l52r4aDBbTuvMYXitWmw+eOW2zcvi37euM/5+ZvDorf3idMn2WckfOuOVi/d1NTEd7/7Xe6++2769+8ftW/37t1MmjSJJ554ggMOaJki/4477uCuu+4K59W6++6746a8r6ysDL8mGAzyk5/8hD/96U8tUsnv27ePr3/969x5551ceOGF/PKXv+Svf/0rK1as4Morr2TChAmt3k+hKKoaSs77UAJ+mBHTvOarSK/PJJnuA+CmiPbb34/x5rzGdHC33norI0aMCKedj3Tttddy+eWXc9JJJ6V0rvnz53PZZZcBiVPeP/jgg5xzzjn069evxevLy8s5+2ynSXvkyJGMHz+esrIyRo4cmdPU814oqhqKqr4MvFxdXX11Ti44+7To5zOyEMgOiBn1NaMbTN8FJUX1XcB0RCnUJIDmEV6T53py2Xnz5vGXv/yFjz76qMW+xx57jLVr1/KnP/3Jk2uFLFiwgHfffZcHH3yQvXv30tjYSNeuXbnjjjsoKysLp7wvKSkJp54vKSkJp57PNI19rtmnUqYa65whvSHiS3xsW/UeGf383ywlizGZ2LVrF5MnT+bxxx+nqqoqat/q1av5+c9/zpNPPhlegjeeqqoq9uzZE36eKOV9pCeffJJ169axZs0a7rrrLq644gruuCPFgApJ0+cXkqKqoeR0YuOvo1dc4xebs3etqfPh98fAzi+c5xpwmtssM7ExaXn44YfZunUrU6dOjdp+yy238NZbb1FXV9dixcX77ruPU045Jfx81KhR+Hw+jj76aCZNmpQw5b2XkqXPLyR5S1+fTdXV1bpo0aLsXeD2vtBU1/z8l9ugtDzx8V6J7a/JRhObMVmUSfp6r5u8TAdOX1+QIoMJ5CaYAJSUOTm+Qv6tF/zsSyjvkpvrG5MPFkjaDetDSdfsb0Q/P2Rs7q49fTv8JGL0SKABfn2IzaI3xhQECyjp2hjRlOargClv5fb68dZIuXuIM1vfGGPyqKgCStbnoQRjlu+dFicRZLaV+Jxhw7H+eFbuy2JMBoqx37Y9yebvv6gCStZzef2qe/PjQ8ZCWWXiY7Mp3hyUDR/AHQMgGMh9eYxJUWVlJTt27LCgkieqyo4dO6Jm8XvJOuVT9cDx0c9z3dQVa/qu6AAHUL8bHj4FrouTXNKYAtCvXz82bNjAtm3b8l2UDquysjLujH0vWEBJxQPHw7aIFCileaqZRCopcVLj3xczHn3rcvjDGXD1G/FfZ0welZWVMWjQoHwXw2RJUTV5ZU1kMOk7Gn5ZIKOqug+Mv/2rD+FXB7cckWaMMVlUVAElJ8khr85zU1ekEh8ccgxx38ZgkzMibccXOS+WMaZjKqqAkpVO+djOw5Is5uzKxJS3YcauxM1wT1wUf7sxxnisqAJKVvwhw0Wycq3XUfG3714DD56Q06IYYzomCyit2dgyxXVBmvJm4lrKnk25LYsxpkOygJLMnHPyXYL0JBossH+Xs9pdKMmeMcZkgQWUZPZ8le8SpC9RLaWhFta+B/95BMw6LadFMsZ0DBZQEqnbCbvWND8/eCj0GZ234qTsl1uSpLVXqNsGm/5hQcUY4zmb2JjIozHNQ9e9X3gjvJI57ARYvyDx/t1rc1cWY0yHYDWURLauaH5c3rV9BROAq/43+QJc/obclcUY0yEUVUDJ2sTGHy/19ny5JAne4qZ9zsqTfzw7t+UxxhStogooWcs2XJml7MW5cOCAxPua6mDj4tyVxRhT1IoqoHjm14dGP/eV5accXqg6BA49FkoSdJcFGuHf+zjLCRtjTBtYp3wsVWjc2/y88sD8lcUL33/V+fnIubBuAWic9VL8+52fc84BEVvD2xiTEQsosdZ/EP2859D8lMNroSARm4o/0rr3nJ+PnAsahIY9cM078Rf0MsaYGPZJEevxf4l+ftXr+SlHtnTu0foxGxY66Vq2fAxfFFB2ZWNMQbOAEitek1AxmTwXBpwMvorExwQaYdeXzuM9G3NTLmNMu2cBJVagMd8lyL7Jc6HfsakdO/93ThOY5QEzxrTCAkqk2DkZnVJoHmqvJs+F3qNaP27nF/DVYqcZrG5n9stljGm3LKBEatwX/fzgIfkpR65UHgBlXVo5SJ1RYIFGeOLCnBTLGNM+FXxAEZF/EZE/iMgzInJmVi+2ZVnz437HOelLitnkuXDIGKdPJRWbllrTlzEmobwEFBGZIyJbReSTmO1ni8hKEVklItMAVPV/VPVq4FpgYo5KCD/4v9xcKt8mz3X+tVpTAVBYOx9mf936VYwxLeSrhvIoENVhISI+4AHgW8Bw4FIRGR5xyC/d/dkR+eFY1TdrlylYh4yBigNSO3b3OiewrJ2f3TIZY9qVvAQUVX0HiO3hHQesUtXVqtoI/Bm4QBz/Abymqtlbj3ft+82PE6UpKWaT58It61M7tm5b82NbCdIY4yqkPpRDgchPtA3uthuAM4DviMi1iV4sIlNEZJGILNq2bVuiw5IINj/siDWUkFRrKSEahKb92SmLMaZdKfiv4qr6e+D3KRw3G5gNUF1drWldJPYb9veeSevlReWW9c7vY8PC1ObkNO51MhbPPAjKu0Afdyiy5QMzpsMppBrKV8BhEc/7udtSlvF6KJuWpHd8sZs8F/qNS+81GnDWrd/+efTSycaYDqOQAsqHwNdEZJCIlAOXAC+lc4KM10OJzC4M6Tf7FKPJc6G0k/tEUn/dvq1Qu8HpWzHGdCj5Gjb8NLAAOFJENojIVarqB64HXgc+BZ5V1eVpnrftKzb6ytvfcr/Z0nMolFZC/xOhrHN6r23Y2/oxxpiikpc+FFW9NMH2V4FX23Del4GXq6urr870HPQdnfFLi055Fzi02qmtPHJumsOEg9F9U/H6VEL7rb/FmKJQ8J3yOZdoDfaOKPKDfvJcmJFmU+K6BU7Nps+o+MEj6Aex2qAxxaKoPj0zavL6r5hsLt/twCO8WpNusNUANO5xFu5aO98JIKEZ9gE/rP97c5p8Y0y7V1QBJaNO+c3LWj/GOA47Hsq7Zv769R80z7B/9Dxnm623YkzRKKqAkhF/zKQ8X3l+ytEeSEkb+5gipgdtiUjjZjPtjSkKRRVQ2jzKq8/RUJ7maKaOJJRIMqwNfz5Bf5uLY4wpLEUVUNJu8tr4j+jn6Q6N7aj6jIZu/WHAiZl3qkfWDLeugBndndn2YJmMjWmnOvYorx1fRD8v9vVPvFLR1fkXqq08ci6oOp3vmdjv5glV91yblzWncLGhxca0G0VVQ2mT8KxwkxER8FW0/TyBRmflzI3/sFqKMe1M2jUUEekOHKaqBTc8SkTOB84/4ogjUnvB/HubH/c+KitlKkqxtYXQ81/3c7IMNNVlfu7Ny5zhxk37muexBBqhYQ9UVGV+XmNM1qUUUERkHjDBPX4xsFVE3lPVm7JYtrSlPVN+Z0ST19VvZKdQHUnfo5109puWZB5U/PXNj1WdYBJogDuPcGbth1gTmDEFJ9Umr26qWgtcBDyuqsfhrFHSvtlII+9JCRxyjEeLlAWdYAJOcNm8zOYNGVPAUg0opSLSF/h/wCtZLE9upbLeh0ldi2HFHgo0OO+XBmHT0ub+FRsRZkzBSDWgzMTJArxKVT8UkcHAP7NXrMykNQ8l8kPooBT7XExqJs/NzhIA/npnqYHGPU7/ijGmoKQaUDap6ihVvQ5AVVcDv81esTKT1jyUyKaT0uSjkybOWsDEWfYBVlDUzWa8aanTHGY1FWPyLtWAcl+K29qPqPU60lhAyhQIt0+lcU/04Ip0A4sFImM8k7TnVEROAE4EeopI5IiuA4B2nnc8Iq/UBcmXrK9r9FNR2s5vN9e69oaGGggGYnaUAEFvrtFQ6/zcuxl6xDRb2oRIY3KutaE45UBX97jISQC1wHeyVajc0NYPAfY2+Pn4q1oO7mpJI9PSuQd0Oh7qd0GgCXasAtTJlRa75LIX1s4HSpwZ/I37nG39T4ieeW/BxZisShpQVPVvwN9E5FFVXZujMhWURr/zbXrnvsxGhIX6Xp655gTPytRuiEAnNz/Xvu1OcOk5DDZ+5Exe9FzQqRGpOtdORaAJgk1ZKIsxHU+qkwUqRGQ2MDDyNar69WwUKlMpz5SPbTPvfHDcwybOWkBTIHnzTIcOGOnocbiTs+sHf3U+7H9zWHOTlZea9gHiVEA3LHTSwcQLLqG/gW0rnFn4xpg2SzWg/DfwMPBfQDa+Wnoi5ZnysZPjug9IeOiSdbsBCKoTPJ655oSoEV8rNtUyvG/bhsgWZVCKTBwJUFoJVYc0f7j3GQVr3yPVpsf0qPMv0Oj8K69yRoXFEwom1udiTJulGlD8qvpQVkuSS6E29las2FQb1X28YlMtE2ctYMWmWuoa/HSuKKWuwc+KTYm/aTcFgq3WcopabGCJ3J6tWkqsxj0tF/RKdca9BRpjUpZqQHlZRK4DXgAaQhtVdWdWSpVtmvhbcWRtIRhzXF1Dy1QtAYW99X6aAkHKfC1HYa/cvId9jZlV6oqy5hKrrHPbkkmmKjJwadCtmSRpCrMAYkzaUg0oV7o/b47YpsBgb4uTKy1rDPEmLtY1RAeCgBJVGwl12Ctw7u/e5f9uGt/iHJkGk5BQrajdB5VEH9AlpTgf7Nlo+orjkXPdhdW05TVDtZbQWizGmLSkFFBUdVC2C5I3CVKETJy1IO5H3J765lpKwN8cmAJxaj2JZteHms0Ahvc9gLpGP1/trmfwtLmIwBe/aW4eaq1JrV27Zb3z89f93KHE7u9QfNkZBRZq6opcLTJk22ctm99iF/uKdz6w2owxrlTT118Rb7uqPu5tcfIg9KEGUX0jbe1oTybyOis21bK/MYA/6HyYtjZ9siibwUScLMWhIFLeJTt9K+sWJG7urNsef3vDHmfwQGQfkAUQY+JKtcnr2IjHlcA3gI+A9h9QYqhCXWOA5RtTSDAZYePu/Yyc8Trg1DqeueaEuDWLibMWUNfgJ6CEg0rsR1yoBrOv3k8QWvTlhI6BIgssSMu8ar5y77JCp1rrWTsf/niW2xyXg6Y4q+mYIpFqk9cNkc9F5EDgz1kpURukNA8l4pumApe4/RP1TQHqGgNO70pQ2duQXpPL/qYgNAXxRfTzRnbih4INOH0xkT/jqWvwh3t69qVZlnanzyg3L9c+J4Dcst55n9YtcIYbBwNZmgiZRMNe6HRg68dtXuaMVgObjW86vExXQdoHFFy/SkrzUBIMF/10Uy2BYNu/jYaCxMRZC6ICRqhWEisY1KjrBhQ++LLl4LnY/pii6leZPLc5gEQq79Lcf7F2fm7LtHUF9IlZFjr2b6d2o1N78pU7wXDzMu9qG1ZrMe1Qqn0oL9Nc9/cBw4Bns1WorIpoPmqkOT+Xlw0bi9bspHNF9K82UW0k1VFgoQCyt96PAF0qvVgRsYCEgkohiQwga98Hgs4wZ3DKumu189jn/h2FgkqfUakHhNY6/mOPTeWcxuRJqp9Kd0U89gNrVXVDFsqTA82f7Dt88VOutPkKWWx2Dw129QeC+Eos7X72xL6JbgNkNubMBP0QsOWoTfuXah/K30SkN82d8wW3WmPKIjLdTvNPQd1Pf3+yDo00BYk/CbKtIs+5vylIRWlJcXXOR37zjv0WnqtZ9WFp/D2oNn+LiKylQPwMAZHbg34nSO1clXkNxGoupkCktMCWiPw/YCFwMc668h+ISDtPXw/zG4eweruThsXrSoWH8Qlw5r/EnjOoyopNtcXVn5JIn1HunKES4v7ZSg7Xq/nNYdH9PY17cL5GBFsOHli3IPlyxaGUMPvjJJ2I7JMxph1ItcnrF8CxqroVQER6Am8Az2WrYNm2RbsRVNi+16MhqXngwRiC9kek5TwV8WVv7kommvY7NY+SNvRz+evd8xT5CD9TVFL9iy8JBRPXDlJfPrhgLP/1yYxwHz8fODW83daLbwdCzTmhIbqR+rvNfbkaCdZa4Ao2wd6tsG+bU2MRX3MHf2RiymTnqdngnGf/ruZtqrB5aWaTLK1ZzORAqgHlf0XkdeBp9/lE4NXsFCl7BjatDj/eoD3DjxetaZ85LgNBjepXie1PCc19+XjGWbkvXEe3a3XLbY373Nn6rdQ6HjkX9m5xn0RUQ3eucvoANy2FvkdnVq50Rp+lcpxX1zNFIWktQ0SOEJGTVPVmYBYwyv23AJidg/KZVqh2sKavPqPiD7GdPJeCrTSrOhmONb18yDsAABllSURBVJh8CGDDXtizOfH+cJCJsX8n1O1oWxmN8UBrNZR7gVsAVPV54HkAERnp7js/q6VzrjUYpw+nm6p6NhDgveCI8GOvO9BzKTSb/uKH3uezLd6tPNja6LG8jS4LfdP9z8Ob82+FtlV0LZx+lCjBiNGFCf7YSiud5iwUqvo0H7f985brt8Su57N1hXdFtUEApg1a+0rXW1U/jt3obhuY6UVFZI6IbBWRT2K2ny0iK0VklYhMc6+1WlWvyvRaiazRvl6fMq/acUzMTM+h+S6Bt4J+4r+Lmni1SWMKTGs1lGTJjDq14bqPAvcTkVxSRHzAA8A3gQ3AhyLykqp68vVr4qwFPKlx0pZnIIerd6Rs0Vqn89YnzX0nkan2Y4USUIYSWcZqt0ORy7q468q3M5EJMGOHGbdIjhlszn0Wb5ACRM/Aj+y/2LzMqcX95rDmnGnQ8hjIfF0Y6zfpsFqroSwSkRZ5sUTkB8DiTC+qqu8AsT3h44BVbo2kESf55AWZXiOer4I9PDmPzU8vUH1GQY+v5bsUueGvT6/mUrsxehlkY7KgtRrKjcALIvI9mgNINVAOXOhxWQ4F1kc83wAcJyI9gNuBMSJyi6r+Jt6LRWQKMAWgf//+cS8woGSbpwU2Bai8S0zK+xLirdDZ7qWb0j806izBgnLGeCFpQFHVLcCJInI6EEq9OldV38p6yZrLsAO4NoXjZuOOPKuurm7RItWWJhxfiRAMariZq1OFr2BTygfUaeqKTKN/1G3/y/7GANUDDwo3b63YVMueej+L1uyMSq0fEhqO3K6WH46cq+LVGir5EDus+KtFiY+NHIQQmUY/1CTW4vi9zT8jO+BDj0OvK61wOvub6qGsMn4TWipNW7vXQUkOsxiYvEo1l9fbwNtZLstXQGSDcD93W8pSWg8F2K1d0i5ce+ZFWv6C5KuAAxP0IQDWOBlPin8L/gY3x9hq6D0888vVrMv8tabdKaQc6B8CXxORQTiB5BLgu+mcINl6KHOYEX7834HxaRWsa4WPqsoyNuxyOvWH9Kpi2YbdBT3cOLJs+5ucJp9QbWR43wPCNZDQypGJXh+7cmWoMz8QVA7v2TU7hU9FsnT3fUY5HdvlXaD3UbDu/dyWLV8aap1Z+ZEZA9bOj+m4j/mjjeyAj5x0KRHdq5HDlkMLiqWTbj+VYyJrOfnq1LfBBG2Wl5lgIvI0zuTII0Vkg4hcpap+4HrgdeBT4FlVXZ7mec8Xkdk1NcmX7/2+77UMS+7oSGnj/QlqN3WNAT5tFyPBJOanMSZb8hJQVPVSVe2rqmWq2k9V/+huf1VVh6jq4ap6ewbnfVlVp3Tr1i3pcT7JvGohwHNTT6R64EEZn6M9STaxO1GwKSgiUF4F5XmsTRnTQRRorgrvbQmmsD54EsX2/bY9xALPiDj/TAQ3FUyijMaZTqZsrINtK+MPCDBFr5D6UNosWaf8B4EhnFfitKX/zn9RaufDaXEe1e9A7vzO0Zz8H29RUdZ+Y3BoBBhEt6S31hd0+C1zCShUVZYyvG/0sNPIFCzppGNJN3VLKsc3H9PcvzKx8ZdM15sZ0bcbbPxH+5z0mI54iScTpaMJ9S0FmxJnan7I/X2HlgiA5kXEQkJ9NLesd0ak1blZlntFdOY/cq7TNxF5bKxHzm2e1FnepWVfTWS/WbaSVobWoMn0/B28H6b9fjrGkazJq4Gy8ON/BJOPAguJ/E7bp1slncp9lNg33farxJfbhbg6olAbqa3j0iEVVUBJJBBUnktzZFcxC6ax6H1joAgnBZr0qTo1nUBTxLagO2O/I7WfmmQ6RJPXqq17OUya1weLrK3E8gl0riilrsFPeWlJeMhtrGeuOYFht77WYn+JOKPAmiLakbpW+NhbQBMh96Wx3r0qhEpe1+BnxabaqGHGoXxgE2ctYNGanXSuaP6TGjnjdeoa/FQPPIhFa3YSUDhukDOY4YMvd4aPyXi9lmTNCqF9kYunuc0ngbXvUeI2+omURPQXFOmsek+4vxe/mw8vshkt0ASo00QVdP82GmpaNqPNcFsOKg5oOQkzNDQ5tCAZtGxag+ZcZDMPchZWi83CHFpsLV4es3Sbo2KPj5z8CS2b7SInf7ZFLpvNvCqzq6hqKImavDRm7P0O9Sb9RKmv5a+vS0UplWXRzSoS0UzmK4AWM/tCGfEm2OgvYzxTVAHFGGNM/nSYgHJuyd9bPWbNHeeG55d0riil/0Gdo/YP73tA1Cin4X0PiDsaNfIYnzjPqypLw+eNlI83YGNNfUavU4W99f52N+T4Vz3ujGo+qJdOBChhRfnI6AN7j8Dk2OaPnX6YSBpwUr/k2u51EQuhmUwUVUBJNlP+x6UvZOWancp8DD64feUGq9nf1PpBcWjEv6JVGlrmp6j+axSuhpr4iTwDeQgoNetsMbM2KqpO+US5vAIPnISUZP6HsmzD7oT7fCVCz6oKVm8v8vkNxA8ke+r94Q720POB06I7E0Md8kDUsbHH+4TwfJfQvkjHDTooak5MaH/s8aHzhQYArNhUy8RZC/jpJufaR6LgzsmpLfFThdOr8uWOffQKlFDhBhMf0Yupxev+Srav0CmtlzuVY8K/hYZaEgbiyA76ePNi4s2D0UD0sTO6NXfYa6C5g158zfNvYvOYzejW4j3aN/MQ1pQNZkTTCqdzv7wL+xoDhL8Wbl7WPIAAnGNu2xk9R6ahtrk8kQMBQiLn1MQOFIi8TkOtM0ghNDghtH/dguh5OPHm4MTLlRY6RyhrdP8TWu/cjzf4IUNFFVC8EDlx7m+fb+PKOQvj7kvVF79xJ9hFrJAY+aFaUZZ4JBlAmS96xJjJ3PfdBKFzmMGRrAFgJQMZw2eUEuSRbtfxs523UUclAFXUtXrOgPsBWmqjw4zpGPV6QSmXwhm2azKzN43hziH+YJBGfzCtuTcdRSq1qoKreXn0PpZpY+LmrUybvVSd5rt89P8UiKIKKIn6UKok+pvmFu3u+bUP7FRGaUQW4u6dyziwU+L5LvGUlgi+Eok6j2m2fGPi7Maa4INmc009Df4gu+sS9xs1UkqAErb4+ra5jCbbvKkJ9vOvBzR+gGran+FZ1RlgEFodswMqqoCSaB6KL6L1/xn/adTg/dyDIb27UhmR5+trvboypHf861QmyAdWWVZCp7KShPtN+hr8zgdQskXG/JTyEUey29cj6bkU5/Mn3mdQvLNbpahwddYkzZnx8qGlIvSGd+DEmB2uD8XL/+ORCRFFhBGHdAv3s0yMnKFNdP/Lmb/9G59vdYYn/vHKY/ndm/8M96+EhPpZDu5awaYMh/ma+L7PDPY0+MNZESIXXwtZyUBK8TOGz6O2x/aZhPpburA/6ouLMR2RfRU2xhjjiaKvoeht3RC87zMxqUt1kFrouHhLEsceEzmkOHZ4ccgHX+6M6lSua/BHlSW8/HFFzAsVAqoI2vpXrvY8bjgrcjjaLTSkOM0mqs4a2SQVhIZaOic82j3/zIPiXye0LXL5ZH999PNEw4VDGvc1D3deO98ZRgzOMWvnO0OTQ0OE1y1whiu7xwQa9lIvnehS7ms+Z+S1I4cfh4YPz3QXB7xtZ3P+NI+a6ayGYjq0bPZzvB4Ym72Tm/yJndkP3v4h5SqD874dzQMQPJrQWVQBJdU15duDQw+s5LDunaxjt0A94T+DnzZdk/SYf/X/MOq5vZdFIt7Mfs/Pn4Oa3ooXnMXVAK96l4uqySvRTPlsStQJn+45Yl/fr7tTCd9sHfI59f0+f2HRGmdARHdqWVRxLUGEz8qP4rDGL1ihAxgua7nV/30AXgg66+wsq/gBXSImQqpCPeW5vwFjUhVIf15Xa4qqhmKMMSZ/LKAYY4zxhAUUY4wxnrCAYowxxhMWUIwxxniiqEZ5JfKZ9qcvu6K29ehSzo59jVHrXSTW+qy12NT2yVLd33z2UK5+fFHS45ONGisR2t2qiekIpDCfMFW+EsEfVLbU1scdthsEvtd0KyJCr9376VVVQUCddzw0ja1By7mx8na+2LO3xe9dnaVVuLjqCV7d8+3w9qU6uMW1bm6awl3lsz26M5Oq2P+97wVGcLJvedS2L4J9OKJkc+KTJBvzXV9Li0+RyHkdgUYI+p0sxIHGludqbQ5IaP/erdGvDTQRTnAZmeE48pgcj1UvqoAiIucD5x9xxBHhbQ1axl46tTi2Z1VFOKCU+iTu731YnyoqSkvoc0Bl1sr8td5Vab+mvLSEhqZgUWeO8moUvt+NAPsaE8+mbgw461Cu37U/nCJfgb105otgXx70X8A/t8ZfGja0guWn2/ZDxJ/JKu2HInwcHMgurWKAbOHN4DHe3JRpk9hgAiQPJq3Z8jEtv5ZGPPfXw77tziqUtRvjzGNp7X+yu39HdF45/PvxAWU0ucElJJjgcfYVVUCJNw9FgRptXqL3uJJP8Ql0c1PLi0BlmS/u+XodUMnoww7MWnkP7FxGz6rYvB+tK/OV0OQPtkhp4nO/itl6XJnzR/zymijlG413Jz1ewz9LGFj/FGsqvxuxVzi/8dcpXXeHVtFD9qRZWpM7yf5TpfAfLlTL0EBqx6dBCuirpfWhGGOM8YQFFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT1hAMcYY44miGjacjngp4wtJaKLjI+99ycyXV4S3dy73hdee/2jdLpoihrl2riilzl0r3YYOt1PlXZ1Fj9JchdCYQmA1FGOMMZ4o+BqKiHQBHgQagXmq+mSei2SMMSaOvNRQRGSOiGwVkU9itp8tIitFZJWITHM3XwQ8p6pXAxNyXlhjjDEpyVeT16PA2ZEbRMQHPAB8CxgOXCoiw4F+wHr3MGtYNsaYApWXgKKq7wA7YzaPA1ap6mpVbQT+DFwAbMAJKpBueQN+OpU0cVnpm83XTiFzsOm4vBrLENTU/84CaRxr8sWjvwwNxkkOmYKa9dHPI5JBltHK2vBNdbDXTX65f7cz4CMLOcWgsDrlD6W5JgJOIDkUeB74tog8BLyc6MUiMkVEFonIom3btgFQv7FlVtHV2tfLMmfdGcN6p3TcET27tH6QaZU/0LbsrB8EhwLwWODMuPtvbZrEm4ExfBgcwv8ETgTg5qZr4hxpQaYo1ddk9rrda6Of+/en/tqgH3asch5/+Ifm7bUbMytLEgXfKa+q+4DJKRw3G5gNUF1d7YTeiHUGDq9/ggA+7ip7KDsFzZLDDupMVWUp++r9CRNRC9CjawVb9jQkOMLkysTG6Un3PxE4kycigs2NTdfThf38gpixJlKS1tDhJvVRJtltEf44OJAH/RfwUPnvsnqdopbj9UlaaEojEGWgkGooXwGHRTzv525LmYicLyKza2oy/BZgjDEmY4UUUD4EviYig0SkHLgEeCmdE6jqy6o6pVu3blkpoDHGmMTyNWz4aWABcKSIbBCRq1TVD1wPvA58Cjyrqi07QYwxxhSkvPShqOqlCba/Crya6XnjLQFsjDEmNwqpyavNrMnLGGPyp+BHeaUjlRrK0YceyBfXngs0J2BMJpVjMlU9oHvKx1aW+6hrDDCwR5dwmX7+wsc89cE6ykud7wUfzzgrfPzAaXO9LWwH0ffATny5fV9OryklvuYnNyyGqj7wyLmweZmzrc8oWDs/p2UyHjt+Krx2M1T1ht3r8l2arLEaijHGGE8UVUAxxhiTP0UVUGweijHG5E9RBRRr8jLGmPwpqoBijDEmfyygGGOM8URRBRTrQzHGmPwpqoBifSjGGJM/RRVQjDHG5I8FFGOMMZ6wgGKMMcYTRRVQrFPeGGPyp6gCinXKG2NM/hRVQDHGGJM/FlCMMcZ4wgKKMcYYTxRVQLFOeWOMyZ+iCijWKW+MMflTVAHFGGNM/lhAMcYY4wkLKMYYYzxhAcUYY4wnLKAYY4zxhAUUY4wxniiqgGLzUIwxJn+KKqDYPBRjjMmfogooxhhj8qc03wXIuQsfyncJ0vbxjLP45Ksazrtvftz9t543nMuOH5DjUhnPlHVuuW3y3OjnM6zWbQqf1VCMMcZ4wgKKMcYYT1hAMcYY4wkLKMYYYzxhAcUYY4wnLKAYY4zxhAUUY4wxnij4gCIig0XkjyLyXL7LYowxJrGsBhQRmSMiW0Xkk5jtZ4vIShFZJSLTkp1DVVer6lXZLKcxxpi2y/ZM+UeB+4HHQxtExAc8AHwT2AB8KCIvAT7gNzGv/76qbs1yGY0xxnhAVDW7FxAZCLyiqke5z08AZqjqWe7zWwBUNTaYxJ7nOVX9TpL9U4Ap7tMjgZVtLnxhOhjYnu9CZJHdX/tm99d+HamqVW05QT5yeR0KrI94vgE4LtHBItIDuB0YIyK3JAo8qjobmO1lQQuRiCxS1ep8lyNb7P7aN7u/9ktEFrX1HAWfHFJVdwDX5rscxhhjksvHKK+vgMMinvdztxljjGnH8hFQPgS+JiKDRKQcuAR4KQ/laK+KvVnP7q99s/trv9p8b1ntlBeRp4HTcDqytgC3qeofReQc4F6ckV1zVPX2rBXCGGNMTmR9lJcxxpiOoeBnyhtjjGkfLKAUOBFZIyIfi8iS0LA+ETlIRP4qIv90f3bPdzlTFS97QqL7Ecfv3YwKy0TkmPyVPDUJ7m+GiHzlvodL3Cbf0L5b3PtbKSJn5afUqRGRw0TkbRFZISLLReTH7vaieP+S3F+xvH+VIrJQRJa69zfT3T5IRD5w7+MZt28bEalwn69y9w9s9SKqav8K+B+wBjg4Ztt/AtPcx9OA/8h3OdO4n1OBY4BPWrsf4BzgNUCA44EP8l3+DO9vBvDTOMcOB5YCFcAg4AvAl+97SHJvfYFj3MdVwOfuPRTF+5fk/orl/ROgq/u4DPjAfV+eBS5xtz8MTHUfXwc87D6+BHimtWtYDaV9ugB4zH38GPAveSxLWlT1HWBnzOZE93MB8Lg6/g4cKCJ9c1PSzCS4v0QuAP6sqg2q+iWwChiXtcK1kapuUtWP3Md7gE9xJioXxfuX5P4SaW/vn6rqXvdpmftPga8DoeS7se9f6H19DviGiEiya1hAKXwK/J+ILHbTywD0VtVN7uPNQO/8FM0zie4nXlaFZP/BC9n1brPPnIgmynZ7f27zxxicb7lF9/7F3B8UyfsnIj4RWQJsBf6KU6varap+95DIewjfn7u/BuiR7PwWUArfyap6DPAt4IcicmrkTnXqo0UzVK/Y7sf1EHA4MBrYBNyd3+K0jYh0Bf4C3KiqtZH7iuH9i3N/RfP+qWpAVUfjTCgfBwz18vwWUAqcqn7l/twKvIDzR7Al1HTg/mzvGZkT3U9RZFVQ1S3uf+Qg8Aeam0Xa3f2JSBnOh+2Tqvq8u7lo3r9491dM71+Iqu4G3gZOwGmKDKXhiryH8P25+7sBO5Kd1wJKARORLiJSFXoMnAl8gpNZ4Er3sCuBF/NTQs8kup+XgCvc0ULHAzURTSvtRky/wYU47yE493eJO5pmEPA1YGGuy5cqt/38j8CnqvrbiF1F8f4lur8iev96isiB7uNOOEuIfIoTWEKZ3GPfv9D7+h3gLbcGmli+Rx7Yv6SjMgbjjCJZCiwHfuFu7wG8CfwTeAM4KN9lTeOensZpNmjCaa+9KtH94IxKeQCnnfdjoDrf5c/w/p5wy7/M/U/aN+L4X7j3txL4Vr7L38q9nYzTnLUMWOL+O6dY3r8k91cs798o4B/ufXwCTHe3D8YJhKuA/wYq3O2V7vNV7v7BrV3DZsobY4zxhDV5GWOM8YQFFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT1hAMcYY4wkLKMZkmYgMFJH9bg6l0LbeIvKUiKx287QtEJELk5zj7dj06CJyo4g8JCKd3LTqjSJycDbvxZhkLKAYkxtfqJNDKTQj+3+Ad1R1sKqOxUkP3i/J6592j4l0CfC0qu53z70xC+U2JmUWUIyJICJPu4sKLRSRtSJybhYu83WgUVUfDm1Q1bWqep9bhsvc6y8RkVki4sNJH35uxOJHA4FDgHezUD5jMmIBxZhoRwOrVXUc8D3gtixcYwTwUbwdIjIMmAic5NY6AsD3VHUnTvqLb7mHXgI8q5bqwhSQ0tYPMaZjEJFKoCcw0920AuguIpOB44CzgNeBf6jqLA+v+wBOHqlGnAWNxgIfumsZdaI5e2+o2etF9+dVXpXBGC9YQDGm2VHAP1W13n1+DLBUVR8RkReBMlW9NvZFIjIcmKCqd4jIfcDP1VnxL5HlwLdDT1T1h25n+iKchIqPqeotcV73InCPuzZ7Z1VdnMlNGpMt1uRlTLOjgf4iUukuFzATuMfdNxZI9AF+LE5mWoBurQQTgLeAShGZGrGts/vzTeA7ItILQEQOEpEBAOos3/o2MAentmJMQbGAYkyzo4HncZZ9/RB4SFXfc/e1FlBWuEGoVW6/x78A40XkSxFZiNPU9TNVXQH8EmfZ52U4y7RGrsfxtFtOCyim4Fj6emNcIvI3YIqqroyz72ng+6q6X0R6A+ep6h/dfXNx1j6pBUaq6tkxrx0IvKKqR2W5/Gtw1hzZns3rGJOI9aEY0+xwnEWiWlDVSyOejgG+hPCSsTtU9Zok5w0A3URkSWguipfc1fcWAGVA0OvzG5Mqq6EYY4zxhPWhGGOM8YQFFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT1hAMcYY4wkLKMYYYzxhAcUYY4wnLKAYY4zxxP8HlcncgrmwNj4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = hist.plot1d(output['pt_lead'], overlay='dataset')\n",
    "ax.set_yscale('log')\n",
    "ax.set_ylim(0.1, 5e3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.1, 20000.0)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAELCAYAAAD+9XA2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deZxU1Znw8d9TvULTICAqAgqIqCwtSIMxqLhEwQ0HEwMaFwiJWzRxnOQNGpcmoxNn1OhERSCJ4i6Mia/7+GqUqJHIkrAjhrBIA7LvDb1UPe8f91Z3VXVVd3X3rb5V1c/38+lP1T13O7cL6ulzz7nPEVXFGGOMaamA3xUwxhiTHSygGGOM8YQFFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT+T6XQEvichlwGXFxcU/7N+/v9/VMcaYjLFo0aIdqtqtJceQbHwOpbS0VBcuXOh3NYwxJmOIyCJVLW3JMeyWlzHGGE9YQDHGGOMJCyjGGGM8kVWd8saY9FZdXU15eTmHDx/2uyptVmFhIT179iQvL8/zY1tAMca0mvLycoqLi+nduzci4nd12hxVZefOnZSXl9OnTx/Pj2+3vIwxrebw4cN07drVgolPRISuXbumrIWYEQFFRIpEZKGIXOp3XYwxLWPBxF+p/P37ElBE5GkR2SYiy2PKx4jIahFZIyJTIlb9HJjTurU0xvglJyeHIUOGMGjQIK688koqKiribvfNb36zlWuW2DnnnENjz7899thjCa+luebOnctnn33m6TGby68WyixgTGSBiOQATwIXAQOAq0RkgIhcAKwEtrV2JY0x/mjXrh2LFy9m+fLl5OfnM3369Kj1NTU1ACn5Ig0fOxUsoKSAqn4M7IopHgGsUdW1qloFvAJcDpwDfAO4GvihiGTEbTpjjDfOOuss1qxZw9y5cznrrLMYO3YsAwYMAKBDhw6A86U6atQoLr/8cvr27cuUKVN48cUXGTFiBIMHD+af//wnAG+++Sann346Q4cO5Vvf+hZbt24FoKysjGuvvZaRI0dy7bXXcvbZZ7N48eLaOpx55pksWbIkql6HDh1iwoQJnHLKKYwbN45Dhw7Vrrv55pspLS1l4MCB3HfffQD85je/YfPmzZx77rmce+65CbcDmDJlCgMGDKCkpISf/vSnAGzfvp1vf/vbDB8+nOHDh/OXv/yF9evXM336dB599FGGDBnCJ5984unvvslU1ZcfoDewPGL5O8DvIpavBZ6IWJ4IXNrA8W4AFgILjzvuODXGpJ+VK1cmtV1RUZGqqlZXV+vYsWN12rRp+tFHH2n79u117dq19bb76KOPtFOnTrp582Y9fPiwHnvssXrvvfeqqupjjz2mP/nJT1RVddeuXRoKhVRV9be//a3ecccdqqp633336WmnnaYVFRWqqjpr1qzafVavXq3Dhg2rV8dHHnlEJ02apKqqS5Ys0ZycHF2wYIGqqu7cuVNVVWtqanTUqFG6ZMkSVVU9/vjjdfv27bXHiLfdjh07tH///rX13L17t6qqXnXVVfrJJ5+oquqGDRv05JNPrq37Qw89lNTvNSze5wAs1BZ+r2fMX/uqOktV32pg/UxVLVXV0m7dWpTfzBjjs0OHDjFkyBBKS0s57rjjmDx5MgAjRoxIONx1+PDhdO/enYKCAk444QQuvPBCAAYPHsz69esBZ9jy6NGjGTx4MA899BArVqyo3X/s2LG0a9cOgCuvvJK33nqL6upqnn76aSZOnFjvfB9//DHXXHMNACUlJZSUlNSumzNnDqeddhpDhw5lxYoVrFy5Mm6d423XqVMnCgsLmTx5Mn/84x9p3749AB988AG33norQ4YMYezYsezbt48DBw404beaeun0HMomoFfEck+3LGnhbMP9+vXzsl7GmFYW7kOJVVRUlHCfgoKC2veBQKB2ORAI1PaL3Hbbbdxxxx2MHTuWuXPnUlZWFvfY7du354ILLuD1119nzpw5LFq0KOm6r1u3jocffpgFCxbQuXNnJk6cGHeYbqLtcnNzmT9/Pn/605949dVXeeKJJ/jwww8JhUL89a9/pbCwMOm6tLZ0aqEsAE4UkT4ikg9MAN5oygFU9U1VvaFTp04pqaAxJrPt3buXHj16APDss882uO0PfvADfvzjHzN8+HA6d+5cb/3ZZ5/NSy+9BMDy5ctZunQpAPv27aOoqIhOnTqxdetW3n333dp9iouL2b9/f4PbHThwgL1793LxxRfz6KOP1vbdXHjhhTz++OO1xwoH3Mhj+s2vYcMvA/OAk0SkXEQmq2oNcCvwHrAKmKOqKxo6TpzjXiYiM/fu3et9pY0xGa+srIwrr7ySYcOGceSRRza47bBhw+jYsSOTJk2Ku/7mm2/mwIEDnHLKKdx7770MGzYMgFNPPZWhQ4dy8sknc/XVVzNy5MjafW644QbGjBnDueeem3C7/fv3c+mll1JSUsKZZ57Jr3/9a8Dp1F+4cCElJSUMGDCgduTbZZddxmuvvZYWnfI2H4oxptWsWrWKU045xe9qJGXz5s2cc845fPHFFwQC6XQzp+XifQ42H4oxxqTAc889x+mnn84DDzyQdcEkldKpU77FrFPeGOOF6667juuuu87vamScrAq91ilvjDH+yaqAYp3yxhjjn6wKKNZCMcYY/2RVQDHGZJ/xM+YxfsY8v6thkmABxRjTpoRT4w8cOJBTTz2VRx55hFAo1OzjhRNUxpo4cSKvvvpqg/tOnDiR9u3bRz2YePvttyMi7Nixo9l18ktWBRTrQzHGNCac1mXFihW8//77vPvuu0ydOtW3+vTr14/XX38dgFAoxIcfflj7NH+myaqAYn0oxpimOOqoo5g5cyZPPPEEqsrhw4eZNGkSgwcPZujQoXz00UcAzJo1i1tvvbV2v0svvZS5c+fWLv/rv/4rAwcO5Pzzz2f79u31zrNo0SJGjRrFsGHDGD16NFu2bKldN2HCBGbPng04afhHjhxJbq7zRMf69esZNGhQ7bYPP/xwbf6xxYsX841vfIOSkhLGjRvH7t27AWeir5///OeMGDGC/v37t+rT81n1HIoxJnNMfXMFKzfvq1e+ckt0WUWlk9hxcNl7UeUDunest++AYzty32UDm1SPvn37EgwG2bZtGy+88AIiwrJly/jiiy+48MIL+fLLLxvc/+DBg5SWlvLoo4/yy1/+kqlTp/LEE0/Urq+urua2227j9ddfp1u3bsyePZtf/OIXPP300wD079+fN954g927d/Pyyy9zzTXXROX/SuS6667j8ccfZ9SoUdx7771MnTqVxx57DHAmCZs/fz7vvPMOU6dO5YMPPmjS76S5LKAYY4zr008/5bbbbgPg5JNP5vjjj280oAQCAcaPHw/ANddcwxVXXBG1fvXq1SxfvpwLLrgAgGAwSPfu3aO2ueKKK3jllVf4/PPPmTFjRqP13Lt3L3v27GHUqFEAXH/99Vx55ZVRxwMnH1k4dX9ryKqAYk/KG5M5km1JhEd4zb7xjJTUY+3ateTk5HDUUUcl3CY3Nzeq4z5eOvowEYlaVlUGDhzIvHmJR6qNHz+eYcOGcf3110elemnKeSOFU/fn5OSkdErjWNaHYoxps7Zv385NN93Erbfeiohw1lln8eKLLwLw5Zdf8tVXX3HSSSfRu3dvFi9eTCgUYuPGjcyfP7/2GKFQqHY010svvcSZZ54ZdY6TTjqJ7du31waU6urqqIm9AI4//ngeeOABbrnllqjyo48+mm3btrFz504qKyt56y1njsFOnTrRuXPn2v6R559/vra14qesaqEYY0xjwrNBVldXk5uby7XXXssdd9wBwC233MLNN9/M4MGDyc3NZdasWRQUFDBy5Ej69OnDgAEDOOWUUzjttNNqj1dUVMT8+fO5//77Oeqoo2o72MPy8/N59dVX+fGPf8zevXupqanh9ttvZ+DA6BbajTfeWK+ueXl53HvvvYwYMYIePXpw8skn16579tlnuemmm6ioqKBv374888wzXv6amsXS1xtjWk1z0ten+pZXW5Sq9PXWQjHGpDULJJkjq/pQ7MFGY4zxT1YFFOuUN8YY/2RVQDHGGOMfCyjGGGM8YQHFGJPenrnE+TFpzwKKMabNeO211xgyZEjUTyAQ4N133+XJJ5+MKh80aBAiwqpVq6KOsX79el566aVmnX/fvn307NkzKtFkNrFhw8aYNmPcuHGMGzeudnnmzJm8+OKLjB49mkAgwI9+9KPadXfddRdDhgyp97xGOKBcffXVTT7/Pffcw9lnn938C0hzWdVCsWHDxphkffnll/zyl7/k+eefj8qfBfDxxx8zZ84cpk2bVm+/KVOm8MknnzBkyBAeffTRhCnvYy1atIitW7dy4YUXRpV36NCBn/3sZwwcOJBvfetbzJ8/n3POOYe+ffvyxhtvAI2nz08XWdVCUdU3gTdLS0t/6HddjDGNeHcKfL2sfvnXS6OXqw46r7/qFV1+TEn9fY8ZDBc92Oipq6urufrqq3nkkUc47rjjotbt2bOHiRMn8vzzz9OxY/0U+Q8++CAPP/xwbV6tRx55JG7K+8LCwtp9QqEQ//Zv/8YLL7xQL5X8wYMHOe+883jooYcYN24cd999N++//z4rV67k+uuvZ+zYsY1eT7rIqhaKMcYk45577mHgwIG1aecj3XTTTVx77bWMHDkyqWN9+umnXHPNNUDilPfTpk3j4osvpmfPnvX2z8/PZ8yYMQAMHjyYUaNGkZeXx+DBg1s19bwXsqqFYozJIEm0JIC6EV6T3vbktHPnzuUPf/gDf/vb3+qte/bZZ9mwYQMvvPCCJ+cKmzdvHp988gnTpk3jwIEDVFVV0aFDBx588EHy8vJqU94HAoHa1POBQKA29Xxz09i3Ngsoxpg2Y/fu3UyaNImXXnqJ4uLiqHVr167lrrvu4pNPPqmdgjee4uJi9u/fX7scTnl/3nnnRaW8jxROiQ9Of8jChQt58MEkAyrQu3dvpk2bRigUYtOmTVHp89OJBRRjTJsxffp0tm3bxs033xxVfuedd/Lhhx9SUVFRb8bFxx9/nLPOOqt2uaSkhJycHE499VQmTpyYMOW9lxpKn59OLH29MabVNCd9vde3vIylrzfGtFUWSDJG2gcUETkF+AlwJPAnVX3KkwM/cTrs+MJ5n9cehl4D6z6GMQ9C597QpY8npzHGmLbCl2HDIvK0iGwTkeUx5WNEZLWIrBGRKQCqukpVbwK+CyQ3jq8x075ZF0wAqitg/kzY/gU8/y/wmyGw+e+wY40npzPGmLbAr+dQZgFjIgtEJAd4ErgIGABcJSID3HVjgbeBd1p85sr9sG1F49vNPAeeGAZlneDwvhaf1hjjyMZ+20ySyt+/LwFFVT8GdsUUjwDWqOpaVa0CXgEud7d/Q1UvAr7X4pP/V9/o5SHXwLgZ0LFH4n0e7AV/fxE2/Q0OW1oXY5qrsLCQnTt3WlDxiaqyc+fOqKf4vZROfSg9gI0Ry+XA6SJyDnAFUEADLRQRuQG4AaiXSiFKsKru/f9ZB+27OO9PnQA1lVB5AB7qW3+/128JnwnK9jR6McaY+nr27El5eTnbt2/3uyptVmFhYdwn9r2QTgElLlWdC8xNYruZwExwhg3H3ejgjrr3UzZCYUyentwC5+fu7RDIgcUvwhu3xZ7JuQ12/ZvQJ3uzhhqTCnl5efTpYwNeslU65fLaBERmf+vpliWt0WzD8550NwzUDyaRcvOdgHLadXDvLugc5z/As5fBrwc1pXrGGJPV0imgLABOFJE+IpIPTADeaMoBVPVNVb2hU6dO8Tf49NfO6w8+iL8+nkAO/GQxlO2FI46PXrdvo9Na2VQ/J5AxxrQ1fg0bfhmYB5wkIuUiMllVa4BbgfeAVcAcVU1iOFbUcRO3UHZvqHvfY1jzKn77UvjZWjj3F9Hlvz0X/vtUsI5GY0wb1nZSr0ztAhqEQB7cuyP+jskKBeHfuznHi3XPDsjJa9nxjTGmlXmReiWdbnmlVrvOzuu/rW75sQI5cPc2yMmvv+7fj4QtS+JPHGSMMVksqwJKg7e8KtxWSUEHb06Wkwv3bIf79gASvW7G2TD9TG/OY4wxGSKrAkrCTvnNi+veBzy+HSXucyl3bY5XIW/PZYwxaSyrAkrCFsrvzndez7sbAim65Pwi5/mVSFOPgPuPSc35jDEmzWRVQEnYQgk502gy9NrUViA33+mUH/SdurKaQ/CffWC6PQRpjMluWRVQ4or8Im/XJfXny8mD8++JLju0C3ZZ5mJjTHbL/oCyNyI9WG6cUVmp0Lk3dO0fXVZ10Embb4wxWSqrAkrcPpRz73Jej27lNCm3LagfVLatsFT4xpislVUBJW4fykf/4bxe9XLrV+i2BXBE7+iyRbNavx7GGNMKsiqgxHXInXaluLs/5799SfTy+/fA2j/7UxdjjEmh7A4o+7fWvZcc/+pxbEzusOfGwn8PherD/tTHGGNSIKsCSr0+lP8+tW5lqp4/ScYNH8LP/hldtnstPHA0PGVP1BtjskNWBZR6fSg1h/ytUKSiI6HbyfXLd69t/boYY0wKZFVASSicGNJvP/ocfrw4uqzqIBzY5k99jDHGQ2k/BXCL5BRAsBKu/h+/a1KnSx8oOgoORgSRh0+EQC70GA6T/9e/uhljTAtkdwslWOm89jjN33rE+tk/6peFaqKDjDHGZJisCigJk0MGfBzhlUjnvvXLdv2zfpkxxmSIrAoo9TrlJQAFHf2tVCI/+Xv88rIjWrcexhjjkawKKPUUFMOQq/2uRWLdh8Ixp8YUqjO9sDHGZJjs7pRPdzfOdV7LOgOhuvJglR+1McaYFsnuFkqmKNsNuYXRZb+70J+6GGNMM1lASRd3b4X2Ebe6yj/3ry7GGNMMFlDSyXWvRS//9nx/6mGMMc2QVQEl4bDhTNH1xOjlTQuhYpc/dTHGmCbKqoCScE75TJFXCGUxwXBHnIcgjTEmDWVVQMlK7/zU7xoYY0xSLKCko9x2de9rbAixMSYzWEBJR92H1L3fs8G/ehhjTBNkd0A5vBdWvu53LZouMuNwOs3pYowxDcjegFJ10HkN1vhbDy88c4nfNTDGmEZlb0BRdV7zi/ytR3MVRiSJ3PApPH2Rf3UxxpgkpH1AEZF/EZHfishsEWl6PpLhk1NQq1YQO13wpoX+1MMYY5LkS0ARkadFZJuILI8pHyMiq0VkjYhMAVDV/6uqPwRuAsY39VwPv7/am0q3tsnvRS8Hq+D+Y+z2lzEmbfnVQpkFjIksEJEc4EngImAAcJWIDIjY5G53fZPUBLX5tfRbICYZdM0hqDzgT12MMaYRTQ4oItJZREpaclJV/RiIzSkyAlijqmtVtQp4BbhcHP8JvKuqf2vquYIhJRjK0KDSo7R+2b7y1q+HMcYkIamAIiJzRaSjiHQB/gb8VkR+7XFdegAbI5bL3bLbgG8B3xGRmxqo4w0islBEFm7fvt3jqvnkmlfrl1XsaP16GGNMEpJtoXRS1X3AFcBzqno6zpd8yqnqb1R1mKrepKrTG9hupqqWqmppt25ZMuNhXlH9215g/SjGmLSUbEDJFZHuwHeBt1JUl01Ar4jlnm5Z0hJlG1bN0FtegQD0+kb98i1LLKgYY9JOsgFlKvAeTh/HAhHpC3idBncBcKKI9BGRfGAC8EZTDpAo2/Bf12ZwCvhJb9cvq9oPh3a3fl2MMaYByQaULapaoqq3AKjqWqDZfSgi8jIwDzhJRMpFZLKq1gC34gSuVcAcVV3RxOPGbaF8ve9wc6uaHvKL65dtWwFTu1hLxRiTNpINKI8nWZYUVb1KVburap6q9lTV37vl76hqf1U9QVUfaMZx47ZQnpq7prlVTQ93lUNOfv1yDUKwuvXrY4wxccTp8a0jImcA3wS6icgdEas6AjmprJiX1m4/6HcVWi630Hm4MVb5504rJd6tMWOMaUWNtVDygQ44gac44mcf8J3UVq3pEnbK+1QfT925ESRBDP9qHvyqV/x1xhjTShpsoajqn4E/i8gsVU37iTlU9U3gzdLS0h/GWYeI+FArD923C8riTG+sQdBQXX+KtVaMMT5oMKBEKBCRmUDvyH1U9bxUVKq5ROQy4LJ+/frVWzd+xjxWfb2fAd07MvvGM1q/cqlWdQC+XgrHlFhgMcb4ItmA8j/AdOB3QDB11WmZhloo89fvJgCEMvWZlGRU7oMNnzkp+48e5HdtjDFtTLIBpUZVn0ppTVpBCFi4fjfjZ8wDyM6WCiGntbJxHtRUQm6BtViMMa0i2WHDb4rILSLSXUS6hH9SWrMUyfj2SfGxSWzkXuWzY+uKvl5qz6wYY1Iq2YByPfAz4DNgkfuTdjM+RY7yOlyT+M7c5+t28fm6XbUtlYzSuU/y225dATVVsOcrp9PeGGNSKKlbXqrahG8x/0T2oSSTsn7lln2pr5TXRCC/g3M7K9TIQ41V++E/j4fqCkAgFBFk7TaYMcZjSQUUEbkuXrmqPudtdVrX/sM1fHvaX/jDLSP9rkryJr3tBIPqCticxPQw1RXuG4XNi1JaNWNM25Zsp/zwiPeFwPk486KkfUBp7MmTypoMvRWUk+/8xHt6PpHIbYPVEMiYZAfGmAyQ7C2v2yKXReQInBkV00q851ACAWlwoPPyzfv47vTPmHPTN1NfQa+Eb1M9cwkc2AY7v0x+3/CtrvL50L6L3foyxnimuXPKHwTSrl8lXnLIXp3bNbpf5FDiNuHrpYBCxU6/a2KMySLJ9qG8Sd2I2xzgFGBOqirlpb7dOsC2hrcJ3/TKuOdTJr0NS2bDazckv8+GT4n6OyIUtFtfxhhPJNuH8nDE+xpgg6qWp6A+nks2fdfKLfsY0L1jaiuTUkLyT9lE9BttnAddToi/md0OM8Y0QVK3vNwkkV/gZBruDDShJ9h/yVzkwcM1LFy/K/OGEue5t/Ry8pp/jF3/rP/Q4zOXuLfGjDEmOUkFFBH5LjAfuBJnXvnPRSRj0tevffASCnIbvtQQEFRnKPF3p3+Wwlp6rP9oOKI3dB8KBR0BgdzG+43qCafAjwws2Zz3zBjjuWQ75X8BDFfV61X1OmAEcE/qqtU8iWZsBMjLST51/cot+zKnkz63ADr1rGuhSAB6DGv6cTToJJcMC1Y7D0bu3+xNPY0xWS/ZgBJQ1ciu7Z1N2DctZPxcKMk4psTJNNwS4Zxf5Quc5V3rWl4vY0ybkGxQ+F8ReU9EJorIROBt4J3UVct7y8pGJ73tgcogC9dnUK6vSW/X/dy5sWWd6JUHnNtftR33akkljTFJaTCgiEg/ERmpqj8DZgAl7s88YGYr1M83qrBw/S5OuPPtzAkskQLN7aQPObe/jDGmiRproTyGM388qvpHVb1DVe8AXnPXZazhx3ducH1kQpaM6lMJ61FK44lnkrThM2ulGGMa1VhAOVpVl8UWumW9U1KjFGqXV3e5z//gdB9r0goCOZBb6NHBQnV9K+GfsGcuqT86zBjTJjUWUI5oYF0zxqamVqJhw2G5OdGXW1zY8HOdQa0bShw2fsa8zGuteKFynxNU7NkUY0wCjQWUhSJSb352EfkBziRbaaWhYcPxDOjesdGgkkhGBJZA7LV5NDAv8pkVCzDGGFdj36a3A6+JyPeoCyClQD4wLpUVSzfz1zkd9O0L6n5laZ2qJTzS61e96p4vyW/vzDffElUHnU77qoNOMKk66JSHb4klGmFmaVyMyXoN/smqqltV9ZvAVGC9+zNVVc9Q1a9TXz1vxQ4dnn3jGSwrG03n9o2PiIp8ZryisoaKypqE24alRSvmmBLIK4KOPeDOFqZfq9wfPQJM1ZlaOPxQZLi1EtvPYoxpE5KdD+Uj4KMU16VVnHxMMbsOVkWlYunaoYDdFY1Mp+uqCYYIKQQy6TnJQI4zF33kw50FHaOfjE9KbCoWjVNmjGmrMuppdy90apdHnyOLop6cLy5Irh8lqHCoOuR8jSos27S3RS0QX1owHY52ZnpsKQ1CVUV0Wfj2F9TdAotlrRdjslbzeqSzTH4jiSPjCQEVVcFmZScOB5GVW/ZRUVnD+BnzUjcHS3gO+rCuJ9a93/BpCw8eM32yanSfzVcRwVJDTp4xY0zWsv/hbZrXH3+C21/lC53gcmi3x+czxqQTa6HgdM73ntL80UfjZ8yrbakM6N6xtrWxeOOeJvW1hI8TeQxPxBtZNent6NaEJzT6tpcGnXPUVDrrDu2Gdg1nKIhiI8OMyShpH1BEpC9O+vxOqtriOVi8vrUUnpgr6P5x/vm6XQwue48B3TtSWePcEhpc9h5QN8w4HDQOHq4hRP3ULuHllE9FfEyJB7e9GlG5j6gUMBqCrcudQHNMSV2wsOBhTMbz5ZaXiDwtIttEZHlM+RgRWS0ia0RkCoCqrlXVyamu0+l9ujC0V0OJAeILNb5JQr6MjwpnJU6VhiblqtoPX33mzLVijMk6frVQZgFPAM+FC0QkB3gSuAAoBxaIyBuqutKXGjZBMM536ML1u2rfV1TW0L4gt7YTPqiwYN2u2oBy8HANn6/bVf8grSJAy8JirHjHcq+0cr/zWnMIgu4s0uGWyddLndbMr3o5KfjDNnwKZZ2gLH46HWNM+vClhaKqHwOx36AjgDVui6QKeAW4PNljisgNIrJQRBZu377dw9qmRijB+1Y16W3oNcKvsxtjskw6jfLqAUT8aUo50ENEuorIdGCoiNyZaGdVnamqpapa2q1btyaffPaNZzDj2mZMnRvH/sM1Ua2W2AST6StAq/2T0KDT+ojtwwm3Uv6jZ92skVDXkrHnWIxJW2nfKa+qO4GbktlWRC4DLuvXr1+zztWxXXMnpcpgEgDEmZs+3LeRDg+/h2ogWOl3LYwxTZBOLZRNQK+I5Z5uWdKamm04VmFeTrP2S1YwpA32WUc6XN1KsyZKAAqKvXl6viWS/cUYY9JWOrVQFgAnikgfnEAyAbi6KQdoaQsl1Sqqkg8SS920LikfOpwOYm97hTMaR4pM5RKZ2Tj29lfkCLZ4Q5HDKfePKYneJ7Lchi4b0yx+DRt+GWde+pNEpFxEJqtqDXAr8B6wCpijqiuactyWtlAActIk6WOr/cE+6e26L9f8IjiuDQQwY0xK+NJCUY9Orl8AABLPSURBVNWrEpS/A7zTytVJW5+v2+VPKyW/KPqJ99YW2zoBp7M+3JIpcOeheeaSujJxb1f+qpdT9/wiJ1CGJwOLHIocKVESy1RqzkOc9uCnyQDp1IfSYo1NAdxUkXPQZ71jSpwv3Ulvw5GnQG7azfBsjElzWfWN6cUtLwPkFUJOGo94U9+e3DHGNCCdOuVbzItO+fYFubXPjPQ7qphlm9rAE9qN3UaRnPi3ofxSXdH4NsaYVpdVAUVV3wTeLC0t/WFLj1WUn0P7/NQOI05GOItxOHlkvP6UeOsa2j55bgM2r72ThytdaKj+yLBwwAtnT67cVzcfS+U+J31LpNiRXuHUL18vdfpcIt25sW4UWKSGRoSFj5Go7yaehkasGZMBsuqWl5cCGTXHb4qI1P1kosaGysVbb8/DGNNsFlBiLCsbzel9ung/J0mmOabEGUKcX+R3TVoglDhABKugfD7UHHaX3SwBVQdap2rGZKGsCihej/Iy2aCBgBKqrst6HH5Ni7wzxmQm60NJc/sP10TNJnnCnW8TVOcBzNLeXeK2ogaXvUdFZQ2lvbtElcf2q3jTz5LuEgSIyJFi8RJURop83iVMcpw+mqld6p7byS+q629pykyYDT0HE5kdINH2kZOURT7tn+jp/4ayAkT2/Xjx7Is9P5M+WuGzyKqAkkpezxqS1iK/oMId1+EvNA3ZbSEvNHXos2rm9mWZNiOrbnmlgv0fNp47sM2ZuXLPV8ltH6xyRtnV3pYzJj1lVUBJRR/KwO4ds+uX1FISnjPFfivNdmCr87prXXLb1+vnMSY9ZdUtL6/6UCL7FIoKchnex+mrCPdftBmx91rD995FnP6COP0EClijzpi2yf7MNMYY4wnRLHyQq7S0VP/80QcUPXI8f+13O9+4ZmqzjxU5EioTWyjhdPztC3KpqKybmvj0Pl1YuH4XQYXiwlwGdO9Y+1R+WLgsnIrmfwrv5yTW057D5BZ0sBZKA8L/TKJ+F62ZwiackTn8GSU6d3i72DloCjo6+0qO8zxS7KiwyJFm4YwE4QzPkSPMIo8BTlnVQWc5cj6byNZw7Eiz8LnD+4bPFXZMCSu27OWXXR9idv79dcdraD6cyFF54fNFZk5o6BixEo3Si90/3vvIbRNJlHUhMntDMvP4xBsFGHFuEVmkqqUNH6Rh1kLJcln490JGECywmrYnqwKKPdjoj3KOiloOBzGLZca0LVkVUFKdvr60dxeKC3M5vU+XxjfOQt+njJHM4u+cHHV7YAtHRm0XJGDBxJg2KKsCijHGGP9YQDHGGOMJCyjGGGM8YQHFGGOMJyygNGL2jWfUPjk/+8YzGNC9o881ar6QNrycrO9TlnDM+4uhC1ioJxEkwH7aNe8ExhvhuV4g8fjxRicha0FK1PC+Nnbdfwe2wY7VLfs8k5BVAaU1hw0n84yBV9t4Jfa/dUv/m2/nCCo0H4CNciz3VE/kqerLWnhU45mkcn819q+gBf9KwpOW2Zg//y15BQ5uj/hMUiOrAkqqhw1H6lyUH7UcLzB0KGw8VVomP/y2nmOpcdPBfZXTi+eDF7JbOySxZyZftTEZqJUmkMuqgGKMMcY/FlCMMcZ4wgKKMcYYT1hAMcYY4wkLKMYYYzxhAaWJIp9LiZSX0/DIpRxx5iA5vY+TYDKs/zHF5Iizvrgwl/b5OZ7XOZVW6vHOcymuCgr5XvU9DOUVLuJJt1SgbK8zLwZAIKsmCjXGuCygGGOM8UTa/6koIkXANKAKmKuqL/pcJWOMMXH40kIRkadFZJuILI8pHyMiq0VkjYhMcYuvAF5V1R8CY1u9ssYYY5LiVwtlFvAE8Fy4QERygCeBC4ByYIGIvAH0BJa5m7XShNzZIwSIQk0wOodPeJ54cNbtqaiiJhgiIEJIlaDC4epg1H41wRAiQjCk1ARDLAn2IaA1BFFWBXvUndM9XzUhyLF55mu1Zk6reik2Epw7FAStql+3yLQt1RUQqoHDESmNaiqdbUJBd1913h/a4ywf2u3sE1Z1EEScfTRUt05DcHiPs28gB0KhurqruuepduamD1bX5aIK1QDi1msPO4IdOFhZA7jniXfs8H7BKtCgUx75uwpWQ/VBp48v3N9XfcjZtiHBKme7/CKo3O/8LgrcbB2H9kCBmz2i6kDdcSv3QU4BSMD53TQkVBM/B1eopq7uh3ZBsAZyGvlKD9XUXffhvZDXvuHtm8iXgKKqH4tI75jiEcAaVV0LICKvAJfjBJeewGIaaFGJyA3ADQDHHXec95VOICDOF2hOQCBY958yN5DcV2hBbvQlSQq+eRU4VB39D3Llln217w9Vh1i99UC9/ZaUR+dEizzGyws2cn/199ylIPfznXrnC5ELObClXX+OjTxQIBdCqc0plJ5Sm5gvSs2hmIIEAaXeduHyiMSSWxY7X2hblzlfnHntYNsKJ9BEqj4I25bD+k+c9ZHn3vL36G13r3de92923n/5v3DyJbD+47o6bVsF21Y6y18vjTlXxLm3LkeCJSzfvA8K/173pb1vM+xZD2s+gP6jnbI9G+qurfogHNxZd74dq50ABBDIc4+9rPGcaDtWO1/OHY6GA1udsiOOhwPbnd9HUbe632OOm7Lp66WQ3wFyC6FiB1TsgvYJZoINX2tswNizIeJ3tRJWvQ6Dvt1wXTXk7Fd9yLm2dp0b3r6J0qlTvgewMWK53C37I/BtEXkKeDPRzqo6U1VLVbW0W7duqa1phPzcAAEgN8d5DQAdCnIY0uuIRvftUJBDYV70qK5AKiJKCmzcVdHoNpXkM+jw71hWODx6RU4e5BfH2aN5174x1Hqfd5sU+ddx+K/b2GAS6eCOxo9Z6f5BU+1+uYe/iCt21m1TdcD50k/CyMCyugWNqeOBbRHn3R+9Y2SwCAcTqPuDJ5kEm+GW2+G6P9KoroBgZfT62ONVHag7Z1JJG2P+KIi9ln1bkjiGe95wC/HQ7uT2SVLad8qr6kFgUjLbishlwGX9+vVLbaXqnTf2VZwWS6P7ZUbwaIkDtI/f7PLw2vfhbbPdGNM86dRC2QT0ilju6ZYlrTWzDRtjjImWTgFlAXCiiPQRkXxgAvBGUw7QmvOhGGOMiebXsOGXgXnASSJSLiKTVbUGuBV4D1gFzFHVFQ0dJ5a1UIwxxj9+jfK6KkH5O8A7rVwdY4wxHkj7TvmmaM1O+TNO6MrbS7fQ84h2dO1QwOwbz2D8jHlRw3GXlY1m/Ix5APXyfyUqBxhc9l7UcyJZ475dzuszl0QPAy06Cg5ug9wCZ1jlMSXOEMuq+kOZTQaSnMaf5TBZIZ36UFrMbnkZY4x/siqgWKe8Mcb4J6sCirVQjDHGP1kVUIwxxvjHAooxxhhPZFVAsT4UY4zxT1YFFOtDMcYY/2RVQDHGGOMfCyjGGGM8kVUBxfpQjDHGP1kVUKwPxRhj/JNVAcUYY4x/LKAYY4zxhAUUY4wxnsiqgGKd8sYY45+sCijWKW+MMf7JqoBijDHGPxZQjDHGeMICijHGGE9YQDHGGOMJCyjGGGM8kVUBxYYNG2OMf7IqoNiwYWOM8U9WBRRjjDH+yfW7ApnuJ9/qz2WnHgvA7BvPqLc+XllD5QDLykZz8wuLeHf5195U0m/n3x29POlt57XMbUmeeAEsfhHG/ApKv++U7d0Ejw6A4mNh/2anrM8oWPfn5M5ZMh6Wzm553U3L9SyFjZ9Hlx1xvPN60hhYNCs15z3xQlj8QnRZ136wZXFqztepJ+z8h/P+hPPqyou6Jd6noBgO7W7e+Zp7LUcPat75kmAtFGOMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT1hAMcYY44m0Dygi0ldEfi8ir/pdF2OMMYmlNKCIyNMisk1ElseUjxGR1SKyRkSmNHQMVV2rqpNTWU9jjDEtl+rnUGYBTwDPhQtEJAd4ErgAKAcWiMgbQA7wq5j9v6+q21JcR2OMMR4QVU3tCUR6A2+p6iB3+QygTFVHu8t3AqhqbDCJPc6rqvqdBtbfANzgLp4ErG5x5dPTkcAOvyuRQnZ9mc2uL3OdpKrFLTmAH0/K9wA2RiyXA6cn2lhEugIPAENF5M5EgUdVZwIzvaxoOhKRhapa6nc9UsWuL7PZ9WUuEVnY0mOkfeoVVd0J3OR3PYwxxjTMj1Fem4BeEcs93TJjjDEZzI+AsgA4UUT6iEg+MAF4w4d6ZKpsv61n15fZ7PoyV4uvLaWd8iLyMnAOTkfWVuA+Vf29iFwMPIYzsutpVX0gZZUwxhjTKlI+yssYY0zbkPZPyhtjjMkMFlDSnIisF5FlIrI4PKxPRLqIyPsi8g/3tbPf9UxWvOwJia5HHL9xMyosFZHT/Kt5chJcX5mIbHI/w8XuLd/wujvd61stIqP9qXVyRKSXiHwkIitFZIWI/MQtz4rPr4Hry5bPr1BE5ovIEvf6prrlfUTkc/c6Zrt924hIgbu8xl3fu9GTqKr9pPEPsB44Mqbsv4Ap7vspwH/6Xc8mXM/ZwGnA8sauB7gYeBcQ4BvA537Xv5nXVwb8NM62A4AlQAHQB/gnkOP3NTRwbd2B09z3xcCX7jVkxefXwPVly+cnQAf3fR7wufu5zAEmuOXTgZvd97cA0933E4DZjZ3DWiiZ6XLgWff9s8C/+FiXJlHVj4FdMcWJrudy4Dl1/BU4QkS6t05NmyfB9SVyOfCKqlaq6jpgDTAiZZVrIVXdoqp/c9/vB1bhPKicFZ9fA9eXSKZ9fqqqB9zFPPdHgfOAcPLd2M8v/Lm+CpwvItLQOSygpD8F/p+ILHLTywAcrapb3PdfA0f7UzXPJLqeeFkVGvoPns5udW/7PB1xizJjr8+9/TEU56/crPv8Yq4PsuTzE5EcEVkMbAPex2lV7VHVGneTyGuovT53/V6ga0PHt4CS/s5U1dOAi4AficjZkSvVaY9mzVC9bLse11PACcAQYAvwiL/VaRkR6QD8AbhdVfdFrsuGzy/O9WXN56eqQVUdgvNA+QjgZC+PbwElzanqJvd1G/Aazj+CreFbB+5rpmdkTnQ9WZFVQVW3uv+RQ8BvqbstknHXJyJ5OF+2L6rqH93irPn84l1fNn1+Yaq6B/gIOAPnVmQ4DVfkNdRen7u+E7CzoeNaQEljIlIkIsXh98CFwHKczALXu5tdD7zuTw09k+h63gCuc0cLfQPYG3FrJWPE9BuMw/kMwbm+Ce5omj7AicD81q5fstz7578HVqnqryNWZcXnl+j6sujz6yYiR7jv2+FMIbIKJ7CEM7nHfn7hz/U7wIduCzQxv0ce2E+DozL64owiWQKsAH7hlncF/gT8A/gA6OJ3XZtwTS/j3DaoxrlfOznR9eCMSnkS5z7vMqDU7/o38/qed+u/1P1P2j1i+1+417cauMjv+jdybWfi3M5aCix2fy7Ols+vgevLls+vBPi7ex3LgXvd8r44gXAN8D9AgVte6C6vcdf3bewc9qS8McYYT9gtL2OMMZ6wgGKMMcYTFlCMMcZ4wgKKMcYYT1hAMcYY4wkLKMYYYzxhAcWYFBOR3iJyyM2hFC47WkReEpG1bp62eSIyroFjfBSbHl1EbheRp0SknZtWvUpEjkzltRjTEAsoxrSOf6qTQyn8RPb/BT5W1b6qOgwnPXjPBvZ/2d0m0gTgZVU95B57cwrqbUzSLKAYE0FEXnYnFZovIhtE5JIUnOY8oEpVp4cLVHWDqj7u1uEa9/yLRWSGiOTgpA+/JGLyo97AscAnKaifMc1iAcWYaKcCa1V1BPA94L4UnGMg8Ld4K0TkFGA8MNJtdQSB76nqLpz0Fxe5m04A5qilujBpJLfxTYxpG0SkEOgGTHWLVgKdRWQScDowGngP+LuqzvDwvE/i5JGqwpnQaBiwwJ3LqB112XvDt71ed18ne1UHY7xgAcWYOoOAf6jqYXf5NGCJqj4jIq8Deap6U+xOIjIAGKuqD4rI48Bd6sz4l8gK4NvhBVX9kduZvhAnoeKzqnpnnP1eBx5152Zvr6qLmnORxqSK3fIyps6pwHEiUuhOFzAVeNRdNwxI9AU+HCczLUCnRoIJwIdAoYjcHFHW3n39E/AdETkKQES6iMjxAOpM3/oR8DROa8WYtGIBxZg6pwJ/xJn2dQHwlKr+xV3XWEBZ6QahRrn9Hv8CjBKRdSIyH+dW189VdSVwN860z0txpmmNnI/jZbeeFlBM2rH09ca4ROTPwA2qujrOupeB76vqIRE5GrhUVX/vrnsbZ+6TfcBgVR0Ts29v4C1VHZTi+q/HmXNkRyrPY0wi1odiTJ0TcCaJqkdVr4pYHAqsg9opY3eq6o0NHDcIdBKRxeFnUbzkzr43D8gDQl4f35hkWQvFGGOMJ6wPxRhjjCcsoBhjjPGEBRRjjDGesIBijDHGExZQjDHGeMICijHGGE9YQDHGGOMJCyjGGGM8YQHFGGOMJ/4/CjS9KtGpGBIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = hist.plot1d(output['pt_trail'], overlay='dataset')\n",
    "ax.set_yscale('log')\n",
    "ax.set_ylim(0.1, 2e4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Events/s: 31703.135042905807\n"
     ]
    }
   ],
   "source": [
    "print(\"Events/s:\", output['cutflow']['all events']/elapsed)"
   ]
  }
 ],
 "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}