chi2_comb_test.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 581,
   "id": "7bf97c44",
   "metadata": {},
   "outputs": [],
   "source": [
    "import uproot\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import hist\n",
    "import mplhep as hep\n",
    "# Conversion to hist\n",
    "def tgasym_to_err(tgasym, restoreNorm=False):\n",
    "    x, y = tgasym.values(\"x\"), tgasym.values(\"y\")\n",
    "    xlo, xhi = tgasym.errors(\"low\", axis=\"x\"), tgasym.errors(\"high\", axis=\"x\")\n",
    "    ylo, yhi = tgasym.errors(\"low\", axis=\"y\"), tgasym.errors(\"high\", axis=\"y\")\n",
    "    bins = np.r_[(x - xlo), (x + xhi)[-1]]\n",
    "    binwidth = xlo + xhi\n",
    "    if restoreNorm:\n",
    "        y = y * binwidth\n",
    "        ylo = ylo * binwidth\n",
    "        yhi = yhi * binwidth\n",
    "    return y, bins, ylo, yhi\n",
    "\n",
    "\n",
    "def tgasym_to_hist(tgasym, restoreNorm=True):\n",
    "    y, bins, ylo, yhi = tgasym_to_err(tgasym, restoreNorm)\n",
    "    h = hist.new.Var(bins, flow=False).Weight()\n",
    "    h.view().value = y\n",
    "    h.view().variance = y\n",
    "    return h"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 567,
   "id": "3a7ad20c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "ColormeshArtists(pcolormesh=<matplotlib.collections.QuadMesh object at 0x7f239d8c1090>, cbar=<matplotlib.colorbar.Colorbar object at 0x7f239d9db880>, text=[])"
      ]
     },
     "execution_count": 567,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAokAAAGwCAYAAADFSv/ZAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABAP0lEQVR4nO3de3hU9Z3H8c+ZyY1bwiWSBIwERUUEEg0So6ViTQ1stVLdLtoLkLJsq4lV46WlViLax+CllLUbzWpF7Xat1G61rbpxa+TyWIOUKOtlAQ1FicWEi0IgSi4zZ/9QpowT8AycH2eG8349z3mezJkz3/Ody5l853c5x7Jt2xYAAABwgIDXCQAAACDxUCQCAAAgBkUiAAAAYlAkAgAAIAZFIgAAAGJQJAIAACAGRSIAAABipHidgGnhcFhbt27VoEGDZFmW1+kAAAAHbNvWnj17NGLECAUC7rRp7du3T93d3a7EOlJpaWnKyMjwOo1DOuaLxK1btyo/P9/rNAAAwGFobW3V8ccff8Rx9u3bp9GjBqptW8iFrI5cbm6uNm/enNCF4jFfJA4aNEiSNDWnQimBNPcChxLjQ+apcNhQXPcvAmTbhnI1wAoGzQR26Zf4Z9k9Pe4H7e11P6YkO4mO28DAge4HNfXZMsVE5w/XGDPDwHvVG+7WiraHI//Hj1R3d7fatoX0bnOBMgd5O9quY09Yo4rfUXd3N0Wil/Z3MacE0twtEu3k+WdjjqnCy+dFYiDJisSAgf8OlqFcreQ5bgNufl9FgiZbkWjgs8WVaM0wOJzL7aFiAwdZGjjI2+FnYSO/gNx3zBeJAAAA+4XssEIe/1YIJUnDBUUiAADwjbBshT0ed+D1/p3iFDgAAACIQUsiAADwjbDCxkbUx5NDMqBIBAAAvhGybYU8nsDk9f6dorsZAAAAMWhJBAAAvsHEFecoEgEAgG+EZStEkegI3c0AAACIQUsiAADwDbqbnaNIBAAAvsHsZufobgYAAEAMWhIBAIBvhD9dvM4hGfinSAyFJDvkdRZwwDZx4fOwmaZ9Kxh0P2jAUAO/obiW5X5c29RrEHL/O8AylasJluV1Bt7jNTDzGoQNfG8b6pENJcDsZq/375R/ikQAAOB7IfuTxesckkES/QQGAADA0UJLIgAA8A3GJDpHkQgAAHwjLEsheTs2Nezx/p2iuxkAAAAxaEkEAAC+EbaNnfAirhySAUUiAADwjVACdDd7vX+n6G4GAABADFoSAQCAb9CS6BxFIgAA8I2wbSlsezy72eP9O0V3MwAAAGLQkggAAHyD7mbnKBIBAIBvhBRQyOOO1JCne3eOIhEAAPiGnQBjEm3GJAIAAMANdXV1KigoUEZGhkpKSrRmzZpDbr9r1y5VVlYqLy9P6enpOuWUU/Tss8/GtU9aEg9XwFB9bRn6dWEbOL172NAlykPJculzSYHk+DUoydxnK5leAwOfLdvEsSXJMvEdY+ozADPHgan3qzdZOjvNSMYxicuWLVN1dbXq6+tVUlKiJUuWqLy8XBs3btTw4cNjtu/u7taXv/xlDR8+XL/97W81cuRIvfvuuxo8eHBc+6VIBAAAvhGyAwrZHo9J/PS3ZUdHR9T69PR0paenx2y/ePFizZs3TxUVFZKk+vp6PfPMM1q6dKl++MMfxmy/dOlSffDBB3rppZeUmpoqSSooKIg7T7qbAQAAPJCfn6+srKzIUltbG7NNd3e3mpubVVZWFlkXCARUVlampqamPuP+4Q9/UGlpqSorK5WTk6Px48frjjvuUCgUXysyLYkAAMA3wrIU9riNLKxPmhJbW1uVmZkZWd9XK+KOHTsUCoWUk5MTtT4nJ0cbNmzoM/5f//pXvfDCC/rmN7+pZ599Vi0tLbrqqqvU09Ojmpoax3lSJAIAAN9IpDGJmZmZUUWiW8LhsIYPH64HHnhAwWBQxcXF+tvf/qa7776bIhEAAOBYkJ2drWAwqPb29qj17e3tys3N7fMxeXl5Sk1NVTAYjKw77bTT1NbWpu7ubqWlpTnaN2MSAQCAb+yfuOL14lRaWpqKi4vV2NgYWRcOh9XY2KjS0tI+H3PuueeqpaVF4QPOQvLWW28pLy/PcYEoUSQCAAAf+WRMovdLPKqrq/Xggw/q0Ucf1fr163XllVeqs7MzMtt51qxZmj9/fmT7K6+8Uh988IGuueYavfXWW3rmmWd0xx13qLKyMq790t0MAACQwGbOnKnt27drwYIFamtrU1FRkRoaGiKTWbZs2aLAAedWzc/P13PPPafrrrtOEydO1MiRI3XNNdfoBz/4QVz7pUgEAAC+EU6Aazfvn90cj6qqKlVVVfV534oVK2LWlZaWavXq1XHv50AUiQAAwDcS42TaZq7U5DaKRAAA4BthBRLmPImJjokrAAAAiOFpkVhbW6uzzjpLgwYN0vDhwzVjxgxt3Lgxapt9+/apsrJSw4YN08CBA3XZZZfFnCsIAADAiZBtJcSSDDwtEleuXKnKykqtXr1af/rTn9TT06MLL7xQnZ2dkW2uu+46/fGPf9QTTzyhlStXauvWrbr00ks9zBoAACSr0KcTV7xekoGnYxIbGhqibj/yyCMaPny4mpub9cUvflG7d+/WQw89pMcee0xf+tKXJEkPP/ywTjvtNK1evVpnn312TMyuri51dXVFbnd0dJh9EgAAAMeghJq4snv3bknS0KFDJUnNzc3q6elRWVlZZJuxY8fqhBNOUFNTU59FYm1trRYuXBgb3LY/WdwSCrkX60AHXELHVQecdd29mGYG3lop7n8sbdvA85dkd/e4HzRo5hdmcnRuIMIy8I4F+BQYY+B7S7297seEwnZAYY9nN4eTZHZzwrR3hsNhXXvttTr33HM1fvx4SVJbW5vS0tI0ePDgqG1zcnLU1tbWZ5z58+dr9+7dkaW1tdV06gAAIEl43c1Md/NhqKys1BtvvKEXX3zxiOKkp6crPT3dpawAAAD8KSGKxKqqKj399NNatWqVjj/++Mj63NxcdXd3a9euXVGtie3t7crNzfUgUwAAkMzCkuezi80MgHKfp+2dtm2rqqpKTz75pF544QWNHj066v7i4mKlpqaqsbExsm7jxo3asmWLSktLj3a6AAAgye0/mbbXSzLwtCWxsrJSjz32mH7/+99r0KBBkXGGWVlZ6tevn7KysjR37lxVV1dr6NChyszM1NVXX63S0tI+J60AAADAHZ4Wiffff78kaerUqVHrH374Yc2ZM0eS9LOf/UyBQECXXXaZurq6VF5ervvuu+8oZwoAAI4FiXHtZloSP5ftYAp4RkaG6urqVFdXdxQyAgAAx7KwLIU9PimY1/t3KiEmrgAAABwNtCQ6lxxZAgAA4KiiJREAAPhGIpzM2uv9O0WRCAAAfCNsWwp7fZ5Ej/fvVHKUsgAAADiqaEkEAAC+EU6A7mZOpg0AAJBgwnZAYY9nF3u9f6eSI0sAAAAcVf5pSbTtTxa3BAzV15ahwazBoPsxe0Pux5RkhwzE7e11P6Zk5P2yUgwdlqbimnhtTR0HAQNxjeVq4DvGVK7JxMR3oWTmOAi7+D8LESFZCnl8Mmuv9++Uf4pEAADge3Q3O5ccWQIAAOCooiURAAD4Rkjed/eaGazlPopEAADgG3Q3O0eRCAAAfCNkBxTyuEjzev9OJUeWAAAAOKpoSQQAAL5hy1LY4zGJNqfAAQAASCx0NzuXHFkCAADgqKIlEQAA+EbYthS2ve3u9Xr/TlEkAgAA3wgpoJDHHale79+p5MgSAAAARxUtiQAAwDfobnaOIhEAAPhGWAGFPe5I9Xr/TiVHlgAAADiqfNOSaPf2yrbcq4mttFTXYiWtgKHm8l7b/ZiWoVxTDRxCKUH3Y0rmXgPbwPtlIqYkhULux3TxeyWKqdcgmQQNHAsmPgMw8/1i6CsrZFsKedzd6/X+nfJNkQgAAMCYROcoEgEAgG/YdkBhj694YnPFFQAAACQrikQAAOAbIVkJscSrrq5OBQUFysjIUElJidasWXPQbR955BFZlhW1ZGRkxL1PupsBAIBvhG3vxwSG45yXtmzZMlVXV6u+vl4lJSVasmSJysvLtXHjRg0fPrzPx2RmZmrjxo2R29ZhTC6iJREAAMADHR0dUUtXV1ef2y1evFjz5s1TRUWFxo0bp/r6evXv319Lly49aGzLspSbmxtZcnJy4s6PIhEAAPhG+NOJK14vkpSfn6+srKzIUltbG5Nvd3e3mpubVVZWFlkXCARUVlampqamgz7PvXv3atSoUcrPz9cll1yiN998M+7Xiu5mAADgG2FZCps6CWMcOUhSa2urMjMzI+vT09Njtt2xY4dCoVBMS2BOTo42bNjQZ/xTTz1VS5cu1cSJE7V7927dc889Ouecc/Tmm2/q+OOPd5wnRSIAAIAHMjMzo4pEt5SWlqq0tDRy+5xzztFpp52mf//3f9ftt9/uOA5FIgAA8I1ku+JKdna2gsGg2tvbo9a3t7crNzfXUYzU1FSdccYZamlpiStPxiQCAADf8Hos4oFjEp1IS0tTcXGxGhsb//4cwmE1NjZGtRYeSigU0uuvv668vLy4XitaEgEAABJYdXW1Zs+erUmTJmny5MlasmSJOjs7VVFRIUmaNWuWRo4cGZn4ctttt+nss8/WmDFjtGvXLt19991699139c///M9x7ZciEQAA+EZYCXDt5jgnzsycOVPbt2/XggUL1NbWpqKiIjU0NEQms2zZskWBwN9bJz/88EPNmzdPbW1tGjJkiIqLi/XSSy9p3Lhxce3Xsm07zlM6JpeOjg5lZWXpgsGzlGKluRbXSkt1LVaUlCSq23t6jIS1uw3EDYfdjylJqe6/X1aqoc+Wqbjd3a6HtLvcjylJ9scfux/UMjNqJ5A91P2gKUH3Y5oUNJBvKOR+TEky8a803jMuO5Uk//Z7w116fuu/a/fu3a5M7thfD3y9cZZSB7hXDxyOns5uPXHBL117bqYkUUUCAABwZMJ2ArQkerx/p5i4AgAAgBj+aUns7XW3W8hAN6NRh3HNxs8VMPQbI2ggrqEuJvvjfUbimmCsG9tEl6CJz4BkZDiHZeo4MHHMGuoaV8BQq4iprmG4z0QXtqFe8XhnF5vKIRkkWaUDAABw+Ohudi45SlkAAAAcVbQkAgAA30ikazcnOopEAADgG3Q3O0d3MwAAAGLQkggAAHyDlkTnKBIBAIBvUCQ6R3czAAAAYtCSCAAAfIOWROcoEgEAgG/Y8v4UNIYuJuM6ikQAAOAbtCQ6x5hEAAAAxKAlEQAA+AYtic5RJAIAAN+gSHSO7mYAAADE8E9LYkqKZLn4dFMMvXQpQTNxTeg1FDfs/rwvOxx2PaYkyXL/16BlGfrtFjD0yzWcTK9BEv0uNvF+BQ09/5Ch4wswgJZE5/xTJAIAAN+zbUu2x0Wa1/t3Kol+VgMAAOBooSURAAD4RliW5yfT9nr/TlEkAgAA32BMonN0NwMAACAGLYkAAMA3mLjiHEUiAADwDbqbnaNIBAAAvkFLonOMSQQAAEAMWhIBAIBv2AnQ3ZwsLYkUiQAAwDdsSbb7V3+NO4dk4Gl386pVq3TxxRdrxIgRsixLTz31VNT9c+bMkWVZUcu0adO8SRYAAMBHPG1J7OzsVGFhob7zne/o0ksv7XObadOm6eGHH47cTk9PP1rpAQCAY0xYliyuuOKIp0Xi9OnTNX369ENuk56ertzc3KOUEQAAOJYxu9m5hJ/dvGLFCg0fPlynnnqqrrzySu3cufOQ23d1damjoyNqAQAAQHwSeuLKtGnTdOmll2r06NHatGmTfvSjH2n69OlqampSMBjs8zG1tbVauHBh7B29vZLlYk3c2+terAMFEr5u/ztTI39DIfdjGsrVSjFwCAUM/cK0DMVNNfAamIgpyeqXYSSuEYG+v+OOSE+P+zElKWzou8DEcWvqOAiH3Y9p6jvWRK4mhA38L9AnJ7K2OJm2IwldJF5++eWRvydMmKCJEyfqpJNO0ooVK3TBBRf0+Zj58+eruro6crujo0P5+fnGcwUAAInPthNgdnOSTG9OomYr6cQTT1R2drZaWloOuk16eroyMzOjFgAAAMQnoVsSP+u9997Tzp07lZeX53UqAAAgCTFxxTlPi8S9e/dGtQpu3rxZ69at09ChQzV06FAtXLhQl112mXJzc7Vp0ybddNNNGjNmjMrLyz3MGgAAJCuKROc8LRLXrl2r888/P3J7/1jC2bNn6/7779drr72mRx99VLt27dKIESN04YUX6vbbb+dciQAA4LAwccU5T4vEqVOnyj7E6M3nnnvuKGYDAACA/ZJqTCIAAMCRYHazc0k1uxkAAOBIfFIkWh4v8eddV1engoICZWRkqKSkRGvWrHH0uMcff1yWZWnGjBlx75MiEQAAIIEtW7ZM1dXVqqmp0SuvvKLCwkKVl5dr27Zth3zcO++8oxtuuEFTpkw5rP1SJAIAAN/wvhXx77OrP3sZ4a6urj5zXrx4sebNm6eKigqNGzdO9fX16t+/v5YuXXrQ5xkKhfTNb35TCxcu1IknnnhYrxVFIgAA8A07QRZJys/PV1ZWVmSpra2Nybe7u1vNzc0qKyuLrAsEAiorK1NTU9NBn+dtt92m4cOHa+7cuXG+Qn/HxBUAAAAPtLa2Rl0Zrq9T/O3YsUOhUEg5OTlR63NycrRhw4Y+47744ot66KGHtG7duiPKjyIRAAD4RiKdTNvE5YP37Nmjb3/723rwwQeVnZ19RLEoEgEAgH8c2N/rZQ4OZWdnKxgMqr29PWp9e3u7cnNzY7bftGmT3nnnHV188cWRdeFwWJKUkpKijRs36qSTTnK0b8YkAgAA/0iESStxtGSmpaWpuLhYjY2NkXXhcFiNjY0qLS2N2X7s2LF6/fXXtW7dusjy1a9+Veeff77WrVun/Px8x/v2TUuiNXCgrECaewFDIfdiHcDu/MhI3KQ5c6cku7fX9ZiB/v1djylJyjBwichPf/G5rrvHSNjwiCPrzujLx3lm3q+PhwXdD2qo1yr7+S2ux7R37XY9pkl2j/vfBQoaahtJou9YK+j+cWAb+N4K292ux0xW1dXVmj17tiZNmqTJkydryZIl6uzsVEVFhSRp1qxZGjlypGpra5WRkaHx48dHPX7w4MGSFLP+88RVJPb09GjatGmqr6/XySefHNeOAAAAvJaMV1yZOXOmtm/frgULFqitrU1FRUVqaGiITGbZsmWLAgH3fwDFVSSmpqbqtddecz0JAACAoyGRJq7Eo6qqSlVVVX3et2LFikM+9pFHHol7f9JhjEn81re+pYceeuiwdgYAAIDkEPeYxN7eXi1dulTPP/+8iouLNWDAgKj7Fy9e7FpyAAAAropz4oixHJJA3EXiG2+8oTPPPFOS9NZbb0XdZ1nJ8aQBAIA/JeOYRK/EXSQuX77cRB4AAABIIIc9FaalpUXPPfecPv74Y0mSnSxlMQAA8C+vL9qcCCfzdijuInHnzp264IILdMopp+gf/uEf9P7770uS5s6dq+uvv971BAEAANzi9Ym0E2F2tVNxF4nXXXedUlNTtWXLFvU/4ATFM2fOVENDg6vJAQAAwBtxj0n8n//5Hz333HM6/vjjo9affPLJevfdd11LDAAAwIgk6e71WtxFYmdnZ1QL4n4ffPCB0tMNXKIMAADAJYnQ3ev1/p2Ku7t5ypQp+uUvfxm5bVmWwuGw7rrrLp1//vmuJgcAAOAqryesJNHElbhbEu+66y5dcMEFWrt2rbq7u3XTTTfpzTff1AcffKA///nPJnIEAADAURZ3S+L48eP11ltv6Qtf+IIuueQSdXZ26tJLL9Wrr76qk046yUSOAAAALrESZEl8cbckSlJWVpZuvvlmt3MBAAAwKxG6e73ev0OHVSR++OGHeuihh7R+/XpJ0rhx41RRUaGhQ4e6mpyrQiHJDnmdxbElHDYU1/2jx7bN5GqZOIl8yNDnNIkum2kHzeQaSnM/brA7Sb7tDbJNfReYOG5NTRgw8L1lim0Zer9wzIm7u3nVqlUqKCjQvffeqw8//FAffvih7r33Xo0ePVqrVq0ykSMAAIA7vJ6wkggtmQ7F3ZJYWVmpmTNn6v7771cwGJQkhUIhXXXVVaqsrNTrr7/uepIAAACusC1zLcrx5JAE4m5JbGlp0fXXXx8pECUpGAyqurpaLS0triYHAAAAb8RdJJ555pmRsYgHWr9+vQoLC11JCgAAwATbTowlGTjqbn7ttdcif3//+9/XNddco5aWFp199tmSpNWrV6uurk6LFi0ykyUAAIAbEmFMoNf7d8hRkVhUVCTLsmQfUPredNNNMdt94xvf0MyZM93LDgAAAJ5wVCRu3rzZdB4AAADmMXHFMUdF4qhRo0znAQAAYJxlf7J4nUMyOKyTaW/dulUvvviitm3bpvBnTqL6/e9/35XEAAAAXMeYRMfiLhIfeeQRffe731VaWpqGDRsm64CrOFiWRZEIAABwDIi7SLzlllu0YMECzZ8/X4FA3GfQAQAA8A5jEh2Lu0j86KOPdPnll1MgAgCA5EN3s2NxV3pz587VE088YSIXAAAAJIi4WxJra2t10UUXqaGhQRMmTFBqamrU/YsXL3YtOQAAAFfRkujYYRWJzz33nE499VRJipm4AgAAkLAoEh2Lu0j86U9/qqVLl2rOnDkG0jEoFJbs8Odv51R6mnuxDmCs0DZwoUi7t9f1mJKkri7XQ9r73I9pimUZGu+bkW4mbq+Lx9WnUj4KuR5TkgZscz+mFTL0bR90/3NgpRv6DLj53XqgUOrnbxOvJPqOTSph9z8Dlh10PSbiE3eRmJ6ernPPPddELgAAAGYxu9mxuH+qXnPNNfr5z39uIhcAAACj9l9xxeslGcTdkrhmzRq98MILevrpp3X66afHTFz53e9+51pyAAAA8EbcReLgwYN16aWXmsgFAADALCauOBZ3kfjwww+byAMAAAAJJO4iEQAAIFlZ8n5MYHJMWzmMInH06NGHPE3LX//61yNKCAAAAN6Lu0i89tpro2739PTo1VdfVUNDg2688Ua38gIAAHAfp8BxLO4i8ZprrulzfV1dndauXXvECQEAABjDxBXHXDul//Tp0/Vf//VfboUDAACAh1ybuPLb3/5WQ4cOdSscAACA+2hJdCzuIvGMM86Imrhi27ba2tq0fft23Xfffa4mBwAA4KZEuOKJ1/t3Ku4iccaMGVG3A4GAjjvuOE2dOlVjx451Ky8AAAB8qq6uTnfffbfa2tpUWFion//855o8eXKf2/7ud7/THXfcoZaWFvX09Ojkk0/W9ddfr29/+9tx7TPuIrGmpibehwAAACSGJOxuXrZsmaqrq1VfX6+SkhItWbJE5eXl2rhxo4YPHx6z/dChQ3XzzTdr7NixSktL09NPP62KigoNHz5c5eXljvd7WGMSw+GwWlpatG3bNoXD4aj7vvjFLx5OSAAAAPMSqEjs6OiIWp2enq709PSYzRcvXqx58+apoqJCklRfX69nnnlGS5cu1Q9/+MOY7adOnRp1+5prrtGjjz6qF1980WyRuHr1an3jG9/Qu+++K9uOfpUty1IoFIo35FFh22HZdvjzN3TI2BmO0tPMxA27f0RYtpmjzERUK+DaRP7ouP36uR80Jeh+TElKMXOBJSvs3nG1X9q2j1yPKUlpKWY+B0nD0PeLqe8Cmfh/coiLQRwRA8eBMQb+H7j5/zUimV7Tw5Sfnx91u6amRrfeemvUuu7ubjU3N2v+/PmRdYFAQGVlZWpqavrcfdi2rRdeeEEbN27UnXfeGVd+cf/X+N73vqdJkybpmWeeUV5e3iGvvgIAAJBIEmniSmtrqzIzMyPr+2pF3LFjh0KhkHJycqLW5+TkaMOGDQfdx+7duzVy5Eh1dXUpGAzqvvvu05e//OW48oy7SHz77bf129/+VmPGjIn3oQAAAN5KoCuuZGZmRhWJbho0aJDWrVunvXv3qrGxUdXV1TrxxBNjuqIPJe4isaSkRC0tLRSJAAAg+STQmEQnsrOzFQwG1d7eHrW+vb1dubm5B31cIBCI1GpFRUVav369amtrzRaJV199ta6//nq1tbVpwoQJSk1Njbp/4sSJ8YYEAABAH9LS0lRcXKzGxsbIaQjD4bAaGxtVVVXlOE44HFZXV1dc+467SLzsssskSd/5znci6yzLkm3bCT1xBQAAIJHGJDpVXV2t2bNna9KkSZo8ebKWLFmizs7OyGznWbNmaeTIkaqtrZUk1dbWatKkSTrppJPU1dWlZ599Vv/xH/+h+++/P679xl0kbt68Od6HAAAAJIYk626WpJkzZ2r79u1asGCB2traVFRUpIaGhshkli1btihwwFk8Ojs7ddVVV+m9995Tv379NHbsWP3qV7/SzJkz49qvZX/2PDbHmI6ODmVlZemCIbOVEnDv9A9WRoZrsaKYOv2JgVMeKM5ma6fCuzs+f6M4WWmGTv2RZWDAcZKdAsfuHzsb74hjGjplkZLoFDiBHbvdD9rT635MSeIUOMl1upYkOQVOb7hbjR8+qt27d7syuWN/PXDigjsUMPU/3KHwvn36620/cu25mWLmvwYAAEAiSoDuZs9bMh2iSAQAAP6RhN3NXkmevhcAAAAcNbQkAgAA/6Al0bHDbkns7u7We++9py1btkQt8Vi1apUuvvhijRgxQpZl6amnnoq637ZtLViwQHl5eerXr5/Kysr09ttvH27KAADA5/afAsfrJRnEXSS+/fbbmjJlivr166dRo0Zp9OjRGj16tAoKCjR69Oi4YnV2dqqwsFB1dXV93n/XXXfp3nvvVX19vV5++WUNGDBA5eXl2rdvX7xpAwAAIA5xdzfPmTNHKSkpevrpp5WXlyfrCE4nMH36dE2fPr3P+2zb1pIlS/TjH/9Yl1xyiSTpl7/8pXJycvTUU0/p8ssvP+z9AgAA4NDiLhLXrVun5uZmjR071kQ+EZs3b1ZbW5vKysoi67KyslRSUqKmpqaDFoldXV1Rl53p6HD/nHsAACBJMSbRsbiLxHHjxmnHjh0mconS1tYmSZGzie+Xk5MTua8vtbW1WrhwYewdPT3unkT1M9esdi+uoblEJk7QbOpk2t3drscMDhzgekxJ6j3hOPeDGjrZb6DLzImUA7s6XY9p7TPz2VKvgdfA1Im/0wx8xyTTyfolKWggX1Mn/jbx2TLx/CUjr4Fl4GTilqHPVSKMCfR6/07F/e1255136qabbtKKFSu0c+dOdXR0RC1emz9/vnbv3h1ZWltbvU4JAAAg6cTdbLW/+/eCCy6IWm/btizLUsilyyjl5uZKktrb25WXlxdZ397erqKiooM+Lj09Xenp7l8mDAAAHCOSpCXPa3EXicuXLzeRR4zRo0crNzdXjY2NkaKwo6NDL7/8sq688sqjkgMAADjGMCbRsbiLxPPOO8+1ne/du1ctLS2R25s3b9a6des0dOhQnXDCCbr22mv1k5/8RCeffLJGjx6tW265RSNGjNCMGTNcywEAAACxHBWJr732msaPH69AIKDXXnvtkNtOnDjR8c7Xrl2r888/P3K7urpakjR79mw98sgjuummm9TZ2al/+Zd/0a5du/SFL3xBDQ0NysjIcLwPAACA/Zi44pyjIrGoqEhtbW0aPny4ioqKZFmW7D5mR8U7JnHq1Kl9xjkw3m233abbbrvNcUwAAICDorvZMUdF4ubNm3XcccdF/gYAAMCxzVGROGrUqD7/BgAASCZ0NzsX98SVnTt3atiwYZKk1tZWPfjgg/r444/11a9+VVOmTHE9QQAAANfQ3eyY45Npv/766yooKNDw4cM1duxYrVu3TmeddZZ+9rOf6YEHHtD555+vp556ymCqAAAAOFocF4k33XSTJkyYoFWrVmnq1Km66KKL9JWvfEW7d+/Whx9+qO9+97tatGiRyVwBAACOjJ0gSxJw3N38l7/8RS+88IImTpyowsJCPfDAA7rqqqsU+PS6pVdffbXOPvtsY4kCAAAcKcYkOue4SPzggw8il8obOHCgBgwYoCFDhkTuHzJkiPbs2eN+hgAAAG5JhJY8r/fvkOPuZumT8xYe6jYAAACODXHNbp4zZ47S09MlSfv27dP3vvc9DRgwQJLU1dXlfnYAAABuoiXRMcdF4uzZs6Nuf+tb34rZZtasWUeekSFW//6yAmnuBTzElWKOyD4zxbbd0+N6zHDHXtdjmtJdeKKRuFsuTHc9ZvoHroeUJI3480dG4gZ6e90Paug4UDCuzhNvffSx+zHDZr63bDtsJK6Vmup+0F7nVwWLh4nXIKn66sIGPgMmYooxifFwXCQ+/PDDJvMAAABAAon7ZNoAAABJi+5mxygSAQCAb9Dd7FwSDdABAADA0UJLIgAA8A+6mx2jSAQAAP5BkegY3c0AAACIQUsiAADwDUven4PS6/07RZEIAAD8g+5mxygSAQCAb3AKHOcYkwgAAIAYtCQCAAD/oLvZMYpEAADgL0lSpHmN7mYAAADEoCURAAD4BhNXnPNPkZg1SAqmuxdvV4d7sQ7UGzIS1rLcbzQO5hznekxJ6jolz/WYrWUuvvcHGDFpq+sxt64d4XpMSZJt6FvJwGfWtsOux5QkKyXNQFAzZzyz93W5HzNk5vvF2GfLhJ5eM3EDBjrmgmaOA4Xdf79MHLN22NTnVd53N3u9f4fobgYAAEAMikQAAOAb+7ubvV7iVVdXp4KCAmVkZKikpERr1qw56LYPPvigpkyZoiFDhmjIkCEqKys75PYHQ5EIAAD8w06QJQ7Lli1TdXW1ampq9Morr6iwsFDl5eXatm1bn9uvWLFCV1xxhZYvX66mpibl5+frwgsv1N/+9re49kuRCAAA4IGOjo6opaur77HIixcv1rx581RRUaFx48apvr5e/fv319KlS/vc/j//8z911VVXqaioSGPHjtUvfvELhcNhNTY2xpUfRSIAAPANr7uZD+xuzs/PV1ZWVmSpra2Nybe7u1vNzc0qKyuLrAsEAiorK1NTU5Oj5/zRRx+pp6dHQ4cOjeu18s/sZgAAgASa3dza2qrMzMzI6vT02DNx7NixQ6FQSDk5OVHrc3JytGHDBke7+8EPfqARI0ZEFZpOUCQCAAD/SKAiMTMzM6pINGHRokV6/PHHtWLFCmVkZMT1WIpEAACABJWdna1gMKj29vao9e3t7crNzT3kY++55x4tWrRIzz//vCZOnBj3vhmTCAAAfMPrsYjxngInLS1NxcXFUZNO9k9CKS0tPejj7rrrLt1+++1qaGjQpEmTDuu1oiURAAD4RwJ1NztVXV2t2bNna9KkSZo8ebKWLFmizs5OVVRUSJJmzZqlkSNHRia+3HnnnVqwYIEee+wxFRQUqK2tTZI0cOBADRw40PF+KRIBAAAS2MyZM7V9+3YtWLBAbW1tKioqUkNDQ2Qyy5YtWxQ44NKQ999/v7q7u/WP//iPUXFqamp06623Ot4vRSIAAPANy7ZleXy98cPZf1VVlaqqqvq8b8WKFVG333nnncPIKhZFIgAA8I8k7G72ChNXAAAAEIOWRAAA4Bvxzi42lUMy8E+R2PmRFAi5Fi68Z69rsQ5kBYNG4tonn+B6zL99Kcv1mJJ09dynXI95Tv9NrseUpEt/db3rMU94vu9rdx6p4MYtRuLa3T3uxwyHXY9pimWZ6ZCxe3vdD2pZ7seUue8t9Rh4DUyNRQu6/zkwcWxJMvM5CJj5bBlBd7NjdDcDAAAghn9aEgEAgO/R3ewcRSIAAPAPupsdo0gEAAC+QUuic4xJBAAAQAxaEgEAgH/Q3ewYRSIAAPCVZOnu9RrdzQAAAIhBSyIAAPAP2zZ3UvV4ckgCFIkAAMA3mN3sHN3NAAAAiEFLIgAA8A9mNztGkQgAAHzDCn+yeJ1DMqC7GQAAADFoSQQAAP5Bd7Nj/ikS09OkQJpr4QKDBroWK0pqqpGwHw/v73rMvaNDrseUpHlZ77se88HdJ7keU5Lyn+9yPWZwebPrMSXJzLslBdLcO65Mszvd7+OxA5brMSXJSjHw9RwMuh9Tkgy9BgoZ6JOzDL1flvsdc7YM9UkaeL9MPH8TMSVmN8fDP0UiAAAA50l0jDGJAAAAiEFLIgAA8A26m52jSAQAAP7BxBXH6G4GAABADFoSAQCAb9Dd7BxFIgAA8A9mNztGdzMAAABi0JIIAAB8g+5m5ygSAQCAfzC72TG6mwEAABCDlkQAAOAbdDc7l9Atibfeeqssy4paxo4d63VaAAAgWYXtxFiSQMK3JJ5++ul6/vnnI7dTUhI+ZQAAkKgYk+hYwldcKSkpys3Ndbx9V1eXurq6Irc7OjpMpAUAAHBMS/gi8e2339aIESOUkZGh0tJS1dbW6oQTTjjo9rW1tVq4cGHMertfuuxgunuJ9Xcx1gHs1KCRuAparoccsMVMric++S+uxxy42Uyu+W+943rMUEqq6zElKTgky0hcE+xQyEzgnl73YwbMjNqx0tMMBHX/e8CoYNj9mEnSzSdJVpqZ74LkYebzasn7MYHJciQm9JjEkpISPfLII2poaND999+vzZs3a8qUKdqzZ89BHzN//nzt3r07srS2th7FjAEAQELbf8UVr5ckkNAtidOnT4/8PXHiRJWUlGjUqFH6zW9+o7lz5/b5mPT0dKWnm2nlAwAA8IuELhI/a/DgwTrllFPU0tLidSoAACAJcQoc5xK6u/mz9u7dq02bNikvL8/rVAAAQDKyE2RJAgldJN5www1auXKl3nnnHb300kv62te+pmAwqCuuuMLr1AAAAI5pCd3d/N577+mKK67Qzp07ddxxx+kLX/iCVq9ereOOO87r1AAAQBKybFuWxxNHvN6/Uwndkvj4449r69at6urq0nvvvafHH39cJ510ktdpAQCAZBVOkCVOdXV1KigoUEZGhkpKSrRmzZqDbvvmm2/qsssuU0FBgSzL0pIlS+LfoRK8SAQAAPC7ZcuWqbq6WjU1NXrllVdUWFio8vJybdu2rc/tP/roI5144olatGhRXBck+SyKRAAA4Bv7u5u9XqRPrgp34HLgFeMOtHjxYs2bN08VFRUaN26c6uvr1b9/fy1durTP7c866yzdfffduvzyy4/otIAUiQAAwD+8ntV8wOzm/Px8ZWVlRZba2tqYdLu7u9Xc3KyysrLIukAgoLKyMjU1Nbn0ovQtoSeuAAAAuCoRrnjy6f5bW1uVmZkZWd1Xq9+OHTsUCoWUk5MTtT4nJ0cbNmwwmiZFIgAAgAcyMzOjisREQ5EIAAB8I9muuJKdna1gMKj29vao9e3t7Uc0KcUJxiQCAAD/2N/d7PXiUFpamoqLi9XY2BhZFw6H1djYqNLSUhOvUAQtiQAAAAmsurpas2fP1qRJkzR58mQtWbJEnZ2dqqiokCTNmjVLI0eOjEx86e7u1v/93/9F/v7b3/6mdevWaeDAgRozZozj/fqmSAwNypCVkuF1Gp/LNtS2a/W637ae/Uav6zElKft/3c81Y9te12NKkvb1fbqCIxHIHOh6TEnSgAFm4va6/zmwurtdjylJdugwzmDrldRU92MGDX3BmHpdvZ5cEI9QyP2YAUPvVzhZjgMzz98Kf7J4Kd79z5w5U9u3b9eCBQvU1tamoqIiNTQ0RCazbNmyRYEDPi9bt27VGWecEbl9zz336J577tF5552nFStWON6vb4pEAACARJrdHI+qqipVVVX1ed9nC7+CggLZLjxHxiQCAAAgBi2JAADAPw44mbWnOSQBikQAAOAbB14Wz8sckgHdzQAAAIhBSyIAAPCPJJ244gWKRAAA4B+2JK/PApQcNSJFIgAA8A/GJDrHmEQAAADEoCURAAD4hy3vxwQmR0MiRSIAAPARJq44RnczAAAAYtCSCAAA/CMsyUqAHJIARSIAAPANZjc7R3czAAAAYtCSCAAA/IOJK475pkjszMtQSmqGa/FSuswMKLB6zXxwUvf0uB4z2NntekxJCmzf7XpMe3eH6zElSWlproe0UlNdjylJsg0NggkYGNxjGRowZCqu3yXJPzyjAgY65vi8mkGR6BjdzQAAAIjhm5ZEAAAAWhKdo0gEAAD+wSlwHKNIBAAAvsEpcJxjTCIAAABi0JIIAAD8gzGJjlEkAgAA/wjbkuVxkRZOjiKR7mYAAADEoCURAAD4B93NjlEkAgAAH0mAIlFe798ZupsBAAAQg5ZEAADgH3Q3O0aRCAAA/CNsy/PuXmY3AwAAIFnRkggAAPzDDn+yeJ1DEvBNkbjrlICC6e41nKZ1mGmEzdrcaySuZeDzaPWE3A8qST09roe0Q2ZytVKC7gc1NVZlX5eZuEEDr4GJmJKstFQjcY3o7nY/Zoqhr/xeM99bsqzkiCmZOW4NHQe+x5hEx3xTJAIAADAm0TnGJAIAACAGLYkAAMA/6G52jCIRAAD4hy3vi7TkqBHpbgYAAEAsWhIBAIB/0N3sGEUiAADwj3BYksfnKQwnx3kS6W4GAABADFoSAQCAf9Dd7BhFIgAA8A+KRMfobgYAAEhwdXV1KigoUEZGhkpKSrRmzZpDbv/EE09o7NixysjI0IQJE/Tss8/GvU+KRAAA4B9hOzGWOCxbtkzV1dWqqanRK6+8osLCQpWXl2vbtm19bv/SSy/piiuu0Ny5c/Xqq69qxowZmjFjht5444249kuRCAAAfMO2wwmxxGPx4sWaN2+eKioqNG7cONXX16t///5aunRpn9v/67/+q6ZNm6Ybb7xRp512mm6//XadeeaZ+rd/+7e49kuRCAAA/MNOgFbET8ckdnR0RC1dXV0x6XZ3d6u5uVllZWWRdYFAQGVlZWpqaurzKTY1NUVtL0nl5eUH3f5gKBIBAAA8kJ+fr6ysrMhSW1sbs82OHTsUCoWUk5MTtT4nJ0dtbW19xm1ra4tr+4PxzezmgZN3KNg/3bV4e5uOcy3WgXrTzdTt7j3zo2BAf9dDBlIMfdQHup+rdnW4H1OS3fmRkbhKS3U9pGUZ+v0aNBA3GHQ/piR7z173g5o6Dnp7zcQ1kK+Vaug1sCz3YybJCZclmTkOLDPH1ieteIkxu7m1tVWZmZmR1enpifXf2jdFIgAAgMJhyfK4AP90TGJmZmZUkdiX7OxsBYNBtbe3R61vb29Xbm5un4/Jzc2Na/uDobsZAAAgQaWlpam4uFiNjY2RdeFwWI2NjSotLe3zMaWlpVHbS9Kf/vSng25/MLQkAgAA/0ig7manqqurNXv2bE2aNEmTJ0/WkiVL1NnZqYqKCknSrFmzNHLkyMiYxmuuuUbnnXeefvrTn+orX/mKHn/8ca1du1YPPPBAXPulSAQAAL5hh8OyPe5ujvcUODNnztT27du1YMECtbW1qaioSA0NDZHJKVu2bFEg8PfO4XPOOUePPfaYfvzjH+tHP/qRTj75ZD311FMaP358XPulSAQAAEhwVVVVqqqq6vO+FStWxKz7+te/rq9//etHtE+KRAAA4B9J2N3sFYpEAADgH2FbsigSnWB2MwAAAGLQkggAAPzDtiV5fZ7E5GhJpEgEAAC+YYdt2R53N9sUiQAAAAnGDsv7lsTkuORiUoxJrKurU0FBgTIyMlRSUqI1a9Z4nRIAAMAxLeGLxGXLlqm6ulo1NTV65ZVXVFhYqPLycm3bts3r1AAAQJKxw3ZCLMkg4YvExYsXa968eaqoqNC4ceNUX1+v/v37a+nSpV6nBgAAko0dTowlCST0mMTu7m41Nzdr/vz5kXWBQEBlZWVqamrq8zFdXV3q6uqK3N69e7ckKfRRV5/bH65Q1z5X4+3X22Pmg9Pb636+dsjd13Q/K2wgbrjb/ZiSFAq6H9NQrrZt6DUw8IvYsgz9fjUR1zLwGZCh98vUPya711Bc9/O1woZeA8tKjpimGDgOej/9LnR7kkevejw/l3averxNwKGELhJ37NihUCgUuTbhfjk5OdqwYUOfj6mtrdXChQtj1r82q85IjgAAwJw9e/YoKyvriOOkpaUpNzdXL7Y960JWRy43N1dpaWlep3FICV0kHo758+eruro6cnvXrl0aNWqUtmzZ4sqHDGZ1dHQoPz9fra2tyszM9DodfA7er+TBe5VceL8+aUHcs2ePRowY4Uq8jIwMbd68Wd3dhnpV4pSWlqaMjAyv0zikhC4Ss7OzFQwG1d7eHrW+vb1dubm5fT4mPT1d6enpMeuzsrJ8e6Alo8zMTN6vJML7lTx4r5KL398vtxt3MjIyEr4wSyQJPXElLS1NxcXFamxsjKwLh8NqbGxUaWmph5kBAAAc2xK6JVGSqqurNXv2bE2aNEmTJ0/WkiVL1NnZqYqKCq9TAwAAOGYlfJE4c+ZMbd++XQsWLFBbW5uKiorU0NAQM5nlYNLT01VTU9NnFzQSD+9XcuH9Sh68V8mF9wuJwLKT5QKCAAAAOGoSekwiAAAAvEGRCAAAgBgUiQAAAIhBkQgAAIAYx3SRWFdXp4KCAmVkZKikpERr1qzxOiX04dZbb5VlWVHL2LFjvU4Ln1q1apUuvvhijRgxQpZl6amnnoq637ZtLViwQHl5eerXr5/Kysr09ttve5MsPvf9mjNnTszxNm3aNG+S9bna2lqdddZZGjRokIYPH64ZM2Zo48aNUdvs27dPlZWVGjZsmAYOHKjLLrss5gITgCnHbJG4bNkyVVdXq6amRq+88ooKCwtVXl6ubdu2eZ0a+nD66afr/fffjywvvvii1ynhU52dnSosLFRdXd/XP7/rrrt07733qr6+Xi+//LIGDBig8vJy7du37yhnCunz3y9JmjZtWtTx9utf//ooZoj9Vq5cqcrKSq1evVp/+tOf1NPTowsvvFCdnZ2Rba677jr98Y9/1BNPPKGVK1dq69atuvTSSz3MGr5iH6MmT55sV1ZWRm6HQiF7xIgRdm1trYdZoS81NTV2YWGh12nAAUn2k08+GbkdDoft3Nxc++67746s27Vrl52enm7/+te/9iBDHOiz75dt2/bs2bPtSy65xJN8cGjbtm2zJdkrV660bfuTYyk1NdV+4oknItusX7/elmQ3NTV5lSZ85JhsSezu7lZzc7PKysoi6wKBgMrKytTU1ORhZjiYt99+WyNGjNCJJ56ob37zm9qyZYvXKcGBzZs3q62tLepYy8rKUklJCcdaAluxYoWGDx+uU089VVdeeaV27tzpdUqQtHv3bknS0KFDJUnNzc3q6emJOr7Gjh2rE044geMLR8UxWSTu2LFDoVAo5qosOTk5amtr8ygrHExJSYkeeeQRNTQ06P7779fmzZs1ZcoU7dmzx+vU8Dn2H08ca8lj2rRp+uUvf6nGxkbdeeedWrlypaZPn65QKOR1ar4WDod17bXX6txzz9X48eMlfXJ8paWlafDgwVHbcnzhaEn4y/Lh2Dd9+vTI3xMnTlRJSYlGjRql3/zmN5o7d66HmQHHnssvvzzy94QJEzRx4kSddNJJWrFihS644AIPM/O3yspKvfHGG4zHRkI5JlsSs7OzFQwGY2aAtbe3Kzc316Os4NTgwYN1yimnqKWlxetU8Dn2H08ca8nrxBNPVHZ2Nsebh6qqqvT0009r+fLlOv744yPrc3Nz1d3drV27dkVtz/GFo+WYLBLT0tJUXFysxsbGyLpwOKzGxkaVlpZ6mBmc2Lt3rzZt2qS8vDyvU8HnGD16tHJzc6OOtY6ODr388ssca0nivffe086dOznePGDbtqqqqvTkk0/qhRde0OjRo6PuLy4uVmpqatTxtXHjRm3ZsoXjC0fFMdvdXF1drdmzZ2vSpEmaPHmylixZos7OTlVUVHidGj7jhhtu0MUXX6xRo0Zp69atqqmpUTAY1BVXXOF1atAnRfuBrUybN2/WunXrNHToUJ1wwgm69tpr9ZOf/EQnn3yyRo8erVtuuUUjRozQjBkzvEvaxw71fg0dOlQLFy7UZZddptzcXG3atEk33XSTxowZo/Lycg+z9qfKyko99thj+v3vf69BgwZFxhlmZWWpX79+ysrK0ty5c1VdXa2hQ4cqMzNTV199tUpLS3X22Wd7nD18wevp1Sb9/Oc/t0844QQ7LS3Nnjx5sr169WqvU0IfZs6caefl5dlpaWn2yJEj7ZkzZ9otLS1ep4VPLV++3JYUs8yePdu27U9Og3PLLbfYOTk5dnp6un3BBRfYGzdu9DZpHzvU+/XRRx/ZF154oX3cccfZqamp9qhRo+x58+bZbW1tXqftS329T5Lshx9+OLLNxx9/bF911VX2kCFD7P79+9tf+9rX7Pfff9+7pOErlm3b9tEvTQEAAJDIjskxiQAAADgyFIkAAACIQZEIAACAGBSJAAAAiEGRCAAAgBgUiQAAAIhBkQgAAIAYFIkAAACIQZEIwBXvvPOOLMvSunXrvE7lsCR7/gDgNopEAJ9rzpw5siwrsgwbNkzTpk3Ta6+9FtkmPz9f77//vsaPH+9hpgAAt1AkAnBk2rRpev/99/X++++rsbFRKSkpuuiiiyL3B4NB5ebmKiUlxcMsE093d7fXKQDAYaFIBOBIenq6cnNzlZubq6KiIv3whz9Ua2urtm/fLim2u3bFihWyLEuNjY2aNGmS+vfvr3POOUcbN2486D72x/jd736n888/X/3791dhYaGampoi29x6660qKiqKetySJUtUUFAQuT1nzhzNmDFDd9xxh3JycjR48GDddttt6u3t1Y033qihQ4fq+OOP18MPPxyTw4YNG3TOOecoIyND48eP18qVK6Puf+ONNzR9+nQNHDhQOTk5+va3v60dO3ZE7p86daqqqqp07bXXKjs7W+Xl5U5fYgBIKBSJAOK2d+9e/epXv9KYMWM0bNiwQ257880366c//anWrl2rlJQUfec73/nc+DfffLNuuOEGrVu3TqeccoquuOIK9fb2xpXjCy+8oK1bt2rVqlVavHixampqdNFFF2nIkCF6+eWX9b3vfU/f/e539d5770U97sYbb9T111+vV199VaWlpbr44ou1c+dOSdKuXbv0pS99SWeccYbWrl2rhoYGtbe365/+6Z+iYjz66KNKS0vTn//8Z9XX18eVNwAkDBsAPsfs2bPtYDBoDxgwwB4wYIAtyc7Ly7Obm5sj22zevNmWZL/66qu2bdv28uXLbUn2888/H9nmmWeesSXZH3/8cZ/72R/jF7/4RWTdm2++aUuy169fb9u2bdfU1NiFhYVRj/vZz35mjxo1KirfUaNG2aFQKLLu1FNPtadMmRK53dvbaw8YMMD+9a9/HbXvRYsWRbbp6emxjz/+ePvOO++0bdu2b7/9dvvCCy+M2ndra6styd64caNt27Z93nnn2WeccUbfLyQAJBFaEgE4cv7552vdunVat26d1qxZo/Lyck2fPl3vvvvuIR83ceLEyN95eXmSpG3btrn+mM86/fTTFQj8/SsuJydHEyZMiNwOBoMaNmxYTNzS0tLI3ykpKZo0aZLWr18vSfrf//1fLV++XAMHDowsY8eOlSRt2rQp8rji4uK4cgWARMQIcwCODBgwQGPGjInc/sUvfqGsrCw9+OCD+slPfnLQx6Wmpkb+tixLkhQOhw+5r0M9JhAIyLbtqO17enoOGWN/nL7WfV4uB9q7d68uvvhi3XnnnTH37S9mpU9eKwBIdrQkAjgslmUpEAjo448/Pqr7Pe6449TW1hZVKLp5bsPVq1dH/u7t7VVzc7NOO+00SdKZZ56pN998UwUFBRozZkzUQmEI4FhDkQjAka6uLrW1tamtrU3r16/X1VdfHWlZO5qmTp2q7du366677tKmTZtUV1en//7v/3Ytfl1dnZ588klt2LBBlZWV+vDDDyOTbSorK/XBBx/oiiuu0F/+8hdt2rRJzz33nCoqKhQKhVzLAQASAUUiAEcaGhqUl5envLw8lZSU6C9/+YueeOIJTZ069ajmcdppp+m+++5TXV2dCgsLtWbNGt1www2uxV+0aJEWLVqkwsJCvfjii/rDH/6g7OxsSdKIESP05z//WaFQSBdeeKEmTJiga6+9VoMHD44a/wgAxwLL/uzgHgAAAPgeP30BAAAQgyIRAAAAMSgSAQAAEIMiEQAAADEoEgEAABCDIhEAAAAxKBIBAAAQgyIRAAAAMSgSAQAAEIMiEQAAADEoEgEAABDj/wEZzb1+T3cEgQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 684.8x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# https://github.com/andrzejnovak/combine_postfits/tree/master/tests/fitDiags\n",
    "f = uproot.open('/home/anovak/cat/combine_postfits/tests/fitDiags/fit_diag_A.root')\n",
    "f.keys();\n",
    "\n",
    "tdir = f['shapes_fit_s/ptbin3pass2016']\n",
    "cov = tdir['total_covar'].to_hist()\n",
    "cov.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 582,
   "id": "c3743422",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "ColormeshArtists(pcolormesh=<matplotlib.collections.QuadMesh object at 0x7f239c4dfd90>, cbar=<matplotlib.colorbar.Colorbar object at 0x7f239cfff4f0>, text=[])"
      ]
     },
     "execution_count": 582,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAApoAAAGwCAYAAAAJ08UyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABOn0lEQVR4nO3de1iUdf7/8dcMCJgCJiaHUrG0NEXwyFpt6soV2tHN3cy1xMPPDoum0ppaedp2QzsYHSi/tp2+12Za+y0r26VVTN0SNUHWLLV0TSwFtBISE3Dm/v3hOjWByujMfQ/cz8d13VfMPZ+53+97sPHt5zQOwzAMAQAAAH7mtDoBAAAANE0UmgAAAAgICk0AAAAEBIUmAAAAAoJCEwAAAAFBoQkAAICAoNAEAABAQIRanUCgud1u7d+/X5GRkXI4HFanAwAAGsAwDH3//fdKSEiQ0+mffrFjx46ppqbGL9c6V2FhYYqIiLA6jYBr8oXm/v371a5dO6vTAAAAZ2Hfvn266KKLzvk6x44dU8cOLVVa7vJDVucuLi5Oe/bsafLFZpMvNCMjIyVJV4f9WqGOZhZnAwAAGuK4Uat1NW95/h4/VzU1NSotd2lvYaKiIq2dOVj5vVsden+pmpoaCs3G7uRweaijmUIdYRZnAwAAfOHvaW8tIx1qGWntVDq37DOVr8kXmgAAACe5DLdchvU52AWrzgEAgG24ZQTF4avc3FwlJiYqIiJCqamp2rRp0ynbfvrppxo+fLgSExPlcDiUk5NTp83J535+ZGZmetoMHDiwzvN33XWXT3lTaAIAAASxZcuWKSsrS3PmzFFRUZGSk5OVnp6u8vLyetsfPXpUF198sebPn6+4uLh623z88cc6cOCA51i5cqUk6be//a1XuwkTJni1e+SRR3zKnaFzAABgG265ZfXAta8ZLFy4UBMmTNDYsWMlSYsWLdJ7772nF198UTNmzKjTvm/fvurbt68k1fu8JF1wwQVej+fPn69LLrlEAwYM8Dp/3nnnnbJYbQh6NAEAgG24DCMoDkmqrKz0Oqqrq+vkW1NTo8LCQqWlpXnOOZ1OpaWlqaCgwC/vSU1Njf76179q3LhxdRZfvfrqq2rTpo26d++umTNn6ujRoz5dmx5NAAAAC/x8n+85c+Zo7ty5XucOHTokl8ul2NhYr/OxsbHasWOHX/JYvny5Dh8+rDFjxnid/93vfqcOHTooISFBW7du1fTp07Vz5069+eabDb42hSYAALCNs12M4+8cpBOb0UdFRXnOh4eHW5LPCy+8oKFDhyohIcHr/B133OH5OSkpSfHx8Ro8eLB2796tSy65pEHXptAEAAC24ZYhV5AUmlFRUV6FZn3atGmjkJAQlZWVeZ0vKys7p7mTJ+3du1erVq1qUC9lamqqJGnXrl0NLjSZowkAABCkwsLC1Lt3b+Xn53vOud1u5efnq3///ud8/Zdeeklt27bVddddd8a2xcXFkqT4+PgGX58eTQAAYBvBNHTeUFlZWcrIyFCfPn3Ur18/5eTkqKqqyrMKffTo0brwwguVnZ0t6cTins8++8zz89dff63i4mK1bNlSnTp1+jEPt1svvfSSMjIyFBrqXRLu3r1bS5Ys0bXXXquYmBht3bpVU6dO1dVXX60ePXo0OHcKTQAAYBs/XfVtZQ6+GDFihA4ePKjZs2ertLRUKSkpysvL8ywQKikpkdP54yD1/v371bNnT8/jxx57TI899pgGDBigNWvWeM6vWrVKJSUlGjduXJ2YYWFhWrVqlaeobdeunYYPH64HH3zQp9wdhmHxux1glZWVio6O1q/Cb+G7zgEAaCSOGzVaXf26KioqzjiPsSFO1gOfb49VZKS1Mwe//96tS7uW+e3eghk9mgAAwDbc/z2szsEuKDQBAIBtuIJg1bnV8c1EoQkAAGzDZZw4rM7BLiydpLBu3TrdcMMNSkhIkMPh0PLly+u02b59u2688UZFR0erRYsW6tu3r0pKSsxPFgAAAD6xtNCsqqpScnKycnNz631+9+7duuqqq9SlSxetWbNGW7du1axZsxQREWFypgAAoClwB8lhF5YOnQ8dOlRDhw495fMPPPCArr32Wj3yyCOecw3diR4AAODn3HLIJYflOdhF0H4zkNvt1nvvvadLL71U6enpatu2rVJTU+sdXv+p6upqVVZWeh0AAAAwX9AWmuXl5Tpy5Ijmz5+vIUOG6J///Kd+/etf6+abb9batWtP+brs7GxFR0d7jnbt2pmYNQAACGZuIzgOuwjaVedu94kZDDfddJOmTp0qSUpJSdH69eu1aNEiDRgwoN7XzZw5U1lZWZ7HlZWVFJsAAECS5AqCoXOr45spaAvNNm3aKDQ0VJdffrnX+a5du+rDDz885evCw8MVHh4e6PQAAABwBkFbaIaFhalv377auXOn1/nPP/9cHTp0sCgrAADQmNGjaS5LC80jR45o165dnsd79uxRcXGxWrdurfbt22vatGkaMWKErr76ag0aNEh5eXl69913vb4QHgAAoKHchkNuw+JV5xbHN5OlhebmzZs1aNAgz+OTcyszMjL08ssv69e//rUWLVqk7Oxs3XPPPbrsssv0f//3f7rqqqusShkAAAANZGmhOXDgQBnG6ZdejRs3TuPGjTMpIwAA0JQxdG6uoJ2jCQAA4G8uOeWyeHdHl6XRzUWhCQAAbMMIgjmaho3maAbthu0AAABo3OjRBAAAtsEcTXNRaAIAANtwGU65DIvnaNroKygZOgcAAEBA0KMJAABswy2H3Bb3s7llny5NCk0AAGAbzNE0F0PnAAAACAh6NAEAgG0Ex2Ighs4BAACanBNzNK0durY6vpkYOgcAAEBA0KMJAABswx0E33XOqnMAAIAmiDma5qLQBAAAtuGWk300TcQcTQAAAAQEPZoAAMA2XIZDLsPiDdstjm8mCk0AAGAbriBYDORi6BwAAAA4N/RoAgAA23AbTrktXnXuZtU5AABA08PQubkYOgcAAEBA0KMJAABswy3rV327LY1uLgpNAABgG8GxYbt9BpTtc6cAAAAwFT2aAADANoLju87t089nnzsFAAC255YjKA5f5ebmKjExUREREUpNTdWmTZtO2fbTTz/V8OHDlZiYKIfDoZycnDpt5s6dK4fD4XV06dLFq82xY8eUmZmpmJgYtWzZUsOHD1dZWZlPeVNoAgAA2zjZo2n14Ytly5YpKytLc+bMUVFRkZKTk5Wenq7y8vJ62x89elQXX3yx5s+fr7i4uFNet1u3bjpw4IDn+PDDD72enzp1qt5991298cYbWrt2rfbv36+bb77Zp9wZOgcAALBAZWWl1+Pw8HCFh4fXabdw4UJNmDBBY8eOlSQtWrRI7733nl588UXNmDGjTvu+ffuqb9++klTv8yeFhoaeshCtqKjQCy+8oCVLluhXv/qVJOmll15S165dtWHDBv3iF79o0D3SowkAAGzj5IbtVh+S1K5dO0VHR3uO7OzsOvnW1NSosLBQaWlpnnNOp1NpaWkqKCg4p/fiiy++UEJCgi6++GKNGjVKJSUlnucKCwtVW1vrFbdLly5q3769T3Hp0QQAALbhNhxyW72P5n/j79u3T1FRUZ7z9fVmHjp0SC6XS7GxsV7nY2NjtWPHjrPOITU1VS+//LIuu+wyHThwQPPmzdMvf/lLbdu2TZGRkSotLVVYWJhatWpVJ25paWmD41BoAgAAWCAqKsqr0DTT0KFDPT/36NFDqamp6tChg15//XWNHz/eb3EoNAEAgG24g+C7zn3ZsL1NmzYKCQmps9q7rKzstAt9fNWqVStdeuml2rVrlyQpLi5ONTU1Onz4sFevpq9xLX2n161bpxtuuEEJCQlyOBxavnz5Kdveddddp1yiDwAA0BBuwxkUR0OFhYWpd+/eys/P//Ee3G7l5+erf//+fntfjhw5ot27dys+Pl6S1Lt3bzVr1swr7s6dO1VSUuJTXEt7NKuqqpScnKxx48addrn8W2+9pQ0bNighIcHE7AAAAKyXlZWljIwM9enTR/369VNOTo6qqqo8q9BHjx6tCy+80LOYqKamRp999pnn56+//lrFxcVq2bKlOnXqJEn6wx/+oBtuuEEdOnTQ/v37NWfOHIWEhGjkyJGSpOjoaI0fP15ZWVlq3bq1oqKiNGnSJPXv37/BK84liwvNoUOHes0RqM/XX3+tSZMm6f3339d11113xmtWV1erurra8/jnWwcAAAD7cskh11lsmO7vHHwxYsQIHTx4ULNnz1ZpaalSUlKUl5fnWSBUUlIip/PHXtL9+/erZ8+ensePPfaYHnvsMQ0YMEBr1qyRJH311VcaOXKkvvnmG11wwQW66qqrtGHDBl1wwQWe1z3xxBNyOp0aPny4qqurlZ6ermeffdan3IN6jqbb7dbtt9+uadOmqVu3bg16TXZ2tubNmxfgzGA5p7UfEk2aw4L31jDMjylJboviArCMr0PXgcrBVxMnTtTEiRPrfe5k8XhSYmKijDN8ri5duvSMMSMiIpSbm6vc3NwG5/lzQb2P5oIFCxQaGqp77rmnwa+ZOXOmKioqPMe+ffsCmCEAAABOJWh7NAsLC/Xkk0+qqKhIDh96WE61qz4AAIBLvg9dByIHuwjaHs1//etfKi8vV/v27RUaGqrQ0FDt3btX9957rxITE61ODwAANEJWrzYPhqF7MwVtj+btt9/u9bVHkpSenq7bb7/ds8oKAADAFy7DKZfFhZ7V8c1kaaF55MgRz8agkrRnzx4VFxerdevWat++vWJiYrzaN2vWTHFxcbrsssvMThUAAAA+srTQ3Lx5swYNGuR5nJWVJUnKyMjQyy+/bFFWAACgqTLkkNviOZqGxfHNZGmhOXDgwDMuv/+pL7/8MnDJAACAJo+hc3PZ504BAABgqqBdDAQAAOBvbsMht2Ht0LXV8c1EoQkAAGzDJadcFg/oWh3fTPa5UwAAAJiKHk0AAGAbDJ2bi0ITAADYhltOuS0e0LU6vpnsc6cAAAAwFT2aAADANlyGQy6Lh66tjm8mCk0AAGAbzNE0F4UmAACwDcNwym3xN/MYfDMQAAAAcG7o0QQAALbhkkMuWTxH0+L4ZqLQBAAAtuE2rJ8j6TYsDW8qhs4BAAAQEPRo4tw5bTIE4DD/Ph1OC/4tGBJiekijttb0mJLs809tO3WfAGfgDoLFQFbHNxOFJgAAsA23HHJbPEfS6vhmsk9JDQAAAFPRowkAAGyDbwYyF4UmAACwDeZomss+dwoAAABT0aMJAABsw60g+K5zGy0GotAEAAC2YQTBqnODQhMAAKDpcRtB0KNpo8VAzNEEAABAQNCjCQAAbINV5+ai0AQAALbB0Lm57FNSAwAAwFT0aAIAANvgu87NRaEJAABsg6FzczF0DgAAgICgRxMAANgGPZrmotAEAAC2QaFpLkuHztetW6cbbrhBCQkJcjgcWr58uee52tpaTZ8+XUlJSWrRooUSEhI0evRo7d+/37qEAQAALJCbm6vExERFREQoNTVVmzZtOmXbTz/9VMOHD1diYqIcDodycnLqtMnOzlbfvn0VGRmptm3batiwYdq5c6dXm4EDB8rhcHgdd911l095W1poVlVVKTk5Wbm5uXWeO3r0qIqKijRr1iwVFRXpzTff1M6dO3XjjTdakCkAAGgKTvZoWn34YtmyZcrKytKcOXNUVFSk5ORkpaenq7y8vN72R48e1cUXX6z58+crLi6u3jZr165VZmamNmzYoJUrV6q2tlbXXHONqqqqvNpNmDBBBw4c8ByPPPKIT7lbOnQ+dOhQDR06tN7noqOjtXLlSq9zzzzzjPr166eSkhK1b9/ejBQBAEATYsj67YUMH9svXLhQEyZM0NixYyVJixYt0nvvvacXX3xRM2bMqNO+b9++6tu3ryTV+7wk5eXleT1++eWX1bZtWxUWFurqq6/2nD/vvPNOWaw2RKNadV5RUSGHw6FWrVqdsk11dbUqKyu9DgAAACm4ejR/Xq9UV1fXybempkaFhYVKS0vznHM6nUpLS1NBQYHf3peKigpJUuvWrb3Ov/rqq2rTpo26d++umTNn6ujRoz5dt9EsBjp27JimT5+ukSNHKioq6pTtsrOzNW/evLpPGCf+DWMah/k1vOFymR5TknTcxPf1vxwhIebHDAszPaYc9pgwbsXvU9J/PxdMDmnF/6chFvQpuMz/XAAam3bt2nk9njNnjubOnet17tChQ3K5XIqNjfU6Hxsbqx07dvglD7fbrSlTpujKK69U9+7dPed/97vfqUOHDkpISNDWrVs1ffp07dy5U2+++WaDr90oCs3a2lrdcsstMgxDzz333Gnbzpw5U1lZWZ7HlZWVdX6RAADAnoJp1fm+ffu8Os/Cw8MtySczM1Pbtm3Thx9+6HX+jjvu8PyclJSk+Ph4DR48WLt379Yll1zSoGsHfaF5ssjcu3evVq9efdreTOnEL8mqXxQAAAhuwVRoRkVFnbGuadOmjUJCQlRWVuZ1vqys7JzmTp40ceJErVixQuvWrdNFF1102rapqamSpF27djW40AzqOZoni8wvvvhCq1atUkxMjNUpAQAAmCYsLEy9e/dWfn6+55zb7VZ+fr769+9/1tc1DEMTJ07UW2+9pdWrV6tjx45nfE1xcbEkKT4+vsFxLO3RPHLkiHbt2uV5vGfPHhUXF6t169aKj4/Xb37zGxUVFWnFihVyuVwqLS2VdGKiapgV8+UAAECjFkw9mg2VlZWljIwM9enTR/369VNOTo6qqqo8q9BHjx6tCy+8UNnZ2ZJOLCD67LPPPD9//fXXKi4uVsuWLdWpUydJJ4bLlyxZorfffluRkZGeGis6OlrNmzfX7t27tWTJEl177bWKiYnR1q1bNXXqVF199dXq0aNHg3O3tNDcvHmzBg0a5Hl8cm5lRkaG5s6dq3feeUeSlJKS4vW6Dz74QAMHDjQrTQAA0EQYhkOGxYWmr/FHjBihgwcPavbs2SotLVVKSory8vI8C4RKSkrkdP44SL1//3717NnT8/ixxx7TY489pgEDBmjNmjWS5Fnz8vN66qWXXtKYMWMUFhamVatWeYradu3aafjw4XrwwQd9yt1hGBYsuzRRZWWloqOj9auw3yrU0cy8wHZadW6w6jxwQS34MHRaMKPGbdEKZbusOrfizxGrznGOjhs1Wl39uioqKs44j7EhTtYDV749UaEtrF3LcbyqWh/d9Izf7i2YBf1iIAAAAH9xy2H5hu1WxzcThSYAALCNxjhHszEL6lXnAAAAaLzo0QQAALbRGBcDNWYUmgAAwDYYOjcXhSYAALANejTNxRxNAAAABAQ9mgAAwDaMIBg6t1OPJoUmAACwDUOWfFdDnRzsgqFzAAAABAQ9mgAAwDbccsjBNwOZhkITAADYBqvOzcXQOQAAAAKCHk0AAGAbbsMhBxu2m4ZCEwAA2IZhBMGqcxstO7dPoelwSI6mPVPA0cyiX6fDJv8yCzX//XWc19z0mGoeYX7MI1Xmx5RkHP3B/KBW/P9ixd9qTgvu022jv72BRsI+hSYAALA9FgOZi0ITAADYBoWmuSg0AQCAbbAYyFxNe9IiAAAALEOPJgAAsA1WnZuLQhMAANjGiULT6jmaloY3FUPnAAAACAh6NAEAgG2w6txcFJoAAMA2jP8eVudgFwydAwAAICDo0QQAALbB0Lm5KDQBAIB9MHZuKgpNAABgH0HQoymr45uIOZoAAAAICJ8KzdraWg0ePFhffPFFoPIBAAAImJPfDGT1YRc+DZ03a9ZMW7duDVQuAAAAAcViIHP5PHR+22236YUXXghELgAAAGhCfF4MdPz4cb344otatWqVevfurRYtWng9v3DhwgZfa926dXr00UdVWFioAwcO6K233tKwYcM8zxuGoTlz5uj555/X4cOHdeWVV+q5555T586dfU0bAADgxEIcq3sUrY5vIp8LzW3btqlXr16SpM8//9zrOYfDtzeuqqpKycnJGjdunG6++eY6zz/yyCN66qmn9Morr6hjx46aNWuW0tPT9dlnnykiIsLX1AEAgM0FwxxJq+ObyedC84MPPvBb8KFDh2ro0KH1PmcYhnJycvTggw/qpptukiT97//+r2JjY7V8+XLdeuut9b6uurpa1dXVnseVlZV+yxcAAAANd9b7aO7atUu7d+/W1VdfrebNm8swDJ97NE9nz549Ki0tVVpamudcdHS0UlNTVVBQcMpCMzs7W/Pmzav7hGFIcvstvzMKCTEv1kkul/kxJRlu8/9p5mhm/hawjvOamx7THdfa9JjVbcy/z4hSa0YoHGXfmB/06A/mx7Sg+8Sw6PMIOCM2bDeVz4uBvvnmGw0ePFiXXnqprr32Wh04cECSNH78eN17771+S6y0tFSSFBsb63U+NjbW81x9Zs6cqYqKCs+xb98+v+UEAAAat5Orzq0+7MLnQnPq1Klq1qyZSkpKdN5553nOjxgxQnl5eX5N7myEh4crKirK6wAAAID5fB5//Oc//6n3339fF110kdf5zp07a+/evX5LLC4uTpJUVlam+Ph4z/mysjKlpKT4LQ4AALAZGw1dW83nHs2qqiqvnsyTvv32W4WHh/slKUnq2LGj4uLilJ+f7zlXWVmpjRs3qn///n6LAwAA7MPqIfOzHTrPzc1VYmKiIiIilJqaqk2bNp2y7aeffqrhw4crMTFRDodDOTk5Z3XNY8eOKTMzUzExMWrZsqWGDx+usrIyn/L2udD85S9/qf/93//1PHY4HHK73XrkkUc0aNAgn6515MgRFRcXq7i4WNKJBUDFxcUqKSmRw+HQlClT9Kc//UnvvPOOPvnkE40ePVoJCQlee20CAAA0mBEkhw+WLVumrKwszZkzR0VFRUpOTlZ6errKy8vrbX/06FFdfPHFmj9/vmeE+GyuOXXqVL377rt64403tHbtWu3fv7/e7ShPx2EYvi1H3LZtmwYPHqxevXpp9erVuvHGG/Xpp5/q22+/1UcffaRLLrmkwddas2ZNvcVpRkaGXn75Zc+G7YsXL9bhw4d11VVX6dlnn9Wll17a4BiVlZWKjo7Wr8J+q1BHswa/7pyx6jygrFh17jy/lekx7bPqvMr0mJI1q84NVp0HjgWfRQic40aNVle/roqKCr+stzhZD7RbNEfO5tbuxe3+4Zj23TWvwfeWmpqqvn376plnnjnxerdb7dq106RJkzRjxozTvjYxMVFTpkzRlClTfLpmRUWFLrjgAi1ZskS/+c1vJEk7duxQ165dVVBQoF/84hcNulefezS7d++uzz//XFdddZVuuukmVVVV6eabb9aWLVt8KjIlaeDAgTIMo87x8ssvSzrRW/rHP/5RpaWlOnbsmFatWuVTkQkAAODNESTHieL3p8dP9wE/qaamRoWFhV7bPTqdTqWlpamgoOCs3oGGXLOwsFC1tbVebbp06aL27dv7FPesuoWio6P1wAMPnM1LAQAArBNE+2i2a9fO6/ScOXM0d+5cr3OHDh2Sy+Wqd7vHHTt2nFX4hlyztLRUYWFhatWqVZ02p9tm8ufOqtD87rvv9MILL2j79u2SpMsvv1xjx45V69bmD/MBAAA0Rvv27fMaOvfnoupg4fPQ+bp165SYmKinnnpK3333nb777js99dRT6tixo9atWxeIHAEAAPzD6kVAP+lR/fm+3/UVmm3atFFISEid1d5lZWWnXOhzJg25ZlxcnGpqanT48OFziutzoZmZmakRI0Zoz549evPNN/Xmm2/qP//5j2699VZlZmb6ejkAAADzGI7gOBooLCxMvXv39tru0e12Kz8//6y3e2zINXv37q1mzZp5tdm5c6dKSkp8iuvz0PmuXbv0t7/9TSE/WVUdEhKirKwsr22PAAAAcO6ysrKUkZGhPn36qF+/fsrJyVFVVZXGjh0rSRo9erQuvPBCZWdnSzqx2Oezzz7z/Pz111+ruLhYLVu2VKdOnRp0zejoaI0fP15ZWVlq3bq1oqKiNGnSJPXv37/BK86lsyg0e/Xqpe3bt+uyyy7zOr99+3YlJyf7ejkAAADTGIYlO37VycEXI0aM0MGDBzV79myVlpYqJSVFeXl5nsU8JSUlcjp/HKTev3+/evbs6Xn82GOP6bHHHtOAAQO0Zs2aBl1Tkp544gk5nU4NHz5c1dXVSk9P17PPPutT7g3aR3Pr1q2en7dv36777rtPkyZN8lS0GzZsUG5urubPn68RI0b4lECgsY9m4LGPZuCwj2ZgsY9mAEOyjybOUaD20bzo6XlBsY/mV5Pm+O3eglmD/rZOSUmRw+HQT2vS++67r0673/3ud0FXaAIAAMAaDSo09+zZE+g8AAAAAs/HxTgBy8EmGlRodujQIdB5AAAABJzDOHFYnYNdnNVEt/379+vDDz9UeXm53G6313P33HOPXxIDAADwuyD6ZiA78LnQfPnll3XnnXcqLCxMMTExcjh+7P51OBwUmgAAAJB0FoXmrFmzNHv2bM2cOdNrKT0AAEDQY46mqXwuNI8ePapbb72VIhMAADQ+DJ2byudCc/z48XrjjTc0Y8aMQOSDc2DJvnUWxXVYsE+pK6GN6TG/HhRpesxjPc3f57HFpvNNjylJCavNj2mbvTsdFvTYWNH/wd6dwGn5XGhmZ2fr+uuvV15enpKSktSsmfcm6AsXLvRbcgAAAH5Fj6apzqrQfP/99z1fQfnzxUAAAABBi0LTVD4Xmo8//rhefPFFjRkzJgDpAAAAoKnwudAMDw/XlVdeGYhcAAAAAotV56byeer05MmT9fTTTwciFwAAgIA6+c1AVh924XOP5qZNm7R69WqtWLFC3bp1q7MY6M033/RbcgAAAGi8fC40W7VqpZtvvjkQuQAAAAQWi4FM5XOh+dJLLwUiDwAAADQxPheaAAAAjZVD1s+RtM9SoLMoNDt27Hja/TL/85//nFNCAAAAaBp8LjSnTJni9bi2tlZbtmxRXl6epk2b5q+8AAAA/I/tjUzlc6E5efLkes/n5uZq8+bN55wQAABAwLAYyFQ+76N5KkOHDtX//d//+etyAAAAaOT8thjob3/7m1q3bu2vywEAAPgfPZqm8rnQ7Nmzp9diIMMwVFpaqoMHD+rZZ5/1a3IAAAD+FAzfzGN1fDP5XGgOGzbM67HT6dQFF1yggQMHqkuXLv7KCwAAAI2cz4XmnDlzApEHAABA4DF0bqqzmqPpdru1a9culZeXy+12ez139dVX+yUxSXK5XJo7d67++te/qrS0VAkJCRozZowefPDB0+7lCQAAUC8KTVP5XGhu2LBBv/vd77R3714Zhvc75XA45HK5/JbcggUL9Nxzz+mVV15Rt27dtHnzZo0dO1bR0dG65557/BYHAAAA/udzoXnXXXepT58+eu+99xQfHx/QnsX169frpptu0nXXXSdJSkxM1GuvvaZNmzYFLCYAAGi6WAxkLp8LzS+++EJ/+9vf1KlTp0Dk4+WKK67Q4sWL9fnnn+vSSy/Vv//9b3344YdauHDhKV9TXV2t6upqz+PKysqA5wkAABoJvhnIVD4Xmqmpqdq1a5cpheaMGTNUWVmpLl26KCQkRC6XS3/+8581atSoU74mOztb8+bNq/uEwyE5/LY/fXCy6P4cIRYEDTH/Xo/FNzc9ZshV35ke84u+S02P2VljTI8pScc+a2F6zOZVx0yPqWPVZ27jZw6n+f+PGrW1pse01WS7poI5mqbyudCcNGmS7r33XpWWliopKUnNmjXzer5Hjx5+S+7111/Xq6++qiVLlqhbt24qLi7WlClTlJCQoIyMjHpfM3PmTGVlZXkeV1ZWql27dn7LCQAAAA3jc6E5fPhwSdK4ceM85xwOhwzD8PtioGnTpmnGjBm69dZbJUlJSUnau3evsrOzT1lohoeHKzw83G85AACApoM5mubyudDcs2dPIPKo19GjR+X82fBLSEhInS2VAAAAGoShc1P5XGh26NAhEHnU64YbbtCf//xntW/fXt26ddOWLVu0cOFCr95UAAAABKez2rDdLE8//bRmzZql3//+9yovL1dCQoLuvPNOzZ492+rUAABAYxQEQ+f0aAaJyMhI5eTkKCcnx+pUAABAU8DQuama+H4/AAAAsEpQ92gCAAD4FT2apjrrQrOmpkbl5eV1VoC3b9/+nJMCAAAIBLY3MtdZfQXluHHjtH79eq/zgdhHEwAAAI2Xz3M0x4wZI6fTqRUrVqiwsFBFRUUqKirSli1bVFRUFIgcAQAAbC03N1eJiYmKiIhQamqqNm3adNr2b7zxhrp06aKIiAglJSXp73//u9fzDoej3uPRRx/1tElMTKzz/Pz5833K2+cezeLiYhUWFqpLly6+vhQAAMBajXCO5rJly5SVlaVFixYpNTVVOTk5Sk9P186dO9W2bds67devX6+RI0cqOztb119/vZYsWaJhw4apqKhI3bt3lyQdOHDA6zX/+Mc/NH78eM83QJ70xz/+URMmTPA8joyM9Cl3n3s0L7/8ch06dMjXlwEAAFju5BxNqw9Jqqys9Dqqq6vrzXnhwoWaMGGCxo4dq8svv1yLFi3SeeedpxdffLHe9k8++aSGDBmiadOmqWvXrnrooYfUq1cvPfPMM542cXFxXsfbb7+tQYMG6eKLL/a6VmRkpFe7Fi1a+PR++1xoLliwQPfdd5/WrFmjb775ps6bBAAAgDNr166doqOjPUd2dnadNjU1NSosLFRaWprnnNPpVFpamgoKCuq9bkFBgVd7SUpPTz9l+7KyMr333nsaP358nefmz5+vmJgY9ezZU48++qiOHz/uyy36PnR+MvHBgwd7nWcxEAAAaBSsHjr/r3379ikqKsrzODw8vE6bQ4cOyeVyKTY21ut8bGysduzYUe91S0tL621fWlpab/tXXnlFkZGRuvnmm73O33PPPerVq5dat26t9evXa+bMmTpw4IAWLlzYoPuTzqLQ/OCDD3x9CQAAQHAIojmaUVFRXoWmVV588UWNGjVKERERXuezsrI8P/fo0UNhYWG68847lZ2dXW9RXB+fC80BAwb4+hIAAACchTZt2igkJERlZWVe58vKyhQXF1fva+Li4hrc/l//+pd27typZcuWnTGX1NRUHT9+XF9++aUuu+yyBuXfoEJz69at6t69u5xOp7Zu3Xratj169GhQYAAAALM1tg3bw8LC1Lt3b+Xn52vYsGGSJLfbrfz8fE2cOLHe1/Tv31/5+fmaMmWK59zKlSvVv3//Om1feOEF9e7dW8nJyWfMpbi4WE6ns96V7qfSoEIzJSVFpaWlatu2rVJSUuRwOGQYdd8l5mgCAICgFkRD5w2VlZWljIwM9enTR/369VNOTo6qqqo0duxYSdLo0aN14YUXehYTTZ48WQMGDNDjjz+u6667TkuXLtXmzZu1ePFir+tWVlbqjTfe0OOPP14nZkFBgTZu3KhBgwYpMjJSBQUFmjp1qm677Tadf/75Dc69QYXmnj17dMEFF3h+BgAAgDlGjBihgwcPavbs2SotLVVKSory8vI8C35KSkrkdP64kdAVV1yhJUuW6MEHH9T999+vzp07a/ny5Z49NE9aunSpDMPQyJEj68QMDw/X0qVLNXfuXFVXV6tjx46aOnWq17zNhnAY9XVNNiGVlZWKjo7Wr8JvUagjzLzATod5sf7LqPVtywH/BXafuY2fORo4Cdmffhjc/cyN/Ozw//ve9JjFfZeaHrPzmjGmx5Skdi+HmB6z+e5vTI9pHPrW9JhWMGprzQ/qMv/zzy6OGzVaXf26Kioq/LJg5mQ9cOkfHlZIeMSZXxBArupj+vyx+/12b8HM58VA33zzjWJiYiSdWJb//PPP64cfftCNN96oX/7yl35PsNFymF9oOsKamR5TkhxOn7djPfeYFhSaYd/VmB6z9uOGD0/4S4puNT1mxJbmpseUpPBDFeYH/eGY+THdFhRDVvRhNO1+E/hLIxw6b8waXCF88sknSkxMVNu2bdWlSxcVFxerb9++euKJJ7R48WINGjRIy5cvD2CqAAAAaEwaXGjed999SkpK0rp16zRw4EBdf/31uu6661RRUaHvvvtOd955p89ftA4AAGAqI0gOm2jw0PnHH3+s1atXq0ePHkpOTtbixYv1+9//3jP5dNKkSfrFL34RsEQBAADOVWPb3qixa3Ch+e2333o2+mzZsqVatGjhtbz9/PPP1/ffm79wAQAAoMGCoUfR6vgm8mkVh+NnC1x+/hgAAAA4yadV52PGjPF8t+WxY8d01113qUWLFpKk6upq/2cHAADgT/RomqrBhWZGRobX49tuu61Om9GjR597RgAAAAHCHE1zNbjQfOmllwKZBwAAAJoYnzdsBwAAaLQYOjcVhSYAALANhs7NZf53BwIAAMAW6NEEAAD2wdC5qSg0AQCAfVBomoqhcwAAAAQEPZoAAMA2HP89rM7BLig0AQCAfTB0bioKTQAAYBtsb2SuoJ+j+fXXX+u2225TTEyMmjdvrqSkJG3evNnqtAAAAHAGQd2j+d133+nKK6/UoEGD9I9//EMXXHCBvvjiC51//vlWpwYAABojhs5NFdSF5oIFC9SuXTuv71nv2LHjaV9TXV2t6upqz+PKysqA5QcAABohGxV6VgvqQvOdd95Renq6fvvb32rt2rW68MIL9fvf/14TJkw45Wuys7M1b968OucdYc3kcDQLZLreQkLMi/VfjrgLTI8pScfPP8/0mO4w899fZ43L9JgXrjlqesyaf0eaHjPiwPemx5QkZ+m3psc0jv5gekwZ5v+tarjdpscEEHyCeo7mf/7zHz333HPq3Lmz3n//fd19992655579Morr5zyNTNnzlRFRYXn2Ldvn4kZAwCAYHZyMZDVh10EdY+m2+1Wnz599PDDD0uSevbsqW3btmnRokXKyMio9zXh4eEKDw83M00AANBYMEfTVEHdoxkfH6/LL7/c61zXrl1VUlJiUUYAAABoqKDu0bzyyiu1c+dOr3Off/65OnToYFFGAACgMQuGoWur45spqHs0p06dqg0bNujhhx/Wrl27tGTJEi1evFiZmZlWpwYAABojI0gOmwjqQrNv375666239Nprr6l79+566KGHlJOTo1GjRlmdGgAAAM4gqIfOJen666/X9ddfb3UaAACgCWDo3FxBX2gCAAD4TTAMXVsd30QUmgAAwD4oNE0V1HM0AQAA0HjRowkAAGyDOZrmotAEAAD2wdC5qRg6BwAAQEDQowkAAGzDYRhyGNZ2KVod30wUmgAAwD4YOjcVQ+cAAAAICHo0AQCAbbDq3Fz0aAIAAPswguTwUW5urhITExUREaHU1FRt2rTptO3feOMNdenSRREREUpKStLf//53r+fHjBkjh8PhdQwZMsSrzbfffqtRo0YpKipKrVq10vjx43XkyBGf8qbQBAAACGLLli1TVlaW5syZo6KiIiUnJys9PV3l5eX1tl+/fr1Gjhyp8ePHa8uWLRo2bJiGDRumbdu2ebUbMmSIDhw44Dlee+01r+dHjRqlTz/9VCtXrtSKFSu0bt063XHHHT7lTqEJAABs4+TQudWHLxYuXKgJEyZo7Nixuvzyy7Vo0SKdd955evHFF+tt/+STT2rIkCGaNm2aunbtqoceeki9evXSM88849UuPDxccXFxnuP888/3PLd9+3bl5eXpL3/5i1JTU3XVVVfp6aef1tKlS7V///4G506hCQAA7MPqIfOfDJ1XVlZ6HdXV1XXSrampUWFhodLS0jznnE6n0tLSVFBQUO8tFhQUeLWXpPT09Drt16xZo7Zt2+qyyy7T3XffrW+++cbrGq1atVKfPn0859LS0uR0OrVx48Z649bHNouBHOHhcjjDrE4joI52bGVJ3MrEZqbHNBymh1TMp8dMjxm2ba/pMUN++MH0mHK5zY8pye204A+Sw4KYx4+bHzMkxPyYQAME02Kgdu3aeZ2fM2eO5s6d63Xu0KFDcrlcio2N9TofGxurHTt21Hv90tLSetuXlpZ6Hg8ZMkQ333yzOnbsqN27d+v+++/X0KFDVVBQoJCQEJWWlqpt27Ze1wgNDVXr1q29rnMmtik0AQAAgsm+ffsUFRXleRweHm5a7FtvvdXzc1JSknr06KFLLrlEa9as0eDBg/0Wh6FzAABgH1YPmf9k6DwqKsrrqK/QbNOmjUJCQlRWVuZ1vqysTHFxcfXeYlxcnE/tJeniiy9WmzZttGvXLs81fr7Y6Pjx4/r2229Pe52fo9AEAAC20pgWAoWFhal3797Kz8/3nHO73crPz1f//v3rfU3//v292kvSypUrT9lekr766it98803io+P91zj8OHDKiws9LRZvXq13G63UlNTG5w/hSYAAEAQy8rK0vPPP69XXnlF27dv1913362qqiqNHTtWkjR69GjNnDnT037y5MnKy8vT448/rh07dmju3LnavHmzJk6cKEk6cuSIpk2bpg0bNujLL79Ufn6+brrpJnXq1Enp6emSpK5du2rIkCGaMGGCNm3apI8++kgTJ07UrbfeqoSEhAbnzhxNAABgH4Zx4rA6Bx+MGDFCBw8e1OzZs1VaWqqUlBTl5eV5FvyUlJTI6fyx7/CKK67QkiVL9OCDD+r+++9X586dtXz5cnXv3l2SFBISoq1bt+qVV17R4cOHlZCQoGuuuUYPPfSQ1/D9q6++qokTJ2rw4MFyOp0aPny4nnrqKZ9ydxiG1e92YFVWVio6OlppbcYptKmvOu/T0ZK4rDoPHCtWnbtttOpcrDoPHLusOnc36b9CLXXcqNHq6tdVUVHhtWDmbJ2sB/r85k8KbRbhhwzP3vHaY9r8twf9dm/BjKFzAAAABARD5wAAwD5+surb0hxsgkITAADYhsN94rA6B7tg6BwAAAABQY8mAACwD4bOTUWhCQAAbCOYvuvcDig0AQCAfTTCfTQbM+ZoAgAAICDo0QQAALbB0Lm5KDQBAIB9sBjIVAydAwAAICAaVaE5f/58ORwOTZkyxepUAABAI3Ry6Nzqwy4azdD5xx9/rP/5n/9Rjx49rE4FAAA0Vqw6N1Wj6NE8cuSIRo0apeeff17nn3++1ekAAACgARpFoZmZmanrrrtOaWlpZ2xbXV2tyspKrwMAAECyfsicofMgs3TpUhUVFenjjz9uUPvs7GzNmzev7hMOx4nDLBZ0ixshJt7fT+Na8M8Vh9v8mM4al+kxjepqe8R0mf/eSpIjtJn5McPMj2nF++sICTE9ptw2+tsbZ49V56YK6h7Nffv2afLkyXr11VcVERHRoNfMnDlTFRUVnmPfvn0BzhIAAAD1CeoezcLCQpWXl6tXr16ecy6XS+vWrdMzzzyj6upqhfzsX83h4eEKDw83O1UAANAIBMPQtdXxzRTUhebgwYP1ySefeJ0bO3asunTpounTp9cpMgEAAE7LbVg/zcLq+CYK6kIzMjJS3bt39zrXokULxcTE1DkPAABwRszRNFVQz9EEAABA4xXUPZr1WbNmjdUpAACARsoh6+dIWrNHjDUaXaEJAABw1vhmIFMxdA4AAICAoEcTAADYBtsbmYtCEwAA2Aerzk3F0DkAAAACgh5NAABgGw7DkMPixThWxzcThSYAALAP938Pq3OwCYbOAQAAEBD0aAIAANtg6NxcFJoAAMA+WHVuKgpNAABgH3wzkKmYowkAAICAoEcTAADYBt8MZC4KTQAAYB8MnZvKNoWmUVsrw+GwOo2Aiij9waLIzU2P6HCZ/z9p6HdHTY/pdluw2ZrD/Bk1jhDTQ56I67TgM8GKv2As+J3KbZ+/SAGcmm0KTQAAAIf7xGF1DnZBoQkAAOyDoXNTseocAAAAAUGhCQAA7MMIksNHubm5SkxMVEREhFJTU7Vp06bTtn/jjTfUpUsXRUREKCkpSX//+989z9XW1mr69OlKSkpSixYtlJCQoNGjR2v//v1e10hMTJTD4fA65s+f71PeFJoAAMA2Tn4FpdWHL5YtW6asrCzNmTNHRUVFSk5OVnp6usrLy+ttv379eo0cOVLjx4/Xli1bNGzYMA0bNkzbtm2TJB09elRFRUWaNWuWioqK9Oabb2rnzp268cYb61zrj3/8ow4cOOA5Jk2a5FPuzNEEAACwQGVlpdfj8PBwhYeH12m3cOFCTZgwQWPHjpUkLVq0SO+9955efPFFzZgxo077J598UkOGDNG0adMkSQ899JBWrlypZ555RosWLVJ0dLRWrlzp9ZpnnnlG/fr1U0lJidq3b+85HxkZqbi4uLO+R3o0AQCAfZxcDGT1Ialdu3aKjo72HNnZ2XXSrampUWFhodLS0jznnE6n0tLSVFBQUO8tFhQUeLWXpPT09FO2l6SKigo5HA61atXK6/z8+fMVExOjnj176tFHH9Xx48cb+k5LokcTAADYiSHJ6u2F/jtyvm/fPkVFRXlO19ebeejQIblcLsXGxnqdj42N1Y4dO+q9fGlpab3tS0tL621/7NgxTZ8+XSNHjvTK55577lGvXr3UunVrrV+/XjNnztSBAwe0cOHCBt2mRKEJAABs5GzmSAYiB0mKioryKuysUFtbq1tuuUWGYei5557zei4rK8vzc48ePRQWFqY777xT2dnZ9RbF9WHoHAAAIEi1adNGISEhKisr8zpfVlZ2yrmTcXFxDWp/ssjcu3evVq5cecaiNzU1VcePH9eXX37Z4PwpNAEAgH0Ysn5+pg8dqmFhYerdu7fy8/M959xut/Lz89W/f/96X9O/f3+v9pK0cuVKr/Yni8wvvvhCq1atUkxMzBlzKS4ultPpVNu2bRucP0PnAADAPhrhNwNlZWUpIyNDffr0Ub9+/ZSTk6OqqirPKvTRo0frwgsv9Cwmmjx5sgYMGKDHH39c1113nZYuXarNmzdr8eLFkk4Umb/5zW9UVFSkFStWyOVyeeZvtm7dWmFhYSooKNDGjRs1aNAgRUZGqqCgQFOnTtVtt92m888/v8G5U2gCAAAEsREjRujgwYOaPXu2SktLlZKSory8PM+Cn5KSEjmdPw5SX3HFFVqyZIkefPBB3X///ercubOWL1+u7t27S5K+/vprvfPOO5KklJQUr1gffPCBBg4cqPDwcC1dulRz585VdXW1OnbsqKlTp3rN22wIh2FYXdYHVmVlpaKjozW41WiFOsKsTieg3J3aWRL3WFxz02M6XOb/sT1vz2HTY7pLvjY9plFTa3pMGdYsAXWEhJgf1IKYRq1v25H4gyXvLZqU40aNVle/roqKCr8smDlZD/wqabpCQxq2kCVQjruqtfqTBX67t2BGjyYAALCNYFp1bgcsBgIAAEBA0KMJAADsoxEuBmrMKDQBAIB9UGiaKuiHzrOzs9W3b19FRkaqbdu2GjZsmHbu3Gl1WgAAADiDoC80165dq8zMTG3YsEErV65UbW2trrnmGlVVVVmdGgAAaGws36w9CHpUTRT0Q+d5eXlej19++WW1bdtWhYWFuvrqq+u0r66uVnV1tedxZWVlwHMEAACNhFuSIwhysImgLzR/rqKiQtKJnevrk52drXnz5tU5b1RXy3CY9y8IK/aQc+7aZ3pMSWqxt9H9MTorVuxFaAVHM3v8Pi3jNr8ngz0tgR+xvZG5gn7o/KfcbremTJmiK6+80rO7/c/NnDlTFRUVnmPfPmuKLwAAALtrVF0XmZmZ2rZtmz788MNTtgkPD1d4uLU7/gMAgCAVDHMkrY5vokZTaE6cOFErVqzQunXrdNFFF1mdDgAAaIzchmTiVLpT5mATQV9oGoahSZMm6a233tKaNWvUsWNHq1MCAABAAwR9oZmZmaklS5bo7bffVmRkpEpLSyVJ0dHRat68ucXZAQCARoWhc1MFfaH53HPPSZIGDhzodf6ll17SmDFjzE8IAAA0YkFQaMrq+OYJ+kLTsPwPAwAAAM5G0BeaAAAAfsPQuakoNAEAgH24DVk+dG2jVeeNasN2AAAANB70aAIAAPsw3CcOq3OwCQpNAABgH8zRNBWFJgAAsA/maJqKOZoAAAAICHo0AQCAfTB0bioKTQAAYB+GrC/07FNnMnQOAACAwKBHEwAA2AdD56ai0AQAAPbhdkuyeB9Lt3320WToHAAAAAFBjyYAALAPhs5NRaEZIIbLZXpMR6hFv06Hw/yYVvxPasHvVC6bDK/Y6OvY5LBgIMmK99eK+wQagkLTVHwSAAAAICDo0QQAAPbBV1CaikITAADYhmG4ZVg8Xcfq+Gai0AQAAPZhGNb3KDJHEwAAADg39GgCAAD7MIJgjqaNejQpNAEAgH243ZLD4jmSNpqjydA5AAAAAoIeTQAAYB8MnZuKHk0AAGAbhtsdFIevcnNzlZiYqIiICKWmpmrTpk2nbf/GG2+oS5cuioiIUFJSkv7+9797vw+GodmzZys+Pl7NmzdXWlqavvjiC6823377rUaNGqWoqCi1atVK48eP15EjR3zKm0ITAAAgiC1btkxZWVmaM2eOioqKlJycrPT0dJWXl9fbfv369Ro5cqTGjx+vLVu2aNiwYRo2bJi2bdvmafPII4/oqaee0qJFi7Rx40a1aNFC6enpOnbsmKfNqFGj9Omnn2rlypVasWKF1q1bpzvuuMOn3B2G0bT7bysrKxUdHa1fNR+hUEeY1ekElCM83Jq4Yc3MD2rBH1vjWLX5MWtqTY9pCRtNjOe7zoGGOW7UaHX166qoqFBUVNQ5Xy+Y6oHjRo1W/7CswfeWmpqqvn376plnnpEkud1utWvXTpMmTdKMGTPqtB8xYoSqqqq0YsUKz7lf/OIXSklJ0aJFi2QYhhISEnTvvffqD3/4gySpoqJCsbGxevnll3Xrrbdq+/btuvzyy/Xxxx+rT58+kqS8vDxde+21+uqrr5SQkNCge+WTAAAA2IfbCI5DJ4rfnx7V1XU7NGpqalRYWKi0tDTPOafTqbS0NBUUFNR7iwUFBV7tJSk9Pd3Tfs+ePSotLfVqEx0drdTUVE+bgoICtWrVylNkSlJaWpqcTqc2btzY4LebQhMAAMAC7dq1U3R0tOfIzs6u0+bQoUNyuVyKjY31Oh8bG6vS0tJ6r1taWnra9if/e6Y2bdu29Xo+NDRUrVu3PmXc+rDqHAAA2IdhSLJ6H80TPZr79u3zGjoPt2gKXCBRaAIAANsw3IYMh7XLU04uj4mKijrjHM02bdooJCREZWVlXufLysoUFxdX72vi4uJO2/7kf8vKyhQfH+/VJiUlxdPm54uNjh8/rm+//faUcevD0DkAALAPwx0cRwOFhYWpd+/eys/P95xzu93Kz89X//79631N//79vdpL0sqVKz3tO3bsqLi4OK82lZWV2rhxo6dN//79dfjwYRUWFnrarF69Wm63W6mpqQ3Ov1EUmr7uHQUAANBUZGVl6fnnn9crr7yi7du36+6771ZVVZXGjh0rSRo9erRmzpzpaT958mTl5eXp8ccf144dOzR37lxt3rxZEydOlCQ5HA5NmTJFf/rTn/TOO+/ok08+0ejRo5WQkKBhw4ZJkrp27aohQ4ZowoQJ2rRpkz766CNNnDhRt956a4NXnEuNYOj85N5RixYtUmpqqnJycpSenq6dO3fWmaQKAABwOsE0dN5QI0aM0MGDBzV79myVlpYqJSVFeXl5nsU8JSUlcjp/7Du84oortGTJEj344IO6//771blzZy1fvlzdu3f3tLnvvvtUVVWlO+64Q4cPH9ZVV12lvLw8RUREeNq8+uqrmjhxogYPHiyn06nhw4frqaee8in3oN9H09e9o34umPbNCjT20QxwSPbRDBz20Qws9tFEIxSofTQH6iaFOiz4e+snjhu1WqO3/XZvwSyoezRP7h310+7gM+0dVV1d7bUPVUVFhaQTv9SmzmE4rInrtuDfKlYUmkaNBTGb/p9bSbb63l/Jgv9PLXl/rfk8QtNx8u9tf/eHHVet5V91flw2+WxXkBeap9s7aseOHfW+Jjs7W/Pmzatzft2xNwOSY1D5weoEAADwr++//17R0dHnfJ2wsDDFxcXpw9K/n7mxCeLi4hQW1rRHWqUgLzTPxsyZM5WVleV5fPjwYXXo0EElJSV++YMarCorK9WuXbs6e3I1RXa5V+6z6bHLvXKfTY8V92oYhr7//nufFp6cTkREhPbs2aOaGvNHp+oTFhbmNR+yqQrqQvNs9o4KDw+vd8PT6OjoJv9BIDVsT66mwi73yn02PXa5V+6z6TH7Xv3dQRQREWGL4i6YBPVs7bPZOwoAAADBIah7NKUTe0dlZGSoT58+6tevn3Jycrz2jgIAAEBwCvpC80x7R51JeHi45syZ0yS/P/Sn7HKfkn3ulftseuxyr9xn02One4V/Bf0+mgAAAGicgnqOJgAAABovCk0AAAAEBIUmAAAAAoJCEwAAAAHR5AvN3NxcJSYmKiIiQqmpqdq0aZPVKflVdna2+vbtq8jISLVt21bDhg3Tzp07rU4r4ObPny+Hw6EpU6ZYnUpAfP3117rtttsUExOj5s2bKykpSZs3b7Y6Lb9yuVyaNWuWOnbsqObNm+uSSy7RQw895PfvNbbCunXrdMMNNyghIUEOh0PLly/3et4wDM2ePVvx8fFq3ry50tLS9MUXX1iT7Dk43X3W1tZq+vTpSkpKUosWLZSQkKDRo0dr//791iV8ls70+/ypu+66Sw6HQzk5Oabl508Nudft27frxhtvVHR0tFq0aKG+ffuqpKTE/GTRKDTpQnPZsmXKysrSnDlzVFRUpOTkZKWnp6u8vNzq1Pxm7dq1yszM1IYNG7Ry5UrV1tbqmmuuUVVVldWpBczHH3+s//mf/1GPHj2sTiUgvvvuO1155ZVq1qyZ/vGPf+izzz7T448/rvPPP9/q1PxqwYIFeu655/TMM89o+/btWrBggR555BE9/fTTVqd2zqqqqpScnKzc3Nx6n3/kkUf01FNPadGiRdq4caNatGih9PR0HTt2zORMz83p7vPo0aMqKirSrFmzVFRUpDfffFM7d+7UjTfeaEGm5+ZMv8+T3nrrLW3YsMFvX5lohTPd6+7du3XVVVepS5cuWrNmjbZu3apZs2bxbTs4NaMJ69evn5GZmel57HK5jISEBCM7O9vCrAKrvLzckGSsXbvW6lQC4vvvvzc6d+5srFy50hgwYIAxefJkq1Pyu+nTpxtXXXWV1WkE3HXXXWeMGzfO69zNN99sjBo1yqKMAkOS8dZbb3keu91uIy4uznj00Uc95w4fPmyEh4cbr732mgUZ+sfP77M+mzZtMiQZe/fuNSepADjVfX711VfGhRdeaGzbts3o0KGD8cQTT5iem7/Vd68jRowwbrvtNmsSQqPUZHs0a2pqVFhYqLS0NM85p9OptLQ0FRQUWJhZYFVUVEiSWrdubXEmgZGZmanrrrvO6/fa1Lzzzjvq06ePfvvb36pt27bq2bOnnn/+eavT8rsrrrhC+fn5+vzzzyVJ//73v/Xhhx9q6NChFmcWWHv27FFpaanXn+Ho6GilpqY26c8m6cTnk8PhUKtWraxOxa/cbrduv/12TZs2Td26dbM6nYBxu9167733dOmllyo9PV1t27ZVamrqaacSAE220Dx06JBcLledbxCKjY1VaWmpRVkFltvt1pQpU3TllVeqe/fuVqfjd0uXLlVRUZGys7OtTiWg/vOf/+i5555T586d9f777+vuu+/WPffco1deecXq1PxqxowZuvXWW9WlSxc1a9ZMPXv21JQpUzRq1CirUwuok58/dvpskqRjx45p+vTpGjlypKKioqxOx68WLFig0NBQ3XPPPVanElDl5eU6cuSI5s+fryFDhuif//ynfv3rX+vmm2/W2rVrrU4PQSrov4ISDZeZmalt27bpww8/tDoVv9u3b58mT56slStXNvm5QG63W3369NHDDz8sSerZs6e2bdumRYsWKSMjw+Ls/Of111/Xq6++qiVLlqhbt24qLi7WlClTlJCQ0KTuEycWBt1yyy0yDEPPPfec1en4VWFhoZ588kkVFRXJ4XBYnU5Aud1uSdJNN92kqVOnSpJSUlK0fv16LVq0SAMGDLAyPQSpJtuj2aZNG4WEhKisrMzrfFlZmeLi4izKKnAmTpyoFStW6IMPPtBFF11kdTp+V1hYqPLycvXq1UuhoaEKDQ3V2rVr9dRTTyk0NFQul8vqFP0mPj5el19+ude5rl27NrlVndOmTfP0aiYlJen222/X1KlTm3yP9cnPH7t8Np0sMvfu3auVK1c2ud7Mf/3rXyovL1f79u09n0179+7Vvffeq8TERKvT86s2bdooNDTUFp9P8J8mW2iGhYWpd+/eys/P95xzu93Kz89X//79LczMvwzD0MSJE/XWW29p9erV6tixo9UpBcTgwYP1ySefqLi42HP06dNHo0aNUnFxsUJCQqxO0W+uvPLKOltUff755+rQoYNFGQXG0aNH5XR6fwSFhIR4ek2aqo4dOyouLs7rs6myslIbN25sUp9N0o9F5hdffKFVq1YpJibG6pT87vbbb9fWrVu9PpsSEhI0bdo0vf/++1an51dhYWHq27evLT6f4D9Neug8KytLGRkZ6tOnj/r166ecnBxVVVVp7NixVqfmN5mZmVqyZInefvttRUZGeuZ4RUdHq3nz5hZn5z+RkZF15p22aNFCMTExTW4+6tSpU3XFFVfo4Ycf1i233KJNmzZp8eLFWrx4sdWp+dUNN9ygP//5z2rfvr26deumLVu2aOHChRo3bpzVqZ2zI0eOaNeuXZ7He/bsUXFxsVq3bq327dtrypQp+tOf/qTOnTurY8eOmjVrlhISEjRs2DDrkj4Lp7vP+Ph4/eY3v1FRUZFWrFghl8vl+Xxq3bq1wsLCrErbZ2f6ff68gG7WrJni4uJ02WWXmZ3qOTvTvU6bNk0jRozQ1VdfrUGDBikvL0/vvvuu1qxZY13SCG5WL3sPtKefftpo3769ERYWZvTr18/YsGGD1Sn5laR6j5deesnq1AKuqW5vZBiG8e677xrdu3c3wsPDjS5duhiLFy+2OiW/q6ysNCZPnmy0b9/eiIiIMC6++GLjgQceMKqrq61O7Zx98MEH9f5/mZGRYRjGiS2OZs2aZcTGxhrh4eHG4MGDjZ07d1qb9Fk43X3u2bPnlJ9PH3zwgdWp++RMv8+fa8zbGzXkXl944QWjU6dORkREhJGcnGwsX77cuoQR9ByG0QS+hgMAAABBp8nO0QQAAIC1KDQBAAAQEBSaAAAACAgKTQAAAAQEhSYAAAACgkITAAAAAUGhCQAAgICg0AQAAEBAUGgC8Isvv/xSDodDxcXFVqdyVhp7/gAQjCg0AZzRmDFj5HA4PEdMTIyGDBmirVu3etq0a9dOBw4caHLfPQ8AOHsUmgAaZMiQITpw4IAOHDig/Px8hYaG6vrrr/c8HxISori4OIWGhlqYZfCpqamxOgUAsAyFJoAGCQ8PV1xcnOLi4pSSkqIZM2Zo3759OnjwoKS6Q89r1qyRw+FQfn6++vTpo/POO09XXHGFdu7cecoYJ6/x5ptvatCgQTrvvPOUnJysgoICT5u5c+cqJSXF63U5OTlKTEz0PB4zZoyGDRumhx9+WLGxsWrVqpX++Mc/6vjx45o2bZpat26tiy66SC+99FKdHHbs2KErrrhCERER6t69u9auXev1/LZt2zR06FC1bNlSsbGxuv3223Xo0CHP8wMHDtTEiRM1ZcoUtWnTRunp6Q19iwGgyaHQBOCzI0eO6K9//as6deqkmJiY07Z94IEH9Pjjj2vz5s0KDQ3VuHHjznj9Bx54QH/4wx9UXFysSy+9VCNHjtTx48d9ynH16tXav3+/1q1bp4ULF2rOnDm6/vrrdf7552vjxo266667dOedd+qrr77yet20adN07733asuWLerfv79uuOEGffPNN5Kkw4cP61e/+pV69uypzZs3Ky8vT2VlZbrlllu8rvHKK68oLCxMH330kRYtWuRT3gDQpBgAcAYZGRlGSEiI0aJFC6NFixaGJCM+Pt4oLCz0tNmzZ48hydiyZYthGIbxwQcfGJKMVatWedq89957hiTjhx9+qDfOyWv85S9/8Zz79NNPDUnG9u3bDcMwjDlz5hjJycler3viiSeMDh06eOXboUMHw+Vyec5ddtllxi9/+UvP4+PHjxstWrQwXnvtNa/Y8+fP97Spra01LrroImPBggWGYRjGQw89ZFxzzTVesfft22dIMnbu3GkYhmEMGDDA6NmzZ/1vJADYDD2aABpk0KBBKi4uVnFxsTZt2qT09HQNHTpUe/fuPe3revTo4fk5Pj5eklReXu731/xct27d5HT++BEXGxurpKQkz+OQkBDFxMTUuW7//v09P4eGhqpPnz7avn27JOnf//63PvjgA7Vs2dJzdOnSRZK0e/duz+t69+7tU64A0FQxax9Ag7Ro0UKdOnXyPP7LX/6i6OhoPf/88/rTn/50ytc1a9bM87PD4ZAkud3u08Y63WucTqcMw/BqX1tbe9prnLxOfefOlMtPHTlyRDfccIMWLFhQ57mTBbF04r0CADBHE8BZcjgccjqd+uGHH0yNe8EFF6i0tNSr2PTn3pcbNmzw/Hz8+HEVFhaqa9eukqRevXrp008/VWJiojp16uR1UFwCQF0UmgAapLq6WqWlpSotLdX27ds1adIkTw+fmQYOHKiDBw/qkUce0e7du5Wbm6t//OMffrt+bm6u3nrrLe3YsUOZmZn67rvvPAuYMjMz9e2332rkyJH6+OOPtXv3br3//vsaO3asXC6X33IAgKaCQhNAg+Tl5Sk+Pl7x8fFKTU3Vxx9/rDfeeEMDBw40NY+uXbvq2WefVW5urpKTk7Vp0yb94Q9/8Nv158+fr/nz5ys5OVkffvih3nnnHbVp00aSlJCQoI8++kgul0vXXHONkpKSNGXKFLVq1cprPigA4ASH8fPJTgAAAIAf8E9wAAAABASFJgAAAAKCQhMAAAABQaEJAACAgKDQBAAAQEBQaAIAACAgKDQBAAAQEBSaAAAACAgKTQAAAAQEhSYAAAACgkITAAAAAfH/ARUMa3cC9yXHAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 684.8x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# https://github.com/andrzejnovak/combine_postfits/tree/master/tests/fitDiags\n",
    "f = uproot.open('/home/anovak/cat/combine_postfits/tests/fitDiags/fit_diag_C.root')\n",
    "f.keys();\n",
    "\n",
    "tdir = f['shapes_fit_s/loosepasshadel2017']\n",
    "# tdir = f['shapes_fit_s/ptbin0pass2016']\n",
    "cov = tdir['total_covar'].to_hist()\n",
    "cov.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 583,
   "id": "6e1997c5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.00098322, 0.00125216, 0.02497574, 0.06100075, 0.05948158,\n",
       "       0.02749619, 0.07300204, 0.17876614, 0.06373522, 0.00813034,\n",
       "       0.00160794, 0.00142784, 0.00028312, 0.00003777, 0.00000102,\n",
       "       0.        , 0.        ])"
      ]
     },
     "execution_count": 583,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "tdir['total_background'].to_hist().variances()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 584,
   "id": "b6e30ee1",
   "metadata": {},
   "outputs": [],
   "source": [
    "hbkg = tdir['total'].to_hist()\n",
    "hbkg *= hbkg.axes[0].widths\n",
    "hdata = tgasym_to_hist(tdir['data'])\n",
    "diff = hdata.values() - hbkg.values()\n",
    "# hdata.view()[10] = (100, 100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 585,
   "id": "30082e84",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_28880/3268034680.py:7: RuntimeWarning: divide by zero encountered in divide\n",
      "  hep.histplot(diff/np.sqrt(hdata.variances()), yerr=1, ax=rax, histtype='errorbar')\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkEAAAGdCAYAAAAVEKdkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABWSElEQVR4nO3deVwTd/oH8M8kkHCjoEJQEBQqtVatZ1G39cBi67GubtvtuipV++uBVMXdqt2tR4/16OVq0V4orq212lZbta1r8Wi79YRlW9dKoQVB5agXp4aQmd8fSEyAICDJJJnP+/XKSzOZTJ6BSXjynWe+jyBJkgQiIiIihVHJHQARERGRHJgEERERkSIxCSIiIiJFYhJEREREisQkiIiIiBSJSRAREREpEpMgIiIiUiQmQURERKRIbnIH4MhEUcT58+fh6+sLQRDkDoeIiIiaQZIklJeXIyQkBCqV9fEeJkFNOH/+PEJDQ+UOg4iIiFqhoKAAXbp0sfo4k6Am+Pr6Aqj9Ifr5+ckcDRERETVHWVkZQkNDTX/HrWES1IS6U2B+fn5MgoiIiJzMzUpZWBhNREREisQkiIiIiBSJSRAREREpEmuCiIiIzEiShJqaGhiNRrlDISvUajXc3NxuefoaJkFERETXVVdXo7CwEFVVVXKHQjfh5eUFnU4HjUbT6m0wCSIiIkLtBLm5ublQq9UICQmBRqPhRLkOSJIkVFdX49dff0Vubi6ioqKanBCxKUyCiIiIUDsKJIoiQkND4eXlJXc41ARPT0+4u7vjzJkzqK6uhoeHR6u2w8JoIiIiM60dVSD7aovfE3/TRGYKCwuxdOlSFBYWyh0KERHZGJMgIjOFhYVYtmwZkyAicjl5eXkQBAGZmZlyh+IwWBNERER0M7t2AQBEScIlUYQBQKBKBY1Z4XS1JOGiKMIdQIBKBZXZY+WiiHJJgq8gwNfsNI7V7Y0f36Lw4uPjsWnTJtP9gIAADBw4EKtWrULv3r1bvLstMXz4cPTt2xerV6+2WJ6amoq5c+fiypUrNn39W8GRICIiomawWwLUSmPGjEFhYSEKCwuRlpYGNzc3jBs37pa26eqYBBEREd2EPROgaklqVYxarRbBwcEIDg5G3759sXDhQhQUFODXX39tdH2j0YgZM2YgOjoa+fn5AIDTp09j2LBh8PDwQM+ePfHVV19BEATs3LmzVTGZi4+Px8SJE/HKK69Ap9MhMDAQCQkJMBgMpnX0ej0WLFiA0NBQaLVaREZGIiUl5ZZf2xqeDiNluj60bf4BpZZUGLvzFwBA1f5DuJCff8tD1ETkGuyZAF0URehuMd6Kigq89957iIyMRGBgYIPH9Xo9HnnkEeTl5eGbb75Bx44dYTQaMXHiRISFheHo0aMoLy/H/PnzbzESSwcOHIBOp8OBAweQk5ODhx9+GH379sVjjz0GAJg2bRoOHz6MNWvWoE+fPsjNzcWFCxfaNAZzTIJIsep/QFWZzZBfJrXdEDUROT97JkDurYxx9+7d8PHxAQBUVlZCp9Nh9+7dDS4lr6iowNixY6HX63HgwAH4+/sDAPbt24eff/4ZBw8eRHBwMADgpZdewujRo1sZUUPt27fHG2+8AbVajejoaIwdOxZpaWl47LHH8NNPP2Hbtm3Yt28fYmNjAQDdunVrs9duDE+HkSJZ+4CqU4O2G6ImIudnzwQooJXz34wYMQKZmZnIzMzEsWPHEBcXh/vvvx9nzpyxWO+RRx5BZWUl/vWvf5kSIADIyspCaGioKQECgEGDBrUqFmvuuOMOqNVq032dToeSkhIAQGZmJtRqNe699942fc2mMAkiRbL2AVXHX2j8A4qIlMmeCZCqlaPP3t7eiIyMRGRkJAYOHIh3330XlZWVeOeddyzWe+CBB/D999/j8OHDrXqd+vz8/FBaWtpg+ZUrVyySLABwd7cc5xIEAeL1z1ZPT882iaclmASRIjX2AXXZLMlxa8MhaiJyHY6aADVGEASoVCpcvXrVYvmTTz6JFStWYMKECTh06JBpeY8ePVBQUIDi4mLTsuPHj9/0dXr06IGMjIwGyzMyMnDbbbc1O94777wToihaxGRrrAkiRbL2AVVfWwxRE5FrcPQESK/Xo6ioCABw+fJlvPHGG6ioqMD4Ri7oSExMhNFoxLhx4/DFF19g2LBhGD16NLp3747p06dj1apVKC8vx9/+9jcAaLKR7JNPPok33ngDTz/9NGbNmgWtVos9e/bggw8+wK7rF6E0R3h4OKZPn44ZM2aYCqPPnDmDkpISPPTQQy38aTQPP9VJ0cw/oOonObb8hkZEzsWeCVB5K0+9f/nll9DpdNDpdBg8eDCOHz+O7du3Y/jw4Y2uP3fuXCxbtgwPPPAAvvvuO6jVauzcuRMVFRUYOHAgZs2ahb/+9a8A0GSD0m7duuHrr7/G6dOnERsbi8GDB2Pbtm3Yvn07xowZ06J9WL9+PX7/+9/jqaeeQnR0NB577DFUVla2aBstIUgSqz2tKSsrg7+/P0pLS+Hn5yd3ONSWdu1q8AFVIwrovvMXFG2ai4Mvv4bbIiMafkDxEnkil3Xt2jXk5uYiIiKi1V3JXc2///1vDBs2DDk5Oejevbvc4Vho6vfV3L/fPB1GitTYN7Qas8dLJY4AEZHy7NixAz4+PoiKikJOTg7mzJmDoUOHOlwC1FaYBJEi3Wyqejc0PkTta78QiYjsrry8HAsWLEB+fj46dOiA2NhYvPrqq3KHZTNMgkiRbjZVvZ9g5Zy/PYMkIrKzadOmYdq0aXKHYTcsjCZFaiwBumRWjChYKXokIiLXwSSIFKm5U9XfbGZpIiJyXjwdRopmngDVT3KYABERuTZ+spNiNTVPR5XEBIiIyNXx050U6WYTIVbdpLcYERE5P6dNgpYuXQpBECxu0dHRpsevXbuGhIQEBAYGwsfHB5MnT7boh0LKdrOZoL2szPxKRNQcVdU1CF+4B+EL96CquubmTyBZOG0SBAB33HEHCgsLTbdvv/3W9Ni8efOwa9cubN++HYcOHcL58+cxadIkGaMlRyBJEqqqayAZRXgbJRhqROgNRugNRlwy+6ByN0qm5Vera1BkqEFVjVj7XI4IEZGTycvLgyAIyMzMlDsUh+LUhdFubm4IDg5usLy0tBQpKSnYsmULRo4cCQDYuHEjbr/9dhw5cgR33323vUMlB3HVYETPxXsR+3MmzHMZLx8NtL4ewPVrxF7d9xM6HS+FIAD+AV5wc1fhysUq7D2hxann4+Clceq3DhG5mPj4eGzatMl0PyAgAAMHDsSqVavQu3dvGSNzbE49EpSdnY2QkBB069YNU6ZMQX5+PgAgPT0dBoMBsbGxpnWjo6MRFhaGw4cPW92eXq9HWVmZxY1cU/0EyNtXi8qKaot16idANQaeDiMixzVmzBjTmZG0tDS4ublh3Lhxcofl0Jz26+zgwYORmpqKHj16oLCwEMuWLcNvfvMbnDx5EkVFRdBoNGjXrp3Fc4KCglBUVGR1m8uXL8eyZctsHDk5itcf6guDWkCFJMFHEKASBezccwYAMC82CtFR3VAjAQGCANEoYd62THkDJiJqglarNZ0dCQ4OxsKFC/Gb3/wGv/76a4N1jUYjHnvsMXz33Xf417/+hbCwMJw+fRqzZs3CiRMn0K1bN6xZswajR4/Gjh07MHHiRDvvjX04bRJ0//33m/7fu3dvDB48GF27dsW2bdvg6enZqm0uWrQISUlJpvtlZWUIDQ295VjJMRnUAvRqAQGCCr4qFaqMNx67phYgqFUIvj6ztN5gtL4hInJZkiTh6k3e/40VPl+tvvGcixV6VGkarnOz0+qe7mqL2etboqKiAu+99x4iIyMRGBiIyspK02N6vR6PPPII8vLy8M0336Bjx44wGo2YOHEiwsLCcPToUZSXl2P+/Pmtem1n4rRJUH3t2rXDbbfdhpycHIwePRrV1dW4cuWKxWhQcXFxozVEdbRaLbRarR2iJbl5+WhQIUmmBKi+GjRsreHm7tRnj4moFerqCG/Fb1YdbNXzWlp/uHv3bvj4+AAAKisrodPpsHv3bqjMPuMqKiowduxY6PV6HDhwAP7+/gCAffv24eeff8bBgwdNfydfeukljB49ulWxOwuX+VSvqKjAzz//DJ1Oh/79+8Pd3R1paWmmx7OyspCfn4+YmBgZoyRH4e2rhU8T8wD5C5YJkEGS0C7Qy64xEhG1xIgRI5CZmYnMzEwcO3YMcXFxuP/++3HmzBnTOo888ggqKyvxr3/9y5QAAbV/I0NDQy0GCgYNGmTX+OXgtCNBf/7znzF+/Hh07doV58+fx5IlS6BWq/HII4/A398fM2fORFJSEgICAuDn54fExETExMTwyjACAFSW6+FjluSIkoTLZvMAudVvripJLIwmUiBPdzVOPR/X5DrWTofVjQB988xweGrUDdZpzumwlvD29kZkZKTp/rvvvgt/f3+88847mDVrFgDggQcewHvvvYfDhw+brp5WMqdNgs6ePYtHHnkEFy9eRMeOHTFs2DAcOXIEHTt2BAC8/vrrUKlUmDx5MvR6PeLi4rBu3TqZoyZHUWV2JVjdRIg/FZ5HefouAMCbn+9GJ58H0VWnw0VRhJsAlF6qAtrLFTERyUEQhJsmK409bp4YBfpoZZlWQxAEqFQqXL161bTsySefRK9evTBhwgTs2bMH9957LwCgR48eKCgoQHFxMYKCggAAx48ft3vM9ua0SdDWrVubfNzDwwPJyclITk62U0TkjOoSoM1pafhzcjJE1I4AfXDwALYc2I9XZ8/Gn0aOhDcEcI5EInJker3edAX05cuX8cYbb6CiogLjx4+3WC8xMRFGoxHjxo3DF198gWHDhmH06NHo3r07pk+fjlWrVqG8vBx/+9vfAKDVxdnOwGVqgohaqi4Byjp/vjYBkiRAqj3lZRRFiJKE+W+8gUvFxY221iAiciRffvkldDoddDodBg8ejOPHj2P79u0YPnx4g3Xnzp2LZcuW4YEHHsB3330HtVqNnTt3oqKiAgMHDsSsWbPw17/+FUDtoIKrYhJEiiQIwGVIMAD4dP9+q990BEHAxn377BscEVELpaamQpIk062srAzHjh3D5MmTAQDh4eGQJAl9+/Y1PScpKQllZWUYMmQIgNpJhb/99lvo9Xr8+OOPaN++9vy/eZ2Rq3Ha02FEt8I/wAs1EhCsUqGgpMRqPzAJQF5JCYDay+qJiFzVjh074OPjg6ioKOTk5GDOnDkYOnQounfvLndoNsMkiBTJzV2FAEGARhAQ3qlT7UhQI4mQACC8UydUSBK8fTmHFBE1j5fGDXkrxsodRouUl5djwYIFyM/PR4cOHRAbG4tXX31V7rBsiqfDSJGuXKyC+/VTYDNGj7Y+EiRJeDg2FhWShMpyvT1DJCKyq2nTpuGnn37CtWvXcPbsWaSmpiIwMFDusGyKSRApkvmcP111Orw6e3Zt8bNQ+5ZQqVRQCQKSZ89Gp+Bg+AiCxWX1RETk/JgEkaJVSxIuiiL+NHIkMt5YD+87RgAA/jh8BNLXrcOEkSPhKwgWEysSEZFrYBJEimW4ngC5AwhQqRAVEgLf/rXzaUy7/wF0Cg6Gb73WGkRE5Dr46U6K5OauwiVJMiVA9ecBqpKkBgkQB4OIiFwLkyBSpHaBXnATGk+AAMCrkeaq/gFsoEpE5Ep4iTwp0oc+d2DJ+DiozPr5eAHYE5aB/puAwJHDgX79TI9J1TX4+DjnCSKi5issLMRbb72Fxx9/HDqdTu5wqBEcCSIiIrKBwsJCLFu2DIWFhXKHQlYwCSIiInJy8fHxEAQBgiDA3d0dQUFBGD16NDZs2ABRFG++getSU1PRrl072wXqYJgEERERtbHs7GysXbsWALB27VpkZ2fb/DXHjBmDwsJC5OXl4YsvvsCIESMwZ84cjBs3DjU1NTZ/fWfEJIiIiKgNbdy4EdHR0di8eTMAYPPmzYiOjkZqaqpNX1er1SI4OBidO3dGv3798Oyzz+LTTz/FF198YXrt1157DXfeeSe8vb0RGhqKp556ChUVFQCAgwcP4tFHH0VpaalpVGnp0qWmfRgwYAB8fX0RHByMP/7xjyi53lfRmTEJIiIiaiPZ2dmYNWsWRFGE0WgEABiNRoiiiJkzZyInJ8eu8YwcORJ9+vTBJ598AqB2Nvw1a9bgf//7HzZt2oT9+/fjmWeeAQAMGTIEq1evhp+fHwoLC1FYWIg///nPAACDwYAXXngB//3vf7Fz507k5eUhPj7ervtiC7w6jIiIqI1s2LChtiFzIwRBQEpKCpYvX27XmKKjo/H9998DAObOnWtaHh4ejhdffBFPPPEE1q1bB41GA39/fwiCgODgYIttzJgxw/T/bt26Yc2aNRg4cCAqKirg4+Njl/2wBY4EERERtZG8vLwmGzLn5eXZN6Drr1uXmH311VcYNWoUOnfuDF9fX0ydOhUXL15EVVVVk9tIT0/H+PHjERYWBl9fX9x7770AgPz8fJvHb0tMgoiIiNpIeHh4kyNB4eHh9g0IwI8//oiIiAjk5eVh3Lhx6N27Nz7++GOkp6cjOTkZAFBdbb1BdGVlJeLi4uDn54f3338fx48fx44dO276PGfAJIiIiKiNzJgxo8mRoJkzZ9o1nv379+OHH37A5MmTkZ6eDlEU8eqrr+Luu+/GbbfdhvPnz1usr9FoTLVMdU6fPo2LFy9ixYoV+M1vfoPo6GiXKIoGmAQRERG1maioKKSkpEClUkGtVgMA1Go1VCoVUlJSEBkZabPX1uv1KCoqwrlz55CRkYG///3v+O1vf4tx48Zh2rRpiIyMhMFgwNq1a/HLL79g8+bNePPNNy22ER4ejoqKCqSlpeHChQuoqqpCWFgYNBqN6XmfffYZXnjhBZvthz0xCSIiImpD8fHxyMrKwtSpUwEAU6dORVZWls2vpvryyy+h0+kQHh6OMWPG4MCBA1izZg0+/fRTqNVq9OnTB6+99hpWrlyJXr164f33329QpD1kyBA88cQTePjhh9GxY0esWrUKHTt2RGpqKrZv346ePXtixYoVeOWVV2y6L/YiSNbG7QhlZWXw9/dHaWkp/Pz85A6H2kBVdQ16Lt4LADj1fBy8NJYXSGZkZKB///5IT09HP7PeYTd7HhE5v2vXriE3NxcRERHw8PC45e1Z+zyhttHU76u5f785EkRERESKxCSIiIjIBnQ6HZYsWcIO8g6MY/pEREQ2oNPpTG0nyDFxJIiIiIgUiUkQERERKRKTICIiIjO8aNo5tMXviUkQERERAHd3dwC4aR8tcgx1v6e631trsDCaiIgItTM7t2vXztQSwsvLy2ofMJKPJEmoqqpCSUkJ2rVrZ5qZuzVcIglasWIFFi1ahDlz5mD16tUAaidRmj9/PrZu3Qq9Xo+4uDisW7cOQUFB8gZLDo2XtBIpW3BwMAC4TG8sV9auXTvT76u1nD4JOn78ON566y307t3bYvm8efOwZ88ebN++Hf7+/pg9ezYmTZqEf//73zJFSs6Al7QSKZsgCNDpdOjUqRMMBoPc4ZAV7u7utzQCVMepk6CKigpMmTIF77zzDl588UXT8tLSUqSkpGDLli0YOXIkAGDjxo24/fbbceTIEdx9991yhUxERE5ArVa3yR9ZcmxOXRidkJCAsWPHIjY21mJ5eno6DAaDxfLo6GiEhYXh8OHDVren1+tRVlZmcSMiIiLX5LQjQVu3bkVGRgaOHz/e4LGioiJoNBq0a9fOYnlQUBCKioqsbnP58uVYtmxZW4dKREREDsgpR4IKCgowZ84cvP/++23S6bfOokWLUFpaaroVFBS02baJiIjIsThlEpSeno6SkhL069cPbm5ucHNzw6FDh7BmzRq4ubkhKCgI1dXVuHLlisXziouLm6wk12q18PPzs7gRERGRa3LK02GjRo3CDz/8YLHs0UcfRXR0NBYsWIDQ0FC4u7sjLS0NkydPBgBkZWUhPz8fMTExcoRMNlL92WfQmM3jUS1JOGsQMfJ47eRZp+4W4XW9trFcFHHJYMT4olPYFXyXHOESEZEDccokyNfXF7169bJY5u3tjcDAQNPymTNnIikpCQEBAfDz80NiYiJiYmJ4ZZiLqZ8AXRRFNDZ3aLkoolyS4CMIqKqotl+ARETksJwyCWqO119/HSqVCpMnT7aYLJFck3kC5KuyPMtblwD5CgI0nPyViIiuc5kk6ODBgxb3PTw8kJycjOTkZHkCIrsxT4ACVCpcE29kOuWiCKOqNgHyVamgNxrlC5SIiByKUxZGE9WpnwCp6vX5qZBuJEB12AqIiIgAJkHk5JpKgADAp14CJEoS/AO87BghERE5KiZB5NQaS4DKRdH0//oJ0GVIcHPnYU9EREyCyMk1lgBVSFKD9URJwiVRRI0EXLlYZc8QiYjIQTEJIqdWPwGquwzeXF0CZAAQIAioMYggIiJymavDSNnML4NX1zsFVpcABapUkEQmQEREVItJEDklSZJw1WCE2mBEhSSh4voIkEYAqoxG1A1yFhtq4K6qHQGSRBH6GiZBRERUi0kQOaWrBiN6Lt6L8UUn4e2rRWW53jQTtFFQAd0HAgDWff0zKi5WWpwCY2E0EREBTILIydVPgIDaeYBqKi6hIvMLnO/YGV5qT9Njbu4qtAvkJfJERMQkiJxcZbkeLz3QE1q32tEdUZJw3ijhy88LUPrvD/D4ilcQE9UdAGCQJFySJLgJwJLxcfB0V8sZOhERyYxJEDm1qopqaN1U0LqrTUXQAm5cHebupobWXY1qSUK5KMILQu1l9Roe+kRESsfiCHIJFpfB12ugerPWGkREpExMgsjp1b8MXmOW5NTYIQEqLCzE0qVLUVhY2ObbJiIi22ESRE5NEIDLkBpNgACgVLL9CFBhYSGWLVvGJIiIyMkwCSKn5h/ghRqp8QQIqC16a6q3GBERKReTIHJqbu4qBAiCRQJUbdY7zE9omACVN9JbjIiIlIdJEDm1Kxer4F4vAbpkNtIjNJIA+bIwmoiIwCSInJz5TNDmV4HVZ54A+ap42BMREecJIhdhngDVT3KYABERUWOYBJFTS4scDOO4OHhp3KAzW76nawb6bwI8R9wD33794CtbhERE5Kj4tZiIiIgUiUkQERERKRKTICIiIlIk1gQRtcConKNQ7zZAdFPhkiiizAiM3XUGAHD1wNfAuXMNe5VNmCBv0ERE1CiOBBG1EJu1EhG5Bo4EEbVAXa8yAQICVSrUiPZt1kpERG2HI0FELdBUrzJ7NGslIqK2wySIqAUa61Vmegxs1kpE5EyYBBG1QGO9yuqwWSsRkXNhEkTUAvV7lbFZKxGR82ISRNQKbNZKROT8nPYTev369ejduzf8/Pzg5+eHmJgYfPHFF6bHr127hoSEBAQGBsLHxweTJ09GcXGxjBFTaxQWFmLp0qUoLCxs9nOys7Oxdu1aAMDatWuRnZ3dpjEZzBKg9mzWSkTktJz2U7pLly5YsWIF0tPTceLECYwcORK//e1v8b///Q8AMG/ePOzatQvbt2/HoUOHcP78eUyaNEnmqKmlCgsLsWzZsmYnQRs3bkR0dDQ2b94MANi8eTOio6ORmpraJvG4uatQbDRCMorwNkow1Nw4HVZmqMElgxFaowSNUYLeYITeYERVdQ0k1gYRETkcp50naPz48Rb3X3rpJaxfvx5HjhxBly5dkJKSgi1btmDkyJEAav843n777Thy5AjuvvtuOUImG8vOzsasWbMgmtXpGI1GAMDMmTMxbNgwREZG3tJrtAv0wrq0HJReqoIkAUZBBXi3BwBs+c9ZtDt9DVUV1ab1BQH46pg7Tj1f2+meiIgch9OOBJkzGo3YunUrKisrERMTg/T0dBgMBsTGxprWiY6ORlhYGA4fPmx1O3q9HmVlZRY3ch4bNmywKE42JwgCUlJSbvk1agyiKQGq7+pVQ4MEyD/A65Zfk4iIbMOpv5r+8MMPiImJwbVr1+Dj44MdO3agZ8+eyMzMhEajQbt27SzWDwoKQlFRkdXtLV++HMuWLbNx1GQreXl5Vk87SZKEvLy8Vm/b012NU8/HQdhVbXEZ/K81EvZ+ng8AmD20G2KiugOoba1RYhSx7uufW/2aRERkW049EtSjRw9kZmbi6NGjePLJJzF9+nScOnWq1dtbtGgRSktLTbeCgoI2jJZsLTw8vMmRoPDw8FZvWxAEeGnc4Klxg9ZdDa27GtVqAVklhShP3wUASPnXF8j/tRjubipUqgUIKgFXLla1+jWJiMi2nDoJ0mg0iIyMRP/+/bF8+XL06dMH//jHPxAcHIzq6mpcuXLFYv3i4mIEBwdb3Z5WqzVdbVZ3I+cxY8aMJkeCZs6c2WavVS6KePerr3Df07NR+b8DAIAPDh5A9JNPIvmrr2qbqwqCxbxCRETkWJw6CapPFEXo9Xr0798f7u7uSEtLMz2WlZWF/Px8xMTEyBgh2VJUVBRSUlKgUqmgVqsBAGq1GiqVCikpKbdcFF2nXBTx33Pn8OfkZIiSBEi1iY5RFCFKEua+8QYuFxdbzCxNRESOx2mToEWLFuHrr79GXl4efvjhByxatAgHDx7ElClT4O/vj5kzZyIpKQkHDhxAeno6Hn30UcTExPDKMBcXHx+PrKwsTJ06FQAwdepUZGVlIT4+vk22XzcP0I79+5s89bZ53742eT0iIrIdpy2MLikpwbRp01BYWAh/f3/07t0be/fuxejRowEAr7/+OlQqFSZPngy9Xo+4uDisW7dO5qjJHiIjI5GYmIjU1FQkJia22QgQANNEiOdLSqyfegOQV1ICoHZeISIickxOmwTd7HJnDw8PJCcnIzk52U4RkRLUzQQd3qlT7UhQI4mQACC8UycYJAntAnmJPBGRo+LXVKIWqGuFMWP06CaLsKeOHo1LksTCaCIiB8YkiKgVuut0WD17du2cQULt20ilUkElCHgrMRHtg4LgJgCll3iJPBGRo2ISRNRCoiThkiji9yNH4vgb6+F9xwgAwB+Hj8AP69dj7IgRtc1VITQ6szQRETkGp60JIpLF+PFQAehw/a4OwNd3D0T//v0x7+VV6Nmvn2nVquoapB3XyBElERE1A0eCiIiISJGYBBEREZEi8XQYOYzCnTsRqFJBc30SwiojMHbnLwCAorQDEM+eNTUvVRuMGF90EruC75ItXiIicm4cCSKHYZ4AAUC1WVWxn6Cy6N5eIUnw9tXaNT4iInItTILIYdRPgC6JN+bYMW9RUS6KqJAkVJbr7RofERG5FiZB5HCqJQkXRRHujTxW17vLRxBQVVFt99iIiMh1MAkih2KeAF0sKkJ5+i4AwJuf78Z/zp419e7yYYd2IiK6RUyCSHaSJKGqugYV1TUoqq6BZBTx8b6v0D/xKVT+7wAAYMvBAxiQkIBPvkqDxihBX8N2FEREdGt4dRjJ7qrBiJ6L9+LhimzUGETknynAx3vXQ4KE2p7sgHi9Pujp5GQcypHg5xMADgYREdGt4EgQOYwag4jSS1X4KS+zthV7YwTgp7z/QBAA/wB2aCciotZjEkQOo/RSFV57sC96dFBZz4EEAT06qPD8H/pizujbcOr5OHi6q+0aJxERuQaeDiOHIUmA1k2F7sFBtZfEN9J9VACgC+oEQa1CsEoFjYaHMBERtQ5HgsjhPDRqFCQr7dclScLDsbENJlasT6fTYcmSJdDpdLYKk4iInByTIHIoFZKEIJ0OybNn184QLdQeoipV7YzRr86ejQGdOzeZAAG1SdDSpUsdOgkqLCzE0qVLUVhYKHcoRESKxCSIHIaXjwYV1+cBemL0aGQmr4f3HSMAAJPvHY5vkpPxVGzsTRMgZ1FYWIhly5YxCSIikgmTIHIY3r5a+AgCfFW1h2V3XQh8+48HAPzp/gcajABVWzllRkRE1BysKiWHsdW7JxaPjwOuFzt7AdgTloH+m4Auo0ZA06+fxfoaGWJsDOuPiIicE5MgoltUV39ERETOhafDiIiISJGYBBEREZEiMQkiIiIiRWISRERERIrEJIiIiIgUiUkQERERKRKTICIiIlIkJkFERESkSEyCiIiISJGcMglavnw5Bg4cCF9fX3Tq1AkTJ05EVlaWxTrXrl1DQkICAgMD4ePjg8mTJ6O4uFimiImIiMjROGUSdOjQISQkJODIkSPYt28fDAYD7rvvPlRWVprWmTdvHnbt2oXt27fj0KFDOH/+PCZNmiRj1ERERORInLJ32JdffmlxPzU1FZ06dUJ6ejruuecelJaWIiUlBVu2bMHIkSMBABs3bsTtt9+OI0eO4O6775YjbCIiInIgTjkSVF9paSkAICAgAACQnp4Og8GA2NhY0zrR0dEICwvD4cOHrW5Hr9ejrKzM4kZERESuySlHgsyJooi5c+di6NCh6NWrFwCgqKgIGo0G7dq1s1g3KCgIRUVFVre1fPlyLFu2zJbhKsOuXab/ipKES6KIMiMw8rg7AODU3SK81EC1JOGiKEIyioj9ORNfdR/cYFM6nQ5LliyBTqezW/j2kJ2djbVr1wIA1q5di2effRZRUVEyR0VEpCxOPxKUkJCAkydPYuvWrbe8rUWLFqG0tNR0KygoaIMIlasuATIACFBZHmp1CZA7gPYQIEmNb0On02Hp0qUulQRt3LgR0dHR2Lx5MwBg8+bNiI6ORmpqqryBEREpjFMnQbNnz8bu3btx4MABdOnSxbQ8ODgY1dXVuHLlisX6xcXFCA4Otro9rVYLPz8/ixu1jnkCFKhSQSMIpsfME6AAlQoqs8dcXXZ2NmbNmgVRFGE0GgEARqMRoihi5syZyMnJkTlCIiLlcMokSJIkzJ49Gzt27MD+/fsRERFh8Xj//v3h7u6OtLQ007KsrCzk5+cjJibG3uEqTlMJEABcUmgCBAAbNmyAYGWfBUFASkqKnSMiIlIup6wJSkhIwJYtW/Dpp5/C19fXVOfj7+8PT09P+Pv7Y+bMmUhKSkJAQAD8/PyQmJiImJgYXhlmB00lQAAaTYC8fDT2C1BGeXl5kKyc+5MkCXl5efYNiIhIwZwyCVq/fj0AYPjw4RbLN27ciPj4eADA66+/DpVKhcmTJ0Ov1yMuLg7r1q2zc6TKZO0UWJ329RKgCkmCt6/WniHKJjw8vMmRoPDwcPsGRESkYE6ZBFn7Jm3Ow8MDycnJSE5OtkNEZK6xBOiSKAJQA4BFAlQuiqiQJFSW6wFve0dqfzNmzMCqVasafUySJMycOdPOERERKZdT1gSRY7NWBF1fuSiiXJLgIwioqqi2X4AyioqKQkpKClQqFdTq2qRQrVZDpVIhJSUFkZGRMkdIRKQcTILIZiwug693iXxdAuQrCPBRWHF0fHw8srKyMHXqVADA1KlTkZWVZTqVS0RE9sEkiNqMJEmoqq6B3mBERXUNiqprIBlFeBslGGpE03qXqmtwyWCE1ihBY5SgN3tMKSIjI5GYmAgASExM5AgQEZEMnLImiBzTVYMRPRfvRdyZTLQL9EKNQUTppSpIEmAUVKgJ6o6KzC/wYkkk3EWNxSkwhQ0GERGRA+BIELW5+glQHWPFJZT++wNcuHCxQQLkH+AlQ6RERKRkHAmiNldjEPHcqCh4uqlNy36tkbD383wAwOyh3RAT1R1A7cSKlyGhRgJemBAHT3d1o9skIiJqa0yCqM2VXqqCp5sa2usJTbkoQm82IuR+/bG6maUFCAhWqaDR8HAkIiL74ekwanPmp8DML4M3d7PWGkRERLbGr94KVlhYiLfeeguPP/54i7q0N/d55pfBq80ukZcUlgDF/nwUV3bpUapWmfa3ygiM3fkLAOBs2gGE5ec37KU2fnyj22vt742IiCxxJEjBCgsLsWzZMhQWFrb58yokyZQA+dabI6hMUk4CBNQWfddI1vfXDQ17qZWL1qcNaO3vjYiILHEkiKwSP/uswYhNUyMYaoMRo3Iy4OWjQYUkIUBQNUiAAKAGjbfWcNUWqm7uKniLEqQaEfrryypqbpwz9Kg3j1KFJKFCkqCuroGnu9pqrzEiIro1TILIqpudsmpsBMPLRwNvXy186o0AiWaFQv5CwwTooijCVU/sXLlYhWe2/dd0381dBZ9Ab+B62vf6V9l4/0QZgBs/v8pyPXad0OLU83HwYsE4EZFN8HQYWXWzbvB+QuPd4CvL9RaF0KIk4bLZ6R23ZvYWcxU1BrN9d1eZ5lGqzzwBUkovNSIiOfErJlllvRt8LcFKN3jzP+DmV4HVZ54ABTRy2szZebqrcer5OKh31+69QZJwSZLgJgBao4Bt1+dNmj/6NtwZ2Q0V16+iczdKmLctU8bIiYiUgUkQWXWr3eDNE6D6SU79BEjlgnUvgiDUnspyV6NaklAuivCCgACVCtfUN/bXoBagVwumGiq9wShj1EREysEkiG7KPGGpX+hsfhm85vrf9bTIwTCOi4NW44YOZuvu6ZqB/psAzxH3QNOvn8vWANV3s4SvqpGr6FwwJyQicjhMgqhJ9f+AXxNv/HWukkSLy+D1Ro5gNEYzYUKDhM8LQMbgQrwVfhmD/vRH+JrN92OsrsFXx1y5SoqIyDEwCSKrOIJhWzqdDkuXLpU7DCIixXK9alRqMzer2fFq5DJ4doMnIiJnwSSIrGosATKfydhLsEyALkOCm3vjh1R2djbWrl0LAFi7di2ys7NtEzQREVEzMQmiBiRJQlV1Dbyvz2SsNxihNxhxsboGl82uXDLU1C6/Wl2DIkMNrhpFXLlY1WB7GzduRHR0NDZv3gwA2Lx5M6Kjo5GammqvXSIiImqANUHUwFWDET0X78WonP+YltVN5FdaYQDgDQB4dd9PCDpRCv8AL7i5q3DlYlWDSQCzs7Mxa9YsiGYjSMbrBdQzZ87EsGHDEBkZafudIiIiqocjQXRT5jMZXzWbB0gQ0GQCBAAbNmyw2vtKEASkpKTYLG4iIqKmcCRIoerX6Dz77LOIiopqsN6LD/VBtUqAj1B7++HseXz0zicAAI0hGxN63Yn+IZ3hbpboGMfFwdNdDQDIy8uDZNZqw5wkScjLy2vjPXNtzf29ERHRzXEkSIGaW6Pj5aNBtUpAgLsagRo3bDl0AIPnzkbl/w4AAHZ8cwgj5zyNj745BK27Glp3NQQ3Fbw0bqbRn/Dw8CZHgsLDw222n66GtVVERG2LSZDCmNfo1NXmGI1GiKKImTNnIicnx7SueTf47PPnMWvt2tpu8FLtaS9RFCFKEmauWYOc8+dN8wqZmzFjRpMjQTNnzrTRnrqWlvzeiIioeZgEKUxLanTMu8Fv2Levyee9vW9fo73FoqKikJKSApVKBbW69hSZWq2GSqVCSkoKi6KbibVVRERtj0mQwrSkRse8GWpuSYn15wH4qbjYajf4+Ph4ZGVlYerUqQCAqVOnIisrC/Hx8a3dDcVhbRURUdtjEqQwranRESUJQR07Wn8egK6dOjXZDT4yMhKJiYkAgMTERI4AtRBrq4iI2h6TIKXYtQvVn32GcaGhkMSGl7IDgCSKeLhrV6h378aonKMAahOgS6KIh0aNanIkIuG++6wmQHTrWFtFRNT2mAQpRF3Rco+QELyTmFibsFxve6G6PoKTPHs2+nbubHqOIACXIcEAYGCXLkh5+ulGn/dOYiJuCwmRY7cUg7VVRERtz2mToK+//hrjx49HSEgIBEHAzp07LR6XJAmLFy+GTqeDp6cnYmNjFd2vyrwZ6ozYWGQmr4f3HSMAAL+7dzjS163DE6NHm9avmwixRgICVSpoBAHxo0ZZPG/yvcPx4/r1mBEba3peuZVRJrp1rK0iImpbTjtZYmVlJfr06YMZM2Zg0qRJDR5ftWoV1qxZg02bNiEiIgLPPfcc4uLicOrUKXh4eMgQsbzqN0PtrguBb//xqDyZhscfGGsxAlTXDf7j9nfihQlx0GhuHCZ3Avh68ED0798fC195Gbf162fxOr722BkFq6utSk1NZW0VEdEtctok6P7778f999/f6GOSJGH16tX429/+ht/+9rcAgH/+858ICgrCzp078Yc//MGeoTqEtuoGT0RE5Cpc8i9dbm4uioqKEGt2msbf3x+DBw/G4cOHrT5Pr9ejrKzM4uYq6idAFY0U2dYVQddIaLQbPBERkStxySSoqKgIABAUFGSxPCgoyPRYY5YvXw5/f3/TLTQ01KZx2oMkSaiqroHeYITeYMTF6hpcMhihNd5Iggw1RlytrkGRoQZVNSK8RanRZqhERESuxGlPh9nCokWLkJSUZLpfVlbm9InQVYMRPRfvxaicDItu8OWVNYB3ewDAa1/9hKhfDE12gyciInI1LjkSFBwcDAAoLi62WF5cXGx6rDFarRZ+fn4WN1dhngCZzwQNAN5+WrskQDqdDkuWLIFOp7PZaxARETWXS44ERUREIDg4GGlpaejbty+A2lGdo0eP4sknn5Q3OJl4+2ox5zfdEOBWO8dMlRHouecMAODBgWG4L6o73M3qhgySBNX4OHi6q9ssBp1Oh6VLl7bZ9oiIiG6F0yZBFRUVFp2zc3NzkZmZiYCAAISFhWHu3Ll48cUXERUVZbpEPiQkBBMnTpQvaBlVlusR4KaG9npSYxBu1AQFuqnhY3YZfLUkoVwUodM47eFBRER0U077V+7EiRMYMWKE6X5dLc/06dORmpqKZ555BpWVlfi///s/XLlyBcOGDcOXX36pyDmCAMtmqKIk4bLZJfJuZiNAdTNL1+8GT0RE5GqcNgkaPny41V5KQG1Tyeeffx7PP/+8HaNyfHWXwRsAqH0C4D/0EQS3ry2QNk+AGusGX4e1PfLiz5+IqG04bRJELWeeAAWoVHDzCUC7YVOgCxAbJEBNNUNlbY+8+PMnImobTIIUoq4ZqgABgSoVakTLU2ClzUyAiIiIXAWTIIX4qvtgrBkfB6/rxc4aAHkTbzzeTo6giIiIZOSS8wQRERER3QxHgogc1Kico6jZqcdFtYAKSYKPIEAlCrgrvfbavfR+1dCrJdRIQIAgmOZ5Mo4bB093NQSe1iQiahKTIAdTWFiIt956C48//jiv/iH89fNTFjN9GwUV0H0gAOD5tJ/g4S40mOk77Zg7Tj1/49QnwOOKiKgxPB3mYAoLC7Fs2TIUFhbKHQrJrKlWJwAabXXi5l77lq6qNqKqusZ0+yW/AMuWLcMv+QUWy81vTU05QUTkijgSRORgPN3VOPV8HIC4m6xp+XhVtREDXvwKky//gIWPHkVdTmMUVNh7vVnu6y9tgr9mt+k5ggD4B3jh4/Z3Nhg9IiJydfzEI3IwgiDcUjJSeqkK5oM6nj4a4Pp9fVVN7aWBuJEA1Y0eEREpDZMgIhdRN4Kk3m0wLauQJPxaI2HnFwUAgPmjb0NMVPfa1imQcNUoYtWOk4CPXFETEcmHSZCLYOErmUaQJk00LdMCCASwJyoD/TcBfqNHQNuvHwDAE0BVdQ32Zni2+LUkScJVg7FN4r5VhYWF2PDuO5gx6zGHOPZ5ZR6R82AS5CLqCqonTJjgEH8IyLVdNRjRc/FexP581HRKrbEi7XaBXqgxiA1O0Vkr+jY/Rdec7dXVOxVtegH5P5WhnWeHW9rerca3t+tA1lYRORG+U4mo1eRMgOq7etUAjVG+BKhdoFeLf35EJC8mQUTUam7uKjxzT3f4uKlNywyShEuSBDcBaA/BohddhSSZJn70MVteV6NUf+LHm23v1xoJez/PBwDMHtoNMVHdb2l7rY2v2GjEurScW/pZEpH9MQlyINnZ2Vi7di0AYO3atXj22WcRFRUlc1Tk7Gx5XH3ocweWTIyD1uz0jxbW66zrapQaY60yqantdQGwJ7JhvVNrt9fa+FTVNfg4Y6+VNYjIUfHaWAexceNGREdHY/PmzQCAzZs3Izo6GqmpqfIGRk6NxxURkXUcCXIA2dnZmDVrFkTxRo2B0Vh75c3MmTMxbNgwdO/eHVcNRlzZ9VmDIfkqIzB2zxkAQMG/9uP23DOmIf6aGhHji05hV/Bddt4rkltzjquQsHCZoiMikh+TIAewYcMGq5fUCoKAlJQUPLfsBfRcvBcPV+Q0KMqExg3Q+AMAUtKy8dmJMtNDdUWepDzNPa6IiJSKp8McQF5entW+TZIkIS8vz3S/0atSAjzN1r/xXPOrXEh5WnJcEREpEUeCHEB4eHiT39jDw8NN92sMIl5/qC+0birTVS7Z587jrXd3AgD8jDlIGt4HwTqdxVUuK8fFwdNd3ehrkGtqyXHlrHgxARHdCo4EOYAZM2Y0+Y195syZpvtpkYPhNnECtJMmwmfy75BWXor7kubg2o8HAQAffnMIveckYld5Kbr+fhICJ/8O2kkT4aVx4yy2CtOS4wpo2Hm+6Zv8s0Wz6JuIbhVHghxAVFQUUlJSMHPmTAiCAKPRCLVaDUmSkJKSgsjISFRV1zR4XnMKXyMjI+22H+RYWnJcjS/6DxbEHzU9VxAA7wBvfOwegorML3BfOx18PbwB3JgYcLJBxMft75Rl33jsW8cWOkTNx5EgBxEfH4+srCxMnToVADB16lRkZWUhPj7e6nOaU/hKytbc48razMjGikso/fcHuHqtAkDDmZblwmPfuroWOoWFhXKHQuTwOBLkQCIjI5GYmIjU1FQkJibe9JssC1+pOZo6ruo6zwNxjT530n8yMGwTsPClmbjrrn6NriNHrRmPfSJqC0yCnJgSCl/Jtkyd563wdHcz/etITUF57BNRW+DpMCfW0sJXIlfBY5+I2gKTIAej0+mwZMmSZhU01hW+qlQqqNW1pyTUajVUKpWp8JUIaNlx5Qx47BNRW2AS5GB0Oh2WLl3a7D9WrSmoJuVp6XEFNJyDJzs721bhtQqP/YYc/Xcml8LCQixdupTF4tfx53GDIFkbUyaUlZXB398fpaWl8PPza5NtSpKEqwYj1Lt3Q5QkXIZk6gVmEAXcle4OADjWrxrXVBLcBKA9BBiMEuZty0Ra5GCcej6uQX1GRkYG+vfvj/T0dPTr13gBK1Fzbdy4EbNmzWr00npHSzIc4divqq5Bz8W1XeQbe3/agzP9zuzNEY4RR6KEn0dz/347TqWjQlw1GNFz8V7E/pxhugy5rhWGUVAB3QcCAFZ+nQMYjCi9VIXmpKmudrqD5ONsc/Dw2He+3xmRo+DpMJnUT4Dqq5uHpX4vMGtac7qDqDHONgcPj33n+50ROQqXHwlKTk7Gyy+/jKKiIvTp0wdr167FoEGD5A4LH7e/Eyf+FgsvjbU5Vhqft2Ul5JmXhZSDc/DcGjlaiuT8ktvk7yznl9xGZ523B093tawte9hfzhJ/HpZcOgn68MMPkZSUhDfffBODBw/G6tWrERcXh6ysLHTq1Enu8OClUTvU3CtEAOfguVUDXvzK7q95Od8I0cppc1EC/pVvNNUs2ZtcNVKAZZ0UUNtf7p///Kdi66T482jIpQujBw8ejIEDB+KNN94AAIiiiNDQUCQmJmLhwoU3fb4tCqMdoYCSqCnZ2dmIjo62qC+po1KpkJWVxfqSeszf13IwXDqH8+8+gUYLCAUBIY+9Bff2IfYPDLjJiLft5GRn467evawex5k//A/dFXQcO+rPw1Z/AxVfGF1dXY309HQsWrTItEylUiE2NhaHDx9u9Dl6vR56vd50//LlywCAs2fPtlkSdNVQg5qyCwCAc+fOmmbkJXIUnp6eePnllzF//nwIggBJkkz/vvzyy/Dw8MDZs2flDtOhSJKEfY/fIWMEd+DjyJVYtOCZBr+z5StWYvKDo+0azdVqEQ/84xsAQN9nttr1teuUHtnW6B98oPYL8ZA/zYf/3Q/ZOSr5OOrP48izo2yy3bKyMgA3LhCwSnJR586dkwBI3333ncXyv/zlL9KgQYMafc6SJUskALzxxhtvvPHGmwvcjh071mSuwGEIM4sWLUJSUpLp/uXLlxEeHo6CgoI2GwkiIiIi2yorK0NoaOhNT927bBLUoUMHqNVqFBcXWywvLi5GcHBwo8/RarXQarUNlvv5+TEJIiIicjJ1bXWscdl5gjQaDfr374+0tDTTMlEUkZaWhpiYGBkjIyIiIkfgsiNBAJCUlITp06djwIABGDRoEFavXo3Kyko8+uijcodGREREMnPpJOjhhx/Gr7/+isWLF6OoqAh9+/bFl19+iaCgILlDIyIiIpm59DxBt8oW8wQRERFQUnYN7x/Nx5TBYejk5yF3OORimvv322VrghxZSdk1vL7vJ5SUXZM7FCIiWZSU6/GPtGyUlOtvvjKRjTAJkgHf/ERERPJjEkRERESK5LJJ0PLlyzFw4ED4+vqiU6dOmDhxIrKysuQOi4iIiByEyyZBhw4dQkJCAo4cOYJ9+/bBYDDgvvvuQ2VlpdyhERERkQNw2Uvkv/zyS4v7qamp6NSpE9LT03HPPffIFBURERE5CpdNguorLS0FAAQEBFhdp34X+boutEREROR6XPZ0mDlRFDF37lwMHToUvXr1srre8uXL4e/vb7qFhobaMUoiIiKyJ0UkQQkJCTh58iS2bt3a5HqLFi1CaWmp6VZQUGCnCImIiMjeXP502OzZs7F79258/fXX6NKlS5PrWusiT0RERK7HZZMgSZKQmJiIHTt24ODBg4iIiJA7JCIiInIgLns6LCEhAe+99x62bNkCX19fFBUVoaioCFevXpU7NCIi2bF9D5ELJ0Hr169HaWkphg8fDp1OZ7p9+OGHcodGRCQ7tu8hcvHTYeS42EGaiIjk5rIjQeTY+C2UiIjkxiSIiIjsKvdCJVK/ywMApH6Xh9wLbGdE8mASZGd88xORkm07UYBRrx7EJxlnAQCfZJzFqFcPYvsJzsumNI5QnM8kyI745iciJcu9UImFH38PUQLE62Wbdf9f8PH3yOOXQkVxhLIIJkF2wjc/ESndthMFEASh0ccEQcCH/EJIdsYkyE745icipTt7+arVK3clScLZy5zHjeyLSZCd8M1PRErXpb1nk18Gu7T3tHNEpHRMguyEb34iUrqHBoQ2+WXw4QGhdo6IlI5JkJ3wzU9EShfRwRsrJ/eGSgBU178T1v1/5eTeCO/gLW+ApDhMguyEb34iIuDBAaHYP384JvXrAgCY1K8L9s8fjgf5RZBkwCTIjvjmJyICwjt4I35IOAAgfkg4vwSSbJgE2Rnf/ERERI6BSRAREREpEpMgIiKFYfseolo2S4LUarWtNk1Et8gRevaQPNi+h+gGmyVB1i4HJ+K3UPk5Qs8esj+27yGy1OIkaMyYMViwYAHef/99/PDDD6ipqWl0PWsTA5Ky8VsokXzYvofIUouToJEjR+L06dP44Ycf8Oyzz8Lf3x/Dhw9HSkqKLeIjF8JvoUTyYvseIkstToI++OADfPrpp1ixYgV27dqFPXv2oH///sjKysJzzz1nixjJRfBbKNENctRlsX0PkaUWJ0Genp7Izs423R8+fDj279+PlStXYs+ePW0aHLkWfgslukGOuiy27yFH4Si1oW4tfcK6devw0EMPYciQIejbty+ys7Ph7u4OQRBgMBhsESO5CNO30EY+hO39LbSk7BreP5qPKYPD0MnPw26vS45DicdAXfueBR9/D6D2VHRdGx+27yF72XaiAAuvH4NAbW3oJxlnsXJyb7t3UGjxSFDfvn1x/PhxDB8+HPn5+dDpdNizZw+qqqrw+9//3hYxkotwpG+hvDqKlHoMsH0PycnRakNbnAS9++67cHNzw4MPPogXXngB8+bNQ8eOHeHl5YUlS5bYIkZyEWwiS+QY2L6H5OJotaEtToJ2796N/fv3m+5XVVXhD3/4Q5sGRa6L30KJiJTL0WpDW5wE/fOf/8TixYtx+vRp/PTTT7jnnnsQFxdni9jIRfFbKBGRMjnaFYrNLoyeN28e+vbtiz59+uDdd9/FlClTIIoiNm7ciL59+9owRGprSiwIJSLH0slXizmjotDJVyt3KGRHDw0IxVuHfm70MTmuUGz2SNCIESNQUFCAv//973jwwQeRm5uLkJAQ7N27l5fGt5Dcb36lFoQSkePo5OeBeaNv4xcxhXG02tBmjwRNmDABsbGx8PLyAgBcu3YNJ0+exPfff499+/Zh7NixNgvS1dS9+YmIiJTmwQGhGBgegDcO5OCj9LOY1K8LZo+IlKU0otkjQYmJiejUqRP69u2LnJwczJ8/HytWrIBer8fq1attGCIRERG5EkepDW12EvT555/jwoULWL9+PYYNG4bOnTtj2rRp+Prrr7F48WJbxnhLkpOTER4eDg8PDwwePBjHjh2TOyQiIiJyAM1Ogvz9/eHh4YGYmBj4+/vj2WefxYQJE/Dee+/h888/t2WMrfbhhx8iKSkJS5YsQUZGBvr06YO4uDiUlJTIHRqRbBxlunoiIrk1Own69ddfsXPnTuTm5sLb+8awlVqttnrNv9xee+01PPbYY3j00UfRs2dPvPnmm/Dy8sKGDRvkDo1IFttOFGDUqwfxScZZALXT1Y969SC2s3ktESlQswujk5KSsGvXLixfvhy//PILhgwZgh49eqBHjx64ePFis7bxzjvv4PDhw+jevTt69+6N8ePHtzrwm6murkZ6ejoWLVpkWqZSqRAbG4vDhw+3aFtV1TVwq65p6xBlc81gNP1bJdN+yR2D3K8vh7yLN6arr1P3/wUff49enf3QNVA5czbJfQwo/fWJbHkMNnd7gtTKYZzc3FycPHnSdHv//fctHler1TAajRbLEhMT0b59e9x333145513sGnTpta8dLOcP38enTt3xnfffYeYmBjT8meeeQaHDh3C0aNHGzxHr9dDr79x2XhZWRlCQ0MROncbVFovm8VKREREbUfUV6Fg9UMoLS2Fn5+f1fVa3EW+TkREBCIiIlo0mtOuXTsAwIABA/DNN9+09qVtZvny5Vi2bJncYRAREZEdtDoJao0XXngB586dw9y5cxEVFWXT1+rQoQPUajWKi4stlhcXFyM4OLjR5yxatAhJSUmm+3UjQcf+OqrJTNLZnDpfht+/eRgfPRGDniHy7JfcMcj9+nJ4bd9P2PhtHoyNDP6qBQGPDgtHkoLmr5LrGPgk4xwWf3oSQO3pyLoJ416Y2Au/u6uz3eJQ4nuAHIstj8GysjLoVt98vVYnQUVFRVaTCWvefvtt/N///R/efPNNnDt3rrUv3SwajQb9+/dHWloaJk6cCAAQRRFpaWmYPXt2o8/RarXQahvO4uylcYOXxq75ok15uKtN/8q1X3LGkHuhEluP1xYCbz1egIQRkYhQQP+yKYO7YsO3uY0+JkHCnwZ3danjvClyHQO5Fyqx+NOTjdZlPbfzJIZ272C3+VIc4XOAlM2Wx2BNM7fX4gaqde67775mr3vlyhU89thj2Lp1K9577z3897//xZIlS1r70s2WlJRkqj368ccf8eSTT6KyshKPPvqozV+bHJOSr45ypOnqS8qu4fV9P6Gk7JrdXrOOnMfAthMFTTaP/NCOx6Hc7XscgZzHITmGVqdeLamnbteuHd555x3s3r0bHTp0wLFjx+zSZuPhhx/Gr7/+isWLF6OoqAh9+/bFl19+iaCgIJu/Njme3AtNXx01MDzAbomAXE1sHWW6+rr+daN7Btl1/+U+Bs5evmr1s1OSJJy9fNVmr10f2/fIdxyS42j1SJC1bzNNKSwsxNtvv42SkhK4udln+HX27Nk4c+YM9Ho9jh49isGDB9vldcnxONK3cDmb2DrKdPVykPsY6NLes8nX79Le06avT0SWWp0Etcb333+PLl264N5778VHH31kz5cmByPHULwjfQsnech9DDw0ILTJ1394QKhNX5+ILNk1CTK/RD46OtqeL03XOUrLhLqheHsOQfNbOMl9DDhSXRYR3UISpFarW/ycF154AY8//jjmzZsHjUbT2pemVlJyUTDAb+HkGMfAgwNCsX/+cEzq1wUAMKlfF+yfPxwP8vgjO2Nh+C0kQf/5z39atH5aWhrGjRuHOXPmoEePHpg3b15rX5pawbwgtK4QtO7/Cz7+HnkKaKLJb+HkKMeAkuuyyHHIWZsIOMYVinY7HbZhwwbs3LkTH330EaKiovDSSy/Z66UJ8heEOgp+CyceA0SOQY6yiPrslgQFBgaarggbO3YsLly4YK+XJshfEOpI+C2ceAwQEdCKeYIMBgMOHjwIDw8P9OzZE4GBgc163r///W8899xz6NevH+66665WXWJPrWcqCG0kEWJRMBERKVGLR4ImTZqE7du343e/+x3uvvtudO7cGWPGjLG6ft2l8Onp6Zg+fTpqamrw9ttvIycnp/VRU4s5QkEoERGRI2nxSFB+fj527dqFY8eOITMzE8nJyThz5ozV9R9++GGsW7cOV65cQb9+/TBp0iQ8+OCDtxQ0tVxdQeiCj78HYNm4kUXBRESkRC0eCfLwqC1g0mg0qK6uRkJCAr799lur60uShKNHj0KtVmPTpk0YPHgw8vLyWh0wtR4LQomIiG5o8UjQ008/jUuXLmHy5Ml44oknMHTo0CaLnG+77TZs2LDBdP/48eN46qmn8Pnnn7cuYroldQWhH6WfZUEoEREpWqtqggICArBgwQLcc889OH36dJMtMAICAnDq1CnT/YEDB+LcuXOti5aIiIiojbR4JGjYsGFIT08HAMTHxwMATp8+bXX91atX43e/+x3uuece9OrVC99//z1CQkJaFy0RERFRG2n2SNCuXbuwcuVKVFRUoKDAcmK9hx9+2OrzBg0ahMzMTIwdOxaVlZUYMGAAPvzww9ZHTERERNQGmj0S1KtXLxQUFODChQuYNm0a8vPz0blzZwQHB8Pd3b3J53p6emLixImYOHHircZL5PTqN7FNGBGJCAXVZil9/8kx8DgkoAVJUEREBJ566in06tUL99xzDwDg3LlzOHPmDHr16mWzAIlcybYTBVh4fZoCoLaJ7ScZZ7Fycm+7XqUnV88eR9l/UjYeh1SnxYXRMTEx2LdvH7755ht4eHhgyJAh8PPzs0VsRC7FkZrYytGzx5H2n5TLkY5DdnGXn81njCaiWkpvYqv0/SfH4EjHodxd3KkVSVB+fj7efvttdOnSBdnZ2Xj22WfRu3dvW8RG5FKU3sRW6ftfn1ynJJWOx2Gt+jVRuQodiW3xJfKNzRg9ZMiQNg+MyNUovYmt0ve/vrpTkmRfPA5ZE2WuxSNBdTNGT5o0CU888QRSUlKanDGayBHJ8S1c6U1slb7/5BiUfhw6Uk2UI2hxEjRlyhS0a9cOCxcubNaM0USOSI7C4LomtirhRvPauv8roYmto+0/T0cpk6Mdh/bmSDVRjqDZSdDly5fxxz/+Ef7+/tBqtYiMjMTp06fx17/+lTVBRM2k9Ca2jrT/ciTC5Bgc6Ti0N9ZEWWp2ErRw4UKEhYXh7NmzqKysNDVAjYmJQVFRkc0CJHI1dU1sASiyia3S958cg1KPQ1NNVCOUUhNlrtlJ0NGjR7FixQr4+vpCo9Hgtttuw4oVK7B48WIsXbrUhiESERFRW1B6TVR9zU6CVKrGV33kkUdw4sSJNguIbI+1EEREyqT0mqj6mp0ElZSU4KOPPsKPP/4Io9Fo8Zi1oTVyTKyFICJSLiXXRNXX7HmC5s+fjy+++AIvv/wysrOzERISgjvuuAM9e/ZESUmJLWMkIiKiNlRXE/VR+llF1UTV1+wkaN68eRb3c3NzcfLkSZw8eRLDhg1rsL61c45EREREjqDFM0bXiYiIQEREBMaPH9/o46IotjooIiIiIltr8WSJziAvLw8zZ85EREQEPD090b17dyxZsgTV1dVyh0ZEREQOotUjQY7s9OnTEEURb731FiIjI3Hy5Ek89thjqKysxCuvvCJ3eEREROQAXDIJGjNmDMaMGWO6361bN2RlZWH9+vVMgoiIiAiAi54Oa0xpaSkCAgKaXEev16OsrMziRkRE1NZyL1Qi9bs8AEDqd3nIVVjjUkehiCQoJycHa9euxeOPP97kesuXL4e/v7/pFhqqvDkTiIjItradKMCoVw/ik4yzAIBPMs5i1KsHsV1hzUsdgVMlQQsXLoQgCE3eTp8+bfGcc+fOYcyYMXjwwQfx2GOPNbn9RYsWobS01HQrKOABSUREbSf3QiUWfvw9RAkQr88kU/f/BR9/jzyOCNmVU9UEzZ8/H/Hx8U2u061bN9P/z58/jxEjRmDIkCF4++23b7p9rVYLrZatJIiIyDa2nSio7bLQyFx6giDgwxMFWDAmWobIlMmpkqCOHTuiY8eOzVr33LlzGDFiBPr374+NGzda7X1GRERkL2cvX22ygenZy1ftHJGyOVUS1Fznzp3D8OHD0bVrV7zyyiv49ddfTY8FBwfLGBlRLaU3sVX6/pNjkOM47NLes8mRoC7tPe0WC7loErRv3z7k5OQgJycHXbp0sXiM7TzIEdQ1sVUqpe8/OQY5jsOHBoTirUM/N/qYJEl4WIFNTOXkkueI4uPjIUlSozciIiK5RHTwxsrJvaESAJVQu6zu/ysn91ZsI1O5uORIEBERkaN6cEAoBoYH4I0DOfgo/Swm9euC2SMimQDJwCVHgoiIiBxZeAdvxA8JBwDEDwmXJQFibR5HgoiIiBSJtXkcCSIiIiKFYhJEREREisQkiIiIiBSJSRAREREpEpMgIiIiUiQmQURERKRITIKIiIhIkZgEERERkSIxCSIiIiJFYhJEREREisQkiIiIiBSJSRAREREpEpMgIiIiUiQmQURERDLo5KvFnFFR6OSrlTsUxXKTOwAiIiIl6uTngXmjb5M7DEXjSBAREREpEpMgIiIiUiSeDmuCJEkAgLKyMpkjISIiouaq+7td93fcGiZBTSgvLwcAhIaGyhwJERERtVR5eTn8/f2tPi5IN0uTFEwURZw/fx6+vr4QBMG0vKysDKGhoSgoKICfn5+MEcpH6T8D7r+y9x/gz0Dp+w/wZ+DI+y9JEsrLyxESEgKVynrlD0eCmqBSqdClSxerj/v5+TncL97elP4z4P4re/8B/gyUvv8AfwaOuv9NjQDVYWE0ERERKRKTICIiIlIkJkGtoNVqsWTJEmi1yp3lU+k/A+6/svcf4M9A6fsP8GfgCvvPwmgiIiJSJI4EERERkSIxCSIiIiJFYhJEREREisQkiIiIiBSJSVArJCcnIzw8HB4eHhg8eDCOHTsmd0h2sXz5cgwcOBC+vr7o1KkTJk6ciKysLLnDks2KFSsgCALmzp0rdyh2de7cOfzpT39CYGAgPD09ceedd+LEiRNyh2UXRqMRzz33HCIiIuDp6Ynu3bvjhRdeuGl/Imf29ddfY/z48QgJCYEgCNi5c6fF45IkYfHixdDpdPD09ERsbCyys7PlCdYGmtp/g8GABQsW4M4774S3tzdCQkIwbdo0nD9/Xr6AbeBmx4C5J554AoIgYPXq1XaL71YwCWqhDz/8EElJSViyZAkyMjLQp08fxMXFoaSkRO7QbO7QoUNISEjAkSNHsG/fPhgMBtx3332orKyUOzS7O378ON566y307t1b7lDs6vLlyxg6dCjc3d3xxRdf4NSpU3j11VfRvn17uUOzi5UrV2L9+vV444038OOPP2LlypVYtWoV1q5dK3doNlNZWYk+ffogOTm50cdXrVqFNWvW4M0338TRo0fh7e2NuLg4XLt2zc6R2kZT+19VVYWMjAw899xzyMjIwCeffIKsrCxMmDBBhkht52bHQJ0dO3bgyJEjCAkJsVNkbUCiFhk0aJCUkJBgum80GqWQkBBp+fLlMkYlj5KSEgmAdOjQIblDsavy8nIpKipK2rdvn3TvvfdKc+bMkTsku1mwYIE0bNgwucOQzdixY6UZM2ZYLJs0aZI0ZcoUmSKyLwDSjh07TPdFUZSCg4Oll19+2bTsypUrklarlT744AMZIrSt+vvfmGPHjkkApDNnztgnKDuz9jM4e/as1LlzZ+nkyZNS165dpddff93usbUGR4JaoLq6Gunp6YiNjTUtU6lUiI2NxeHDh2WMTB6lpaUAgICAAJkjsa+EhASMHTvW4jhQis8++wwDBgzAgw8+iE6dOuGuu+7CO++8I3dYdjNkyBCkpaXhp59+AgD897//xbfffov7779f5sjkkZubi6KiIov3gr+/PwYPHqzIz0Sg9nNREAS0a9dO7lDsRhRFTJ06FX/5y19wxx13yB1Oi7CBagtcuHABRqMRQUFBFsuDgoJw+vRpmaKShyiKmDt3LoYOHYpevXrJHY7dbN26FRkZGTh+/Ljcocjil19+wfr165GUlIRnn30Wx48fx9NPPw2NRoPp06fLHZ7NLVy4EGVlZYiOjoZarYbRaMRLL72EKVOmyB2aLIqKigCg0c/EuseU5Nq1a1iwYAEeeeQRh2woaisrV66Em5sbnn76ablDaTEmQdQqCQkJOHnyJL799lu5Q7GbgoICzJkzB/v27YOHh4fc4chCFEUMGDAAf//73wEAd911F06ePIk333xTEUnQtm3b8P7772PLli244447kJmZiblz5yIkJEQR+0/WGQwGPPTQQ5AkCevXr5c7HLtJT0/HP/7xD2RkZEAQBLnDaTGeDmuBDh06QK1Wo7i42GJ5cXExgoODZYrK/mbPno3du3fjwIED6NKli9zh2E16ejpKSkrQr18/uLm5wc3NDYcOHcKaNWvg5uYGo9Eod4g2p9Pp0LNnT4tlt99+O/Lz82WKyL7+8pe/YOHChfjDH/6AO++8E1OnTsW8efOwfPlyuUOTRd3nntI/E+sSoDNnzmDfvn2KGgX65ptvUFJSgrCwMNPn4pkzZzB//nyEh4fLHd5NMQlqAY1Gg/79+yMtLc20TBRFpKWlISYmRsbI7EOSJMyePRs7duzA/v37ERERIXdIdjVq1Cj88MMPyMzMNN0GDBiAKVOmIDMzE2q1Wu4QbW7o0KENpkX46aef0LVrV5kisq+qqiqoVJYfm2q1GqIoyhSRvCIiIhAcHGzxmVhWVoajR48q4jMRuJEAZWdn46uvvkJgYKDcIdnV1KlT8f3331t8LoaEhOAvf/kL9u7dK3d4N8XTYS2UlJSE6dOnY8CAARg0aBBWr16NyspKPProo3KHZnMJCQnYsmULPv30U/j6+prO+fv7+8PT01Pm6GzP19e3Qf2Tt7c3AgMDFVMXNW/ePAwZMgR///vf8dBDD+HYsWN4++238fbbb8sdml2MHz8eL730EsLCwnDHHXfgP//5D1577TXMmDFD7tBspqKiAjk5Oab7ubm5yMzMREBAAMLCwjB37ly8+OKLiIqKQkREBJ577jmEhIRg4sSJ8gXdhpraf51Oh9///vfIyMjA7t27YTQaTZ+LAQEB0Gg0coXdpm52DNRP/Nzd3REcHIwePXrYO9SWk/vyNGe0du1aKSwsTNJoNNKgQYOkI0eOyB2SXQBo9LZx40a5Q5ON0i6RlyRJ2rVrl9SrVy9Jq9VK0dHR0ttvvy13SHZTVlYmzZkzRwoLC5M8PDykbt26SX/9618lvV4vd2g2c+DAgUbf99OnT5ckqfYy+eeee04KCgqStFqtNGrUKCkrK0veoNtQU/ufm5tr9XPxwIEDcofeZm52DNTnTJfIC5LkwlOdEhEREVnBmiAiIiJSJCZBREREpEhMgoiIiEiRmAQRERGRIjEJIiIiIkViEkRERESKxCSIiIiIFIlJEBERESkSkyAiIiJSJCZBREREpEhMgoiIiEiRmAQRERGRIv0/YtSkSQp7Sq4AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, (ax, rax) = plt.subplots(2, 1, height_ratios=[2,1])\n",
    "fig.subplots_adjust(wspace=0, hspace=0)\n",
    "hbkg.plot(label='Bkg', alpha=1, ax=ax)\n",
    "hbkg.plot(histtype='band', alpha=0.3, label='Bkg Unc', ax=ax, facecolor='red')\n",
    "hdata.plot(label='Data', alpha=1, ax=ax, histtype='errorbar', c='k')\n",
    "rax.axhline(0)\n",
    "hep.histplot(diff/np.sqrt(hdata.variances()), yerr=1, ax=rax, histtype='errorbar')\n",
    "ax.legend()\n",
    "rax.set_ylabel(r\"$\\frac{Data-Bkg}{\\sigma_{Data}}$\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 586,
   "id": "ff334709",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "111.9553304407635"
      ]
     },
     "execution_count": 586,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# chi2 no-corr\n",
    "np.sum(np.nan_to_num((hbkg.values() - hdata.values())**2/(hbkg.variances()+hdata.variances()), posinf=0, neginf=0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 587,
   "id": "faf6ef4c",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_28880/1569011114.py:2: RuntimeWarning: divide by zero encountered in divide\n",
      "  np.sum(np.nan_to_num((hbkg.values() - hdata.values())**2/(hdata.variances()), posinf=0, neginf=0))\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "16.83030987844864"
      ]
     },
     "execution_count": 587,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# chi2 no-corr - data-only - should match \"pull\" plot above\n",
    "np.sum(np.nan_to_num((hbkg.values() - hdata.values())**2/(hdata.variances()), posinf=0, neginf=0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 588,
   "id": "aac77083",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0.10068152,  0.20043927,  2.47086453,  6.03295708,  5.94815826,\n",
       "        2.74485588,  7.3138895 , 17.87661362,  6.3735218 ,  0.81303436,\n",
       "        0.16079427,  0.14357626,  0.02831182,  0.09442087,  0.00255251,\n",
       "        0.00000513,  0.00000003])"
      ]
     },
     "execution_count": 588,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hbkg.variances()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 589,
   "id": "3ee0dbf7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 4.,  3., 12., 12., 13., 25., 30., 48., 29., 19., 12.,  4.,  1.,\n",
       "        3.,  1.,  0.,  0.])"
      ]
     },
     "execution_count": 589,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hdata.variances()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 590,
   "id": "bdaa0416",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.00100682,  0.00096771,  0.00468452,  0.00736983,  0.00675476,\n",
       "         0.0028593 ,  0.00021518, -0.00019035, -0.000652  , -0.00023233,\n",
       "         0.00017659, -0.00032066,  0.00009207,  0.00002338,  0.00002052,\n",
       "        -0.00000025,  0.        ],\n",
       "       [ 0.00096771,  0.00200439,  0.00416386,  0.00641046,  0.00660311,\n",
       "         0.00425275,  0.00160866,  0.00117741,  0.00099808,  0.00055142,\n",
       "         0.00025924,  0.00030966,  0.00011217,  0.00010249,  0.00001627,\n",
       "         0.0000002 , -0.00000001],\n",
       "       [ 0.00468452,  0.00416386,  0.02470865,  0.0385591 ,  0.03749459,\n",
       "         0.01758856,  0.00174191, -0.00126097, -0.00222593, -0.00051703,\n",
       "         0.00089978, -0.00105717,  0.00045223,  0.00022432,  0.00010223,\n",
       "        -0.00000066, -0.00000003],\n",
       "       [ 0.00736983,  0.00641046,  0.0385591 ,  0.06032957,  0.05847818,\n",
       "         0.02715217,  0.0024799 , -0.00214338, -0.00377468, -0.00096147,\n",
       "         0.00140732, -0.0017724 ,  0.00070965,  0.00033617,  0.0001609 ,\n",
       "        -0.00000113, -0.00000004],\n",
       "       [ 0.00675476,  0.00660311,  0.03749459,  0.05847818,  0.05948158,\n",
       "         0.03114725,  0.00604244,  0.00161539,  0.00041527,  0.00105811,\n",
       "         0.0017713 , -0.00040402,  0.00081929,  0.00053551,  0.00015448,\n",
       "        -0.00000002, -0.00000006],\n",
       "       [ 0.0028593 ,  0.00425275,  0.01758856,  0.02715217,  0.03114725,\n",
       "         0.02744856,  0.02872228,  0.03909862,  0.02462115,  0.00965889,\n",
       "         0.00421004,  0.00223054,  0.00137368,  0.00079486,  0.00012538,\n",
       "         0.00000253,  0.00000021],\n",
       "       [ 0.00021518,  0.00160866,  0.00174191,  0.0024799 ,  0.00604244,\n",
       "         0.02872228,  0.07313889,  0.11355486,  0.06754183,  0.02305502,\n",
       "         0.00930174,  0.00297292,  0.00245437,  0.00110938,  0.00018951,\n",
       "         0.00000518,  0.00000087],\n",
       "       [-0.00019035,  0.00117741, -0.00126097, -0.00214338,  0.00161539,\n",
       "         0.03909862,  0.11355486,  0.17876613,  0.10524943,  0.03495969,\n",
       "         0.01412612,  0.00344552,  0.00354776,  0.00154263,  0.00028983,\n",
       "         0.0000076 ,  0.00000142],\n",
       "       [-0.000652  ,  0.00099808, -0.00222593, -0.00377468,  0.00041527,\n",
       "         0.02462115,  0.06754183,  0.10524943,  0.06373522,  0.02117023,\n",
       "         0.00802692,  0.00265663,  0.00188043,  0.00099803,  0.00015922,\n",
       "         0.00000531,  0.0000008 ],\n",
       "       [-0.00023233,  0.00055142, -0.00051703, -0.00096147,  0.00105811,\n",
       "         0.00965889,  0.02305502,  0.03495969,  0.02117023,  0.00813034,\n",
       "         0.00332893,  0.0019866 ,  0.00107631,  0.00043158,  0.00004885,\n",
       "         0.0000018 ,  0.00000024],\n",
       "       [ 0.00017659,  0.00025924,  0.00089978,  0.00140732,  0.0017713 ,\n",
       "         0.00421004,  0.00930174,  0.01412612,  0.00802692,  0.00332893,\n",
       "         0.00160794,  0.00080232,  0.00059472,  0.00017412,  0.00002516,\n",
       "         0.00000049,  0.0000001 ],\n",
       "       [-0.00032066,  0.00030966, -0.00105717, -0.0017724 , -0.00040402,\n",
       "         0.00223054,  0.00297292,  0.00344552,  0.00265663,  0.0019866 ,\n",
       "         0.00080232,  0.00143576,  0.00044125,  0.00015434, -0.00000562,\n",
       "         0.00000047, -0.00000001],\n",
       "       [ 0.00009207,  0.00011217,  0.00045223,  0.00070965,  0.00081929,\n",
       "         0.00137368,  0.00245437,  0.00354776,  0.00188043,  0.00107631,\n",
       "         0.00059472,  0.00044125,  0.00028312,  0.00006403,  0.00000636,\n",
       "         0.00000007,  0.00000002],\n",
       "       [ 0.00002338,  0.00010249,  0.00022432,  0.00033617,  0.00053551,\n",
       "         0.00079486,  0.00110938,  0.00154263,  0.00099803,  0.00043158,\n",
       "         0.00017412,  0.00015434,  0.00006403,  0.00003777,  0.00000274,\n",
       "         0.00000011,  0.00000001],\n",
       "       [ 0.00002052,  0.00001627,  0.00010223,  0.0001609 ,  0.00015448,\n",
       "         0.00012538,  0.00018951,  0.00028983,  0.00015922,  0.00004885,\n",
       "         0.00002516, -0.00000562,  0.00000636,  0.00000274,  0.00000102,\n",
       "         0.00000001,  0.        ],\n",
       "       [-0.00000025,  0.0000002 , -0.00000066, -0.00000113, -0.00000002,\n",
       "         0.00000253,  0.00000518,  0.0000076 ,  0.00000531,  0.0000018 ,\n",
       "         0.00000049,  0.00000047,  0.00000007,  0.00000011,  0.00000001,\n",
       "         0.        ,  0.        ],\n",
       "       [ 0.        , -0.00000001, -0.00000003, -0.00000004, -0.00000006,\n",
       "         0.00000021,  0.00000087,  0.00000142,  0.0000008 ,  0.00000024,\n",
       "         0.0000001 , -0.00000001,  0.00000002,  0.00000001,  0.        ,\n",
       "         0.        ,  0.        ]])"
      ]
     },
     "execution_count": 590,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.set_printoptions(suppress=True)\n",
    "acov = cov.values()\n",
    "acov"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d0db87c1",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 591,
   "id": "268f9c30",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "111.9553304407635"
      ]
     },
     "execution_count": 591,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Baseline chi2 no cov\n",
    "np.sum(diff**2 / (hbkg.variances() + hdata.variances()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 592,
   "id": "2a88c59b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "256445370828.56155"
      ]
     },
     "execution_count": 592,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Chi2 MC only cov according to prescription - seem to high?\n",
    "diff @ np.linalg.inv(acov) @ diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 593,
   "id": "6ed90a83",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "8.213088800632335"
      ]
     },
     "execution_count": 593,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Chi2 if cov matrix from combine is already inverted? - seems too low?\n",
    "diff @ acov @ diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 594,
   "id": "68c60fa1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2532.5493883569693"
      ]
     },
     "execution_count": 594,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Chi2 according to prescription - seem still to high?\n",
    "diff @ (acov+np.diag(hdata.variances())) @ diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 595,
   "id": "93a89062",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.8149768104471574"
      ]
     },
     "execution_count": 595,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Chi2 = undo C-1 from combine, add add variances, redo invert - maybe?\n",
    "diff @ np.linalg.inv(np.linalg.inv(acov)+np.diag(hdata.variances())) @ diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5593a869",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3fbc0030",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e5989a7c",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:combine]",
   "language": "python",
   "name": "conda-env-combine-py"
  },
  "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.10.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}