cms_final_lbt_mart_ads_r1.0_0-10.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "88c93f3a",
   "metadata": {},
   "source": [
    "## Jet Analysis(Nuclear Modification Factor)( $R_{\\mathrm{AA}}$)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b06e4546",
   "metadata": {},
   "source": [
    "# Preparations\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1dfb1c9a",
   "metadata": {},
   "source": [
    "# Loading all the files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "dcfa051a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set File Paths\n",
    "# Please set the path for pp jet data file\n",
    "file_pp = '../../data/music_cms/cms_pp_vir0.5_R1.0_trkpT4.0_maxT250_pt200-1001_m1-7_merged_389140.dat'\n",
    "\n",
    "# Please set the value of jet cone size you used in the jet reconstruction\n",
    "jetR = 1.0\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "d70a0423",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Please set the path for pbpb jet data file\n",
    "file_pbpb = '../../data/music_cms/LBT/R1.0/cms_music_0-10_R1.0_s1-7_35k.dat'\n",
    "file_pbpb1 = '../../data/music_cms/MARTINI/R1.0/cms_martini_0-10_R1.0_s1-6_12790.dat'\n",
    "file_pbpb2 = '../../data/music_cms/ADS/R1.0/cms_ads_0-10_R1.0_s1-4_6461.dat'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "fcee61ca",
   "metadata": {},
   "outputs": [],
   "source": [
    "# preperations\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# pi and 2pi \n",
    "pi = 3.1415926\n",
    "twopi = 2.0*pi\n",
    "\n",
    "# define plot style\n",
    "width = 0.05\n",
    "plotMarkerSize = 8\n",
    "labelfontsize = 15\n",
    "import matplotlib as mpl\n",
    "mpl.rcParams['figure.figsize'] = [8, 6]\n",
    "mpl.rcParams['lines.linewidth'] = 2\n",
    "mpl.rcParams['xtick.top'] = True\n",
    "mpl.rcParams['xtick.labelsize'] = 15\n",
    "mpl.rcParams['xtick.major.width'] = 1.0\n",
    "mpl.rcParams['xtick.minor.width'] = 0.8\n",
    "mpl.rcParams['xtick.minor.visible'] = True\n",
    "mpl.rcParams['xtick.direction'] = \"in\"\n",
    "mpl.rcParams['ytick.right'] = True\n",
    "mpl.rcParams['ytick.labelsize'] = 15\n",
    "mpl.rcParams['ytick.major.width'] = 1.0\n",
    "mpl.rcParams['ytick.minor.width'] = 0.8\n",
    "mpl.rcParams['ytick.minor.visible'] = True\n",
    "mpl.rcParams['ytick.direction'] = \"in\"\n",
    "mpl.rcParams['legend.fontsize'] = 15\n",
    "mpl.rcParams['legend.numpoints'] = 1\n",
    "mpl.rcParams['font.size'] = 15\n",
    "mpl.rcParams['savefig.format'] = \"pdf\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "707d8005",
   "metadata": {},
   "outputs": [],
   "source": [
    "def ratio_error(v1,e1,v2,e2):\n",
    "  #v1, e1: numerator value and error\n",
    "  #v2, e2: denominator value and error  \n",
    "  error1 = e1/v2\n",
    "  error2 = (e2/v2)*(v1/v2)\n",
    "  error = np.sqrt(error1*error1+error2*error2)\n",
    "  return error"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "fe7826ff",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load files\n",
    "data_pp = np.loadtxt(file_pp, delimiter=',')\n",
    "data_pbpb = np.loadtxt(file_pbpb, delimiter=',')\n",
    "data_pbpb1 = np.loadtxt(file_pbpb1, delimiter=',')\n",
    "data_pbpb2 = np.loadtxt(file_pbpb2, delimiter=',')\n",
    "\n",
    "# Indices of the data array\n",
    "i_pp = data_pp[:,0] \n",
    "i_pbpb = data_pbpb[:,0] \n",
    "i_pbpb1 = data_pbpb1[:,0]\n",
    "i_pbpb2 = data_pbpb2[:,0]\n",
    "\n",
    "\n",
    "# Get Indices of jets in the data array\n",
    "jet_id_pp = np.where(i_pp < 0.1)\n",
    "jet_id_pbpb = np.where(i_pbpb < 0.1)\n",
    "jet_id_pbpb1 = np.where(i_pbpb1 < 0.1)\n",
    "jet_id_pbpb2 = np.where(i_pbpb2 < 0.1)\n",
    "\n",
    "# Extract jets\n",
    "jets_pp = data_pp[jet_id_pp]\n",
    "jets_pbpb = data_pbpb[jet_id_pbpb]\n",
    "jets_pbpb1 = data_pbpb1[jet_id_pbpb1]\n",
    "jets_pbpb2 = data_pbpb2[jet_id_pbpb2]\n",
    "\n",
    "# Extract associated charged particles for pp\n",
    "assoc_pp= []\n",
    "for i in range(len(jet_id_pp[0])-1):\n",
    "  chunck = data_pp[jet_id_pp[0][i]+1:jet_id_pp[0][i+1]]\n",
    "  assoc_pp.append(chunck)\n",
    "chunck = data_pp[jet_id_pp[0][-1]+1:]\n",
    "assoc_pp.append(chunck)\n",
    "\n",
    "# Extract associated charged particles for pbpb\n",
    "assoc_pbpb= []\n",
    "for i in range(len(jet_id_pbpb[0])-1):\n",
    "  chunck = data_pbpb[jet_id_pbpb[0][i]+1:jet_id_pbpb[0][i+1]]\n",
    "  assoc_pbpb.append(chunck)\n",
    "chunck = data_pbpb[jet_id_pbpb[0][-1]+1:]\n",
    "assoc_pbpb.append(chunck)\n",
    "\n",
    "assoc_pbpb1= []\n",
    "for i in range(len(jet_id_pbpb1[0])-1):\n",
    "  chunck = data_pbpb1[jet_id_pbpb1[0][i]+1:jet_id_pbpb1[0][i+1]]\n",
    "  assoc_pbpb1.append(chunck)\n",
    "chunck = data_pbpb1[jet_id_pbpb1[0][-1]+1:]\n",
    "assoc_pbpb1.append(chunck)\n",
    "\n",
    "assoc_pbpb2= []\n",
    "for i in range(len(jet_id_pbpb2[0])-1):\n",
    "  chunck = data_pbpb2[jet_id_pbpb2[0][i]+1:jet_id_pbpb2[0][i+1]]\n",
    "  assoc_pbpb2.append(chunck)\n",
    "chunck = data_pbpb2[jet_id_pbpb2[0][-1]+1:]\n",
    "assoc_pbpb2.append(chunck)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a38aafd6",
   "metadata": {},
   "source": [
    "## Jet-$R_{\\mathrm{AA}}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "760626b2",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi4AAAGdCAYAAAA1/PiZAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAACCeklEQVR4nO3dd3gUVdvA4d8hBFJoCeWF0DsRpUakWEA6UhQCBGygAhY+lWKj5gXEhmBFiC8KKiZAQCnSpQiCCqiANEGIIjVIACUESPJ8f2x23c1ukk3YEDZ57uuaK9kzZ86cmQnkyWljRASllFJKKW9QKK8roJRSSinlLg1clFJKKeU1NHBRSimllNfQwEUppZRSXkMDF6WUUkp5jcJ5XYH8qFixYlKvXr28robKBfHx8ZQtWzavq6Hs6DPJPXpvs0/vmaMdO3YkikigJ8s0Oh3a8wIDA+XixYt5XQ2VC8LCwti+fXteV0PZ0WeSe/TeZp/eM0fGGETEeLJM7SrKBWXKlMnrKqhcMnjw4LyugkpHn0nu0XubfXrPcp+2uOSCsLAw0YhbKaVUQactLkoppZTyJic8XaAGLkoppZTKLcc9XaAGLkoppZTyGhq4KKWUUspr6DoueeTChQucPn2aq1ev5nVVlFJu8PX1pVy5cpQoUSKvq6JUgaaBSy6Ij48nLCyMwYMHu5wad+HCBU6dOkXFihXx9/fHGI8OuFZKeZiIcOnSJY4dOwagwYtSWYiKiiIqKgrA4+uD6HToXJDVdOhDhw4REhJCQEDAdayVUupaJSYmcvz4cWrVqpXXVVHKKxhjdohImCfL1DEueeDq1av4+/vndTWUUtnk7++v3btK5TENXPKIdg8p5X30361SeU8Dl3zAGKP/oSqllCoQ8k3gYoypZYyZaYzZZYxJMcZscOOYasYYcbHFuMjbwxiz2xiTZIzZa4zpmysXopRSSqkM5ZvABagPdAEOAL9m89iRQAu7bYz9TmPM7cBCYD3QGfgKiDbGdLjGOueItYUlfUtLRum5aeHChdx9992UKlWKokWLUqdOHYYPH87x4/8ulmitz2effeZ0/Geffeayvr///jsPPvggVapUwc/Pj8qVK9OjRw+++eabTOtjf/0///yzLb1atWqMHDnS9jkyMjLDl2EOGDCAsLB/x5Jt2LDB6d5at8cee8yWr3Xr1rb0woULU7NmTYYNG8aFCxcyrbOIMHnyZCpXroy/vz933nmnQ90BZs+e7fL8M2bMcCpv0aJF1KpVi5SUFAC2bdvGwIEDqVWrFgEBAdStW5f//ve/JCUlOR377bffctttt+Hn50f16tV55513HPYfOHCAp556itDQUAICAqhRowbPPPMM586dc8h36NAhhgwZQoMGDfDx8aF169ZO5zp27BjFihXj8OHDOToHQExMDE2aNKFYsWJUrFiRhx56yOFnL71hw4ZhjHH4WQDLz4P1noaHh2d4vFIq7+Wn6dBLRWQxgDEmluxNwTogIt9lsn8s8I2IPJ32eb0xpj4wDlido9rmAyNGjOCtt95i4MCBDBs2jBIlSrB3715mzJjBkSNH+OKLL2x5ixUrRkxMDA888IBDGdHR0RQrVox//vnHlpaQkEDz5s2pUKECr7zyCiEhIcTFxbFkyRK2bt3KnXfemWW9wsPDqVOnjkevd+7cudSoUcMhrVy5cg6f27Rpw+TJk0lOTmbbtm2MHTuWo0ePEhsbm2G5r776KhMnTuSNN96gXr16TJ06lXbt2vHLL79Qvnx5h7zr1q1zGNidvj6pqamMGzeO5557Dh8fHwDmzZvHb7/9xgsvvEDt2rXZtWsXY8eOZdeuXSxcuNB27KFDh+jYsSNdu3bllVde4YcffmD48OEEBATYArQ1a9bw7bff8sQTT9CgQQMOHz7MmDFj2Lp1K9999x2FCln+FtqzZw/Lly+nefPmGQ5mrVixIn379mXChAnMnj3blu7uOZYsWUK/fv146qmneOONNzhx4gRjxozhnnvuYceOHbZ8Vnv37mXWrFkupzI/9thjdOrUiSeffDLD56SUukGISL7bgFhggxv5qgECdM0kT1HgCvB4uvSHgBSgZPpjmjZtKpnZu3dvpvuzklZnAclgs+zPTUuWLBFAZs2a5bQvOTlZli9f7lDffv36ia+vr5w9e9aW/tdff4mvr6/079/fob5RUVFijJFTp045lZ2ampppvQB59913ndKrVq0qI0aMsH0eP368lC5d2mUZDz/8sNg/w/Xr1wsgu3fvzvTcd911l/Tq1cshbeLEiQLI6dOnXR5z6dIlKVGihPz3v/+1pf3zzz9SpkwZGT16tC3t448/FkD+/vvvTOuwevVq8fPzkwsXLtjS4uPjnfLNnDlTAImLi7OlDR48WGrXri1Xr161pT3xxBNSqVIl230/c+aM0zNYtWqVALJhwwZbWkpKiu37Xr16yV133eWyvps2bZIiRYrImTNnbGnunqNv377SpEkTh3yLFy8WwOW/sbvvvlvGjBnj9LNgz9UzTO9a//0qVZAA28XDv+PzU1fRtfg4bVzMCWPMVGOM/VzlmoAvsD/dMfuwdLV59s96LzFt2jSaNGnCI4884rTPx8eHzp07O6S1aNGCkJAQh7/wFy5cSEhICC1atHDIe+7cOYoUKUJwcLBT2d42CLlp06YAxMXFudy/ZcsWLly4QJ8+fWxpgYGBdOvWjRUrVmT7fHPmzKFDhw4UL17cluaqS6xx48YADt0qK1asoGfPnhQu/G9DbEREBH/++Se//PILAKVLl3Z6Bq7KSt/akZFWrVoRHBxMTMy/w8rcPcfVq1cpWbKkQ75SpUoBWP+4sImNjWX//v28+OKLbtVLKXXjKuiBy2XgfeBRoC0wE3gCsB+cG5T29Vy6YxPS7bc5fvy4wziEyMhIT9Y5z129epUtW7bQqVMnt48xxtC3b1+io6NtadHR0URERDjlbdKkCZcvX+bBBx9kx44dpKameqTeriQnJztt6X/pWaWkpLid18oasKTv8rHav38/Pj4+1K5d2yE9NDSU/fvTx8pQs2ZNChcuTN26dZk5c6bT/nXr1tGyZctM6wSwdetWChUqRM2aNQG4ePEiR48epV69ek71sNYzs7KAHHXNGWNo3rw5a9euzbK+6c/xyCOPsGnTJj755BMuXLjAr7/+ypgxY7j77ru56aabbPkuXbrEiBEjePXVVwkMDMx2HZVS2WM/ZgwI8XT5BTpwEZETIjJURJaIyAYRiQSGA92NMQ1zWm5ISIhDs1Z+C1z++usvLl++TJUqVbJ1XEREBBs3buTUqVOcPHmSjRs3ugxc2rZty7Bhw5g3bx5hYWGUKlWKXr16ZfnLLbv++usvfH19nbZPPvnEZf5GjRo55Z0zZ45DHhEhOTmZy5cvs2nTJl5++WXCwsKoVKmSyzITEhIoVqyYbTyKVVBQEImJiVy5cgWAChUqMHHiRD799FOWLl1K8+bNefzxx5k2bZrtmOPHj3PixAluvvnmTK/75MmTTJo0iQcffNA2Rsc68NXaYmFfD2s9XUlMTOSFF17grrvusrUuZVfDhg354YcfMtyf0TnuueceZs+ezeDBgylZsiR169YlJSXFoVUP4JVXXqFChQpO46uUUrkjMjLSflhFxqPlcyg/Dc71lFhgOtAU2Mm/LSsl0+WztrS4/h+9AMhut03jxo2pVasW8+fPR0SoU6cOjRo1YvPmzU55p06dypNPPsmXX37JN998w8qVK/niiy+YPn06jz/+uEfqX7JkSZfB0H//+19OnDjhlB4TE2NrobCqXr26w+dFixbh6+tr+9yqVSs++uija+7i6tixIx07drR97ty5M0lJSUyaNIlnnnmGQoUKcfLkScB115DVlStX6NOnD8WKFXMIenJCRHj00Uc5ffo0X331VY7LKVOmDKdPn0ZEnO5TZudYv349jz/+OM888wydO3fm1KlTREZGct9997F27Vp8fHw4cuQIU6ZMYf369V7XzaiUck0DF2eS7utvwFWgHrDRLl89IJXsT732eqVLl6Zo0aL88ccf2T62b9++xMTEICL07Zv5Uji1atVi5MiRjBw5kjNnztChQwdGjRrFkCFDPPJLqHDhwg7Tnq1Kly7tMnCpX79+lq0Zd999N6+99hqFCxematWqthaLjAQFBfHPP/+QkpLi0OqSkJBAQEAARYoUyfDY8PBw5s+fT1xcHDVq1LBNby5atKjL/CLCQw89xJ49e/j2228d6mZtaTl//rzDMdaWFlfX8cILL/DFF1+wZs0ap9lN2VG0aFFb15t90JfVOUaMGEH37t157bXXbGmNGjWiXr16LF68mJ49e/Liiy/SuXNn6tata2tVSk1N5fLly5w7d46SJUtqQKOUlynQXUUZsC7isANARC5jWb+ld7p8fYGtInKePGMy2HKXr68vrVq1YtWqVdk+NiIiwjat1VU3UUbKlCnDwIEDSUhI4PTp09k+7/USFBREWFgYjRo1yjJoAahXrx4pKSkcOnTIIX3//v1O403Ss/7CtX61DmZ2td4JwLPPPsvixYtZvHixU9mBgYFUrlzZaSyL9XP6/NOmTWPKlCl88skn3HHHHZnWMyvnzp2jWLFiTkFLVufYv38/jRo1ckirW7cu/v7+/Pbbb4BlTZhFixYRFBRk244ePcp7771HUFCQ7W3PSinvkW8CF2NMgDEm3BgTDlQEylo/G2MC0vIcMsbMsjsm0hjzpjGmpzGmnTFmAjANWCQiu+yKnwi0Nsa8ZYxpbYx5HctidxOu3xXeWJ599lm2b9/uNMYDLH/Rrly50uVxoaGhDB48mMGDB2f4izk+Pt5l+sGDBylatKjTTBJv1rJlS0qUKMGCBQtsaYmJiSxdutRpZlZ6sbGxlClThqpVqwKWRfaKFCnCkSNHnPK+8sorvPfee3z22WfcfvvtLsvr3LkzX3zxhW3hOrCsAVO5cmWHlqa5c+cyYsQIpk6d6jAbKqfi4uKcBva6c46qVavy448/OqTt27ePS5cuUa1aNQD+97//sX79eoftP//5D3369GH9+vWULVv2muuvlLq+8lNXUTlgQbo06+fqQByW67UfBbkfy6q5jwH+wB/AG8DL9oWIyOa0gGgSlllHR4D+IpIni8+ln8li/Ys7qxkuntStWzeGDx/Oo48+yrfffkuPHj0oVqwY+/fvZ8aMGVSrVi3DWUeuVnu1N2fOHObOnctDDz1Ew4YNuXr1KmvXrmX69Ok88cQT+Pn55ajO19olsGvXLoeF8sAyTsY68yYn/Pz8ePHFF5k4cSJBQUG2BehSU1P5v//7P1u+Xr160axZMxo0aEBKSgrz5s1j3rx5vPPOO7apx35+fjRt2pQdO3YwcOBA27Gff/45o0aNYsCAAVSsWJHvvvt3rcWaNWvafnk/99xzzJ07lwcffJBBgwaxbds2Zs6cyQcffGC7dxs3bmTgwIF06NCB5s2bO5RVqVIl2yDkxMREli9fDlhWyL1w4YJtEb4uXboQEBBgO2779u20atXK9tndczz++OMMGzaMkJAQ2xiXCRMmUK1aNbp06QLgsivQuhKzq9V8lVJewNMLw+iW+wvQpcd1WHAuI7GxsdK6dWspUaKE+Pr6Su3atWXEiBFy4sQJh/q5WhTO6t1333Wo/549e+TJJ5+U0NBQKVasmJQsWVKaNGki06dPd1gczZWMzlWuXDkZNWqU7XNOFqBztbVt29aWz53Fy1xJTU2VSZMmScWKFcXPz09uv/12+fHHHx3yvPTSS1KnTh3x9/cXPz8/adKkiXzyySdOZb3xxhtSs2ZNp+vJqP4ff/yxQ95NmzbJrbfeKkWLFpWqVavK22+/7bB//PjxGZY1fvx4W74jR45kmO/IkSO2fKdPnxYfHx+HheXcPUdqaqpMnz5dbrnlFgkICJCQkBDp06eP/Pbbb5neb12ATqnrh1xYgM7IdfwrvaAICwuT7du3Z7h/37591/RXenp50eJyozLG8Pbbb/Pkk0/aFlI7f/48pUuX5v3332fIkCF5XMPcderUKapUqcLmzZu59dZb87o6WZo5cyZTpkzh119/zdNBsqmpqaSmptK2bVvKli2b6SsaPP3vV6n8zBizQ0Scmz6vQb4Z43IjiY+PJywsjKioqOtyPmsUqiyeeeYZfH19+f7771m/fj1DhgzBx8fH1n2Qn/3nP//hscce4+23387rqmRJRHj77bcZPXp0ns/smTBhAr6+vlm+xFMp5Z6oqChrV2123hvolvw0xuWGUbZsWTJrcVG5Z9u2bbbvAwMDufPOOwkNDWXp0qVUrlw5D2t2/YwdO5ZZs2Y5TbG+0Zw8eZL777+fBx98MK+rwuDBg+natSuAy1dNKKWyxzoJwxhzxtNla1dRLrjeXUVKqetH//0q5T7tKlJKKaVUgaaBi1JKKaW8hgYuSimllPIaGrgopZRSymto4JIfGGPZlFJKqXxOAxellFJKeQ0NXHLB9V6ALi9ERkZijKF27dou99euXRtjDJGRkU77qlevjjHG4Y3I1apVwxiT6eZOntmzZxMXF5fh/nbt2tnOOWDAAFt6oUKFqFSpEv369SMuLu6a78+AAQNcvifHasOGDQ718vX1pUaNGjz//PNcvHgRgNmzZ7t1T5RS6kajC9B5meu5AF1CQgKxwEmg/IcfEh4eTlBQ0HU5t5+fH0eOHGH79u0Ov6S3bdtGXFycy5chbt261RYYREdHM3bsWAC++OILLl++bMvXqVMnwsPDeeyxx2xpPj4+Dm8udpWnZs2atl/8U6ZMcXh5H+D0Zul69erx8ccfk5qayv79+xk9ejRdunTh559/pkiRItm9Jdk2d+5catSowdWrV/nxxx8ZPXo0586dIyoqinvuuYetW7fa8sbGxvLmm286pBUtWjTX66iUUtmVmwvQaeDipUSESePGMXXKFDoCtYC1w4bx/NNPM3zkSMZMmJDry6gHBgbSpEkTYmJiHAKXmJgY7r77bnbs2OF0THR0NIGBgdx8880OgUvjxo0d8hUuXJhKlSrRvHnzDM+fUR5r4FK3bt1Mj7degzVPy5YtCQgIoF+/fmzfvp2WLVtmeJwxhvXr11/zG4YbNGjAzTffDMAdd9zBsWPHmD17NlFRUZQtW9b25mbAFgxndU1KKZWfaVeRl5o0bhxfTp3KnqQkYoBJwLyLF9mTlMSXU6cyady461KPiIgI5s+fb3tXkogwf/58IiIinPKmpKQwf/58unfvziOPPMK+ffvYuXPndamnuxo2bAjA0aNH8+T8xYsX5+rVq3lybqWU8gYauHihhIQEpk6ZwtLERELS7QsBliYmMnXKFM6dO5frdenZsyenTp1i8+bNAGzatIn4+Hh69uzplHf9+vWcOnWKiIgIwsPD8fX1JTo6OtfqlpqaSnJyssOWmpqa6TF//PEHYBmHcz2kpKSQnJzMpUuX+Pbbb5k5cyb33XffdTm3Ukp5Iw1cvIV1yrMxxAYH0yEpySlosQoBOiQlsSAoKNenSZcqVYpOnToRExMDWLqJOnXq5DSWBCzdRNb8wcHBdOjQgZiYmFx7s3WPHj3w9fV12CZMmOCULzk5matXr7Jr1y5efPFFOnXqRLNmzZzy2G/wb9Bh3XJyHY0aNcLX15eAgABuv/12qlev7hVvdlZKqbyigYsXOgm4nsvzr1pp+a6HiIgIYmNjuXz5MrGxsS67ia5cucKiRYu47777bINeIyIi+P333x0Gm3rStGnT2LZtm8M2ePBghzw7duzA19eXIkWK0LBhQy5cuODUChQXF+cUAAG0a9fOIW3jxo3ZrmNMTAzbtm3j+++/Jzo6mjNnznDvvfdm2TKklFIFlQ7O9RZ2f82X//BD1g4bBmmDUF05FBhIu2nTYNCgXK9a9+7deeyxxxg9ejQXL16kW7duTnlWrFjBuXPn6NKli60Lq3Xr1hQtWpTo6OhMB8LmVK1atTKdkgwQGhrKJ598wtWrV9m8eTOjRo1iyJAhzJs3z5YnJCSEbdu2ORx36623MmPGDJo2bWpLq1u3brbrWL9+fdvg3GbNmlGnTh2aNm3K8uXL6dq1a7bLU0qp/E4DFy8UHh7O808/zXFw2V10DFidksLM3r2vS30CAwPp2rUr06ZNo3fv3gQGBjrlsbZi9HZRpwULFvDWW2/h4+OT63VNLyAgwBbctGjRgqSkJMaNG8fw4cO57bbbAChSpIjLAKhu3bpZBkbZFRoaCsC+ffs0cFFKKRe0qygX5PYCdEFBQQwfOZJuAQEcT7fvONA9IIDhI0dSqlSpXDm/K0888QTdunXj8ccfd9p38eJFli5dSr9+/Vi/fr3DNnXqVE6dOsW6deuuW10zM2LECMqUKcNrr72WJ+f/5ZdfAKhcuXKenF8ppTxBF6DzMtdjAboxaYNM60+ZQoekJGph6R5anZLC8OHDbfuvl9atW2e4psnixYtJTEzkmWeesbViWLVq1YqXX36Z6Oho2rdv79E6HThwgDJlHP/N+Pn50ahRowyPCQgIYNiwYYwdO5aDBw9muDKwOxISEoiNjXVK79Kli+37Xbt28c8//5Camsrhw4eZOHEiVapUoXPnzjk+r1JK5TVdgE45McYwduJEhg4fTmxwMCeBdtOmMbN37+va0uKO6Ohoateu7RS0APj6+tKnTx8+//xzPvjgA4+uBDty5EintJo1azq8asCVoUOH8sYbbzBlyhRmzpyZ4/MfPnzYZdfYkSNHbN/ff//9gOV5hoSEcMcddzBp0iSXs7KUUkqBya2pqAVZWFiYZNbism/fPttYBo+wTnnWZ6lUrvP4v1+l8jFjzA4R8ehgQG1xyQ80YFFKKVVA6OBcpZRSSnkNDVyUUkop5TU0cFFKKaWU19DARSmllFJeQwOXXJDbC9AppZRSNzJdgM7LXI8F6JRSSqkbVW4uQKctLkoppZTyGhq45AsmbVNKKaXyt3wTuBhjahljZhpjdhljUowxG9w45lZjzMfGmEPGmERjzAFjzHhjjF+6fJHGGHGxdcq1C1JKKaWUk3wTuAD1gS7AAeBXN4/pC9QEXks79n1gODDXRd7zQIt029Zrq/K1S0hI4MMPYeJE+PDDD0lISLgu542MjMQYk+FLCGvXro0xhsjISKd91atXxxjj8M6gatWqYYzJdHMnz+zZs4mLi8twf7t27WznHDBggC29UKFCVKpUiX79+hEXF3fN98datv35rC5dukSJEiVs9U3v4sWLBAYGEhAQwN9//+20f/bs2Q7XVKZMGdq0acM333zjdF0ZbQMGDLDdd/t3Olmfa8eOHZ3OGx4e7vAizQ0bNmCMsb3RWimlrof8NDh3qYgsBjDGxOLeSOZXRcR+4NAGY0wSMNMYU1VEfrfblywi33mwvtdERJg0aRxTp06hY0eoVQvWrh3G888/zfDhIxkzZgLG5G73kZ+fH0eOHGH79u3W0eMAbNu2jbi4OPz8/JyO2bp1qy0wiI6OZuzYsQB88cUXXL582ZavU6dOhIeH89hjj9nSfHx8SElJyTRPzZo1uXjxIgBTpkyhVatWDudP//LCevXq8fHHH5Oamsr+/fsZPXo0Xbp04eeff6ZIkSLZvSUOihUrxoYNGzh16hT/+c9/bOnLli0js3eELVmyhMTERAC+/PJLHnzwQZf51q1bh7+/PydPnmTy5Ml07NiRXbt2MXbsWB5//HFbvlGjRnHu3DmmT59uSytbtmymdV+9ejXbtm3j1ltvdetalVLqesk3gYuIpObgGFejnX9K+xoC/O5i/w1h0qRxfPnlVPbsSSIkxJp6kePHoVu3qQCMHTsxV+sQGBhIkyZNiImJcQhcYmJiuPvuu9mxY4fTMdHR0QQGBnLzzTc7BC6NGzd2yFe4cGEqVapE8+bNMzx/RnmsgUvdunUzPd56DdY8LVu2JCAggH79+rF9+3ZatmyZ4XHGGNavX+/QApFe3bp1+fvvv1mwYAFDhw61pcfExNC9e3c+//xzl8dFR0dTo0YNRITo6OgMA5dbb72VYsWKAdCsWTMqVarEggULGDVqFDVr1rTlCw4OJjU1Nct7YZ+/YsWKvPzyy3z55ZduHaOUUtdLfuoq8pQWQCrwW7r0UsaYM8aYq8aYn4wxPfOgboCle2jq1CksXZpoF7RYhITA0qWJTJ06hXPnzuV6XSIiIpg/f76tBUFEmD9/PhEREU55U1JSmD9/Pt27d+eRRx5h37597Ny5M9frmB0NGzYE4OjRox4pr2/fvsTExNg+//333yxfvtzl/QHLs121ahV9+/YlIiKCNWvWcOZM1rMJQ0JCKFu2rEfqbYxh9OjRLFmyhN27d19zeUop5UkauNgxxpQHxgCfishpu12HgOeB3kAv4DiwMKPg5fjx4w7jCVyN88hB7WxbbGwwHTokOQUtViEh0KFDEgsWBJHbs4169uzJqVOn2Lx5MwCbNm0iPj6enj2db8369es5deoUERERhIeH4+vrS3R0dK7VLTU1leTkZIctNTXzhrk//vgDsIzD8YR+/fqxZcsWW7lffPEFQUFB3HXXXS7zL1y4kCtXrhAREUG/fv1ITk4mNjY2y/NcvHiRs2fPeqzevXv3pnbt2rz88sseKU8pVXBYx8qlDVfI4DdVzmngksYYUwSYD/wDDLPfJyKfichUEVkvIkuArsB3wDhXZYWEhCAits0zgcu/Tp6EDMbE2tSqZcmX20qVKkWnTp1srQoxMTF06tTJaSwJWLpArPmDg4Pp0KEDMTExmY73uBY9evTA19fXYZswYYJTvuTkZK5evcquXbt48cUX6dSpE82aNXPKY7+BpQXJPs3VdYSGhnLLLbcwb948wHJ/+vTpQ6FCrv/pRUdHExoaSoMGDbjllluoX79+hsGd9fzHjh1jyJAhVKhQgYEDB2brHmWkUKFCvPTSSyxYsIBff3V3rLtSSlkCF+vvPyx/6HuUBi6AsYSFn5A2M0lEMp2aI5ansQhoYIzxuQ5VBMS2lS8fxcGDgZnmPnQokPLlo9KOyV0RERHExsZy+fJlYmNjXXaDXLlyhUWLFnHffffZBr1GRETw+++/s3Vr7kzOmjZtGtu2bXPYBg8e7JBnx44d+Pr6UqRIERo2bMiFCxecAoW4uDinAAigXbt2DmkbN250WY+IiAhiYmL466+/WLt2bYbdRCdOnGDDhg0O+yMiIti0aRN//vmnU/5SpUrh6+tLpUqVWLhwIbGxsVkOus2OBx54gCpVqvDKK694rEyllLpWGrhYvAX0AHqIyH43j7FGEtddeHg4q1encDyDOPbYMVi9OoXevXtfl/p0796df/75h9GjR3Px4kW6devmlGfFihWcO3eOLl26cO7cOc6dO0fr1q0pWrRornUX1apVi7CwMIctJF3/WmhoKNu2bWPLli28/vrr/PHHHwwZMsQhT0hIiFMABDBjxgyHtKZNm7qsR0REBD/++COvvPIKFStWzHCQ7Pz580lNTaVTp062e9S5c2dExNZiY++bb77hhx9+4LPPPiM4OJiIiAjbwGRPKFy4MM8//zyfffYZv/9+w45TV0oVMPlmVlFOGWNeAoYCfURks5vHGCxjXXaKSEpW+T0tKCiI4cNH0q3bVKcBusePQ/fuAQwfPpxSpUpdl/oEBgbStWtXpk2bRu/evQkMdG4NsgYnroKpBQsW8NZbb+Hjc50ar+wEBATYZkS1aNGCpKQkxo0bx/Dhw7ntttsAKFKkiMOsKau6deu6TE+vevXqNGvWjGnTpvHcc89lmM96j6znTb9vxIgRDmmNGzemWLFi3HrrrVSuXJm77rqL9957jxdeeCHLOrnrkUceYdKkSbz22mseK1Mppa5FvglcjDEBWBaRA6gIlDDGhKd9Xi4iicaYQ8BGEXk07Zj+wGRgNnDMGGP/p/BvIhKflm8jsBDYDwQCg4DbgHtz9aIyMWaMZaxG/fpT6NAhiVq1LN1Dq1enMHz4cNv+6+WJJ57g8uXLDuuHWF28eJGlS5fSr18/p66an376ieHDh7Nu3Trat29/vaqboREjRvDOO+/w2muvsWjRIo+W+/nnn/PQQw+53H/48GG+//57hg0bRvfu3R32rVixgtdff52DBw9muODfnXfeSZcuXXjrrbd49tlnKVq0qEfqXbRoUUaOHMlLL71E06ZNbd1kSimVV/JN4AKUAxakS7N+rg7EYble+z/rO6R9HZC22RuIJaABy6yiZ4EKWKZK/wjcIyIrrrXSOWWMYezYiQwdOpzY2GBOnoR27aYxc2bv69bSYq9169YZrmmyePFiEhMTeeaZZ5xaE1q1asXLL79MdHS0xwOXAwcOUKaM4zqEfn5+NGrUKMNjAgICGDZsGGPHjs00UMiuPn360KdPnwz3x8TEUKhQIUaOHOnUnXXTTTcxdepUoqOjGTfO5XhwAEaPHk2rVq347LPPePTRRz1Sb4AhQ4YwefJktmzZkuFsKKWUul7yzRgXEYkTEZPBFpeWp5qIDLA7ZkAmx8y2y/eoiNQQEX8RCRSRO/IyaLEXFBTEoEEwdiwMGjQoT4KWrERHR1O7dm2XXSC+vr706dOHRYsWOayc6wkjR46kRYsWDlt4eHiWxw0dOpQSJUowZcoUj9YnM9HR0bRt29YpaAEoV64c7du3z3IsUMuWLWnTpg1Tpkzx6EwtazCnlFI3ApNbU1ELsrCwMNm+fXuG+/ft20doaKgHz2hdq0WfpVK5zfP/fpXKv4wxO0Qk68GA2ZCfuooKMA1YlFJKFQz5pqvoRhIfH09YWBhRUVF5XRWllFLquouKirLOunTnhcfZoi0uuaBs2bJk1lWklFJK5WeDBw9m8ODBGGOyftlaNmmLi1JKKaW8hgYuSimllPIaGrgopZRSymto4KKUUkopr6GBSz5g+HclF6WUUio/08BFKaWUUl5DAxeVI5GRkRhjMnyXT+3atTHGEBkZ6bSvevXqGGM4dOiQLa1atWoYYzLd3Mkze/Zs4uLiMtzfrl072zkHDBhgSy9UqBCVKlWiX79+xMXFeew+Xbx4kcDAQAICAvj777/dOmbZsmUYYxzq8fvvv/Pggw9SpUoV/Pz8qFy5Mj169OCbb75xOn7+/PmUL1/eYdn/hQsXcvfdd1OqVCmKFi1KnTp1GD58OMePH7flyeiebd68Ocv7nr6+SimVW3Qdl1xgXYDOOo89NyUkJEBsLJw8yYflyxMeHk5QUFCuntPKz8+PI0eOsH37dutCQwBs27aNuLg4/Pz8nI7ZunWr7RdcdHQ0Y8eOBeCLL75weFdRp06dCA8P57HHHrOl+fj4kJKSkmmemjVrcvHiRQCmTJlCq1atHM5fsmRJh8/16tXj448/JjU1lf379zN69Gi6dOnCzz//TJEiRbJ7S5wsWbKExMREAL788ksefPDBbJeRkJBA8+bNqVChAq+88gohISHExcWxZMkStm7dyp133umQ/6uvvqJLly4YY+lAHDFiBG+99RYDBw5k2LBhlChRgr179zJjxgyOHDnCF198YTt2xIgRTu9zCg0NZevWrbbPhw8f5v777+f999+nSZMmtvQKFSpk+9qUUvlTVFSUdRFWjy9Ah4jo5uGtadOmkpm9e/dmut8dqampMmbCBPErVUro21cYPVoC+/QRv1KlZMyECZKamnrN58jM+PHjpXTp0tK+fXsZMWKEw77hw4dLhw4dpHTp0jJ+/HiHff/3f/8ngYGBctttt0loaGiG5bs61t08R44cEUCWLl2a6fEPP/ywpH9W0dHRAsi3336b6bGArF+/PtM8IiLdunWTGjVqSPXq1aVz585Z5hcRWbp0qQBy5MgRERGJiooSY4ycOnXKKW/655ySkiJly5aV2NhYERFZsmSJADJr1iynY5OTk2X58uUO1/Tuu+9mWb/du3e7ff35kSf+/SpVUADbxcO/Y7WryEuNmzSJqV9+SdKePRATA5MmcXHePJL27GHql18ybtKk61KPiIgI5s+fb+uWEBHmz59PRESEU96UlBTmz59P9+7deeSRR9i3bx87d+68LvV0V8OGDQE4evToNZeVkJDAqlWr6Nu3LxEREaxZs4YzZxwXkRQRIiMjKVeuHMWLF+ehhx7iwoULDnnOnTtHkSJFCA4OdjqHtVXFatu2bZw7d4727dsDMG3aNJo0acIjjzzidKyPjw+dO3e+1stUSqnrSgMXL5SQkMCUqVNJXLoUQkIcd4aEkLh0KVOmTuXcuXO5XpeePXty6tQpNm/eDMCmTZuIj4+nZ8+eTnnXr1/PqVOniIiIIDw8HF9fX6Kjo3OtbqmpqSQnJztsqampmR7zxx9/AJZxONdq4cKFXLlyhYiICPr160dycjKxsbEOed555x0mTJjA4MGDiY2Nxd/fn+eff94hT5MmTbh8+TIPPvggO3bsyPQavvrqK+644w5KlCjB1atX2bJlC506dXK7zunvmX3XnFJK3Qg0cPESxm4Ljo0lqUMH56DFKiSEpA4dCFqwINenSZcqVYpOnToRExMDQExMDJ06dXIaSwKWMS3W/MHBwXTo0IGYmBiHQaSe1KNHD3x9fR22CRMmOOVLTk7m6tWr7Nq1ixdffJFOnTrRrFkzpzz2G1hakOzT0l9HdHQ0oaGhNGjQgFtuuYX69es7BGopKSm89tprDBkyhEmTJtGxY0dmzpxJaGioQzlt27Zl2LBhzJs3j7CwMEqVKkWvXr1Yu3at07V89dVX3HPPPQD89ddfXL58mSpVqrh9z5555hmH+3XXXXe5faxSSl0PGrh4o5MnIYPZPDa1alnyXQcRERHExsZy+fJlYmNjXXYTXblyhUWLFnHffffZBr1GRETw+++/Owz89KRp06axbds2hy39YOkdO3bg6+tLkSJFaNiwIRcuXHBqBYqLi3MKgADatWvnkLZx40bbMSdOnGDDhg0O9yIiIoJNmzbx559/ApbuqBMnTtCjRw+H87lqrZo6dSq//vorb7zxBq1bt2blypV06NCBGTNmOJzzp59+sgUuVum7kzLz3HPPOdyvWbNmuX2sUkpdDxq4eAmx26LKlyfw4MFM8wceOkRU+fLkTluGo+7du/PPP/8wevRoLl68SLdu3ZzyrFixgnPnztGlSxfOnTvHuXPnaN26NUWLFs217qJatWoRFhbmsIWka6UKDQ1l27ZtbNmyhddff50//viDIUOGOOQJCQlxCoAAZsyY4ZDWtGlT2zHz588nNTWVTp062a63c+fOiAjz5s0D4GRaYFmuXDmH86X/bH89I0eOZMmSJfz+++80atSIUaNG2Vp6li9fTo0aNahbty4ApUuXpmjRorbuL3dUqVLF4X5Zy1JKqRuFBi5eKDw8nJTVq8FuDQ4Hx46Rsno1vXv3vi71CQwMpGvXrkybNo1u3boRGBjolMcanPTu3ZugoCCCgoKoXLkyly9fZsGCBXk2liIgIICwsDBatGjBc889x7hx45g/fz7ff/+9LU+RIkWcAiCAunXrOqQVL17cdoz1em+77Tbb9VqPs+4rX748AKdPn3aoU/rPrpQpU4aBAweSkJBgy2/fTQTg6+tLq1atWLVqVbbvi1JK3ag0cPFCQUFBjBw+nIBu3ZyDl+PHCejenZHDh1OqVKnrVqcnnniCbt268fjjjzvtu3jxIkuXLqVfv36sX7/eYZs6dSqnTp1i3bp1162umRkxYgRlypThtddey3EZhw8f5vvvv2fYsGFO1/v888+zY8cODh48SOXKlSlfvjyLFy92OH7RokUOn+Pj412e5+DBgxQtWpSSJUty5coV1qxZ49RN9Oyzz7J9+3bmzJnjdHxqaiorV67M8XUqpVRe0AXocsH1WIBuwpgxAEypX98yULdWLQIPHSJl9WqGDx9u23+9tG7dmtatW7vct3jxYhITE3nmmWe47bbbHPa1atWKl19+mejoaNsUXk85cOAAZco4rn3k5+dHo0aNMjwmICCAYcOGMXbsWA4ePJjhysCZiYmJoVChQowcOdKpa+qmm25i6tSpREdHM27cOJ5//nlGjhxJmTJluOOOO1i4cCH79u1zOGbOnDnMnTuXhx56iIYNG3L16lXWrl3L9OnTeeKJJ/Dz82Pt2rWIiNNg2m7dujF8+HAeffRRvv32W3r06EGxYsXYv38/M2bMoFq1atmadaSUUu7QBei8bLseC9BZnT17VoiKEiZMkKioKElISPBY2ZmxLkCXGesCcV27dpXatWtnmO+JJ56QkiVLSlJSktOx7pSfnnUBOldbzZo1bflcLUAnInL+/HkpVaqUDB48OMNzk8kCbDfffLO0b98+w2M7d+4s9erVE5G0hQTHjJEyZcpIsWLFpH///jJ37lyHBej27NkjTz75pISGhkqxYsWkZMmS0qRJE5k+fbpcvXpVRESeffZZ6dGjR4bnjI2NldatW0uJEiXE19dXateuLSNGjJATJ044XJMuQJc1XYBOKfeRCwvQGcmlqagFWVhYmGzfvj3D/fv27XOa8notrHNG9EkWXHXq1OG5555j0KBBeV2VfM/T/36Vys+MMTtEJCzrnO7TrqJ8QAMW9euvv+Z1FZRS6rrQwblKKaWU8hoauCillFLKa2jgopRSSimvoYFLHtFB0Up5H/13q1Te08AlD/j6+nLp0qW8roZSKpsuXbpke1eVUipvaOCSC6wL0KUtvuOkXLlyHDt2jMTERP0LTikvICIkJiZy7NixDN8lpZT6V1RUlPU1Jx5fgE7XcckFWa3jAnDhwgVOnz7N1atXr1OtlFLXwtfXl3LlylGiRIm8ropSXkPXcclHSpQoof8BKqWUUtmkXUVKKaWU8hr5JnAxxtQyxsw0xuwyxqQYYza4eVxJY8zHxpgEY8x5Y8xcY0xpF/l6GGN2G2OSjDF7jTF9PX4RSimllMpUvglcgPpAF+AAkJ31z+cDrYHHgAHArcCX9hmMMbcDC4H1QGfgKyDaGNPhGuuslFJKqWzIT2NclorIYgBjTCxujGQ2xrQAOgB3icg3aWnHgO+NMe1EZG1a1rHANyLydNrn9caY+sA4YLWHr0MppZRSGcg3LS4ikpqDwzoDp6xBS1o5PwBH0vZhjCkKtMHSMmMvBmhhjCmZsxorpZRSKrvyU4tLTtQD9rtI35e2D6Am4Osi3z4sgV8dYJv9jgMHoHVrx8x9+sCTT0JiInTp4nzCAQMs25kzEB7uvP+JJ6BvXzh6FB580Hn/iBHQrZvl3EOGOO8fMwbatYOff4Znn3XeP3kytGwJW7bAqFHO+996Cxo1grVrYdIk5/0zZ0LdurB0Kbz5pvP+Tz+FypVh3jz44APn/bGxUKYMzJ5t2dJbvhwCAmD6dJifPoQENmywfJ0yBZYtc9zn7w8rVli+nzgRvv7acX/p0rBwoeX7l16CrVsd91eqBJ99Zvn+2Wct99BenTpgXbJn8GBI/6LmRo0s9w/ggQfgzz8d97doAa+8Yvm+Vy/46y/H/W3bwtixlu87d4b0axd27QojR1q+T/9zB/qzpz97lu/1Z895v/7sWb7P7Z89TyvogUsQcM5FegJQwy4PLvIlpNtvc+XKFTZu3GD7XLVqNaBaTuuolFJKeY24uDg2bvySt98eBhDi6fLz5QJ01jEuItI6i3xrgIsicm+69M+AGiLS0hjTCtgMNBaRn+3y1AIOAh1FxGGcizsL0CmllFL5XW4sQJdvxrjkUALgaoxKEP+2qFi/ps8XlG6/UkoppXJZQQ9c9vPvWBZ79mNffgOuushXD0gle1OvlVJKKXUNCnrgsgIon7ZOCwDGmDAs41tWAIjIZSzrt/ROd2xfYKuInL9OdVVKKaUKvHwzONcYE4BlATqAikAJY4x1nPpyEUk0xhwCNorIowAistUYsxr4xBgzEksLymvAZrs1XAAmAhuMMW9hWZyuS9rWKZcvSymllFJ28k3gApQDFqRLs36uDsRhuV6fdHn6AtOAj7C0QC0DnrbPICKb04KgScATWNZ56Z9+UK5SSimlcle+nFWU13RWkVJKKaWzipRSyusZYzDG5HU1lPJaGrjkgvj4eMLCwoiyLiuolFJKFSBRUVGEhYWBG+8NzC7tKsoF2lWklMqItbVF/+9VBYF2FSmllFKqQNPARSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cFFKKaWU19DARSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cMkFunKuUkqpgkxXzvUyunKuUiojunKuKkh05VyllFJKFWgauCillFLKa2jgopRSSimvoYGLUkoppbyGBi5KKaWU8hoauCillFLKa2jgopRSSimvoYFLLtAF6JRSShVkugCdl9EF6JRSVtYF57Ki/xer/EgXoFNKKaVUgVY4ryuglFIFQ0YtKu61yCilLLTFRSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cFFKKaWU19DARSmllFJeQwOXXKAL0CmlnJkMNqXyH12AzsvoAnRKKStdgE4VZLmxAJ2u46KUUrkofUBiDWQ0UFEqZ7SrSCmllFJeI18FLsaYm4wxXxtjEo0xx40xE4wxPlkcE2mMkQy2l+zyzc4gT73cvzKlVL4hYtmUUjmSb7qKjDFBwFpgL9ADqAm8iSU4G5PJof8DVqZLuxd4AViRLn0/MDBdWlyOKqyUUkqpbMs3gQvwOOAP9BSRC8AaY0wJINIY83pamhMR+RP40z7NGDMW2C8iP6fLflFEvvN81ZVSSinljvzUVdQZWJUuQInBEszc5W4hxpjSQHsg2rPVU0oppdS1yk+BSz0sXTk2IvIHkJi2z129AF9cBy43GWMuGGMuG2M2G2PcDoiUUkopde3yU+ASBJxzkZ6Qts9dEcCPInIwXfpPwAigG3A/4IOlO6pZ+gKOHz+OMca2RUZGZuP0SimllPeKjIy0/f4DQjxdfr5ZgM4YcxV4TkTeSpf+J/CJiIxyo4wKWMa7vCAiU7LIGwDsAXaKyL32+3QBOqVURqzL0eWP/3mVylxuLECXn1pcEoCSLtKD0va5ow+W/1fmZZVRRBKB5UATdyuolFJKqWuTnwKX/aQby2KMqQwEkG7sSyYigM0ictTN/IL+4aSUUkpdN/kpcFkBdDTGFLdL6wtcAjZmdbAxphrQHDdnExlj/IF7gB3ZrqlSSimlciQ/BS4zgMvAImNMO2PMYCASmGo/RdoYc8gYM8vF8RFAMrAg/Q5jTEljzCZjzBBjTFtjTF9gPZZBR5Nz4VqUUkop5UK+WYBORBKMMW2B94ClWGYYTcMSvNgrjGVGUHoRwNcicsbFvstAPJYVeMsBScBW4C4R0VG4Siml1HWSbwIXABHZC9ydRZ5qGaQ3yuSYJKDntdRNKaWUUtfuunYVGev73JVSSimlcsDtwMUY4zTA1VVaBsfWN8a8DvyRjbp5rfj4eMLCwoiKisrrqiillFLXXVRUFGFhYQBlPF222wvQGWN+FJEmdp+LAdtEJDSD/OWwrDD7EBAKTABiReTXa671DU4XoFNKuZKQkEBwbCycPElU+fKEh4cTFJSdhb2V8i55sgCdMWa8MSYVaGSMSU3bUrCsMDvfRf4IY8xy4EegCjAIOCEikwtC0KKUUumJCBMnjqVGjRDu+XoYL10ez9q1w6hRI4SJE8eSX1YwV+p6yHJwroj8F/ivMWaqiAx3o8zPgQ3A7SISB2CM0X+VSqkCa9KkcXz55VT27EkixPbmloscPw7duk0FYOzYiXlWP6W8SXa6ikKBWUAFEalujGkM3Csi49Plawg8gGV68UHgUyBSRKp6tOY3MO0qUkpZJSQkUKNGCHv2JOHvD2k9RZQvD+HhcOkS1K/vx5EjJyhVqlReV1cpj8qNrqLsBC4bgFHA+yLSOG2G0C8iUj+D/IWAdljGuHQHVgPzRMRpgbf8JizMiMYtSimADz+ENWug7i2+TJnqg0/HtiTWakDAwZ2krF7HyOEpHNh9lfbtYdCgvK6tUp5lDB4PXLKzjkugiGyxzmgWETHGXMkos4ikYglWVhtjAoFw4HFcrEyrlFL51cmTcDzel6++rEHSnnVY+4ouAhw/ztRud9Oo+GFOnryap/VUyltkJ3C5ZIwpSdpLBY0xjYBEN4+tBtQH6mSnct6rKaBNLkop8PGZzNYfJpJ68N+gxSYkhMSl6/iudk26dvwv8FKe1FGp3OP55duyE7iMAlYBVY0xS7D8du6dUeYMpkO3zXlVvccOcuNRKaW80s6d0OZu56DFKiSE1DZ3M2rnTkZd35op5ZXcXoBORDYDnYEHgSjgFhHZkj6fTodWSik7Z85Ag4aZ57n5FoiPvz71UcrLZWvJfxFJEJHlIrJMRM5mkO1zwA/LdOhhaS8hLFDToavEx9M0LIyZUVEI6KabbgV4e7BiRXx27SIzPr/8woMVK+Z5XXXTzVPbzKgomubVyrnWFXPTFqFzlfks8LKIvJWWvyE6HVqnQyulADhy5Ag1broJfvvNdXfRsWNQqxZx+/dTtWqB+W9SFRB5snKudZl/ESkkIj7pN+AmYKhd/p0i8hxQFXgVy7iWIGNMrDEmwzExSimVH1WvXp27WrTAtG0Lx4877jx+HNO2LXe1aKFBi1Juuua3Q4tIPNDFRXqqiKwWkQeACsAyYPC1nk8ppbzNurVrubNCBahZE5+uXeHFFy1fa9bkzpAQ1q1dm9dVVMprXHPgApDVoFsRuSgis0WkvSfOp5RS3qRQoUJsWLeOw3v3khIcDNu20T84mLj9+9mwbh2FCnnkv2KlCoTsTIcGwBgTDLwHtAdSgCXACyKSYJenGdAJ6CQiLT1UV6WU8mrVq1eHTz4B4JM8rotS3irbgQswDegDfAWcBloCy4wx9wFj0vaVxTKQ9zsP1VMppZRSKkeBS0fgRRGZApD2zqIoYBtQGViJZUr0ChH5y1MVVUoppZTKScdqOWCz9YNY5lNPxBK0jBORLiLymQYtSimllPK0nI4IS0n3+Vja1zXXUJd8Iz4+nrCwMKKiovK6KkoppdR1FxUVRVheLUDndIBlIbqtWF7Jsy9t+w34HWgqIj95upLeRhegU0plxPoes+z9z6uUd8qNBehyMsblbeBmLINwy+H47+89Y8wW4CfgZ2C/iKReayWVUkoppSAHgYuIDLN+n/YG6FvStgZpX58E/NOyJALFrr2aSimllFI5a3GxEZHTwNdpG2CbZVQbSyBz8zXVTimllFLKzjUFLq6kzTL6NW2L9XT5SimllCq4dJ1ppZRSSnkNDVyUUkop5TU0cFFKKaWU19DAJRfoAnRKKaUKshtqATqVNV2ATimVEV2AThUkubEAnba4KKWUUspraOCilFJKKa+RrwIXY8xNxpivjTGJxpjjxpgJxhifLI6pZowRF1uMi7w9jDG7jTFJxpi9xpi+uXc1Sqn8SNBuIqWuhccXoMsrxpggYC2wF+gB1ATexBKcjXGjiJHAt3afz6Qr/3ZgITAdeBroAkQbYxJEZPU1X4BSSimlspRvAhfgcSzvSOopIheANcaYEkCkMeb1tLTMHBCR7zLZPxb4RkSeTvu83hhTHxgHaOCilFJKXQf5qauoM7AqXYASgyWYuetaCjbGFAXaAPPT7YoBWhhjSl5L+UoppZRyT34KXOoB++0TROQPLG+orufG8R8bY1KMMSeMMVONMf52+2oCvunLB/ZhuYd1cl5tpZRSSrkrPwUuQcA5F+kJafsychl4H3gUaAvMBJ7A0ppiXzYuyk9Itx+A48ePY4yxbZGRkW5UXymllPJ+kZGRtt9/QIiny89PY1xyREROAEPtkjYYY04B040xDUVkZ3bLDAkJ4fjx4x6ro1JKKeUtIiMjbX+wG2M8/sswP7W4JACuxpoE8W/LiLti0742tSsbF+UHpduvlFJKqVyUnwKX/aQby2KMqQwE4Dw2JSuS7utvwNX05ad9TgV+zWb5SimllMqB/BS4rAA6GmOK26X1BS4BG7NZVnja1x0AInIZWA/0TpevL7BVRM5nv7pKKaWUyq78NMZlBpaF4RYZY14DagCRwFT7KdLGmEPARhF5NO1zJFAcy+JzF4A7geeARSKyy678iVjGv7wFfIllAbouQKfcvCillFJK/SvftLiISAKWWUE+wFLgv8A0YHy6rIXT8ljtx7LOy8fAcqA/8EbaV/vyN2NpiWkHrAK6A/111VyllFLq+jEi+tYMTwsLC5Pt27fndTWUUkqpPGWM2SEiYZ4sM9+0uCillFIq/9PAJRfEx8cTFhZGVFRUXldFKaWUuu6ioqIICwsDKOPpsrWrKBdoV5FSSimlXUVKKaWUKuA0cFFKKaWU19DARSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cFFKKaWU19DARSmllFJeQwOXXKAL0CmllCrIdAE6L6ML0CmllFK6AJ1SSimlCjgNXJRSSinlNTRwUUoppZTX0MBFKaWUUl5DAxellFJKeQ0NXJRSSinlNTRwUUoppZTX0MBFKaWUUl5DA5dcoCvnKqWUKsh05VwvoyvnKqWUUrpyrlJKKaUKOA1clFJKKeU1NHBRSimllNfQwEUppZRSXkMDF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwCUX6AJ0SimlCjJdgM7L6AJ0SimllC5Ap5RSSqkCTgMXpZRSSnmNfBW4GGNuMsZ8bYxJNMYcN8ZMMMb4ZHHMrcaYj40xh9KOO2CMGW+M8UuXL9IYIy62Trl7VUoppZSyKpzXFfAUY0wQsBbYC/QAagJvYgnOxmRyaN+0vK8BB4EGwMS0r73S5T0PpA9U9l1r3ZVSSinlnnwTuACPA/5ATxG5AKwxxpQAIo0xr6elufKqiJyx+7zBGJMEzDTGVBWR3+32JYvId7lTfaWUUkplJT91FXUGVqULUGKwBDN3ZXRQuqDF6qe0ryGeq55SSimlrlV+ClzqAfvtE0TkDyAxbV92tABSgd/SpZcyxpwxxlw1xvxkjOmZ49oqpZRSKtvyU+ASBJxzkZ6Qts8txpjyWMbEfCoip+12HQKeB3pjGftyHFjoKng5fvw4xhjbFhkZ6fZFKKWUUt4sMjLS9vuPXOi5yDcL0BljrgLPichb6dL/BD4RkVFulFEEywDfSkBTEUnIJK8BtgD+ItLIfp8uQKeUUkrpAnRZSQBKukgPStuXqbRA5BOgPtAls6AFQCwR3yKgQVZTrpVSSinlGflpVtF+0o1lMcZUBgJIN/YlA29hmUbdXkTcyQ8gaZtSSimlroP81OKyAuhojClul9YXuARszOxAY8xLwFDgARHZ7M7J0lpoegE7RSQlZ1VWSimlVHbkpxaXGcDTwCJjzGtADSASmGo/RdoYcwjYKCKPpn3uD0wGZgPHjDHN7cr8TUTi0/JtBBZiab0JBAYBtwH35upVKaWUUsom3wQuIpJgjGkLvAcsxTLDaBqW4MVeYcB+TEqHtK8D0jZ7A7EENGCZVfQsUAHLVOkfgXtEZMW1114ppZRS7sg3s4puJDqrSCmllNJZRUoppZQq4DRwyQXx8fGEhYURFRWV11VRSimlrruoqCjCwsIAyni6bO0qygXaVaSUUkppV5FSSimlCjgNXJRSSinlNTRwUUoppZTX0MBFKaWUUl5DAxellFJKeQ0NXJRSSinlNTRwUUoppZTX0MAlF+gCdEoppQoyXYDOy+gCdEoppZQuQKeUUkqpAk4DF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwEUppZRSXkMDF6WUUkp5DQ1clFJKKeU1NHDJBboAnVJKqYJMF6DzMroAnVJKKaUL0CmllFKqgNPARSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cFFKKaWU19DARSmllFJeQwMXpZRSSnkNDVyUUkop5TU0cMkFunKuUkqpgkxXzvUyunKuUkoppSvnKqWUUqqA08BFKaWUUl4jXwUuxpibjDFfG2MSjTHHjTETjDE+bhxX0hjzsTEmwRhz3hgz1xhT2kW+HsaY3caYJGPMXmNM39y5EpXX4uLiMMZc900ppVTmCud1BTzFGBMErAX2Aj2AmsCbWIKzMVkcPh+oAzwGpAKvAV8Cd9iVfzuwEJgOPA10AaKNMQkistqT16Ly3vz587ly5Qq+vr55XRWllFJ28lOLy+OAP9BTRNaIyAzgv8BwY0yJjA4yxrQAOgAPi8hCEfkCeAC43RjTzi7rWOAbEXlaRNaLyHPASmBcTiobGRnp8Jd2SEgIvXr14rfffrPlGTBggHVUdq6XYyUiNGrUiDlz5tjS9u7dS9u2bQkICCAkJIRx48aRkpKSZVmHDh1iyJAhNGjQAB8fH1q3bp1hXnfOsX//fm677TZKlixJREQE//zzj8P+b775hooVKzqlDx06lEcffdSNq/+Xtwcts2fPdtmiM2PGjCyPdedZLFiwgO7du1OxYkWKFStG06ZNiY6Otu3/4YcfsmxdKlSoEH///bdb19OmTRu6deuW4f67776b2rVrk5yc7FZ5SinvlW9aXIDOwCoRuWCXFoOl9eQuYGkmx50SkW+sCSLygzHmSNq+tcaYokAbLC0t9mKAj40xJUXkfHYrXLJkSVauXAnA4cOHGTt2LG3btmXPnj0EBgZe93LA0tJw9uxZ+vfvD0BCQgLt2rXjpptuYvHixfz222+MGDGC1NRUJk2alGlZe/bsYfny5TRv3pyrV69mmM/dcwwYMIBatWoxYcIEnnvuOSZPnszkyZMBSE1N5ZlnnuGVV16hWLFiDuWPHDmSevXq8dJLL1GrVq0s78GuXbto0KBBlvm8wbp16/D397d9rlGjRqb53X0WU6dOpXr16kybNo0yZcqwfPly+vfvz5kzZ/i///s/atWqxdatW23533nnHZYtW8bq1f82Tvr6+lK8eHG3riM0NJQ1a9a43Ld+/XrWr1/P559/TuHC+em/NKWUSyKSLzbgNBDpIv0i8Fwmx80HNrhI/wr4Ku37mwABWqfLc2ta+q326U2bNpWsjB8/XkqXLu2QtmnTJgFk/vz5IiLy8MMPS1Zleaocq5YtW8qoUaNsnydPniylSpWS8+fP29Jee+018ff3d0hzJSUlxfZ9r1695K677nKZz51z/P333wLI6dOnRUQkJiZGwsLCbPk//PBDufXWWyU1NdXlOdq2bSvDhw/PtL5Wb7zxhly6dMkhbcWKFdKqVSspWbKkFC9eXJo2bSqrV692q7y88PHHHwsgf//9d7aOc/d5x8fHOx3br18/qVatmstyu3TpIi1atMhWXey988474uPjI5cvX3ba16pVK2nYsGGGz14plXeA7eLh3/f5qasoCDjnIj0hbd+1HGf9mj5fQrr9ABw/ftyhSTwyMjKT0/+radOmgGVgqL0vv/ySevXq4efnx+23387evXtzpZxDhw6xZcsWwsPDbWkrVqygY8eOlCjxb29bREQEly5dYuPGjZnWo1Ah93683DnHlStXAGytBwEBAba0CxcuMGbMGN5+++0MB7j26tWLuXPnkpqammV9kpKS8PPzs33esGED9957L+3bt2fRokV8/vnntGvXLl/+de/u8y5TxnlNqcaNG3P8+HGX5e7evZubb77Z5b6kpCTGjx9PrVq18Pf3JywsjM2bNzvkCQ0NJSUlhUOHDjmkr1y5km+//ZaXX35ZBzcrdYOwH8IAhHi6/PwUuNwwQkJCHKJDdwMXa6BRvnx5W9rvv//O8OHDGTt2LJ9//jnnz5+nY8eOJCUlebycr7/+msDAQBo2bGhL279/P/Xq1XMov0qVKgQEBLB//363risr7pwjODiYatWq8e6773L27Fn7VRmZOHEi7dq1o0WLFhmeo2XLlpw6dYrdu3dnWpeDBw9Sp04dh7S5c+dyzz33MH78eO6++266du3Kq6++Sps2bXJyuW4REZKTk7PcslKzZk0KFy5M3bp1mTlzZpb5r+V5b9261eneAZw/f56jR49yyy23OO1LTk7mnnvu4cMPP+S5555jyZIl1KxZky5dunDq1ClbvtDQUAAOHDjgcPzYsWNp1aoV99xzT5bXppS6PiIjI+17Jlz/NXMN8tOfjAlASRfpQfzbMpLRcWWzOM76NX35Qen2Z5v1l8/hw4d58sknKV68OO3a/Tsm+MyZMyxevJiWLVsCltaUmjVrMnv2bB5//HGPlrNjxw5CQ0MdWkoSEhIoVaqUU72DgoJISMjxZTtw9xzTp0+nd+/ejBo1itq1a/P+++9z6NAh/ve//2UZkNSvXx8fHx9++OEHh8Dsn3/+cRgT89VXX/HYY485HOvv78+qVat46623iIiIcAgIc8ucOXMYOHBglvkkg5WvK1SowMSJE2nWrBkpKSnExMTw+OOPk5iYyLBhwzIsL6fP++uvv+bLL7/ko48+ctr3yy+/ALhscXnzzTf54Ycf2LFjhy3oad26NVWqVCE2NpannnoKgIoVK1KiRAmH4Gnx4sVs376db775xqlcpVT+lZ9aXPYDDn8qGmMqAwFp+9w+Lk09u+N+A666yFcPy/TpX3NQX/766y98fX3x9fWlbt26HD58mHnz5lGhQgVbnnLlytmCDYCqVavStGlTfvjhB4+Xc/LkSZddAFlJ3zrgzoyjnOjcuTOnT5/mwIED7Nu3jypVqjB8+HCGDRtGpUqVeP/996lSpQpVqlRh+vTpDscWLlyYUqVKcfLkSVvaW2+9RZcuXRzypQ9kwNKiM2DAACZMmEDFihXp2LEje/bsse03xjBgwADb5w0bNjjMoMpqvyvdunVj27ZtWW4Z6dixI2PGjKFDhw507tyZOXPm0KdPHyZNmuRWd1l2xMXF0b9/f3r06OFwnVbWwCV9i0tqaipvvvkmDz30EDVq1LD9/BhjqFWrFn/88YdD/nr16tlaXESEcePG0blzZ+644w6UUgVHfmpxWQE8Z4wpLiLWOZZ9gUtAZoMxVgBjjTG3i8hmAGNMGFAjbR8ictkYsx7oDdi3t/cFtkoOZhSBZTbQ2rVrMcZQvnx5QkJCnPrpy5Ur53RcuXLlOHHihMfLSUpKIiAgwCFPUFAQ5887X15CQgJBQZYGp40bNzp0m9x1111s2LAhkyt35M45rAICAmx/ma9Zs4adO3cyb948du7cydixY9myZQsALVq04Pbbb3eYHVS0aFGHrrFnn32WTz/9lOPHjxMSEsIff/xBlSpVnOpRsmRJ3nvvPd5++22+/vprHn/8cR599FG+++47W55NmzZx8OBBateu7fIas9qfXnBwMCVLumpAzLnw8HDmz59PXFxchrOLsvMsAM6ePUvnzp2pWrUqc+fOdVnm7t27+c9//uMUFO/evZv4+HimT5/uFGgCdOrUyeFzaGiorcVlwYIF7N69m9mzZ7s8p1Iq/8pPgcsMLNOVFxljXsMSeEQCU8VuirQx5hCwUUQeBRCRrcaY1cAnxpiR/LsA3WYRWWtX/kRggzHmLSyL03VJ2xz/d82GwoULZ7m+yunTp12m1a9f3+PlBAcHO7RIgOWv3PRjG44ePUpiYqJtLETTpk0d/vp3d4prds6RXkpKCsOGDeP111/H39+fDRs2cPfdd9vyt23blo0bNzoELufOnSM4ONihnPDwcBYtWsTQoUNZtmwZ/fr1y7CePj4+dOjQgY4dO5L+JZrPP/88EydO5JNPPnF5bFb707vWriJXrMFsZoNYs/MsEhMT6dq1K1euXGHZsmVOQa/VL7/84nJ8izVoXrt2rcsgrWrVqg6fQ0NDWbx4MSkpKURGRtKnTx8aN26c4bUopfKnfBO4iEiCMaYt8B6WNVvOAdOwBC/2CgPpXwPQNy3vR1i6z5aRbs0WEdlsjAkHJgFPAEeA/pLLq+aePn2aLVu22Lp5/vjjD3788Ue3fqllt5y6des6rL0Blu6ZN954g7///tsWkMybNw9/f3/uuusuwBKouLvAnSvunCO9Dz74gKCgIPr2/fetC4mJibbvL1686PBLPT4+nsTERKfBo+Hh4QwePJihQ4dy9uxZp1aFIUOG4OvrS+vWrSldujSbNm1i1qxZREVFOeR75JFHeOONN/j1V9e9hlntT8/aVeRJsbGxlClTxikgsOfus0hOTqZ3794cPHiQLVu2uGzRs/rll1948MEHndKtXZmlSpWyzYTLTGhoKOfOnWPq1KkcPHiQxYsXZ3mMUiof8vT8at1yvo5Leg8//LCUKVNGqlevLnPnzpVFixbJLbfcIiEhIbZ1RjxVjojIqlWrHNZKERE5e/aslC9fXtq1aydr1qyRmTNnSmBgoIwePTrLa7x48aIsWLBAFixYIM2bN5ebbrrJ9vnixYs5Psdff/0lZcqUkR07dtjSfvzxRylSpIjMmjVL/ve//0mRIkXk559/tu1fvny5GGPkzJkzTuU1adJEdu3aJVFRUU773nzzTbn11ltt67c0b95cvvzyS4c8ln9GlrVk7r//flm/fr3DmjVZ7c8NPXv2lFdffVWWL18uS5culQceeEAAeeedd2x55syZIz4+PhIXF2dLc/dZDBo0SAB5++23ZevWrQ5bUlKSLd/x48cFkFmzZjnV8erVq9KoUSOpU6eOzJkzRzZu3CixsbHy/PPPy+zZs53y//rrrwKIn5+fDBo0yBO3SSmVy8iFdVzy/Jd8ftw8Gbg0bdpUFi5cKLVr15YiRYpIy5YtZffu3R4vR0Tk8uXLEhwcLJ988olD+p49e6RNmzbi5+cn5cuXlzFjxkhycnKW13jkyBHBskCf03bkyJEcn+P//u//5JFHHnFKf+utt6R8+fJSvnx5h1/QIiJPP/20tG7d2mV5EyZMkGbNmsnJkyezvCZXrIHJlStXpFatWvLBBx+4DFwy2p8bXnrpJalTp474+/uLn5+fNGnSxOm5Whepy8mzqFq1qlvP1hoMf//99y7refToUenfv79UqFBB/Pz8pHr16nL//ffLgQMHnPImJydLkSJFxM/PT44ePZqzG6OUuq40cPGSrUqVKtK0aVOZOXNmhg/zRvX0009Lly5d8roaHpWcnCwVK1aUTz/91OX+vXv3yh133JHj8q2BiYjIrFmzpEqVKi4Dl4z2K6VUfjNz5kxp2rSpAHHi4d+xRsT9wX3KPWFhYZJ+8Ka3+PPPP6lTpw4///yzy8XEvFFMTAxjx45l3759Ga52u2fPHoeBytlhjMH67yg5OZl69epRqVIl28yqrPYrpVR+ZYzZISI5HwTpQr4ZnKs8o1KlSnz00UecOHEi3wQuIsKsWbMyXaI/p0GLtXyrwoULOy1Ln9V+pZRS7tMWl1zgzS0uSimllKfkRotLflo5VymllFL5nAYuSimllPIaGrgopZRSymto4JIP/Prrr0RGRnLlyhWH9I8++oiSJUva3hztSTExMRQtWpSrV696rMxhw4ZhjLFt//nPf3jooYf466+/PHYOq0OHDjFkyBAaNGiAj49Pli89TG/v3r20bduWgIAAQkJCGDduXLZfLikiNGrUiDlz5tjS4uLiMMYQFxeXrbJyKjfPt2DBArp3707FihUpVqwYTZs2JTo62inf0KFDefTRRz1+fqVU/qSBSz6wYsUKpk+fTpEiRRzS27Vrx7fffpvpbJqc2rlzJzfddBO+vr4eK3P37t00b96crVu38u233zJixAjmzp3LU0895bFzWO3Zs4fly5dTt27dbM+eSkhIoF27dhhjWLx4MePGjePNN99k/Pjx2Spn/vz5nD17lv79+2frOG8xdepUihUrxrRp01iyZAlt2rShf//+vPvuuw75Ro4cydy5c3W2lVLKPZ5eGEa3678A3SOPPCJt2rS5Luey6ty5swwcONCjZZYrV06eeuoph7Ru3bpJcHCwR88jIpKSkmL7vlevXtlaEG7y5MlSqlQpOX/+vC3ttddeE39/f4e0rLRs2VJGjRrlkGZdbTj9ara5xRPnO3TokMv0+Ph4p7R+/fpJtWrVnNLbtm0rw4cPz3EdlFI3ltxcgE5bXHJB2bJl2b59O4MHD84039y5c7nzzjspXbo0xYoVo0WLFg4vOUxNTaVYsWLMmDGDMWPGUKlSJUqWLMmgQYNITU0FoFy5cnz00UesX7/e1sWybt0627Effvhhrlzjzz//TKNGjWyf161bR+nSpXn22Wez3WUCcOrUKU6fPk1oaKhDerly5fDxSf9OzGtXqFDOf/RXrFhBx44dKVGihC0tIiKCS5cusXHjRrfKOHToEFu2bCE8PDzLvI0aNXL5s/TAAw9wxx13uF9xD0pKSmLu3Lncfffd1K5d22WeMmXKOKU1btyY48ePO6X36tWLuXPn2n6ulVLebfDgwaQtC3LG02Vr4JKH9uzZw8MPP8yCBQtYsGABFSpU4N5777WNVTl8+DAXL15kypQpXLp0iY8//piRI0fyv//9j8WLFyMiLFmyBD8/P0aPHs3WrVvZunUrLVu2tB3boEEDj9c7Pj6eEydO2AKXd999l65du/Lqq6/y1ltv5SjQ2L17NwD16tWzpaWmprJ161a6devm8hgRITk5OcvN0/bv3+9QT4AqVaoQEBDA/v373Srj66+/JjAwkIYNG2aZ97bbbiP9ukA//PAD0dHRTJ061f2Ke8BPP/3E0KFDqVChAo888gilS5fmq6++cvv4rVu3uuyaa9myJadOnbL9HCilVEZ05dw8NHnyZNv3ycnJ1KhRg3r16vHbb78RGhpq+0/8ySefZPjw4QC0b9+eDz74gIMHD2KMoXz58iQlJdG1a1eaN29uK2/37t0UKlSIm2++2eP13rlzJ8YYbrrpJgYNGsSiRYtYvnx5tge42rNea61atUhOTub48eNMnDiREiVK8Prrr7s8Zs6cOQwcODDLssXDiywmJCRQqlQpp/SgoCASEhLcKmPHjh2Ehoa61fJz22238fHHH5OUlISfnx9gGcjcr18/br311mzVPSfOnz/P3LlzmTVrFj/++CONGjXiv//9Lw888ADBwcFul/P111/z5Zdf8tFHHzntq1+/Pj4+Pvzwww9uBXNKqYJLA5c8cunSJd59910+/fRTDh8+TGJiom1fQEAAALt27aJEiRIMHTrUtk9EOHfunK0ZfteuXRhjuOWWWxzK37VrFzVq1CAwMNDp3OfPn+fEiRNZ1jF9q4LVzz//THBwMD179iQ+Pp7vv/+eWrVq2fZHRkYSGRmZZfnp6wtQrVo1W1qlSpXYtm0bpUuXdnlMt27d2LZtW7bOc6M4efKky64UV5o3b87Vq1f5+eefad68OfPmzePHH390mKHTqFEjfv7550zLyclzWblyJT179sTf35/777+fWbNmOXQRuisuLo7+/fvTo0cPBgwY4LS/cOHClCpVipMnT2a7bKVUwaKBSx4QEbp06cLBgwcZPnw4t9xyC0FBQXz00UfMmTOHypUrA5ZWiDvuuMNhttBvv/3GpUuXbO/WyShA2b17d4bdRAsWLGDQoEFu1dOVnTt3IiJs2rSJtWvXOgQtAP/973+z/QvSOqPo3Xff5erVq6xZs4bx48czadIk3nvvPZfHBAcHU7JkyWydxxOCgoI4f/68U3pCQgJBQUFulZGUlGQLULNSr149SpQowbZt22jUqBEvvPACw4cPp0qVKrY8WQUtkLPnUrRoUQICArh06RLnz5/n3LlzlrezGuN2GWfPnqVz585UrVqVuXPnZnqupKSkbNVPKVXw6BiXPPDNN9+wYcMGFi5cyPDhw2nfvj1hYWHs27ePW265xdZ9sHv3bqdm8507dzp0Ae3atctl0/ru3budWmGsHnvsMbdGbmdk586dPP300zzwwAMMHDiQ+Ph4275hw4YBlhYAdweOpqamsnfvXlq2bElYWBgtWrRg3LhxdOrUifnz52c4YHPOnDn4+vpmuXlavXr1nMayHD16lMTExAxbqdILDg7m3LlzbuUtVKgQzZo1Y9u2bUybNo2kpCRefPFFhzz2gcSmTZu44447aNq0KXfeeSd79+7N0XMBaNOmDceOHWPWrFkcO3aMu+++m5o1azJhwgR+//33LI9PTEyka9euXLlyhWXLlmUarJ07dy5bXU9KqYJJA5c88OeffwI4DFJcuXIlGzZssAUhly5d4tChQ07N8ulbWFwNFLUemxsDcy9fvsz+/ftp0KABUVFRlClTht69e9sGwU6bNg2wtABs2rTJrTIPHjzIpUuXaNy4sUN6//79iY+P54cffnB5nLWrKKvN0zp37syqVav4+++/bWnz5s3D39+fu+66y60y6taty5EjR9w+52233cb69et55ZVXmDBhAsWLF3eZ7+zZs4wcOZIlS5awY8cOXnnlFR599NEcPRerokWLEhERwdq1a/ntt9+4//77+fDDD6levTrt2rXjs88+c3lccnIyvXv35uDBg6xcuZJy5cpleI74+HgSExPzzRvJlVK5yNPzq3UTmjZtKpnZu3evFCpUSB588EFZu3atvPLKK1KvXj0xxsj7778vIiLbtm0TQA4cOOBw7L333is9e/a0fW7VqpXcfvvtsnHjRtm6daukpqbajj148GCm9ciJH3/80aHsuLg4KVOmjAwdOtSWx/Jj5Wj9+vUCyPr16532LViwQADZvXu3Q3p8fLwUKlRIxo8f79FrEBG5ePGiLFiwQBYsWCDNmzeXm266yfb54sWLtnxz5swRHx8fiYuLs6WdPXtWypcvL+3atZM1a9bIzJkzJTAwUEaPHu3W9YqIrFq1SgA5ffq0Q3pG66osWbJEALn55pslOTnZqTzrPV+2bJkEBQVJw4YNbVvlypUd8rhzvqwkJyfLsmXL5N5775XChQu7zDNo0CAB5O2335atW7c6bElJSQ55ly9fLsYYOXPmTLbqoZS6sQHbxcO/Y/P8l3x+3NxZgO6DDz6Q8uXLS4kSJaRPnz6yadMmAWTz5s0iIvLRRx9JQECAw0JpIiI1atSQyMhI2+dNmzbJTTfdJIULF5b//Oc/tmMDAwMlNTU1w/Pn1Mcff+xU9tdffy0+Pj7y8ccfi4jrX5BfffWVALJnzx6nfePGjRM/Pz+5evWq074WLVrIbbfd5rkLSGP9he1qs/8l/vHHH7v8xb5nzx5p06aN+Pn5Sfny5WXMmDEOAUVm1ysicvnyZQkODpZPPvnEZb3Sn2/r1q0CyMqVK12WZ73nS5YskXvuuSfTPO6cLztOnjzpMr1q1apu3WMRkaefflpat26d4zoopW4subkAXZ7/ks+PW1YtLvld8eLFnVaQHTduXIH6xeTO9T799NPSpUsXh7SMAolu3bpJ586dMyzLGpScOXNGypcvLz/++KOIWFYI3rFjh4i4fi7Xe6VeV5KTk6VixYry6aef5lkdlFK5IzdaXHRWkfK4Z555hmbNmlG2bFnbeIotW7bY1qIpCNy53ueee446derw66+/uhzbkZSUxK5du4iNjeXrr79m586dWZ63dOnSLFiwgKeeeoqLFy9y5coVwsPDadKkicvnciNYsGAB/v7+RERE5HVVlFJewFgCIuVJYWFhkn6lU6VciYmJoUKFCrZBvXFxcVSvXp0jR45w4MABOnfuTPXq1Xn77bfp2rWryzLi4+OpV69ejt6ibX8++zV0rqfo6GgqVqzInXfemSfnV0rlHmPMDhEJ82SZOqvoBjBkyBD8/f2Ji4vLl+fLLUePHqVt27aEhoZSv359nn/+ea5nIP7ee+9Rr1492rRpA5CjhdkiIiIynInUsWNHUlNT+e233zIMWg4cOEDz5s0ZNWpUts99o+jXr58GLUopt2mLSy7ITovL9u3bCQ8P5+GHH2bnzp18+eWXuVo3d8+XnJxM4cI3dk/iiRMnOHbsGGFhYVy5coX27dvz9NNP06tXr+ty/tDQUL744gu3125RSqmCRltc8hkRYejQoURFRTF+/HjOnDnDqlWr8ux8xhjGjx9PkyZNePnllx26DjZs2ODwLiJjDJMnT6ZJkyaEhoby/fffZ3rub775hmbNmhEYGIgxhvj4eBo1asScOXOc8l69epVp06bRtGlTihUrRokSJejYsaPTOSpUqEBYmOXfQ5EiRWjQoAFHjx617R86dCiPPvpodm6R2x566CF+++03evbsyXPPPQf8uwjc7t27qVWrFmfPngUs683MmDEjV+qxd+9e2rZtS0BAACEhIYwbNy7Lt3PHxsbSsmVLSpcujZ+fH3Xr1mXSpEm2l3tazZ492/bGcfstJ9fi7jO9ESxYsIDu3btTsWJFihUrRtOmTR1er2AvJSWFDz74gObNm1O8eHH8/PyoU6cOQ4cO5ddff83Wedu0aZPhC0UB25u4c+PFoUp5FU+P9tXNe2cVAbYp3EeOHJGqVava9q1fv17uuusuh7zWqbyff/65tG/fPsNyL1++LGXKlJH7779fvvnmG/n+++8lJiZGKleuLFeuXHHIm5CQILfeequULl1aJkyYIKtXr5bo6Ghp1qyZFClSRDZu3OjyHGfOnJFKlSrJ3r17bWlHjhyRokWL5sp6NiKW6b72s3Gwm248bdo0CQ8Pl/fee0+6d++eK+c/e/asVKhQQdq2bSurV6+WDz74QAICAhzWk3FlxowZMnr0aFm0aJGsW7dOXn31VfHz85OnnnrKIZ91Kvi6desc1mA5depUtuqZ02eaV5o3by79+vWTefPmyddffy0jRowQQN555x2HfP/884/ceeed4u/vLyNGjJCvvvpKvv76a3n99delSpUqEhYWlq3zPvHEE1KrVi2X+9atWyeAfP755zm+LqXyArqOi3ds3hy4/PPPPyLiXuBiXUTs0KFDUrNmzQzLtS62Zr/eR8uWLWXUqFEO+VJTU6V169ZSuXJlOXz4sMO+S5cuSeXKlcXVvU1KSpLWrVvLlClTnPa1bdtWhg8fnvFFX4PMApfU1FRp27athISESHx8fK6cf/LkyVKqVCmHKc6vvfaa+Pv7O017zsqoUaOkZMmSDuvzWAOXv//+O8d1zOkzzUuunle/fv2kWrVqDmldu3aVoKAg2blzp1P+f/75R+bMmZOt877zzjvi4+Mjly9fdtrXqlUradiwYa6szaRUbsqNwEW7inJBfHw8YWFhREVFZZrv/PnzDBkyhODgYEqVKsXw4cN55ZVXnF5a6Anunsv6KoHChQs7vCPI1cvvihYtCoCPj0+Gzde33347HTt2BKB8+fIYY/jmm2/YsmUL4eHhDnlnzZrFhg0bmDVrFtWrV3fY5+fnxwMPPMCOHTsc3o2UkpLC/fffT+PGjRkxYoTT+Xv16sXcuXMzfN+RK2fPnmXIkCGUKVOG4OBgXnjhBdu1vPvuu26V8ffff3P06FEKFy7s9juJsmvFihV07NiREiVK2NIiIiK4dOkSGzduzFZZpUuXduoq8oScPNO85uqt3Y0bN+b48eO2zwsXLmTZsmVERUW5fLVGYGAgDz30kENaUlIS48ePp1atWvj7+xMWFsbmzZtt+0NDQ0lJSeHQoUMOx61cuZJvv/2Wl19+OVsvt1QqL0VFRVm78p3/QV0rT0dCurnX4nLlyhVp3ry53HLLLTJv3jxZtmyZNGrUSCpUqOCwpL9VamqqXL16NcvtWs6FXYvBlStXpGTJknL06FFJTU2VgQMHOrW4WKVvnbH3008/Sfv27eX222+XrVu3ynfffSczZsyQwMBAp1WB69atK02aNMnwnk2fPl0Ah79wH330URkwYECGf4n+/PPPAsjPP/+cYbn2EhISpG7dutK0aVNZvHixvPzyywLI9OnTpXz58nLp0iVb3sxaXB544AGZNGmSLFmyRJo1a+by2VzLMxURKVu2rMvXIQQEBMjrr7+e5bUmJyfLxYsXZdOmTVK1alWnlilri0u5cuXEx8dH6tSpIzNmzMiyXHs5eabX6lrvqyv33Xef3HzzzbbPLVq0cPiclatXr8rdd98tFSpUkBkzZsjq1aulT58+Urx4cVtL5J9//imALFq0yOHYsLAwadWqVbbqq9SNAu0q8o7NncBl8uTJUqZMGfnrr79sad9++60ADkv6W1l/iWS1Xcu50h8/c+ZMqV69utx2223y3HPP5ShwERFp2LChPP/887bPgwYNcur/P3DggADy5ptvZljOhAkTHN6TtHnzZtv7e6zv5Xn77bcdjrl69ar4+PhIVFRUhuXaGzlypJQsWdKhu6BMmTJSsmRJp66ojAKX6Ohouf32222vAHj88cedusVEru2ZiogULlxYpk2b5pResWJFeemll7K81qJFi9rO8dBDDzkFkitXrpSJEyfKqlWrZPny5fLQQw8JIFOnTs2ybJGcPVNPuNb7mt7atWvFGGN7pcXJkyfFGOMUNKYPmOyD6VdffVWKFSvm8O6xK1euSPny5eW9996zpZUoUUImT55s+/zll18KIN98800274JSN4bcCFxu7Pmu+VRqaipvvfUWTz31FMHBwbb0GjVqANjeEG3P+ibk3DyX5WfsX4MHD2bw4MEuy7XPW61atQzXhElOTmbfvn22mTcAJ0+edGqO/+WXXwC46aabMryW/fv34+fnZ5vt1KpVK6c6p1e4cGFKlSrFyZMnM80HlmuaM2cOgwYNcqhfyZIlOX/+PI8//rhD/vTXbK1LRESEwyqwH3zwgcvz5fSZesqWLVtITEzkhx9+YMKECQwdOpTp06fb9nfs2NHWzQeWt2InJSUxadIknnnmGQoVyrynOSfP1BM8eV/j4uLo378/PXr0YMCAAQD8+uuviIjTuj0vvvgir7/+uu3zvHnz6NOnD6mpqbz55ps89NBD1KhRw9ataoyhVq1a/PHHH7Zj6tWrx4EDBwDLz9O4cePo3Lkzd9xxh0euR6n8QAOXPPDLL79w+vRp2rVr55D+559/Aq4Dl+DgYEqWLHldzuVJ+/fv58qVKw7jAJKSkggICHDI9/fffwOuxxcAXLx4kSVLltCpU6dsry9TtGhRl2N00jtw4ADx8fG0b9/eIT0lJYVnnnnGNv7HU3L6TK2CgoI4f/68U3pCQgJBQUFZHt+kSRPAMnanTJkyPPzww4wYMYKaNWtmeEx4eDjz588nLi7OFvxmJDefaWau9b5anT17ls6dO1O1alXmzp3rkA44jC0CePzxx+nduzcbNmzgueeeo3HjxoBlenx8fDzTp093CAytOnXqZPs+NDSU/fv3A5Zp2bt372b27NnXfC1K5Sf5anCuMWaQMeagMSbJGLPDGNPWjWOGGGPWGGNOGWPOG2O+NcZ0cJEvzhgj6bas/4x34dixYwBUqlTJIX3t2rWUKFHC5V+fc+bMwdfXN8vNE+fypF27duHr6+uwSFtwcLDTgNUKFSoAOKzDYm/69On8888/LgfgZuXcuXMOrU0ZOXLkCABVq1a1pW3ZsoW4uLgcrYqblZw+U6t69erZfslZHT16lMTExGwvimcNYqz3ICPWwaHuDBLNzWeamWu9rwCJiYl07dqVK1eusGzZModAu2zZsoDzdVWvXp2wsDASExMpXry4beD7iRMnAMu/uW3btjlt9q2aoaGhHDhwgJSUFCIjI+nTp48tAFJKWeSbFhdjTD9gBhAJbAYGAsuMMbeKyC+ZHDoaWAm8D1wEHgBWGmPuFZEl6fJ+DthPK8nRNIzSpUsDOLwf5syZM0ydOpVbbrnF5S+FnDZ/Z+dcxhiHrpe4uDhat27Njh07aNvWEgOePXuWf/75hypVqgDw2muv0alTpwy7bHbu3Em9evUcflHUrVuXrVu3OuRr2bIlJUqU4KOPPuK+++5z2Ldx40bGjRvHoEGDuP3227N1/fHx8SQmJrp8iWF6Pj4+tmsES1P9888/b/ve0661S6Nz58688cYb/P333xQvXhywdE/4+/tn+BqBjHz77bcATjN/0ouNjaVMmTIOwV1GsvtMExMTuf/++zl8+DApKSn06NGDsmXLcubMGSZNmsT//vc/hg0bRkJCAj4+PlSuXJlff/3VqfXuWu9rcnIyvXv35uDBg2zZsoVy5co57G/atCnBwcF89NFHPPzww07H79ixg0aNGtn+bVkDuFKlStG0adNMzx0aGsq5c+eYOnUqBw8eZPHixTm+DqXyLU8PmsmrDTgAfGT3uRCwG/gsi+PKuEjbAqxPlxYHTHGnLlkNzk1KSpKQkBBp0KCBLFu2TGJiYqRx48YSHBwsTz75ZKbHZld2zkW6AYuuBt1+/PHH8vDDD2d6nL1OnTrJ/fff75BmXdfl9OnTDukffPCBANKnTx9ZsmSJLFu2TJ5++mkpWrSo9OnTx2mxOncsX75cjDFy5swZW9r69esFkPXr1zvkPXPmjPj5+cntt98uK1askAEDBkhoaKjUrl1b+vfv7zAQ90Zw9uxZKV++vLRr107WrFkjM2fOlMDAQIcF6ObMmSM+Pj4SFxdnS+vYsaO88cYbsnz5clm1apWMGzdOAgMDpW/fvg7l9+zZU1599VVZvny5LF26VB544AGnhdgyupdW2XmmCxculMGDBztc308//WSbUXP//fdLs2bNZOvWrbJ79+5cm2kzaNAgAeTtt992WHhv69attrWL/ve//wkgbdu2lejoaNmwYYPExMTIwIEDxc/PT0aMGGEr7+rVq9KoUSOpU6eOzJkzRzZu3CixsbHy/PPPy+zZsx3O/euvvwogfn5+MmjQoFy5PqWuJ3RWUYbBRw0sMwU6pUsfB/yVg/LeBw6kS/NY4CIismXLFmnQoIH4+/vLbbfdJhs3bpSSJUtme7qpO9w9V24ELiEhIfLaa685pF2+fFmCg4NtK+/aW7hwoTRv3lwCAgIEkCJFisjcuXPduErXnn76aWndurVD2ldffSWA7Nmzxyn/vHnzpFq1auLr6ytNmjSRw4cPy+zZsyUwMNB2HRcvXpR7771XGjRoIPXr13c5Y+h62bNnj7Rp00b8/PykfPnyMmbMGNtsJpF/Z9jYB11jxoyR+vXrS2BgoJQsWVIaN24s77zzjlNg+NJLL0mdOnXE399f/Pz8pEmTJk7PLLN7aeXuMz106JDUqlVLRo4cKV999ZUkJydLSkqKlC9fXi5evCgNGjSQzz77TCZPnizvvvtulisE51TVqlUznIlkfx8XL14srVq1koCAAClatKhUr15d+vXrJ1988YXTDK2jR49K//79pUKFCuLn5yfVq1eX+++/32GWkYhlinqRIkXEz89Pjh49mivXp9T1pIFLxoFGl7T/WKqlS++dll42m+X9CCxNlxYHJGDpHjoPxAJVXR2fk5VADx8+LIBs3bo128d66lyAbVpxw4YNJTQ09JoDl4w8/fTT0qVLl0zzrFy5UrB7tUB2JScnS8WKFeXTTz91SB83bpxTMJMdrloGCqrs3susnum5c+dkwYIF0rdvX+nQoYOIiPTo0UM++OAD6dGjhxw/flzat28vvXr1kjVr1njkGpRSuUcDl4wDjfvTApRS6dLbpaXXyUZZj6Qd0yZd+ttAP+AOYDBwDPgDKJm+jAoVKjj8leZqkbD0Fi1aJIUKFbItuZ+bMjpXbrS4ZOTo0aPi7+/v9Bdnev3795fSpUs7dSu5Izo6WmrVquW02Fi7du1kyZIl2S7PylXLQEGVk3uZ0TM9evSoXLx4UUREjh07JqVKlRIRy3ufqlSpYluzpnHjxhISEmLLq5S6sYwfP97+d+Bx8fDv/Bt2cK4xpiRQIat8IrI/qzzZOGdTLINv3xaR9enO84zdx03GmC3Az1gGAb9lnzckJMRheXB37Ny5kxo1anh8ym1enysjlSpV4qOPPuLEiROZDpy1n4aaXSLCrFmznKbarlmzJsdlAtSsWZPt27ezZs0aPvnkE95+++1cfav3jSwn9zKjZ7p7925efPFF2yDx9957D4DWrVszbNgw29vJmzVrxi+//OI0KFcpdWOIjIwkMjISAGNM9n4ZusGIeH62hCcYYx4DPswqn4gYY0wX4CssXUW/25XRG5gPlBORTF+GYoypgWVQ7g/AfSKS4kYd9wA7RMThpSRhYWGyffv2rA6/4WQ0q8h+obXZs2ezYcMGh7Ul0h+X3/35558EBwcTEBDA8ePHqV+/PgkJCXldLaWUuuEYY3aISJgny7xhW1xE5H/A/9zMbm11qQf8bpdeDzjrRtBSDliVdmyEO0GLtZppW4Fnv05MmzZt+PTTT/OwNrkro5YBpZRSue+GbXHJLmPMAWCTiDyW9rkQlq6cXSLyQCbHFQM2ACWBllkFOXbH3ZxW/nARecd+n7e2uCillFKeVKBaXHIgEvjMGBMHfAs8DNQG+lszGGPuAr4G2orIxrTkRUADYABQ0xhjW+9cRL5LO+4eLAvTLQOOY2nJGYNlcO7s3LskpZRSStnLN4GLiESntZ68AIwF9gBdxXHVXAP4pH21sr6YxtWIQWu+o0A5LINwSwF/YVltd5SIXPDQJSillFL5TYinC8xX7yoSkQ9FpJaIFBWRJiLydbr9G0TEiMgGuzST0WaXZ5eItBWRsiLiKyLlRWSAiLgcLR0f71Zvk/JCUVFReV0FlY4+k9yj9zb79J45yXJ2cHblq8DlRnHmzJm8roLKJfqf0o1Hn0nu0XubfXrPcl++GZx7IzHG/I3l3Ukq/ykDaGR6Y9Fnknv03maf3jNHoSLi0UXDNHBRSimllNfQriKllFJKeQ0NXJRSSinlNTRwUUoppZTX0MDFSxhjKhtjvjbG7DPG7DHGvG6MMVkfqW50+mxvLPo8co/e2+zL7/fMGLPRGLPTGLPLGBNrjCmR1TEauHiPZOAFEQkFGgO3AT3ztkrKQ/TZ3lj0eeQevbfZl9/vWXcRaSgiDbCsRv9cVgcUuMDFGDPAGCMutsdz4Vy1jDEz0yLJFGPMhkzy3pQWVScaY44bYyYYY3ys+0XkhIhsT/v+CrALqOzpOns7Y0xhY8yLxpiDxpjLxpg/jTHTcuE8+myzYIzZkMG/NTHGtPDwuQrc8zDGRBhjfjTG/GOMOWaM+cQY4/FVSvPTvTXG3Jt2HZeNMUeMMcNz6Tz55p654u71ZXVtACJyPi1vISAQN15cnG+W/M+Bu4FLdp8P58I56gNdgO8A34wyGWOCgLXAXqAHUBN4E0tgOcZF/tLAvUAHj9fY+83G8mz/i+Wt4ZWBm3LhPPpss/YkkL7ZdwKWvxq3efhcBep5GGO6A9HA+1j+Qq0ATAK+MsY0FZFUD54uX9xbY0wrLO+m+wgYiaXl4jVjTKqIvOXh0+WLe5aJLK8vO9dmjFkO3IrlVT0jsjy7iBSoDcvLFAUolo1jGgOVMtjXBfDJYF8hu+9jgQ0Z5HsJSABK2KU9DyTap6WlFwXWAyPy+l7eaBvQCbgK3KTP9sbbgCLAWeADfR7XfC9jgB3p0rqn/d8WqvfWZf1WAZvSpb2Z9jNZRO9Ztu5llteXnWtL2+cDvA48n9X5C1xXUQ5FAmuMMWXsE40xDwFLsfyF70Tc/6unM7BKHF/YGAP4A3fZnc8Hy8sgfxKRN92ufcHxCLBORPZm45hI9NleL52AICwtBRmJRJ+HO3yB8+nSzqV9zWjgZiQF+942AtakS1uN5Wcyo67LSAr2PXPJzetz69rsykwB5gAPZVVwQQ5cfjPGJBtjDhhjhmSRdwCWv+RXmrQRz8aYHliaHEeJSPp/DNlVD0u3ho2I/IElMq1nlzwT+Bt3mtIKptuAX40x7xljLqT1qy7Kot9/APpsr5cI4E9gUyZ5BqDPwx0fAXcYYx4yxpQwxtTB0lWUWeA+gIJ9b/2AK+nSrJ9DMzhmAAX7nl2LLK/NGBNkjPmPXZZewC9ZFVwQA5cTwFjgQaAblj66GcaYYRkdICIJQEcskfkSY0wXYB7wpoi85oE6BfHvX0v2EtL2WftnHwXCgJ+MMT8bY572wLnzk/JY/qNphOWX5ECgKfCFMa6nD+qzvT6MMQFYujLmS1q7sCv6PNwjIl9h+VmPwtLycgBLU3uvTI4p6Pf2EJZxFPaapX0NdnWA3rNrkuW1pX1dljbIdzdwM/BMVgUXuMG5IrIKS1+n1QpjjB8wxhjzdkZNYCJywhjTDvgW+Ar4n4i8kPs1tp3/WzJuAlYWJm3rISJ/ARhjTgAbsTTpfu3qIH2210U3LDMGMusmAvR5uMMY0waYAbwNrAD+g6Vb4wtjTLu0ZncnBfzezsDyR+ogLOMymgHWWUUZdn0U8HuWq0TkMM7BZJYKYouLK7FYIu5qWeTzxzLAMBUoljZ9yxMSgJIu0oPS9in3JAC7rUFLms1YmoOzmlmkzzZ3RQCHJG1apxv0eWTuTWCJiLwgIhtEZB6W2SatsczgyExBvbcfAR+kbWexzDCamLbvZBbHFtR7di1y7do0cLGQdF+dGGOqYRnI9R2WgUVdsUxF9IT9OPZnYoypDASQro9QZWofrv8yMWTyF5U+29xljCmJZaBelq0tafmroc8jK/WAn+0TROQAliUeamZ0UEG+tyKSIiJDgbJAAyytVN+l7f4uo+MK8j27Rrl2bRq4WIQDZ4DfXe1MGzy0BvgN6C0im7H8VTPQGDPZA+dfAXQ0xhS3S+uL5T+hjR4ov6BYBtySbgbAnVhmYOx0dYA+2+viPixTObMMXPR5uO13oIl9gjEmFEvLQJyrA/TeWohIgojsFpF/sKw1tEVEXP4i1Xt2TXLv2jw9v/tG34CFwAtY/gLsCnyKpaXl/zI5ZhmwHed59T2wjDjvmMFxAViConBgK5bFdayfA+zyBWEZNLwGaAcMBv4BJuX1/fKmDctiZ3+k3etuQH/gKLBGn22ePpeVwM9u5tXn4d59egZLK+KbaddxP5YBukeAQL23Lq+jOZaF59phWTJ/AXABaKA/j9m+l1leX25eW57fgDy44ZPT/oEnYon8dgAPZnFMHaBMBvtaAiaDfdWwBEWutmrp8t4ErEur0wksfa8uFzfSLdNnVQtYDlzE0o86GwjSZ5tnz6NM2n/wL7qZX5+He/fJAE9gWf79InAMy2yXGnpvM7z+plhWbP4HS8DyFXCL/jzm6F66dX25dW0mrXCllFJKqRuejnFRSimllNfQwEUppZRSXkMDF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwEUppZRSXkMDF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwEUplSeMMT+7kScy92uilPImuuS/UuqGZYwRETF5XQ+l1I1DW1yUUnnCGCN2399hjNlkjNlhjPnGGHOTMWZa2r6fjTGb8q6mSqkbiba4KKXyhLU1xRgTDKwAOolIgjGmFTBFRFpoi4tSKj1tcVFK5bUWQG1gfdq4l/eBitkpwBizwRgjaduz6fbda4xZbYz5yxhzxRhzzBgTa4zplM1zLDXG7M5k/3vGmHPGmKLGmEi7+sRm5zxKqcxp4KKUymuFgC0i0shuq5KDctZjCYJirAlp3U0LgWPAY0A74EXAH1hhjKmZjfKjgZuNMTel32GM8QHCgUUichn4X1pdfsrBdSilMqGBi1Iqr20BmhpjGgMYYwoZY5qk7fvbGFPCzXLOish3InIyrZwewLPAoyIyUES+EJFvRORTEbkH6A5cykY9FwOJQD8X+9oA/8ES3CAif4rId8CFbJSvlHKDBi5KKY8zxpw1xkw0xsxN66I5a4yZaYwpmj6viPwF9AbeN8bsBPYA96Xtfhv4IYeDc58FtonIbFc7RWSpiBxPV+87jDEbjTGJafX+0BhTPC3/RWAp0NdFcRHAaWBdDuqplMoGDVyUUh5ljKkCBAHDgLNYgpK3gMHAU2l5yqbtA0BENotISxFpKCKhIjI2LX2siNQTkTuyWYfCWLpqVmfjmFbAWuAklm6fZ4EuwMd22aKB2saYpnbH+QI9gfkikpKdeiqlsq9wXldAKZXvNEz7+p6IvJj2/TpjzJ1AJ2PMV8ByYHIu1qE0UBQ4ap9ojDGAj11Sivw7tfJVLGNt+trlPwZ8bYy5WUR+wTL76RyWFpYdadk6YgnUonPhOpRS6WiLi1LK0xpiGQvyRrr0X4FgETkgIjVF5M3rUJf06z2MAK7abdYWoAAsLTTzjTGFrRuwOS1fUwARuQIsAvqkBUFg6Tr6Hdiay9eilEIDF6WU5zUEvkkbu2KvPHDcRf7c8BdwGaiULv1T4Na0zV4QlpaY6TgGNpcBX6CyXd5ooArQwhjjB/QAYuxabpRSuUi7ipRSntYQcBhMmzbA9W5gmt07ioKBYsAfaZ9fEJFVnqiAiCQbY7YCHYBxdumngFNpdbI/5ByW1plILN1Y6dkHXOvTyogAKgDF0W4ipa4bDVyUUh5jjAkEauLcsvJ02te3ReS/aXkHAK1FZEAuVect4EtjzIMi8mlmGUXkojHmO6CuiEzIIm+KMWY+lkHHFYF9IrLTU5VWSmVOAxellCfdgqULuqYxZhyWlpcOWMaWPCAi565XRURksTHmLWC2MaYNlqnMZ7AM3O2Qlu0fu0OexzIQNxWIBf7G0iV0DzBaRH61yxsN/B+Wadvjc/M6lFKONHBRSnlSQyyLrt2DZRrxaOAA0FtEFl/vyojIMGPMN8CTwCws3TrxWAbSdhGRFXZ5N6fNfPovlrEwPlgG3a4krXvJLu9WY0wcUA3tJlLqutLARSnlSQ2B3SKyi7SZONeRSZsJZD/FGRH5AvjCnQJE5HvArXcYiUj1TCpSCEvLk74gUikP01lFSilPagjk1XiPnlhmAj2TR+e3Nw5LXe7M64oold9oi4tSyiPS1jW5BfgkD04/BEs3EPw7SykvRQHL0r4/m1lGpVT2aOCilPKItO4Zd1+I6OlzH8iL82Yk7R1I12vNGqUKFKNrJimllFLKW+gYF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwEUppZRSXkMDF6WUUkp5DQ1clFJKKeU1NHBRSimllNfQwEUppZRSXuP/AcYRwNV84UFjAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Prepare arrays of Jet-pT\n",
    "jet_pt_pp_in = jets_pp[:,1]\n",
    "jet_pt_pbpb_in = jets_pbpb[:,1]\n",
    "jet_pt_pbpb1_in = jets_pbpb1[:,1]\n",
    "jet_pt_pbpb2_in = jets_pbpb2[:,1]\n",
    "\n",
    "bins_13 = [500,1000]\n",
    "y_axis13= [0.815]\n",
    "yerror_13 =[0.142]\n",
    "\n",
    "pt_min = 500\n",
    "pt_max = 1000\n",
    "\n",
    "n_ev_pp = 389140\n",
    "n_ev_pbpb = 35000\n",
    "n_ev_pbpb1 = 12790\n",
    "n_ev_pbpb2 = 6461\n",
    "\n",
    "sigma = 5.15\n",
    "d_eta = 2*2\n",
    "\n",
    "# Fill Histogram\n",
    "n_pp, pt = np.histogram(jet_pt_pp_in, bins= bins_13 )\n",
    "n_pbpb, pt = np.histogram(jet_pt_pbpb_in, bins= bins_13 )\n",
    "n_pbpb1, pt = np.histogram(jet_pt_pbpb1_in, bins= bins_13 )\n",
    "n_pbpb2, pt = np.histogram(jet_pt_pbpb2_in, bins= bins_13 )\n",
    "\n",
    "# Statistical Errors\n",
    "err_n_pp = np.sqrt(n_pp)\n",
    "err_n_pbpb = np.sqrt(n_pbpb)\n",
    "err_n_pbpb1 = np.sqrt(n_pbpb1)\n",
    "err_n_pbpb2 = np.sqrt(n_pbpb2)\n",
    "\n",
    "# bin width\n",
    "dpt = (pt[1:]-pt[:-1])\n",
    "# bin center\n",
    "pt = pt[0:-1] + 0.5*dpt\n",
    "\n",
    "# Jet Spectrum\n",
    "dn_dpt_pp = (n_pp*sigma)/n_ev_pp/dpt/d_eta\n",
    "dn_dpt_pbpb = (n_pbpb*sigma)/n_ev_pbpb/dpt/d_eta\n",
    "dn_dpt_pbpb1 = (n_pbpb1*sigma)/n_ev_pbpb1/dpt/d_eta\n",
    "dn_dpt_pbpb2 = (n_pbpb2*sigma)/n_ev_pbpb2/dpt/d_eta\n",
    "\n",
    "# Errors \n",
    "err_dn_dpt_pp = (err_n_pp*sigma)/n_ev_pp/dpt/d_eta\n",
    "err_dn_dpt_pbpb = (err_n_pbpb*sigma)/n_ev_pbpb/dpt/d_eta\n",
    "err_dn_dpt_pbpb1 = (err_n_pbpb1*sigma)/n_ev_pbpb1/dpt/d_eta\n",
    "err_dn_dpt_pbpb2 = (err_n_pbpb2*sigma)/n_ev_pbpb2/dpt/d_eta\n",
    "\n",
    "# Generate Plots\n",
    "fig = plt.figure()\n",
    "\n",
    "# Calculate RAA and error\n",
    "raa =  dn_dpt_pbpb/dn_dpt_pp\n",
    "error_raa = ratio_error(dn_dpt_pbpb,err_dn_dpt_pbpb,dn_dpt_pp,err_dn_dpt_pp)\n",
    "\n",
    "raa1 =  dn_dpt_pbpb1/dn_dpt_pp\n",
    "error_raa1 = ratio_error(dn_dpt_pbpb1,err_dn_dpt_pbpb1,dn_dpt_pp,err_dn_dpt_pp)\n",
    "\n",
    "raa2 =  dn_dpt_pbpb2/dn_dpt_pp\n",
    "error_raa2 = ratio_error(dn_dpt_pbpb2,err_dn_dpt_pbpb2,dn_dpt_pp,err_dn_dpt_pp)\n",
    "\n",
    "\n",
    "plt.errorbar(pt, y_axis13,ms=8,mfc='blue',mec='black',mew=2,fmt= 's',label=\"CMS [JHEP 05(2021)284]\",\n",
    "            xerr=None, yerr=yerror_13, color='black')\n",
    "\n",
    "plt.errorbar(pt, raa, ms=9,mfc='red',mec='black',fmt='o', label=r\"MATTER+LBT\",\n",
    "             xerr=0.5*dpt, yerr=error_raa, color='red')\n",
    "\n",
    "plt.errorbar(pt, raa1, ms=9,mfc='yellow',mec='black',fmt='o', label=r\"MATTER+MARTINI\",\n",
    "             xerr=0.5*dpt, yerr=error_raa1, color='yellow')\n",
    "\n",
    "plt.errorbar(pt, raa2, ms=9,mfc='cyan',mec='black',fmt='o', label=r\"MATTER+AdS/CFT\",\n",
    "             xerr=0.5*dpt, yerr=error_raa2, color='cyan')\n",
    "\n",
    "\n",
    "\n",
    "#axes setting\n",
    "plt.legend(loc='upper left')\n",
    "plt.xscale('log')\n",
    "plt.xlabel(r\"$p^{\\mathrm{jet}}_{\\mathrm{T}}$ [GeV]\")\n",
    "plt.ylabel(r\"$R^{\\mathrm{jet}}_{\\mathrm{AA}}$\")\n",
    "plt.xlim(pt_min,pt_max)\n",
    "plt.ylim(-0.3,1.6)\n",
    "plt.axhline(1, color = \"blue\", linestyle=\"dashed\", linewidth=1.5) \n",
    "plt.text(pt_min+5,-0.05, r\"PbPb(0-10%)\"r\"$\\sqrt{s_{\\mathrm{NN}}}=5.02 TeV$\")\n",
    "plt.text(pt_min+5,-0.15,r'$anti-k_{\\mathrm{t}},R=1.0,(|y_{\\mathrm{jet}}|<2)$')\n",
    "plt.text(pt_min+5,-0.25,r\"$\\hat{q}=\\hat{q}^{\\mathrm{run}}_{\\mathrm{HTL}}f(Q^{\\mathrm{2}}),\\alpha^{\\mathrm{fix}}_{\\mathrm{s}}=0.35,Q_{\\mathrm{sw}}=2GeV$\")\n",
    "\n",
    "# save plot to the Desktop\n",
    "plt.tight_layout()\n",
    "plt.savefig('Raa(0-10)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1e74ec30",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "009adc34a3a8452e02639bd458ec0b5a4261dfbc6c9492726e6a0ec46c9767d3"
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}