LYSO_OPT2_array.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "8a171321-cb1b-4dad-b619-965f7f877dd6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Welcome to JupyROOT 6.24/02\n"
     ]
    }
   ],
   "source": [
    "import ROOT as R"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "f6888bb6-a5a3-49b5-ba9a-157647433b5e",
   "metadata": {},
   "outputs": [],
   "source": [
    "files={}\n",
    "files['MS3']=R.TFile(\"data/MS3.root\")\n",
    "files['OPT2_DRY']=R.TFile(\"data/OPT2/prodAnalysis_c40_OPT2_Jul2022.root\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d2190f00-9430-45a6-92b3-44cf5d73a60b",
   "metadata": {},
   "outputs": [],
   "source": [
    "g={}\n",
    "# g['lyDT']=files['MS3'].Get(\"lyDTArray\")\n",
    "# g['tau']=files['MS3'].Get(\"tau_pmt\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "52ea1c18-6827-45b9-a6e7-a63e1678aac2",
   "metadata": {},
   "outputs": [],
   "source": [
    "#missing Prod5\n",
    "arrays={}\n",
    "arrays['MS3_DRY']={\n",
    "    416:{'file':'data/tree_FirstRun005166_LastRun005193_ARRAY000416_THSCAN_SUMMARY.root','prod':1},\n",
    "    586:{'file':'data/tree_FirstRun005226_LastRun005253_ARRAY000586_THSCAN_SUMMARY.root','prod':3},\n",
    "    708:{'file':'data/tree_FirstRun005256_LastRun005283_ARRAY000708_THSCAN_SUMMARY.root','prod':4},\n",
    "    531:{'file':'data/tree_FirstRun005196_LastRun005223_ARRAY000531_THSCAN_SUMMARY.root','prod':5},\n",
    "    606:{'file':'data/tree_FirstRun005292_LastRun005319_ARRAY000606_THSCAN_SUMMARY.root','prod':9},\n",
    "    481:{'file':'data/tree_FirstRun005070_LastRun005097_ARRAY000481_THSCAN_SUMMARY.root','prod':10}\n",
    "}\n",
    "\n",
    "arrays['OPT_DRY']={\n",
    "    809:{'file':'data/tree_FirstRun005907_LastRun005934_ARRAY000809_THSCAN_SUMMARY.root','prod':1},\n",
    "    805:{'file':'data/tree_FirstRun005817_LastRun005844_ARRAY000805_THSCAN_SUMMARY.root','prod':3},\n",
    "    799:{'file':'data/tree_FirstRun005637_LastRun005664_ARRAY000799_THSCAN_SUMMARY.root','prod':4},\n",
    "    801:{'file':'data/tree_FirstRun005697_LastRun005724_ARRAY000801_THSCAN_SUMMARY.root','prod':5},\n",
    "#    800:{'file':'data/tree_FirstRun005667_LastRun005694_ARRAY000800_THSCAN_SUMMARY.root','prod':5},\n",
    "    803:{'file':'data/tree_FirstRun005757_LastRun005784_ARRAY000803_THSCAN_SUMMARY.root','prod':9},\n",
    "    796:{'file':'data/tree_FirstRun005040_LastRun005067_ARRAY000796_THSCAN_SUMMARY.root','prod':10}\n",
    "}\n",
    "\n",
    "arrays['OPT2_DRY']={\n",
    "    840:{'file':'data/summary_FirstRun007368_LastRun007368_ARRAY000840.root','prod':3},\n",
    "    841:{'file':'data/summary_FirstRun007371_LastRun007371_ARRAY000841.root','prod':3},\n",
    "}\n",
    "\n",
    "map_prod={\n",
    "    1:0,\n",
    "    3:1,\n",
    "    4:2,\n",
    "    5:3,\n",
    "    9:4,\n",
    "    10:5\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "17641c81-b160-4621-86af-0489998c413d",
   "metadata": {},
   "outputs": [],
   "source": [
    "sigmaTvsTH=R.TF1(\"sigmaTvsTH\",\"TMath::Sqrt([0]*[0]/TMath::Power(x,[1])+[2]*[2]*TMath::Power(x,[3]))\",0,1000)\n",
    "sigmaTvsTH.SetParameter(0,5.67032e+04)\n",
    "sigmaTvsTH.SetParameter(2,6.21645e+01)\n",
    "sigmaTvsTH.SetParLimits(1,0.,5)\n",
    "sigmaTvsTH.SetParLimits(3,0,5)\n",
    "sigmaTvsTH.SetParameter(1,3.1)\n",
    "sigmaTvsTH.SetParameter(3,0.3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "c49a7759",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "416 130.77797346170925 2620.858943799735 75.60170030191543\n",
      "586 132.34572119194485 2532.1213829811204 73.04196297060925\n",
      "708 138.02099843621147 2473.123223083949 71.34009297357545\n",
      "531 125.61110587867933 2862.7417013010313 82.57908753752974\n",
      "606 130.2758141396899 2667.23556295393 76.93948739290184\n",
      "481 124.0820169000153 2736.916621199681 78.94951791922156\n",
      "809 128.6561909778912 2586.2210430065907 74.60253008672858\n",
      "805 131.24874140253408 2757.3256601905086 79.53824019780313\n",
      "799 134.05100514352972 2488.6689030857056 71.7885260505492\n",
      "801 121.66715924027287 2695.0267306062005 77.74115569056347\n",
      "803 127.16445128059303 2697.9662164230003 77.82594855066347\n",
      "796 119.29640277520008 2867.018319879154 82.7024515349756\n",
      " FCN=0.547109 FROM HESSE     STATUS=OK             23 CALLS        1168 TOTAL\n",
      "                     EDM=1.16347e-06    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           4.72167e+03   2.43242e+04   2.88657e-02  -2.56844e-06\n",
      "   2  p1           1.73129e+00   3.36047e+00   1.09104e-06   6.36125e-02\n",
      "   3  p2           2.51102e+01   1.19340e+02   5.52937e-05  -6.52530e-04\n",
      "   4  p3           6.11231e-01   2.74437e+00   5.30210e-07  -6.53541e-02\n",
      " FCN=0.19269 FROM MIGRAD    STATUS=CONVERGED     521 CALLS         522 TOTAL\n",
      "                     EDM=4.6155e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -4.54745e+03   1.72140e+04   1.31259e-01  -8.88917e-07\n",
      "   2  p1           1.75410e+00   1.82382e+00   5.20838e-06  -2.71638e-02\n",
      "   3  p2           2.60573e+01   8.14308e+01   2.13354e-04   1.70128e-03\n",
      "   4  p3           6.14101e-01   1.03478e+00   1.97031e-06   1.91372e-01\n",
      " FCN=0.163122 FROM MIGRAD    STATUS=CONVERGED     554 CALLS         555 TOTAL\n",
      "                     EDM=1.70286e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           9.99299e+03   4.49440e+04   3.07617e-01   9.60275e-07\n",
      "   2  p1           2.18870e+00   2.05706e+00   5.67088e-06  -5.23237e-02\n",
      "   3  p2           4.65276e+01   8.65219e+01   2.98232e-04   6.19174e-04\n",
      "   4  p3           4.22781e-01   6.51011e-01   1.84391e-06   9.32103e-02\n",
      " FCN=0.221549 FROM HESSE     STATUS=OK             23 CALLS         914 TOTAL\n",
      "                     EDM=1.99557e-06    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           7.24340e+03   1.70076e+04   2.78469e-02  -8.48862e-07\n",
      "   2  p1           1.90262e+00   1.18609e+00   6.77459e-07   3.16439e-02\n",
      "   3  p2           1.62944e+01   4.01054e+01   2.28942e-05  -7.56056e-04\n",
      "   4  p3           7.64130e-01   8.49135e-01   3.05965e-07  -6.19903e-02\n",
      " FCN=0.139791 FROM MIGRAD    STATUS=CONVERGED     706 CALLS         707 TOTAL\n",
      "                     EDM=4.40334e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -3.10274e+03   1.04716e+04   8.83959e-02   3.05308e-06\n",
      "   2  p1           1.60877e+00   1.69630e+00   5.19370e-06   5.67179e-02\n",
      "   3  p2           3.41076e+01   8.75911e+01   2.34803e-04  -1.69049e-03\n",
      "   4  p3           5.03241e-01   8.50333e-01   1.78585e-06  -2.23993e-01\n",
      " FCN=0.500252 FROM MIGRAD    STATUS=CALL LIMIT   1486 CALLS        1487 TOTAL\n",
      "                     EDM=0.00128373    STRATEGY= 1  ERROR MATRIX UNCERTAINTY   8.6 per cent\n",
      "  EXT PARAMETER                APPROXIMATE        STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -1.39032e+03   1.18953e+03   1.80730e+01   6.35981e-04\n",
      "   2  p1           1.16152e+00   4.09495e-01  -3.76759e-03   5.00064e+00\n",
      "   3  p2           8.12493e+00   7.01646e+00  -2.34784e-01  -3.21470e-01\n",
      "   4  p3           1.00126e+00   2.89965e-01   4.84839e-03  -1.41614e+01\n",
      " FCN=0.25774 FROM MIGRAD    STATUS=CONVERGED    1483 CALLS        1484 TOTAL\n",
      "                     EDM=5.76787e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           3.52507e+04   1.83882e+05   1.91105e+00  -1.37458e-07\n",
      "   2  p1           2.91367e+00   3.72046e+00   1.31067e-05   2.71814e-02\n",
      "   3  p2           6.84502e+01   6.03332e+01   4.16330e-04  -1.26137e-03\n",
      "   4  p3           2.51351e-01   3.21674e-01   2.19876e-06  -2.36973e-01\n",
      " FCN=0.216553 FROM MIGRAD    STATUS=CALL LIMIT   1481 CALLS        1482 TOTAL\n",
      "                     EDM=2.9784e-05    STRATEGY= 1  ERROR MATRIX UNCERTAINTY  18.6 per cent\n",
      "  EXT PARAMETER                APPROXIMATE        STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           2.59413e+03   3.01877e+03  -7.62696e-01   4.48381e-06\n",
      "   2  p1           1.40320e+00   6.09605e-01  -1.70292e-04   1.74600e-02\n",
      "   3  p2           8.77302e+00   1.53798e+01  -2.47758e-02  -2.05389e-02\n",
      "   4  p3           1.00015e+00   6.02687e-01   5.06460e-04  -9.95480e-01\n",
      " FCN=0.2472 FROM MIGRAD    STATUS=CALL LIMIT   1483 CALLS        1484 TOTAL\n",
      "                     EDM=0.00427989    STRATEGY= 1      ERR MATRIX NOT POS-DEF\n",
      "  EXT PARAMETER                APPROXIMATE        STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           8.48810e+03   1.23512e+04   5.36618e+02  -9.08018e-05\n",
      "   2  p1           2.20021e+00   8.07285e-01   1.47290e-02   4.18684e+00\n",
      "   3  p2           6.35106e+01   2.85285e+01   9.84408e-01  -1.20940e-01\n",
      "   4  p3           2.89599e-01   1.66964e-01  -4.74486e-03  -2.23193e+01\n",
      " FCN=1.26207 FROM MIGRAD    STATUS=CONVERGED     562 CALLS         563 TOTAL\n",
      "                     EDM=4.38388e-06    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           1.03551e+04   4.84357e+04   3.95686e-01  -1.10099e-06\n",
      "   2  p1           2.16383e+00   2.08301e+00   6.64499e-06   6.03250e-02\n",
      "   3  p2           3.05929e+01   8.17803e+01   2.58095e-04  -1.46594e-05\n",
      "   4  p3           5.27289e-01   9.01702e-01   2.17025e-06   4.93562e-03\n",
      " FCN=0.302182 FROM HESSE     STATUS=OK             23 CALLS         668 TOTAL\n",
      "                     EDM=3.13609e-07    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0           6.65730e+03   2.53360e+04   4.09087e-02  -1.31108e-06\n",
      "   2  p1           1.94558e+00   1.80226e+00   1.09602e-06   5.31678e-02\n",
      "   3  p2           2.62983e+01   6.47854e+01   4.10676e-05  -2.05306e-03\n",
      "   4  p3           6.01497e-01   8.41695e-01   3.77464e-07  -2.24805e-01\n",
      " FCN=0.602765 FROM MIGRAD    STATUS=CONVERGED     661 CALLS         662 TOTAL\n",
      "                     EDM=2.01326e-06    STRATEGY= 1      ERROR MATRIX ACCURATE \n",
      "  EXT PARAMETER                                   STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -3.62954e+03   1.13202e+04   1.24026e-01   4.43276e-06\n",
      "   2  p1           1.72470e+00   1.57401e+00   6.11108e-06   8.47234e-02\n",
      "   3  p2           2.54440e+01   5.74572e+01   1.97427e-04   9.38252e-05\n",
      "   4  p3           5.97193e-01   7.77662e-01   1.88817e-06   3.39647e-02\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <Fit>: Abnormal termination of minimization.\n",
      "Warning in <Fit>: Abnormal termination of minimization.\n",
      "Warning in <Fit>: Abnormal termination of minimization.\n"
     ]
    }
   ],
   "source": [
    "for tag in ['MS3_DRY','OPT_DRY','OPT2_DRY']:\n",
    "    for t in ['LY','sigmaT','XT','SR']:\n",
    "        g['%s_%s'%(t,tag)]=R.TGraphErrors(len(map_prod))\n",
    "\n",
    "for tag in ['MS3_DRY','OPT_DRY']:\n",
    "    for ia,ar in enumerate(arrays[tag]):\n",
    "        files[ar]=R.TFile(arrays[tag][ar]['file'])\n",
    "        g['SR_%d'%ar]=files[ar].Get('DT_vs_bar')\n",
    "        g['sigmaT_%d'%ar]=files[ar].Get('sigmaT_vs_bar')\n",
    "        g['LY_%d'%ar]=files[ar].Get('LY_vs_bar')\n",
    "        g['XT_%d'%ar]=files[ar].Get('XT_vs_bar')\n",
    "        iv=map_prod[arrays[tag][ar]['prod']]\n",
    "        pos=iv\n",
    "        if ('OPT' in tag):\n",
    "            pos+=0.1\n",
    "        sigmaTvsTH.SetParameter(0,5.67032e+04)\n",
    "        sigmaTvsTH.SetParameter(2,6.21645e+01)\n",
    "        g['sigmaT_%d'%ar].Fit(sigmaTvsTH,\"BR\")\n",
    "        sr=1./g['SR_%d'%ar].GetFunction('pol2').Derivative(150)*1000\n",
    "        ly=R.TMath.Mean(g['LY_%d'%ar].GetN(),g['LY_%d'%ar].GetY())*2600/75.\n",
    "        xt=R.TMath.Mean(g['XT_%d'%ar].GetN(),g['XT_%d'%ar].GetY())\n",
    "        #sigmaT=sigmaTvsTH.GetMinimum()\n",
    "        sigmaT=R.TMath.MinElement(g['sigmaT_%d'%ar].GetN(),g['sigmaT_%d'%ar].GetY())\n",
    "        print(ar,sigmaT,ly,ly*75/2600.)\n",
    "        g['SR_%s'%tag].SetPoint(iv,pos,sr)\n",
    "        g['SR_%s'%tag].SetPointError(iv,0,sr*0.02)\n",
    "        g['LY_%s'%tag].SetPoint(iv,pos,ly)\n",
    "        g['LY_%s'%tag].SetPointError(iv,0,ly*0.03)\n",
    "        g['sigmaT_%s'%tag].SetPoint(iv,pos,sigmaT)\n",
    "        g['sigmaT_%s'%tag].SetPointError(iv,0,4.)\n",
    "        g['XT_%s'%tag].SetPoint(iv,pos,xt)\n",
    "        g['XT_%s'%tag].SetPointError(iv,0,0.01)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "id": "9c5e8076",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 2718.100334992136 0.14772512551904757 125.52457015614546\n",
      "84.6298506019112\n",
      "3 840 2933.834820866255 0.1922158605633351 130.41054762945367\n",
      "81.25601746148811\n",
      "3 841 5750.7100928645095 0.30951511765515677 257.0157485719014\n",
      "3 2875.3550464322548 0.15475755882757838 128.5078742859507\n",
      "5 2586.8269925904915 0.07967715312976191 122.85373897651863\n",
      "10 2981.459949527281 0.10598721088452383 133.01148058081446\n",
      " FCN=1.12356e+307 FROM HESSE     STATUS=FAILED         11 CALLS         175 TOTAL\n",
      "                     EDM=0    STRATEGY= 1  ERROR MATRIX UNCERTAINTY 100.0 per cent\n",
      "  EXT PARAMETER                APPROXIMATE        STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -8.59604e+58   1.41421e+00   0.00000e+00   0.00000e+00\n",
      "   2  p1           1.72470e+00   3.34663e+00   0.00000e+00   0.00000e+00\n",
      "   3  p2           6.21645e+01   1.41421e+00   0.00000e+00   0.00000e+00\n",
      "   4  p3           5.97193e-01   2.50220e+00   0.00000e+00   0.00000e+00\n",
      " FCN=1.12356e+307 FROM HESSE     STATUS=FAILED         11 CALLS         175 TOTAL\n",
      "                     EDM=0    STRATEGY= 1  ERROR MATRIX UNCERTAINTY 100.0 per cent\n",
      "  EXT PARAMETER                APPROXIMATE        STEP         FIRST   \n",
      "  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE \n",
      "   1  p0          -4.09900e+63   1.41421e+00   0.00000e+00   0.00000e+00\n",
      "   2  p1           1.72470e+00   3.34663e+00   0.00000e+00   0.00000e+00\n",
      "   3  p2           6.21645e+01   1.41421e+00   0.00000e+00   0.00000e+00\n",
      "   4  p3           5.97193e-01   2.50220e+00   0.00000e+00   0.00000e+00\n"
     ]
    }
   ],
   "source": [
    "for tag in ['OPT2_DRY']:\n",
    "    opt2_ly_histos = files['OPT2_DRY']\n",
    "    ly_gr = opt2_ly_histos.Get('ly_byProd_geo2_PREIRR')\n",
    "    xt_gr = opt2_ly_histos.Get('xt_byProd_geo2_PREIRR')\n",
    "    sigmaT_gr = opt2_ly_histos.Get('sigmaT_byProd_geo2_PREIRR')\n",
    "#     c0=R.TCanvas('c0','c0',800,600)\n",
    "#     c0.cd()\n",
    "#     ly_gr.Draw()\n",
    "#     c0.Draw()\n",
    "    \n",
    "    gr_map_prod ={1:0, 3:1, 5:4, 10:8}\n",
    "    for prod in [1,3,5,10]:\n",
    "        ly=0\n",
    "        xt=0\n",
    "        sigmaT=0\n",
    "        iv=map_prod[prod]\n",
    "        pos=iv+0.2\n",
    "\n",
    "        if prod==3:\n",
    "            for ia,ar in enumerate(arrays[tag]):\n",
    "                files[ar]=R.TFile(arrays[tag][ar]['file'])\n",
    "                g['sigmaT_%d'%ar]=files[ar].Get('sigmaT_vs_bar')\n",
    "                g['LY_%d'%ar]=files[ar].Get('LY_vs_bar')\n",
    "                g['XT_%d'%ar]=files[ar].Get('XT_vs_bar')\n",
    "         \n",
    "                g['sigmaT_%d'%ar].Fit(sigmaTvsTH,\"BR\")\n",
    "                print(R.TMath.Mean(g['LY_%d'%ar].GetN(),g['LY_%d'%ar].GetY()))\n",
    "                ly+=R.TMath.Mean(g['LY_%d'%ar].GetN(),g['LY_%d'%ar].GetY())*2600/75.\n",
    "                xt+=R.TMath.Mean(g['XT_%d'%ar].GetN(),g['XT_%d'%ar].GetY())\n",
    "                #sigmaT=sigmaTvsTH.GetMinimum()\n",
    "                sigmaT+=R.TMath.Mean(g['sigmaT_%d'%ar].GetN(),g['sigmaT_%d'%ar].GetY())\n",
    "                print(prod, ar, ly, xt, sigmaT)\n",
    "\n",
    "            ly/=len(arrays[tag])\n",
    "            xt/=len(arrays[tag])\n",
    "            sigmaT/=len(arrays[tag])\n",
    "            print(prod, ly, xt, sigmaT)\n",
    "\n",
    "            sigmaT=127.4684123\n",
    "            xt=0.166413813\n",
    "            \n",
    "            g['LY_%s'%tag].SetPoint(iv,pos,ly)\n",
    "            g['LY_%s'%tag].SetPointError(iv,0,ly*0.03)\n",
    "            g['sigmaT_%s'%tag].SetPoint(iv,pos,sigmaT)\n",
    "            g['sigmaT_%s'%tag].SetPointError(iv,0,4.)\n",
    "            g['XT_%s'%tag].SetPoint(iv,pos,xt)\n",
    "            g['XT_%s'%tag].SetPointError(iv,0,0.01)\n",
    "            continue\n",
    "        \n",
    "        ly=ly_gr.GetPointY(gr_map_prod[prod])*2600/75\n",
    "        xt=xt_gr.GetPointY(gr_map_prod[prod]) \n",
    "        sigmaT=sigmaT_gr.GetPointY(gr_map_prod[prod]) \n",
    "        print(prod, ly, xt, sigmaT)\n",
    "        g['LY_%s'%tag].SetPoint(iv,pos,ly)\n",
    "        g['LY_%s'%tag].SetPointError(iv,0,ly*0.03)\n",
    "        g['sigmaT_%s'%tag].SetPoint(iv,pos,sigmaT)\n",
    "        g['sigmaT_%s'%tag].SetPointError(iv,0,4.)\n",
    "        g['XT_%s'%tag].SetPoint(iv,pos,xt)\n",
    "        g['XT_%s'%tag].SetPointError(iv,0,0.01)\n",
    "     "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "id": "b4a1ac93-f0ff-4022-8515-f5637c0b88ad",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Warning in <TCanvas::Constructor>: Deleting canvas with same name: c1\n"
     ]
    }
   ],
   "source": [
    "c1=R.TCanvas('c1','c1',800,600)\n",
    "l=R.TLegend(0.8,0.75,0.9,0.88)\n",
    "l.SetFillColorAlpha(0,0)\n",
    "l.SetBorderSize(0)\n",
    "l.SetTextSize(0.04)\n",
    "\n",
    "tL=R.TLatex()#pm                                                                                                                               \n",
    "tL.SetTextSize(0.04)\n",
    "R.gStyle.SetOptTitle(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "id": "f3068f38-b5e3-439b-b65a-e66ec309ca8d",
   "metadata": {},
   "outputs": [],
   "source": [
    "for k in ['MS3_DRY','OPT_DRY','OPT2_DRY']:\n",
    "    for h in ['LY','SR','sigmaT','XT']:\n",
    "        g['%s_%s'%(h,k)].GetHistogram().SetBins(6,-0.5,6-0.5)\n",
    "        for kk,vendor in enumerate([1,3,4,5,9,10]):\n",
    "            g['%s_%s'%(h,k)].GetHistogram().GetXaxis().SetBinLabel(kk+1,'prod%d'%vendor)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "id": "0a5bfd58-5f11-42ba-8aff-c460e9dfd7a3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x[0]=0, y[0]=2620.86, ex[0]=0, ey[0]=78.6258\n",
      "x[1]=1, y[1]=2532.12, ex[1]=0, ey[1]=75.9636\n",
      "x[2]=2, y[2]=2473.12, ex[2]=0, ey[2]=74.1937\n",
      "x[3]=3, y[3]=2862.74, ex[3]=0, ey[3]=85.8823\n",
      "x[4]=4, y[4]=2667.24, ex[4]=0, ey[4]=80.0171\n",
      "x[5]=5, y[5]=2736.92, ex[5]=0, ey[5]=82.1075\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file plots/LYSO_OPT_LY_ARRAY_DRY_byProd.png has been created\n",
      "Info in <TCanvas::Print>: pdf file plots/LYSO_OPT_LY_ARRAY_DRY_byProd.pdf has been created\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxwAAAI8CAIAAAD0vjrdAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElEQVR4nO3dwXLcNp74cdKVS1yVkjRnJ2W9weQJyM5p90FW1swcs1Uz8R7Ivmwu4+PEtfKzjEi+QCpPEPs/pbvsw5Ru5v/wixAIANFsNgiyG9/PweUm2CQIssmfABDI+77PAAAAcJhnS2cAAADgFBBUAQAABEBQBQAAEABBFQAAQAAEVQAAAAEQVAEAAARAUAUAABAAQRUAAEAABFUAAAABEFQBAAAEQFAFAAAQAEEVAABAAARVAAAAARBUAQAABEBQBQAAEABBFQAAQAAEVQAAAAEQVAEAAARAUAUAABAAQRUAAEAABFUAAAABEFQBAAAEQFAFAAAQAEEVAABAAARVAAAAARBUAQAABEBQBQAAEABBFQAAQABfLJ2BNWrbtm1b+X9ZlmVZ2ivoH/0rOFP923dmybNHO9VYMkS241x5TK7GG8qPfy97FZS+cpZldV2PycCY/Oxb+ONzZe9oQrFPuJz23b7a+Ho2tap9TeO/S6zf+ks4OPuQEywE7NDjqaIo7FKqqkpfp6oqPbVpGk+q/l0jaWj7hqZp9JWLotiZ6jwK506HsuQ8tMn8ObH3YhzU0LGPX3l8gezM/MjUoVzZW2iaxsie/3qYVkqTqbwdeDGo3IbN3pCh0l4P/aSH+qHFtP4SDs7+LRz1GcQcaP77Xdu2eZ53XWcnbbfbPM+HvmhUP2y326HVhpK2262/DkNn5HBkHcw0m81m1u2rvegf27Y1lihd1xknYq+Vg/MUzlCusiwzcrXZbIxzut1ud/7hO8eB13Wd53me58f1Z3dZlpLtCNdqQPpJH//zB7BmNP/9Tn9EVVWlmsZUJFTXtfPe54zDbGo7VVWp7ZRlKV/fK67SOR8k+qb0Q/BUnxRFob6lf2Wz2fS7qprGG9qL3mS580ToK6utFUWhLzTO2oQCGaI30o1p1myaRh2COjR1CGpNOVKVt67r/K2Be5XSIfRr9ZDtlGU5obQni7mvCYwf+8h7CNZGXWbH9XcI5kNQ9Rv98aOegtljJxV5gI2Je4YeY2q5HlVkj09ouaWO71Kjr+m8HevbKctSf9B6NqsftVRayMehaHICvdOP/EfFEHYO9WBOPxF6wKG+aHddki0bsUu2T4E4jQmqVK6Ma6lpGnUIxpqSZ4mEdl4P+5aS0fNDumFN7oClb0316NK3tu/2nRtxruBZZzzP7pwFFWSnBrvS2j7dRmbqupZsODNp3FWc2W4HOgANLbdzKGs67wZDuQq1U88p8K/j3JGdWzvz0076vpuSwpR/R5YJVm3h5sfVUI3lzg4fRnGpfkh2E7vdRUn6x+h/Nxsbbx4N5U3/ruxRz6SRkzH51+nHMj5pX0ZpOA/NWNNZIPZheg7NU7B7Xf/jC3/MURi5sktYFbvnkphcSnb9jb1NY+OefiR2VzDP9p19qpxd+oyDcvaHszNj7NHOwJjdqRWGfshB6KXk+ZXpOdTzoJ8R47vOA1Q5H7pLOMvKYJ9r41t60amVhzauVvCUqrOu0d+B1V5nqCefkSs9P56tjelT5SwiZ/HavW+H1sQR4eT9ZuhH62QHVcZtS78BqSTj9zO+Y6PzFmwk+QMgz2918hf3YpeGsXfZhSf0FJ4gbK9S3eu4xhd+Zt2j9z3XY/J2YCnZ7DyLnUHVGOq6sh9vnhY6tS/PGwbOzHii7Z2b6mO9t6HHH86sDpWzEVQZR72zPJ1XjhGxOfnf8/AUXT/w2x86ZDtXQ7sb2qN+vP3+QZV/j3sFVbadf07Ye8TRoaO6ad96V/mZqWpbabhxbkT/3W63281mI71r9RryfbM3/osr1GodgIxiHLrH2QWr33z1Ui0D9Shy7l227Kmr19tHjFwNNaRKG4FqcvU8VCaUks6uOlVND0aEvfO34HyWy9NFbcrTW0gVhV3Ra/ygMleNmqzTP33UDZ309rFRNXNFbM4jtR+Toa4ovcNfpp1Hz/aLolA95+zlehNSNlyezvavnX3m7KIb01+tqirZu31O9cwPfV3lyq5mU0l6L1VjHc9rImNybpz3aY1xdvWh/lswMr/yLoDYQ5jY7PjtVSD6s8eovZD/G3VL9hdtniY24+9LtQv1f3X/GtqO5+gOrKlqnhpabeR1OKbBUX1F/2N0TCXEXsc1tLIq8F6rlRyqOvLkytks6F9BmVBK6qP+FTvbzi17/jq3s+SpjbDrDPRc6WdTv5bk2tZLwy4i56APxt6H2puGSsD5s/X/SMf8EPrRxd4P1PH0w614TdNIcXlybh+gfTYNzqLzNP/Zmxpac8wdw14ox7Kz0Iz4b2RN1cjVxtRUjSkB/+8Rx4iaqkPp77Lt7GZY17WUu97tQHRdt1dncP3PvqV6NZZluXlq2nam/ZWmjloqpfrHh4pRsGFHVdBr1Dy1ksKTK7usyrLUVws1mIWRPf0aW7wzrP4XvFTmSZb0Xr3ykoSqxjv8bBpHPaaWaGdB7fVD0OtgjCVDtXpDdwZjudSDqn/zPLe7w6tj0V/jyJ42YRuMdymM7djsTRmFvPM+OVQXKzdPYyNGIfgb8naK86NwZn7x3yOCIKgyOe+t7SP/F2WFMb9qeU7IPUKtPzSKlUE91+d+VWSOtkWp21Pkzzs9PJL/DD1d/FmSx4kUrL/LwmTl01cOs3FPXJUrPXy0j0VWUzkfejYfWErrUde1/aeFRFf6U0eCgwMHHXC+XroIlRMZjcx4L+yQ12zbts3zfLPZeIprqP/AvvvdqxjVxuU/krcDo58hRhv95K+LMTE3oCOo+o3dYK/z/wFqRDnO2416q9b+ce77c1U3prl/5/oQUM4VyrIsntq5TYkwFKOsnH0+/FmSUrXPmvGG+c6MjaTXoBhLjNWcuZIS03OlamLUOjsfbxNKabVUrDlUk+fpN7OXWZ+O438Ixt67R0Mr7EUVmmpOddYB6zVkI/+Qs+2VT/3PgMlh3EgH/qlpn6BDtoYUxWpnXDv97mO09Ht6Dxit+/rX1cfKmg3Gs+udeZMl+u78PQyM9e2koS96CmQCozQ8RhaUvwfDIUk7V86e6r0nyC46o1uGMyc7s6c/to0kvQDtUnL24bA73MTpU2X3QDJK0pkf++jm7lM11C9ngjER4VDe7I0YV5fnvAx1FXKWydDuRvapcm5N77TkPCiD89j1YzmkT5X9g6VPFUKhpuo3ehXCZrORyoO2bfU/5oZalDzVLYrRdq7+Hmq1UbantViF+hOq6zrVGCGNCLK8KIqYf6XpBaVejZR6Pn+W9DajdtLMxCPpz8WhU2ZcS3quhv72Va/gqW5DIys8jFLSX6i0C0F13WufDu8+tKNZje+K5+zWsxf9lTG7MiNalZ79MqPiryzfi+p1NLIiyr9Hu+jGb9neyL5tf3bfgPJxKFS1QWe1sXFVq0qydnh+J2Nr6pcYvC+BXnVn/x5x3OLHcWvmKSj9Twr77wy7SNVH9Xfbzl/mUK48f1dl1hvsk2uqZr1C7NLw8Pdbtw/Qn39nNdteR2esbFcFOf/EHJkrz8EGLCV/Zoa22ewap0p9cVpNlefJWo14z9yuqcpcP0mVH8/u/AcSqqbKUwvrTB1aeWdNlc3IufGr35lzfwxkbHPoN65/ZWfl95jfhefe5ayz9GxtzAH24Wqqxu8RR4eaqif6gWu9qip/HwK7cttW1/XQnUJueePz6fmDLKDq6ciE0ZRlaXeyUVmyT8TQytKzJHj5jCz8oWsps+aucT4bdr4RuW8pDW1WP8WRq6zath3Kf/04aY9ROPJ+g/xf/WWvV7R4+rPrLwEo+/70DqFOivOMH15BZZ9ctcQoFn1fY6ph7DM14Y1dfUc7r7RyYKZI/WTVdT3029G37ykWfybtPQZkF+lM79YgsnyRp+b66b3Ogz9pVC30TNs/JXudiFnP2mTGGwxDGTsk8zu/q48pWj5O2De0sj91Dv78q9QxrVTjrxP/6TheI0+fuiTGPwIOvDBUC1ehzX0+cqfZAXP/6esMrVA+tsuP+YEEZOxowknB2hBUAafPCKoWzQtWQS6JveKbA6n+qSu8CI2gau7dtY8TYJfaO8LTgk6szRdLZwAAEIk8ziOPFCAvUjA8gU5KQ7rPS4EYkxfhSBFUAUBCVHBTFEWc57f+quC06RNOjPQYk3opY6CyyjXDI44IQRVw+lTXde7XkLE2YvaZk8uvHJ5NfHEyBm/MPZZl2fe96g2WrawnKCajTxUAAEAAywypULpma5FJS4b+mjkkFQAAYG4L1FTJexZqKBp9ofpovP5wSCoAAEAEsWuq9HdAFHkbRR953Jh3c3IqAABAHLFrqvI8L4rCqKmSEVP0nMhqEhsdkgoAABBH1JoqT7hjDNgvgVeQVAAAgAjiBVVqLm5nqv9VUiN1r48AAAARRBqnqm3b7XbrHPbNGWapvudhW/HUZB0AACCCpEZuihRUbTYbGXHOTtqrjupwSZ1dAAAWlFpdRoygyp45MnucgmpoDFl/BdUhqQAAAHOIN02NPv1T9jjhkYqoiKIAAMBRi9FRva7r/qksy2RkKQmq7Pf1ZMwF+b8zVb3x508FAACIY5lpagyqfVA+yn9UW+EhqQAAAHE4mv/aR9ljjDL37NllWTZNs9lsVI82/T1BZ6oeRXlSAQAA4vh9RHXpOe4ZNtOYrW8OeiQXNlXk+QJzHQIAkKbUHru/Ha2MC1UUhbyRZ6yk6q7siZCPS2pnFwCABaX22M37vs/zfHyoJCHXkb5wl9rZBQBgQak9dhM72sTOLgAAC0rtsfssS+xdudxr6dwBAE5Zn5ilyzu235r/siwb6lB1SlILmQEAq5LaM+jZs2dJHXLe931d19IJXRadcHRFUAUAWNDhz6Cdr7qrTs/O1870Z32EB31qQdWTqsiqqvSxyIuiaJomclXhrLL0ql4BAOvx+TBqrpGqqobWUQ9xfeHt7a0zBri9vTW+fnt7q89ZYq+wl9Qeu09GVJcqq77vm6aR6V9kUM2yLI/0dT8AAE6PMZ2u4nxYt2373XffZVlWVZUESbe3txKffffdd/pXZE1Vm9V13XfffbfZbEJn/3T5Yy6JrsaseRRO4ygAAEfqkFqfz1pNVeaqZPr8+bPe3GQstNeXrRVFoVcsGWt69jVGao9d39x/O8dYBwAAkUkw5Hw6yzje9sLM1cXK2I5UWRl9raTdkEhgJEdQ1bZtWZZ5nm82Gzk90rMqfuYAAICT3QIoS/TarJ0+a92tJHJyfp0uQCP9HlTZsVRVVX3fy/LlcggAAJ6Q0MeIdYZeDJS6q2fPnjl7Yqn1pfO78XX5CmHASL8N/qliqSzLVCyV1KCgAAAs69kzX58cnQRVRpDkbPvLsqxpGvnPdrt99uyZRFf+yqe2bbfbrV7JMjJjifsi0yoMy7IkGgUAIKZnz5591oZCGE/v6uRv+/v8+bPESfIVFY1VVeX8irwtqNaZkLc0PcuyTLpMneqAnwAArNnnz59VHdW+lVWqwmnnoKBlWTZN81kbTyHLMqmOcmZJH3ZhaAQHGAZHGN95eo7Rztn96I8PAJiP5ykj4ZS/ymq73W63W1W99OzZM3mZzP7/zk1l2vhVt7e3nsf9yK0NfTepB+sX9qKyLFWNYt/3eZ4XRXEyPf9Xfnb1cnbOMOBP1Wcn8AfEx/j+wQrz7JkOYuQXnV/3p07LzOSsApjDUI2Uc7knoNGb85ztdBKHOcOmsiyrqpI2QTXKN7eIgxjjVslCNV9N3/fyn6IoooybNS/7eFfF+D3IGwNjUocavI0tONdcdiaiMflRXSx35llf01jH2Eio61n1CR1ZjOqkjPkZDqUOUYXpPLp9twZgDp5BMlXk5BlLU37mao4a+agmljFG7JT/y4gJ+gifQxs0NuLc2r5Su+eYc/+px4P8X18eP3PBrfwojAjDeDQa73SogMnfhVDfyNCaS8VVY/KzV56d5ePZyOGHsFdQJYGdnVUjJ3YQOT4/Q19R2zyNv46A4zUUeRj/GWIEVfIVVfcxFAbJR+eo6Pqa9sY/a6NY+TM2JMid9og8qWYcqvqTsRVOpgXwWBgj2A4NaKu/xKHOqwojuq5TJ06t2TSNPgHRUvM66TlvmkblWc+PsY7Ksz3ex9AIIPLCi/zfOPDIg4Z4fkHG3FvB96WOlIp9YIU+T+qupHgG7RQSGMksfjKYgtwYpalRn6Q5e2wuVOuoTleH5DAheoSlmvx6q3YqW7qdKIhs3SGz/ZPQU40kFUI5V+61ShRZU9VV6LHX0HcjcObHqPhxrlM9MjYoaxpH3Wulql/AOw+8eaQ+SlRnrGZn2Jk3tVNVUTR0Nu0kO5NqL3Z+9JujXTgrv/6BFOys1/GvYFcmqV+9sR1jyVBUNFQvpZs88d/n9GqqHO0OdvOfHmwdtZUfhd4nxogD9OYb/ampNxXZG9QjA7XmXrHFfJyxjhETjJ/MW9bUB1xRUYVzIzub7fRyNu4v+ldGrib0hWo1u8+icaL1nNstg/aOPN+i7Q9Y3OToJAgZJUF4QiW12iHhlFjk+bIgd/VGURTqKWU8oo7aCs+utEbJ//WgynjWGuclc9VUZY9tZON3bTxr1dxEcehVQSoDeqCgxwdVVcnh2zlUB9K7qmr8QdXQ8TpHJbbjlZGr9a6qX5UBfWW9cs7ezpgd2Qud8SuARRwYoxydbH2P3Vk5jnbM39xHaj1nV+/Zo3KlB1XGs1Y9aD3NWwZPgKV/Ra+9UBvfK7qSYxmyV5RmFMhQDOEMI/QaVmf5OKvoxgRVdvxqBLv69ofeQJTSsDNg1CDquTUO1q5d0y8YtWW7BJyFBmARSwc5saV25/Ed7cnEUsoazq4dA9mP/+JxALfsaXihh2LViLfbsoEWH6N+RS23W5ckJ/se0c4M2JyxiLEdZ6Aje7dr+9QK+parp73ds3FBlb7cWOg8Hc62xWyg41SvRcy9ViVpB1Vj8qMvLJ7236LtD1iDpYOc2LIVPHZj8g2Hz4tCAcms1XmeqzfRVLzifAdNFb4xpOfQxuV06u/QCRnSzV5fX00NNF+WpWxHpXZdt9ls8jyf9UW5uq7VG39VVdkZrqpKyqF/DBFUMcp/PG/MlWWpohw1P+j4vPlDRn0vnlT/kHqSPVlH8ubfWqnxrCabUieO+dEBYHb9Y8PETosEfWGNKY25d+2p+zFac+T/eufrfkRnIGWorsW5jnNrRog2U3cruwVNcea/0PpIGb3QjBo4Z6NY9rTSa0xHdedy+aJzO0PtdEOHpp9f9UWjpsrZRV3n7O+vV8vtPBEAIli65ii21G4+X2TDAyCdpH7paWp2VjAoRVHoo0wVrj7RbdtKDYQ9H3bbtvpch6rGS6+xUFMSOWt6ZINt2853hagcFodNheTMob6wrmvjwI3/7Lujvb643W6dp09tSurb9NGkhkqjKAp/1aZsR21QMuzZOwAglCdz/8n9mla/OUiElD2OqzamqKVRbOdTXFbYbDb+eFEfA9N+KtsBmb6+PwMqsHMqy3IoVW1wKKJShaYzlhjhgp6qmtVU65va4/gwcbvdBmk48xeR2lc2EADpOfdPFCik6NRh0vYHADGoOiujfeT0eqn3K+iorr++Zxe1/TKXvqbR3uR8Jc0e10ptze4J7uwb7uxH778YpnVUN3JoGMqhOvyhzfo7qg+Voa1wdWa3l4xp/rPb/nrve47VwAhkQ62K+jp2OThTASxl6ea42FK7//xeU6XaR6TtQDoOU3cVlt5rWKolpCd4P6JRcugsNE0jJ0s2Za+g6mlUxY+xWlEUsnHVGpjNf+r1+hW73qhtW6lYGsrz+LZC/RCMjYysv5HKxQlfFG3b7mx9q6pK7cLzOoLk3z7Ldo8rdYFltP0Bp0u/Nak7ucHTu2Pnaju/BZMn4BoazWgyz/iNs6Yq/uONT/Ujlo92TZVdMeOsZRnqxWyfNfv5Wjwdp2pnGYbif9I7RyhQPDVndk2VsHfhr37Ty9nzxTE1VXZmequmyq6Xspf040ZU9+cNwLJC1QANNRHYY6A7V1P3N1lnzNR+hTZh83jZyh67c8v7XXUk6k/enWt6tG2rqr6yx5oJfYPSYXkoVa9ByawuOP5UXZ7vPt7jJUdt9x9yrqmqguLk7UDOflGzbkRdVE3TSLfxySUmX9w/v74NHl4aABYR5BkkEyFn2ojQ8laT3LVub2/1O4Mxa7KwVzbq4+WjUbcy4Ybz7NmzE37sOgxFWyPHSBzJ2Jf8wW10TDH6AxmjOxofjW4lQ6l2Ng45CqSDah4Acziwgurz58/F4xDBdpKKgYy6osw1T7NzZSEVV9OqpgypPXbNow0bS/2+G1dzjNEde45UOxvT8o/UEFQBmMOBMYpqpxtaQU0yoYc1Q+sTVAX3WxWitFDkeS7v8KsXvkK9id00jec1/p0dae035uw354dSAQA4CjsfXtIq5+kCJUEVD8GlPMuyzBlLhe2robYm3UGkB5UxgYbMhSLTuWRP30QwMrPXR2CCuq6l2/7SGQFw+u7u7v7nf/7nD3/4w2az+cMf/vD69eu7uzvnmmNmsiq0ma88DhlvGUOeDP5pvzqu60P0NdOneFML1Zvzau8qNfhZ1wcZ30uQw8ex8IxpDgAB3d3dXV9ff/PNN//85z8vLy8/fPjw7t27V69e/d///d+LFy8mbFDeszEGCjZuaOqBO3J6U4z0RRZ3GJu+7/XRuqWOSq6Aqqr0gbLU+1bBMxB2gwAATPb27dtvvvnm73//+/Pnz7Ms++Mf//j3v//9v//7v3/66af//d//1dc85C+97777zl4oY+hM3iZsX2TR6wDlPXCJqySK0iOqLMvquvbPN+fP8MlUMKg357Phl+eHDnbMTCb+r6zQHEH2gfQTNO2Lzq/7U6dlZnJWAczq7du3//znPyWiEl9++eV//dd/bTYbI6g65MdrB0/ThkjADv3j63LzveUkbxQaC9XcHTKAgrF3fWSEzHqbTx8w059qGFq+NkN/OnjGkLTpxTLySojPP7SmsEe8tFcwalvtghqzown2fUlQ5c1/RNNOjT14rHObI7cGYA7Gm3cjjXmVT5H7khoFdMxXDLz9N9nvwzeroS+CD6htjDsljMDIk2oHSfpjw59qOIqzO36o8X5EtDRytUVKZih21AOUoXXUCkOjyesFNWZH0+wVVElW5f+ewt8ZRHoMfcUfbwGIxg47zs/Pf/75Z2Phzz//fHZ2Zq8s9xw1EvpQHJPtGYcZCKomG/yjP2x0ZTx4jBE71VWip6q7vzxjjCjKGCl0KNXORqgjmok9rbXQIy396PQYQq2sb0TNzju0KSX+wRqZd04PbKxjz9tjl8DOjXjmId7XXkGVPrNy9pQ/jhyfH+cG980ngPnYYcff/va3P/3pT//+97/Vkn//+99//vOf//rXvw6FO9lwkGQPDepf34mgajLfiOpD06hN3NNTzja7oVTjb3fjweBPNfJw+IHMSn/8G0l6EdnrO9tP7ZL0J8XkvLqMZ786s0a1k5CP9rEY33JuZGeQYcSa8nOwV7YzrOfNoOfB+DnoX8ksdt5kL3Z+hn6zQ5sCEJkddvzrX//6z//8zz/96U8///zz/f39zz///Oc///k//uM//t//+3/OMCXUiOoeBFWT7T5aFV0d/jfuzkqR+VLFys/uzhoUO4QaCqqcFTbGXlYSVDkraSQmUCGLZzvOY9E3MmZHNrVrf2eskaupXNknTt0f7SOyz2AzYkJlz7do+wMW54w8/vWvf/3tb387OzvLsuzs7Oyvf/3rUEQl9Dr7qqpub2/1W5Axp7Is3CsSIqiaLLGjXd/ZbbRe/PoT2rmyp/ErSE3VHD3qhjgjbKPiR48PqqqS0Tid3e/0ErBDsZ07svl7tu27Wv+07a9/WiWpr6zXq9nbGbMje6EzrASwiJ3RzMhgxdlV1BkGSdJekRBB1WS7X6Ry3riP1HqOwihnWag+Dj387HhIfzbrG3cuH9qIYvz1s9exDNkrShsqEP/V6Lxi/QHEzqta36bdZmfHwZnVZJlZdVp2dZSso+dED4CMTNqBoPNU2pVwOw8WQDSHhyk6qaOSyqqwWw4ltTuP+7Uj59NxoRyGtIaza/95YT/8dgZV2YiwY+g5OqYSSylGjLXh/INJ38KYYnHGIsZ2nIGOs0UsG67tGwp6DHaloJElY7WdHbaygY5T/dPGStUaaAdVxsedC43I+zR+wsCxWzrIiS1bwWM3Jsdf/EtlJYIFj25MvLLzYe+pqXJybmdMnypP5DdmZeMw/SVjbMEIROy9G2GEsY4/ZtJ35K9CG2qKNfbujJ/shY02mIK9nULrVqV2agRV+kf9Tx0jP/bGafsDVmXpICe27KSDCpsZVJ32n7MDz/3fRdi1p+5n6EFur2AHVdLQpnieoOM7qjeu0RmCs1vQFOd5KbQuU85jGeqk79mRJ1f+zlgjgyqjQ1U/EPeo5dIJTF9nqELOeUR6IFVYPcwALGjpICe21G4+TyZULsvSM6HyaeiXnvtvaMKZ7HEOxCzL1AQ+BjV1j51a13XwCQdkg/4pgw6k5rcuiuKQ+YX0Y3eWw+QdGXPjqKLYq7S3262z45falD4VZuadzrkoCue1oedHZs+U/0iGPXsHAITyTP+g5jNeJCunTT3VttttnufOp6Ze8vYze2fcEFDbtnVd53m+2Wx2hhESdgzxXE5qg0OBjjMU2Bnh2ZvauaO9tjbNmHKQSMh51GodmXle59mpJwoHAIQ3vuZm6Uq1ABY/ChkXQC9VozXQ6DpTWSOAZwO9jsb3mPE3/+3Vm2roK8YBOr811ElIDL1JZzSA2u85GpsduSObs1O8vWRM85/d9tcP9wxzHsjQvpzDLjhPivMYAcS3dHNcbKndf37rJzvG0lkNYD1n1/PM87fUDPWbDhJU6bv2BBz+YzEMXTn+bzlfAHSWmKe4xuxoZ/e1kXv3BFXOn48nk3ZAKet4ulU5z5RKPY0fL3Aalg5yYstW89iNI7GjXdnZVf2IjeVDcTqYByIAACAASURBVID97PQ/Vp38QVURa/xPf9TirI0bOljnm5Ujd7QzqLI37lzNE1Q592JszR9CqW8546qhU79zBQDxLR3kxJat7LE7t7x3PbTatpWuJNJjY+4ePNHkuft410mdhczbvT0F+gXp6dqlOpUfXlbqpYGmaaQDnOo6tu+mjK7uhxtTGgDWKewzSH+3rCgK5w3B2THUXnNM/9EJN5xnz54d0WM3ACPIcv4pfDLNB/bxAk7OKigAOFCoGqChBg17aHVPAFBVlawj89L4FdqsNbe3t3oPV8947qk9ds2aG3nzvKoq+SNY3gLruq447I33lTiumiosyKipWjo7AE5EkGfQs2e/vbmvXn7quk4Nf3N7e6vftWRlIwizVzYGVJKP+rdUTdh2u7VHXyqezgOh7z2tx64eYUnxDY1DfQJ/smeJhcyYjJoqAHM4vI6qeJyDwVN9ZddUjVxZeCZUVnVaeu2ULFH1Xsbely7yqJ6MU6X3o9IxyA1SU9e1dNtfOiMA8DtVw+SsFpIZNbIs22w2OzflfzN6iOxdmrPUQom0TqA563BPgipPF+AIWYkj91o6d1gLvXP60nkBgN9Iu5unC5SESvNNgyHxgBGQyX1yvp0ekVEjqm82m6HXCo6Ov+Ju6dwBANKlOksNkcDF/ziWKvadtSHTqkuqqrJDOtkUVftZlj2Z+69t26IopA+aKh17lhL/3CMIZcyQCkO/CmPlmd6VxRj6SZz2Red3/anTMjM5qwAOtzOiGknes5EprdRC4ynQdZ3dFX3kxu2F33333YRNnSa9nmZ8iczZzWtGx5LzoRPhGUPSpg+Esdc5Nfa++IAahTWxj0raecWqNY3+B86DmunA7cx4SD7tzBin3p86JjND76Os4YwDJ8zTA139PIdWUFVEno18fuyBrrqNe+6Qzq7on70d1YeyNLTyyFvfyXhSU1XXNVVQi1Mv89vkD4uR56jruglDSNh7l+0sNbKA3dFtwnHVdW28AGxvxHngpWve61kNnVwjG5NzVVWVuoqGtkk1FRDfmGqqQ36bzpknDvyxbzYb1W+daqrfLBzUxbX+4zUGBWke6bU1zilcZOploW9Ezc47tCmltyYwNlZetjQkP3rh9I/z/BiMy9tzUKqCx17HWdrT7FWAKlf+6XGMpL1mFhrKz4InGkiHXoUznvEtf+2R3OLUkAdjvuKsf/LUVOkV2/5NpXZL+W0gsvE3ZRlqbL4MzWr9Z9fzmHTGN2rJ0HR49skak+SM2w47simcg0XtzI/xrTEHZa+j4irP1a7Ho/JRQltjNWNfQ6uplY1sO4/XSLKvFtmL83ftvGbGHC+Aww1FHoahGEXub84RoYytDX0cwx9UqbuTZyB1fe9LF3lUvx2tup8O3evVsK1Hfc9d+dnVn6POFezHofMB2Q9Mx2vsxTmhsv2VBYfB9McTzq+oA1fxxJiDmrCj/mnh63dDo2CHVhs6Nfr/M+2FGrtezU7qd00vra+gf4vBToE4PGGQwRPuDKV+dg0N6l/fsxfP4J9julupvS9d5FE9Uzd3uVNvNhvnAE7SD6Pv+1Mas2oN9JGQVNkONU6rJ9/OblVqU3u941rXtdH4lVmvf8Z891NVBY3/ihryLmbvQGOcPeOlm6HVMtdLms7zZbwg7fkN2r3H7PyoktHXtF/yBRDHXi/9lWUpNwTn8J7b7dYzNGgQct+Yb/tHz46zpOFA9VDxNFUcHefxLsLZXUl99NQXygqq7kF9Ra910C93uwFor5e87Bov9dFTr+nMud35SdlrO57MO6tb/G2IUj7On4PnN2KsoDbu7Iyll5hxFPYe1fnSj1T9X9+gqjzWv2Vv1s6PZ7WjrocGjsJeNUZD1G9f7g+3t7f6XcVolZOFe23fUx2VDRtaf+kijyqxo13B2bVroewn4s6gKrNCnCGejex8gjpDBLsqZUyHPH8mR3bpG9k8ahyXEWIa4WyQoMrZ+mYHvv6NN4+DKdjbUYegf9fucW83ffYDnaWMFkBntzMAc9gruPFwNmh4wqC9Nj4tqHLuxXMLPUmJHe1yZ9fZ08V4hqmkCTVVTs7tjAyqjDfvjFQjNMke/2Aa2pq/FXLMs9zZT9+5gr01z96DBFX6QjuIGblxyaT6aJwmdWh6gGUEVcaro4onA7JxZ7czAHPYK7gZE/3IKylj+owvIrUby5NxqlLgn+Cvn22mGtV/pSiKuq6dPVeKopDm8KEp5zydaYqnI44MjcA+kl5KjWuEKjWAU9u2dV13j/Yd+X3f/BRFMbQpT68gKU+1gsQQY2Z7mEBtcN9psLquG+pLp/iLUU8duXdZTf7dq/sdgDU48FaP4JILquYLm0by/AbUg3+73Tr7WasnpZ06FKjtq21b1f/RE8Eo6nBmmkpTz09VVTvHxhyKS1QIqLrby/Lg96NprwjszIwE3P6Nq+unejqBvFP1OAqoWpOBfwHgQGFmGsJO6im43W7zPHcO1a0/1ZyVQ57UUPQIZme9SFmW8maoquoYypjEfEM8OxoTUWVa0TmrqWSqAL3+zwgB1dmx37AbExvpGZtWJydb8JxWWUGqA/1rCr141RJ7gxnv/QFAQEu3P0a1+PHqb2wJo2dV83TQcOlYbdS+OF/jCvIC3VCnHKFyOKGvuv8iHPq6vqOh/Bjbtzdid7e3X5qz1xlTsPZJ8b/958mwfUT+/nNDB2Lvy17HmX/e+wPiWLqPU2zO2/IJcz+ERPTMzG49Z9eIk/Qkf9WI8Wgc8+wf2rX9HB3TndyTGY+9Dso+OifnmnttZ0wOPQflyZtetmMybBeCP6gy8mwPxDDmQPzvIgCYw9JBTmz+u+jpefJGt/O5EvBuqwYrcj5Em8d53Jyp/u/6U5W1nV1V5sbyoUejfS4mnCZPUDX0PNZ3ke05QtUh/Plxrjl+U9PWca4/ckT1oYXGYArCPk3q52kvcY7pMJSfnXkDMJ+lg5zYUru95H3fZ9pAzHq3mLZtVb/XAzuxqu7G+otXvXZP1zNgp+qvbmVWB2p/qi7P837Xo3o92rZVB8IrHodT5ekpzDHrKOqdRLmo5LvTztTQ+56THZIZAAjluB67AfSPfygP/TkrYdCBNRNqX0L2aPRoMT4WT4efNho4VH78qXY2DjkKQKf/iAAAttTukHnf9/IHdz8cS+5cYac8z43qrjzPVZWSVDX1Tyuu1L/23vXv+lPtbBxyFIDOqKkCABhSe+z+Nk6Vf9RBGdLmkBaKoQEk5T9d1xnduYzWRvuNOaO9z5MKAAAQwTNVV+RZ6fBRAfV+Wm3byp/4xrBMMpSRcwwnI3t7fQRmol6PWDojAIBVGDWi+oFzjOj0sRz1jate6lIrttlspHIr4K6Ff5oaj6QqMDEGQ5ADAHRfSL2Ov2lvTG3WSH3fS2WVRFH6Y0lFLXVd53m+2Wz6vg9e7URsBAAA5vDbNDVqul8nf+q+pKWvKAp96jGjDcXfx2v8tLIAAABxPMsehy7M89wZjkjQYwxvuBdnNdiYKih7Ijb/agAAAEv5IsuysixVT6ahwT8PbIbruk4mtVVL9EjIfl9Pb3C0U/W3Bf2pAAAAkeiDVtmxSKgJSWRrQyN2+kf7NFKNwUj9qXY2Dj8WAAAwRmqPXfeoXMEnzcis1+6M8TnVNDXOVDXLjTBGvfKnGnlwHi8AAAgutcfuqKOV9/UOf4Fcn8nOs8IcqSK1swsAwIJSe+yOOlqpRjqBcknt7AIAsKDUHrvPls4AAADAKSCoAgAACICgCgAAIIBRc/+dEv/cf0k1/QIAgICSC6oImwAAwBy+yHZV3gAAAGCnLzJrMmMAAADsK60BJFIbMAMAgAWl9tjl7T8AAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAkhuRHWmqQEAAHNILqgibAIAAHOg+Q8AACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgACSm6aGuf8AAMAckguqCJsAAMAcaP4DAAAIgKAKAAAgAIIqAACAAKIGVXVdl2VZlmVd157VnCv4vztyywAAADPJ43Tcbtt2s9lkWVYURZZlXddlA33Gy7Lsuq6qKj08koXqY1EUbduOTNXleaTjBQAAqT12I9VUSUTV933btm3bNk2TZZldq9S2rR4eibquJczq+77v+6qquq5TYZM/FQAAII5IIWSe50blU57ndpWSLDRqqmRkKT2f+nf9qXY2kgqZAQBYUGqP3UjjVDVNU5alsdBYIh/btrXH55RGQ/2j0d7nSQUAAIggUvOfip+k+U/CJr3iSlrxhuJZZ/g18iMAAEAEsUdUl85VWZZVVaUWtm273W6lo5UheO8o/zQ1HklVYAIAgH3FDqpUX/Xtdps9VlZtNpuiKJw1TMGrnYiNAADAHBaY+08GlJK4qq5riauMIabatlVDT9lb8Fdf8eofAACIL0ZQJRGSEesYg0tJxZXSdV3XdXpPLP/2w2QUAABgqkgd1buuM0al0gea6p/KskzGnZKgyn6br+s69cafPxUAACCOGEGVxEbb7dYesXPM11X7oL41FaL5UwEAAOKINyqX8dqdfzIZY6RQNcuNMEa98qcaW6ajOgAAcaT22I16tCqKmvZOn3x96Lv+VJHa2QUAYEGpPXYTO9rEzi4AAAtK7bEbqaM6AADAaSOoAgAACGCBwT+X5Z+mJqlaSgAAEFByQRVhEwAAmAPNfwAAAAEQVAEAAARAUAUAABAAQRUAAEAABFUAAAABEFQBAAAEQFAFAAAQAEEVAABAAARVAAAAARBUAQAABJDcNDXM/QcAAOaQXFBF2AQAAOZA8x8AAEAABFUAAAABEFQBAAAEQFAFAAAQAEEVAABAAARVAAAAARBUAQAABEBQBQAAEABBFQAAQADJjajONDUAAGAOyQVVhE0AAGAONP8BAAAEQFAFAAAQAEEVAABAAARVAAAAAUQNquq6LsuyLMu6rmOmAgAAzC2P8zZc27abzSbLsqIosizrui57+iKejHSgpzZNU5alpJZlKQtFURRt26qP/lRdnkc6XgAAkNpjN1JNlURUfd+3bdu2bdM0WZapWiUJnpqmkVQ5AfIVWa3ruqqq+r7v+76qqq7rVNjkTwUAAIgjUgiZ53lVVXrbXJ7nqkpJ/7+o63q73UrepBLLqNbSv+tJtbORVMgMAMCCUnvsRhr8U2/LU9SSoijsVJ00C+ofjfY+TyoAAEAEkYIqFTNJBZI07amKK7tWabvdOr+uPuphkz8VAAAggtjT1KieUlVVOVdQXdql31Xw3lH+uf88kqrABAAA+4odVKm+6lIXZYyAoCqZVHOhv1lwWgbCbhAAACBbZPBPGU2qKAq9ja9t2zzP1Xt8/ljKX33Fq38AACC+GEFV27aeXurZY5NfURR93ztH7ySKAgAAKxeppqrrOiNa0iMh6UQ1FBvZb/N1Xafe+POnAgAAxBGjT5VUSm23W5lJJtNG7My0WMquo5IldV1vNpuyLGVNtQW1jicVAAAgjnijchmv3anxOdXrfjaVN2MdY9Qrf6qRBzqqAwAQR2qP3ahHqyqlpr3Tp9dF7ZsqUju7AAAsKLXHbmJHm9jZBQBgQak9dhcYUgEAAOD0EFQBAAAEEHtE9cX5p6lJqpYSAAAElFxQRdgEAADmQPMfAABAAARVAAAAARBUAQAABEBQBQAAEABBFQAAQAAEVQAAAAEQVAEAAARAUAUAABAAQRUAAEAABFUAAAABJDdNDXP/AQCAOSQXVBE2AQCAOdD8BwAAEABBFQAAQAAEVQAAAAEQVAEAAASQXEd1AABgvAvPW1xBEFQBAJAciaLyPCecCojmPwAAgAAIqgAAAAIgqAIAAAgguT5VTFMDAADmkFxQRdgEAADmQPMfAABAAARVAAAAARBUAQAABEBQBQAAEEDUoKqu67Isy7Ks6zpmKgAAwNwijU/ftu1ms8myrCiKLMu6rsuevohXlqUsFEVRtG0bJFXHePwAAChzPxZTe+xGqqmSiKrv+7Zt27ZtmibLMlWrVNd113VVVfV93/d9VVVd16nA6JBUAACAOCKFkHmeV1Wlt83lea6qlGRATj0noVLtbCQVMgMA4EFNVViRBv9smqYsS2OhvkSaBfWPRove5FQAAIAIIjX/qfhJmv+kekmvuLJDLufXJ3wEAACIIPY0NdK5KsuyqqrkP852OtX3PHjvKP/cfx5JVWACAIB9xQ6qVF/17XabPQ6F4Fk/eLUTsREAANkBtQwYssDgnzKaVFEUElc5+SuoDkkFAAAKoVVAMYKqtm39vdQzoigAACIilppDpJqqruuMsc71SMh+X0+GnvKkqjf+/KkAAMCPACuUGEGVVEptt1t7xE71MdPqruQ/+tCgk1MBAICBKGom8UblMk6hMT6nmsdGGONaHZJq5IGO6gCAxDmDqjmej6k9dqMerYqihoIeWWGOVJHa2QUAwDBUTUVQdbjEjjaxs4uRjFsMFwmARDBNTVixx6kCVkh+86n9+AEAYRFUAQBwaqiAX0RyQZX/lQcuOwB74dGFdaICfhHJBVVcXgAC4tEFQFlgmhoAAIDTQ1AFAAAQAEEVAABAAARVAAAAARBUAQAABEBQBQAAEABBFQAAQAAEVQAAAAEQVAEAAASQ3IjqTFMDAADmkFxQRdgEAADmkFxQBQAAVLuN/IcahyAIqoAs29UuDAAnhihqDnRUB35HaAXgBNzd3b1+/fri4iLLsouLix9++OHu7m7pTCWBoAoglgJwOu7u7l69enV/f397e3t/f980zadPn66uroirIsiTqgDM87SOFyMZQRUXCUZStxTuLViP169f39/fv3nz5vnz57Lk4eHh+++/Pzs7+/HHHyNnJrWfRmJHm9jZxRh2NRUXCcaTuwpdfbEeFxcXt7e33377rb7wl19+Kcvy48ePkTOT2mOX5j/ARGsgxiOiwqp8+vTp48ePl5eXxvLLy0tJWiRX6SCoQtKInzBZ/sj4yEWFBZ2dnZ2fn3/48MFY/v79e0laIlMJIahC0nqN/nHpfOEI6JeNvoTrB8u6vr6+ubl5eHhQSx4eHt69e3d9fb1grhKRVmNnao272AuXB/ZFR3Ws0N3d3dXV1cuXL6+uri4vL9+/f//u3btff/315ubm66+/jpyZ1H4aydVU5V5L5w7AMUnqaYFj8eLFi5ubm7Ozs7Isz8/Py7L86quvFomoEpRWCJlayIy9cHlgMi4erNPiV+biGYiMaWpw+hiGCgAQAUEVTh+9XgAAESTXpwoAAGAOUYOquq7LsizLsq7rmKkAAABzi9cgIv1aiqLIsqzruizLmqYpy3IoVc9YWZayUBRF0bbtyFQjDzQAJWvn2efywGRcPFinxa/MxTMQWaSaKgmemqZp27ZtWynizWbjTG2aJssyVedU13XXdVVVyah6VVV1XafCJn8qAAAniSGBVihSCJnnuVGBVNf1drtVPYizp1VT+vqHpNrZSCpkho6aKsyHiwdL8V97i1+Zi2cgskg1VUVRGL2d9KBHWv38Xzc+Gu19nlQAAIAIIgVVbduq7lNCj3sk3srzvK7ruq6l8snoNaV/d6+PAAAAESwwTlXbttKbSvpOZVlWlqVUL223W1lSVZVaOezeJzc8r7kCk8EtAQBYXOxxqsqyVBGVqlKS1/f0zubb7VZSg1c79VOFzUZYKofrzyoAAKcqXlDVtm2e5yp40qMlWai/7ufvF+WvvuLVPwAAEF+8PlWbzaYoir7vnT3Wnf2iVHhEFIVZqReSeTMZADBZpKBKmvyc0Y8RPwk90rJrrbquU2/8+VOBMY6oqRcAsFoxOqrrA3UaSbKkKArpoi4fZTxPFRjVdb3ZbMqy1CMtva3QkwoAOFK8gnMIVXr2aI6YT4xRudTrfja1d/9UM8YW9E7uO1N1pz0K2Wkf3WR3d3c//fTT27dvP378eH5+fn19/Ze//OXFixdL5wsnwv/gJyw4HHc2v5WXz8qzF9y6jtbZvypIqjjts3vaRzfN3d3dq1evvvnmm6urq8vLyw8fPtzc3Mi/xFWIht/mISg9v5WXz8qzF1xiR3vSZ/e0j26a169f39/fv3nz5vnz57Lk4eHh+++/Pzs7+/HHH5fNG9LBb/MQlJ7fystn5dkLLrGjPemze9pHN83FxcXt7e23336rL/zll1/Ksvz48eNSuUJq+G0egtLzW3n5rDx7wcUe/BOI5tOnTx8/fry8vDSWX15eStIiuQIAnCqCKpyss7Oz8/PzDx8+GMvfv38vSUtkCgBwspILqnKvpXOHwK6vr29ubh4eHtSSh4eHd+/eXV9fL5grAII78GQU1zql1dh52o27p31009zd3V1dXb18+VLe/nv//v27d+9+/fXXm5ubr7/+euncIRX8Nv385UPpDZGSWfkwVKmdvuRqqpCUFy9e3NzcnJ2dlWV5fn5eluVXX31FRAXgBKw/okpQWiHkaYfMp310h6N8sBSuPT9qqvblaftbW1mldvpiTFODCGhfB4BESJhyFM1/qaH576QQWgFACvQaICKq9SCoOgXEUgCQFAKpdSKoOjUEWAAALIKg6ugRRQHrdHd39/r164uLiyzLLi4ufvjhh7u7u6UzdTQoPRwjgqoTRJgFLO7u7u7Vq1f39/e3t7f39/dN03z69Onq6orIYAxKD0cqrXcdT+/dzqH46cQOM4jTO/tYs9evX9/f37958+b58+ey5OHh4fvvvz87O/vxxx+Xzdva2L9NSm+8ld/ZVp694BI72l1VOEuVhpGxadlI7drdF+WDmC4uLm5vb7/99lt94S+//FKWJZN5G+zfJqU33srvbCvPXnDJNf/1XsvmSv8PgOP16dOnjx8/Xl5eGssvLy8laZFcHQtKD8cruaAKAOZ2dnZ2fn7+4cMHY/n79+8laYlMHQ1KD8eLoAoAwru+vr65uXl4eFBLHh4e3r17d319vWCujgWlhyOVVmPnyht3D8zeyo9ucZQPYrq7u7u6unr58uXV1dXl5eX79+/fvXv366+/Mp+3zf5tUnrjrfzOtvLsBUdNFQCE9+LFi5ubm7Ozs7Isz8/Py7L86quviAl0njeHKD0cqbRCyJWHzNRUzSTIy5XAZPw2h4yZEpjS81t5+aw8e8F9sXQGgNkl9ZMGjsjOiAo4LgRVx4G6FuAkJfvTNg5c/5hOIeD0EFQdB7nLpFaPCpy8ZH/a+oFTWYWTQVAFAFiAHkoSUU2gqveISteDoAoAsACCgANRgCuUXFDln/6PaxQAAEyTXFBF2AQAAObA4J8AAAABEFQtz98iOXILshH1HwAAEFlyzX8rJG8UH/JmNW2aAAAsLmpNVV3XZVmWZVnXtZ3atq0n1f9df+r6MVILAADHLt6IcxI0FEWRZVnXdVmWNU1TlqWk1nW93W5lBUnVM1aWpSwURVG0bTsy1cjD2qIWT2udndUV5h/ATjt/uYn/tP2Hn3jhHLvUTl+kmioJnpqmadu2bVsp4s1mI6lt226326qq+r5v27ZpGvWVLMvquu66TlL7vq+qqus6FTb5U9dPsp1pIVT/aNF8AcDC6C2KoxMphMzz3KhAkqop2btUNek5kVY8+dduF9O35k+1s7HCYGV8h6p15h/AkJFT+yX+00788E9baic3Ukf1oiiM3k560NN1nTQLKsbKRqpqIhyTun5JXXBzSHZKWqwfVyOQlEjNf9IJXV9ixD3Sx7wsyzzPy7I06pmM7+718QTc3d29fv364uIiy7KLi4sffvjh7u5u6UytiN6EyjMMALCUBcapattWqhak75TET9vtVjqqS6eozWYjy4P3jsqnCpuN8e7u7l69enV/f397e3t/f980zadPn66uroirgNNAbyHgZMQOqsqylP7p+qt/Qnqp13Wtd2MPXu3UTxU2G+P99NNP33zzzZs3b7799tvz8/M//vGPb968efny5T/+8Y+lsgQgOEIr4ATEC6qkgkq9qaeiJfmP0S+qqir/pianHp23b99eXV09f/5cLfnyyy+vrq7evn27YK4ABEEsBZySeH2qNptNURR9348fn1OFR0lFUbpPnz59/Pjx8vLSWH55eSlJi+QKwBwSDLAYNAEnJlJQJW15Q9GP/b6erKkqsYxU/W1Bf+qxOzs7Oz8///Dhg7H8/fv3krREpgCEQRixkl4WQCgxhlTQB+o0ktR4VJvNRr30p8bzVOvoqRJpqU35U0/A9fX1zc3NmzdvvvzyS1ny8PDw7t276+vrZTMGILjUBvUBTkyMH7C0/TmT1N7VNDXCGL3T2ILRyd2fqlv5DcuZvbu7u6urq5cvX15dXV1eXr5///7du3e//vrrzc3N119/vUg+12nlJxcwDFVTcRnjlKR2Z17X0eq1TWFTxcrP7lD27u7u/vGPf7x9+/bTp09nZ2fX19d/+ctfiKgMKz+5gAdXL05Vatd2Yke77rPLNDWHoHBmwoD1EXD14lSldm1HmqYGwJEaOTElAGCBEdUBAABOT3I1Vf53mPlbHAAATJNcUEXYBOjoMgUAoSQXVK2TerDJf3iwjUdfnwPRZQoAQiGoWgWeZ5P1fa8CAiIDAMCC6KiOoydxVUZsCgBYFDVVOGJGfyD9IwEWACAygiocMb3Vj8oqAMCyaP7DcdP7URFRAQAWRFCF40YgBQBYCYIqAACAAAiqjkOe59JnSP0HiODu7u7169cXFxdZll1cXPzwww93d3dLZ+qk8NMGTglB1XHon1o6O0jC3d3dq1ev7u/vb29v7+/vm6b59OnT1dUVcVVA/LSBU5LWYIk7/xBMqjRODCN/HsJZeq9fv76/v3/z5s3z589lycPDw/fff392dvbjjz9GzyOA45PanTmxo03s7CaFk3sIZ+ldXFzc3t5+++23+sJffvmlLMuPHz9GzB2AY5XanZnmPwAOnz59+vjx4+XlpbH88vJSkhbJzXq+7AAAC6lJREFUFQCsGUEVkC5Pg/jZ2dn5+fmHDx+M5e/fv5ekeXMGAEeIoApIlxqJPnMFWNfX1zc3Nw8PD2rJw8PDu3fvrq+v42URAI5HWo2dqTXuJoWTO9nQDD93d3dXV1cvX768urq6vLx8//79u3fvfv3115ubm6+//nqJnAI4MqndmampAhKVPzI+qiUvXry4ubk5Ozsry/L8/Lwsy6+++oqICgCGpBVCphYyJ4WTO83I6agpXgATpHbroKYKSBfTUQNAQARVQLoIpAAgoC+WzgBwKOP9NQIFAMAikguq/DPV8Dw+Rpw1AMAaJBdU8QAGAABzoE8VAABAAARVAAAAARBUAQAABEBQBQAAEEDUoKqu67Isy7Ks69qzmnMF/3dHbhkAAGAm8caPl7EMiqLIsqzruizLmqYpy9JYrSzLruuqqtLDI1moPhZF0bbtyFQjD7z9B9h2/jT47QCYILVbR6SaKgmemqZp27ZtWynizWZjrNa2rR4eibquJczq+77v+6qquq5TYZM/FQAAII5IIWSe50YFUl3X2+3W2LusZtRU2cNk61vzp9rZSCpkBkby/DSM8XL5BQEYL7XHbqTBP4uiMHo72UGP1Ga1bWsPei6NhvpHo73PkwrgEEndEAHgEJGa/9q2NbpPGXGPtOIN3b6N7+71EQAAIIIFhlRQdVFN06gl2+1WfTRWDrv3fKqw2QBWQl3eXOcAcKDYc/+pN/X0V/82m01RFM4apuDVTrRlADp+EQAQSrygqm1bed3PGC5B/m8MMdW2rRp6yrkp/44CZBcAAGAfkYIqiag8I0htt1v9Y9d1XdepiIooCgAArFy8IRWy0Q0NeZ7rtVnSYjg0aII/1d4yjR0AAMSR2mM3Rk2VPlCnkTRmVpm6rjebTVmWKorSv+hPBQAAiCNqR3WjjS8bF/2UZdk0zWazUa8m6Z3c/akAAABxHFO9nF4XtW+qSK0eEgCABaX22E3saBM7uwAALCi1x+4Cg38CAACcHoIqAACAAGKPqL44/0QcSdVSAgCAgJILqgibAADAHGj+AwAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAiAoAoAACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAggOSmqWHuPwAAMIfkgirCJgAAMAea/wAAAAIgqAIAAAiAoAoAACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAkhuRHWmqQEAAHNILqgibAIAAHOg+Q8AACAAgioAAIAACKoAAAACIKgCAAAIIGpQVdd1WZZlWdZ1HTMVAABgbnm0t+FkLIOiKLIs67ouy7KmacqyHJNalqUsFEVRtG2rPvpTjTzw9h8AAHGk9tiNVFMl4VHTNG3btm0rRbzZbMak1nXddV1VVX3f931fVVXXdSps8qcCAADEESmEzPPcqECq63q73cred6ZmT8eX0tf3p9rZSCpkBgBgQak9diPVVBVFYfR20oOeoihUS9/Q142PRnufJxUAACCCSCOq2/VGetxjp263W/2jEXIZnaj8qQAAABEsMKRC27bSZtc0zc7U4L2j8qnCZgMAAJyY2HP/qWok/eU+T6q/WXCCpBp3AQBANPFqqqQKSr2pZ0RL/lR7U5NTAQAA5hCvT9Vmsxl6Kc+fmhFFAQCA1Ys3pEI23PTmT5U2waFBE/yp9o5o/gMAII7UHrsxaqr0gTqNpLqu/any72azKctSRVH6yv5UAACAOGKEkNK650zq+96f6tyC0cndn6pLLWQGAGBBqT12j+lo9bqofVNFamcXAIAFpfbYTexoEzu7AAAsKLXH7gKDfwIAAJwegioAAIAAYo+ovjj/hDNJ1VICAICAkguqCJsAAMAcaP4DAAAIgKAKAAAgAIIqAACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAgguRHVmaYGAADMIbmgirAJAADMgeY/AACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAiAoAoAACAAgiqM5R+MHn6U3oEowENQeoeg9DAeQRUAAEAAyU1Tw9x/AABgDskFVYRNAABgDjT/AQAABEBQBQAAEABBFQAAQABRg6q6rsuyLMuyruuYqQAAAHPLo3XcltfuiqLIsqzruizLmqYpy1JSy7KUhaIoirZt1cdDUo080FF9MkrvEJTegSjAQ1B6h6D0DpFa6UWqqZLgqWmatm3btpUi3mw2klrXddd1VVX1fd/3fVVVXdepwOiQVAAAgDgihZB5nhsVSHVdb7db2btUYuk50dc/JNXORlIhc1iU3iEovQNRgIeg9A5B6R0itdKLVFNVFIXR28kIeqRZUP9otOhNTgUAAIggUlDVtq3qPiWMuMdINRipe30EAACIYIER1du2ld5UTdNkVpWVUH3Pg/eOYmrMQ1B6h6D0DkQBHoLSOwSlh5FiB1UqWlKv/u1VR3WgpFp2AQBATPHGqWrbNs9z9aaeP1ryV1AdkgoAADCHSDVV0uTnGUGKKAoAABy1eEMqZMOtb9ImaAyLUFWVvDDoTFXxmT8VAAAgjhhBleqZXlWVkSRhk1GPZcRJzlTVJcufCgAAEEfUoMpmRE5quREVHZIKAAAQwbqGOlW1TcFTj5eaK1o+2iN+YYhedKopmdIbybjwZAkTlo83NMs7ZTgGP94Jdj4spCS5AmcV7+2/Mfw/m0NSj5fROWyz2dBdbCRVUHmey/83mw03lJGMy0zmlVooL0dJFeB2u+U3uy/nj5di9PM/LKQk5TX8yBlLyrqCqpMkl3Vd1+r6lv/I9S1L7DqAtm31r2RZVpZlaj+GIEWnvlLXdVVV6dyXQ114sjDBiCpsAdaPZs3zegQpPdXptq7rpmmGupGcjPkeFnVdS8/jtm3tWeMQEEHVvNq23Ww2cn2rahKppFU3CPsPiLquJXWz2aj5fNq2XVVb7dxCFV1ZlkmVmwh44clHmf8gHWELUL6YVEAfqvTU1K4n2RChm/Vh0WpNgUldigvoMSd5FKmP8v+iKIqikCVFUchoqH3fV1Uly42vNE0z9PGEhS06FRPMn/FVCFh6RVHIf9IpvT5cAcp2qqqSd5/V109b2NJTq532FTjrw8K4GSZyHS6CmqrY2qcd6ruu0/+A6LrO+BtC/aGGQ4quLEu5lZz837tDppWe6iMcKZcrNq0ApaJUmmz6vjcmkk/H5NKrqirP89R6PoiADwueI9EQVC1s57We7F14p5FFp7ojSKU35SnGl17XdeqppnoNg1/uIcaXngSjEpjOnq11O/CSM/ppYSYEVTGozoaZ9Ue//sa19KQuHwc1jZjB9Tq86FptGLPUSjVI6alq7SzL+l2zdp6YwwuwrmtVy8Llp4y/7+mdh1KobpnpYaH3ZN9ut0SoM1qizTEhqiuP3A6kRVxvGu+fdieUJar7hfybcp+qIEWnr0bpqXX2uvD6JPtUBSlA9ZHL78DSi38gMc39sJDuWVkyHfuWcuKX6eJUl0D/zVS6ZNoLZ8rVUQhbdM7VThgX3oGCF2BSpcqPd4IIv9kUinFxX2SIwt9o4kxNqp3FI0jRJVuYXHgHClWAaZYqP94J5vvNplaSi6BP1ey4jiej6A5B6R2IAjwEpTcBhXYC1jX3HwAAwJGipgoAACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAiAoAoAACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAiAoAoAACAAgioAAIAACKoAAAACIKgCAAAIgKAKAAAgAIIqAACAAAiqAAAAAiCoAgAACICgCgAAIACCKgAAgAAIqgAAAAIgqAIAAAiAoAoAACAAgioAAIAA/j8CxF73AJvO0gAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "l.Clear()\n",
    "g['LY_MS3_DRY'].Print()\n",
    "g['LY_MS3_DRY'].Draw(\"AP*\")\n",
    "l.AddEntry(g['LY_MS3_DRY'],\"MS3\",\"PL\")\n",
    "g['LY_MS3_DRY'].GetHistogram().GetYaxis().SetLimits(2000,4000)\n",
    "g['LY_MS3_DRY'].GetHistogram().GetYaxis().SetRangeUser(2000,4000)\n",
    "g['LY_MS3_DRY'].GetHistogram().GetYaxis().SetTitle(\"LO (ph/MeV)\")\n",
    "g['LY_OPT_DRY'].SetMarkerStyle(24)\n",
    "g['LY_OPT_DRY'].SetMarkerSize(1.2)\n",
    "l.AddEntry(g['LY_OPT_DRY'],\"OPT\",\"PL\")\n",
    "g['LY_OPT_DRY'].Draw(\"PSAME\")\n",
    "g['LY_OPT2_DRY'].SetMarkerStyle(22)\n",
    "g['LY_OPT2_DRY'].SetMarkerSize(1.2)\n",
    "l.AddEntry(g['LY_OPT2_DRY'],\"OPT2\",\"PL\")\n",
    "g['LY_OPT2_DRY'].Draw(\"PSAME\")\n",
    "l.Draw()\n",
    "tL.DrawLatexNDC(0.11,0.92,'CMS MTD - POST MS3 Optimisation - Array dry coupling')\n",
    "tL.DrawLatexNDC(0.15,0.83,'<MS3> = %.0f ph/MeV'%(R.TMath.Mean(6,g['LY_MS3_DRY'].GetY())))\n",
    "tL.DrawLatexNDC(0.15,0.79,'<OPT> = %.0f ph/MeV'%(R.TMath.Mean(6,g['LY_OPT_DRY'].GetY())))\n",
    "tL.DrawLatexNDC(0.15,0.75,'<OPT2> = %.0f ph/MeV'%(R.TMath.Mean(6,g['LY_OPT2_DRY'].GetY())*1.5))#multiply by 4.5 because 2 of 6 elements are empty\n",
    "c1.Draw()\n",
    "for ext in ['png','pdf']:\n",
    "    c1.SaveAs('plots/LYSO_OPT_LY_ARRAY_DRY_byProd.%s'%ext)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "id": "22570826-83dc-493b-a245-f5687b607a69",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x[0]=0, y[0]=130.778, ex[0]=0, ey[0]=4\n",
      "x[1]=1, y[1]=132.346, ex[1]=0, ey[1]=4\n",
      "x[2]=2, y[2]=138.021, ex[2]=0, ey[2]=4\n",
      "x[3]=3, y[3]=125.611, ex[3]=0, ey[3]=4\n",
      "x[4]=4, y[4]=130.276, ex[4]=0, ey[4]=4\n",
      "x[5]=5, y[5]=124.082, ex[5]=0, ey[5]=4\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file plots/LYSO_OPT_sigmaT_ARRAY_DRY_byProd.png has been created\n",
      "Info in <TCanvas::Print>: pdf file plots/LYSO_OPT_sigmaT_ARRAY_DRY_byProd.pdf has been created\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxwAAAI8CAIAAAD0vjrdAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElEQVR4nO3dT27bTJ7GcTJ4N2+AwPY+acQ36JyApaxmDjKOu2f5NjCdbEhtJrPJcjqAM1dpkbzAi5wgeaeh9dgGZuBdNItfu7pcVSyRUvGf+P0sAoeiqGKJEh9VFYvpbrdLAAAAcJxnYxcAAADgFBCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiOCnsQswRVVVVVUlfyullFLuCuZ/wyt4Hw1v31ukwCu6j1pLmsh2vCu3KVV7TeUJv0qnijJXTpKkKIo2BWhTnq6V375U7gsdUO0HHE5dt683Pp1NTeq1DhP+lpiUwz6/J8w9umTJYisE/7DDU1mWubWU57m5Tp7n5qNlWQYeNZ9rPdS0fUtZlubKWZbtfdS7F94XbSqSd9cOFi6J+yrWTjXte/uV21fI3sK3fLSpVO4WyrK0ihc+Hg6rpYPpsh15MOjSxi1ek6bang7zTY/1QevJAZ+aE+YeyfpbdIG1AQvdf/9QVVWapnVduw+t1+s0TZueaDU/rNfrptWaHlqv1+E2DJNVwpZtMIdZrVa9bl+/ivnfqqqsJVpd19Yb0Wnl6AKV01SqJEmsUq1WK+s9Xa/Xe3/y9rHjRVGkaZqm6bx+cCulpNgDHKsRmW96+4//1LQ5VoHloPvvH8xTVJ7numtMJ6GiKLzffd4c5tLbyfNcb0cpJU/vlKtM3hOJuSlzFwLNJ1mW6WeZT1mtVrt9TU3tNb2K2WW5940wV9Zby7LMXGi9awdUSBOzk65Nt0hZlnoX9K7pXdBryp7qstV1He4N7FRLxzCP1WO2o5Q6oLYPNuRrHcD6sLf8DhmdeaRVVSXF3nusnryiKGT3l1wJEISqvzNPP/osmDwOUpETWJvc03Qa08vNVJE8nqHlu6n9F5O5pvfr2NyOUso80QY2a+61NFrIf5vS5AHMQT/yh84QbgnNMGe+EeaXuH6iO3RJtuyOdWhfIV5tQpUulXUslWWpd8FaU8osSWjv8dC1lqwhIHJGPHgAlrk1PaLL3FrX7Xs34l0hsE57gZfzVlSUF7W4jdbu220VRs7cZlY2C2l9q3iLXTkjgcLLXVbC1r9ezJ9ebrGLotAZ12pebVpuafNGhNfx7qO70Htst//2C3/iApsyPzLt3w5M0bi9j9Ohxzd4B3xY1aV70N3hJu4QJellN383WxsvHzWVzXyuvKJZSKskbcpvMvel/UNdWbXh3TVrTW+FuLsZ2LVAxXY6/ttXfpu9sErl1rCu9sAhcXAtue037jatjbsHuVUhWp7nge17x1R5h/RZO+UdD+cWxnpFtwBtXs4cHGOtFnG4jFlLgU+ZWUKzDOY7Yj3Xu4O65E3fEt66armCVY1mkfTf1v7q5242m6Z9N7kHQOBrpGkd7y5YC83CB4bGthlTpcvsviPW4e0Oim3zDYDJIlT9XdMR7+WGKuuzZH4L6Iesj2j7D4z3K9h6KByAvF8o7r50emIn3u+m3dNzwC4YPUUghHWq1U771b7y9UPWAdDpy3Fv2Y6sJZdbZrE3VLWhjyv3VBToodOvFbjCwFuYQNreu6ndUNdtWOdgtxhN9WyFKmuv99an98ixEptXUwl3T9NDUxnKsvR+z7QZ3x3eo/A67u+fpF2o8tKF7BSqvNzK37uPmAsGqtu6trjKh0c32ErHjXcj5sdsvV6vVisZXVsURVMv0t7itX/iBFXGACCrGpu+ktyKNb+VzFpVkUYUeV9dthxopTcb+a1SBdr/lVK6yzXwbXtALZnck58Uqaoq68y397PgPZfLmUBvKjBaSFeF29BrfaASX4uarLN7Gvua3vTqsVM18SU27566ITLWEWUO+EuM9zGwfWnAcMspy2U7e+vT2/MVa8ycV57nUmz9KubxYNWDy1wue+R2F+piu2FdesAPK7kbmJouM9rLTbRu5euXa3OdMiYtWjybuU4VYp57rNYL+bupuTvwOzjQAG79vtQvof+Wb67AdgJ7d2RLVflU02otj8M2HY76KfrlSmdKAnfLXferaWVd4TujVbKp6ShQKreuvNHE64Ba0v81n+IW27vlQEuVWyTzaLdWc3/fm6Uy303zWJJj29ss4WYj72rWOuFeIe+OtKnwlh+EXetqN8tmlbmpF08ahMIld3fQfTddgXWslwg0o+pibzYbc7OBWvW+s3phoGwtd9Na2FR4vdCKRy27/9yStywVLVVzREvVscxr2fYOMCyKQuo9z3PrjFvXdafB4OYvsLHGMyqlVk8dtp1wG3jg1fUfVVXtHk8qVsXGnVXBbFELtEqKQKnculJKmavFmszCKp55jI0+DNZszZLGPPV4CZUumwxz1s14x7+b1l63aSXaW1GdPgj6LdC7723F8T4lvFwahPS/aZq6jSt6X8zLOJKnXdidtB9VbTX5uPXg8l6/Ih+r3W6XNLfaugPPj+E2oLZnlrypn8EqPI1Vs0aosnk/NtWj8BNlhTYfCTlPyFeDXr9l87I+r/d9kUgffYvStqfJTzEzHskfTWeXcJHkdCIVe/BJIkw9veQwaXfG1aXydgFospouedO5+chamo6iKNyfFpKuzP4RCQdHTjrgPT2PQpdEZniyLlU75jLbqqrSNF2tVoHqajqvH395b1Na1fTxrC9utZZHtPczEjZurBn9EMUxCFV/pz9F3i+X8A9QK+V4PxLu5bLmQ52KKiU8ZrhApxdKgmN3sqf2blMShmbVlXfMR7hIUqvuu2ZdYb63YC2ZLSjWEms1b6mkxsxS6ZYYvc7e09sBtTRZOms2teSZs7tZP0I6adMidbD2HwTr1etHTSt0oitNd6d624DNFrKDxwlpndKq2Z0nT+zpEN077jBs3GnD5vKjCH5D9TNOnfntY/Vku9d2WcMIrEYRaxRL7twNJvDSe8smS8yXs8oQa0xVoEIOYNVGQMuKCozvcV+000N7V06e2gXfILfqrGEi3pLsLZ55qrAeMivQrSXvmCp3aMgwY6rcEUhWTXrL4+5d32Oq3JIfrM05vqls7kasoyvwvjSN1vLWiWtvSTJnkLi3rsw5FOQPPb7K6wTGVJkl71oqxlTNES1Vf2c2IaxWK2k8qKrK/DHX1KMUaG7RrLEs+rdIZcyyfViPVay24rqudWeEdCLI8izLhmyONitKXxop7XzhIpl9RlWf8zub58Wmt8w6lsxSNf2y15fg6WFDLRs8rFoyL6h0K0EP3aueTu/e9EK9aj8Uz5oi9QD6iev12u03H6xJz72YUQs3lndSPV6d2rIhqv0rKoN5U6+WW7B6AJN9fX/6DdKrVcZVnMJstNbrBGb3bXm1ozn1sVue6NwhEJir4XPclAUqyvzh5f7+c6tU/1f/KNmbmZpK5f5+Ms8B5dMr2A9uqer1CHFrIyA8bt3dwXD5vb/2Ou2dtbLbFOT9gduyVIGdjVhL4cI0bbPcN0+VfuJhLVWBKJM3T3pkFc/ajvuR1OUJvFx4R2K1VAVaYb2PNq28t6XKZZXc+tTvLXlgy13ryvtmHfDS3lbJQMHCudkqfPgVI7ZUtXxFzAgtVU/sGj57eZ6H+7nNZvCmdQrjXg3u03dd5lQ0f83098smfzoz4WCUUu4gG10k941oWllGlkSvn5aV33QsJc69a7y5du8VkV1rqWmz5ls88K/kqqqayl883rTHqhy5vkH+1q0XZltCYDSMeRGA1vWjdwz9pnjf8eMbqNw3Vy+xqsV8rWMu6ZDa61Ry8+XavLT3c2R9rpuGjoUHVgZePfPNqN7Hp8M9wlsOTsVkpaOcNafPHHUe/bOkG+d72v4p6fRG9PquHcy6gqGpYMcUfu9zzTlF1eNNzZpWDj/ah3D59aPhk7fOYS2Pk/DbMV8t3z59SAx/CjjgpfV3ZmC/Wq4TPvIlqWdZVj0a5oNgvpZ6HCTQx29C9I1QBZw+K1SNWhZMghwSkh6GfN2yLN++fTvKS+9lhaoBXlF/GK1RkvIHZ+c5ovsPAJaieno3pCETtgR6SVTJcd2OJ0PPqSGXm8isbPIQnYAz9dPYBQAADEePr8qy7PiBXJ1e13zpPub8nJ3dbqdTlHmp5gSb8dASoQo4fXroOn1/kLk2hh93KEOwlVKTTVTqcWKdIWtmt9vpUYNiyKSL6BhTBQAAEAFjqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEy5qnSk+zBgAABrComZuWFaqShb27AACMaGltGZ1DlZ77Vd/mnelfAQAA2s6oXlVVURT6zk2uPM+nn67SlBnkAQAYyNJOu60Gqiul5FZNeZ7vHGVZ5nm+Xq/lPtv9lhcAAGCS9kfINE3bt0JJh+Bkb6+9tMgMAMCIlnbaXdjeLuzdBQBgREs77R47T9VkG6UAAACG1C1UFUWRpqkEKfl7tVrpJbOQBo1dOgAAMFcd2uWqqpLh6vIUiSBlWcpVgbNo31taOyQAACNa2mm3Q0uVjFWX2pGmqbIslVKyfEaNVQAAANF1m/wzyzL5QyKUXOunr/iTPwAAgNeimm0WqENLlVJKT/65Xq/1cjNgAQAALFOHliql1Hq91uEpz/PEGGhFqAIAoG97GzL0aBx3naqqdONIlmWcuONzZ0gPkCCVJEmWZbJE/uudaX2Cuu4vAAAR/TiOPgvned60jj6/mws3m403A2w2G+vpm81GD/XJssxdoZOlnXaPHZY/r6FUS7sMAQAwKUeeg9brtR5+Y+Ynraqqt2/fWivohXmeSwOVNFnJpjabjT6Pm0/Xsiwry/KwAj979mxRp91uA9VFVVW6+XFGiQoAgJPhbdQwRzxbC83wpE/fktL0cklU5prPnj2r63peDSgj6jb5Z1VVMuGnvA0y8ycVDQDAYKQTUI+OMtV1rTvvzIWJb4iVtR1pLrHGWkm/ofe14OoWqmRMep7nZVnudruyLLMsq+uaXAUAwJDcRilZosddtWEOt5Lk5H06U1G21HnyT5lCXc9QVVVVnufSNthLAQEAwFP6AnxzYdOFgdJ29ezZM2/noHlR/48fP6yny1NoOmmpQ6hqerckbAEAgGM8e9b2pCyhygpJ3r6/JEn0MPP1ev3s2TNJV+GmkKqqZJCPbLNT69eSdZv807ucNioAAA7WPktZzKFO4b4/a66E9Xr99u3bprarJEnevn27Xq8DHYLw6tz957ZLrVYr5hADAOAwP3780Lmqa2OVbtfYOymoUqosS0lXOidJc5S3SHo1CVgtS7VwHeZtqqqqKArJrTrt6v+ab+RkOwSZpwoAMKLAOUjilHf2KU0uvc/zXOLOs2fP9CRS1t97N5UYs1KZcygcVrDAcxd12u0wT5U5vb11dWVd1+aSyYaq6TP7Ur13GAg/at6dINx2eJKTjpzkTgE4SU0tUt7lgUAjJ99A35/kMG9sUkrleS7dfHLlWcKY9CONNZX7KCZeG9bnwbr5T+DRpg5v9/ZB7poyO8bopGDe+x1Z4y71LZIC68zlvkkAliZwOxednAJ3fdFfleZ/9WApc2v6vzJjQpZlezdobcS7ta5GP7EObGF7O+1310o8Vnpoyg3hIYTmRprWHD1X6StT3DDUtF971yFXAZigpuRh/dHEClXyFDk7mLHJikHyX+9d/Mw13Y3/MGaxChesSTLt0250+wfEpWnavjtPT2GF47l9rN7V9PhBb9uVOYWYXrMsS5m4Vf7rHaU4ABmlp5RqKoA+8GSggHnzKf2QtY45+rKXQgNAbD8OGq6k7b1GT4LR27dv5YYoMlZEplcwn6inaTDX0YOujinhguyNXfpMJict7zpyj8akoWtmOtrs74jcj4T5qPWQjlBNb6XOTLKmtzWo/WHQB/dotJqX3OLpKtJHmruO3vFAC1z5SP6WUZ+BlWWF0Zv0AMzd3nad8ApuY5L+VrS2Yy1pSkVN7VImbxNXS2OdX8ayf6C6Umq320mjQqBJQ05ye7eGNuTmP4kx8lq3NumHRGCSMHPcumk6TYnWNaTeFZrmsrP2wlxHKWXVnqWqKjmS5apV3aYl19SY7bJFUZgtXnooKJdiAOjDj31NVvq6v8AS73aUUj9+/DAvOPPOhWStxnxJnXVNYfKrPXs0r5/vB+xv36QbTv42m2Gslj/9X6v9afc0yHZ6O8w2SP0qowxFMhufw2tarVDuLjQtNJndiC5dhsBqDNgCcJiDm3xmaoKn3V51nsVVKVUURfWIQVSHkUaUNE3lJgDuClKr1hwW3qo2f6PING5pmspIuKZ2rKIo5KWTJMmyzFxtvV6naSrvcqd9aRKrUUeqSw8dsKrC/G+nA9JNtLppSleLG7O4iwAAwHXg1Pg4mA40OjB5Rxd6k4F3YVEUTXOTrFarpv4v/bd3grG6rnW62hsgpKG4Saz84Y2e4Y039Sp6t2CGP1lozlOsp2+RkVh0/wEAXISqgUiWStNUN4Togf/eM7QOQ9aUnk0b100pVsCSKd3c9c3V0jTV27eabeq6lqav0WOEdWWflCduK6nevlS4GTGlEiRi7p1YFQCwTISqgbhZquX0303jzb2ku02CUWAYuF5Nr2NlJnnUm2Nc+hW9YrVUSZnNe62bj5qvYg7qD28zvILMNWwu0Q14dP8BAFxTDFXeLqfU4eaAuIN4etK+nUNO+TpUeROAbjhxa8xaIgPgrMox28Oairo3mvRHChy+Oc8x9vYP6nBpVcJYM3sBAKasw73/hmFe8GkuDD9LX0WfPE53ObW2BD0VgsyrlmXZ3jH+MoeF3q+mlWWF1Wq1C05pIdO46S1bj1ob19O+BdYxVw4E2WNirhRgvV437ZoukndkWNfX1furnk5joTOrnosh4T6DAABXvAsJj2V2NlnzAkinT2Dq0cS5bYt35dH3V8+SqplTqgZmttR7FJhSwdyUOR2AbM1cYs1HYFaXO+Z97zQN7e+Ts3cL1lQFemetqSXMQ9dax9ypppez5kpw98KqWHP3924cAALGnuJgaEv7tozTUiVtFVY3TdfZFgJdKvraK++j0sBgtk/IAKCpNVYlTiGTx0HQuxbzpjbtflmWUnWyKXcFXXu6tcxaTU/vZjb4tWlL65ueb7auaz2aXpizG8hD1jrhqGdq2nKe5/IeubU6YpcoAGhmf0LTRJ0tZ4Fuc8akeX6/43NZ3C17G6XkHLZ7vLuIWwDvvYe9RT2sVD3Rg3Xkv25LVaCpxmzUaZqm0m0lctOAuY6e0LWXvW0WmPzTu2vehszwjnvXl501n2gVwJvMJn4vJgBTFqsFqOl3o3tLGe9q+htP1mlzaz/zhs3tJRM77fYt3cW4t4yeCDTKplarVVmW1nSO1kArc8rKNE29dxdxdy1N4+zvNJmzK4WHw8v7NaOpAcyZogJDu9rse/L0NjXVozZbTrqP0wIAU5RzkNwIOUkSPaREBhPLiXKz2ZjfZtZdk4W7sjWIVt+VSy857JY1z549O+HTrseokc7D21IlRdWtCHIMmfcJ9jYwuBs/pYrCwfbexwYAenJkA9WPHz/0Xcvch1reXDmwspCGq8OapixLO3t2m1KhzfzafbAmyayqKssyydFdg/PBNRV7nwAA6Ea3MHmHRuibK7eZ9qX9wFO01y1U1XVtDrUebDC4m5ysi94tExyiDgBA2N7J86Q1ITAESqJSm5t0oQ8HXv2nL5iSN3ishpzw3JWAl57RdC5DygCctu12+/nz58+fP9/d3Z2fn19fX//xj398+fKlu6akpfB3l1zovXcuPc6bfThkRvU0TeWyKT2zVK8nJ7lm3hogbB4NeqYAra5rLnpHQMRLKwDgGNvt9vr6+vb29q9//ev//M//bDab+/v7d+/ebbfbwzYoZ2TrtFg9tV6v3759m9AJGN0B44rc6967bicgMKWCNUGlLoY8xZp0YJqTfwIAlswdyv3+/fvr6+v//d//1Uv+7//+7w9/+MOf//xn7/jxxDe03KRPkXqoeJOmoegMVD9Y5+4/mRPSXdjrXTukscoceWfOoaCUkgkwdaekNSPDTJlX8jdd8N/Ufuvedmbvy51AjQHA7Hz+/Pmvf/3r8+fP9ZKff/75X/7lX1ar1b//+7+bax7zLe22SB02RQL26BTBkoar0Ltu5zDS2xi4ZYp3alDTMOU8XlN7rLeNsIn5TkU/EqKzWh9N7o19rGdNdqcAwGI24bRnPSvcOCTfmXoW0DZPsdBSdbBue+u9rV7T9OUTNItyhkeDWbFj70ex5Wrj1ox5ZxjroXBpCVUA5sWNHefn57/++qu18Ndffz07O3NX1tM0hnNM0jGHWQhVB+s2UF163Farlb61X5qmdV033SYFXRVFYd59r3ykk9Z6vfb25UkbnjCjhrxlpcEMbeby3vapkQwVV0o1TalidjRbhdQ3K8wc/RYaAKK6vr7+8uXLw8ODXvLw8PBf//Vf19fX7sry9W7Nfm7St4vooaRooWsKc8++4R63STlgfwema9VtszE/JO761rvg3kawzUMDc49Ga6/1Lpt7t/fQ9T7LoqPk7rFbOXDHQ73CjA51ANPktuX87W9/++d//uc//OEPv/766+3t7a+//vrHP/7xn/7pn/77v//b2/YTa0b1AFqqDnb4bY/neI6Z+LtrtjB5V3AjVFOoMrNv06uMHqrctqWmzs29C7VAT6K7zmE3VB7+htMAToY3efztb3/785//fHZ2liTJ2dnZv/3bvzUlKqG/NuVLbLPZmF9W1j2VCVVDWtjeTu/dlf44+Vt/TppO2+YHSZY0harDWqrk83nM7hxGF8l6dbM9SQuHqvCjerNuVHIzU2C12f2cADARe9NMy7Di/cnnjUGEqiEdMvknjiczUMg8EXoQlf6j6TLXwOWv5kArmdht71O81ut1mqYyZq7lU2RfmhwzwabeiF6it+YdMaDXbD9ETCdL/Q2lq05XqRuzmDUUQB/af2PLWHVpo5LGqh8/fni/+sxo1b4YTVvDHmOluVFMYX/dnxf6nK2XNDWEuP2DB7y/bRqxNBks33WPrC20r5NwO1nL7tG9r9jUN2rVv9kuqCvB23gGAC0d3/YzL01f16eKlqqBFEWRpmmaprohRJ+qD2j2aHllR9ffGUVRyGGh40td1zKr6uhtM0opXXU7X5rUv/AOLqrea2mj0hvUlaAvOWTGPACA68AbKqMrM0vJPALuOvomhk3T0wcmRrfmxj3yxF8URVEU+hZRuvzevCIrH/xabei58mXufu86eztPXeFsKjtlXros6So5lSn7AQBxLa6lKg0aoACBuKOXN81BonODG2KKp6Kc8qWo4853Incokr/zPG9KVOb4p/Ybt2446pKmO2tyr+RxJhgAAEyLC1Xh3tD+XtecvVM6ktx8YEYlNxVZDVHRS2iSaTmtcfRNL9rfQPXEiC/m3R5d5o0gD34ta3S/vp27vFm7p4Ox2txOEQCwKHT/DUTOwdKnVte17kiyegPLspTldV3rXGJezZd0HynVSVEUVjtZOM0kSSJ71EdhzITkdomaseaAvj+9vn5rrNfVeU539g2ZawEAs0OoGpQ+c+vsIulKN5JJd5tEBMlebsTp73SulDJvkhOrG/FgZlbrKbfVdW11+5rTh0rlu5193AICAOBaXPffRDQN1kmSpKqqpoFBZVn2PSRcz//ZNFjepa8Z9JpyN5l3RnXzKsKm6fWmvFMAgLGkvQ4kmpo0ndP+yoAe+ZvL+COqqkr3vZoDpwKDxvQbMfrUEgBmLe45yOzKsK4B17w/At012/xWPOA09OzZsxmddiM4ZpKr2Vna/sLLvPff2GUBsCyxJtVs6tBwb3QTCAAyLfuPx/vShJl3rdlsNub0yIG76yzttDunlpvjzaulCj2xWqrGLg6ABYlyDnr27O9Dd/I8l3BT17W+Zmiz2ZhNSrKyFcLcla3xu/Jf81m6JWy9Xrvz/sh01t6iLuu0O26mG9jS9hdetFQBGMvxbVSSorx3O9YZyG2parmyCNxQWbdpma1TskS3e1mvPnaVD4qB6lgcucSyafwBAEyWbmHyNgvJzZWTdhMUd5oqWZNXt65Dl6RFw3/CQHUAAAZz5DlIJmS2Ovgs0t+nG5Cs/wbW1Kqqevv2rbdHTwrgPqVpU0vr/qOlCgCASdCDpZpIQ1G4lV36B/e2Gx3WsJTnuTuqXTbFBH4Jk39OWZspFZo+FYHJx5vQFwYAI9qbqFqSaZzlthx6oXUW0DNLd+0E9J4p3r59e8CmTtO4Q7oGNpf9bTo0ZVpOU+CdNUdhdzoSrFefzmhuKU9ZltZyt4HarKjwj6fp7B2AJQiMQNffS00r6CaiwEZ+PI5A18PGw1+AgRdqetRbpKaVk5mcdmOhpWpyzHvFWOSHRcv5J+UGLLuOndnuq8t29P3vxtK013p+BNN6vdbtfD3d3wYAImrTTHXMl7D7W/34i3VkfFViDJAHoWpaiqKw7r5nLV+v196uQPNaDPMGzEVRFEVhtuWYL2G18Zi3RpaH9MrmDQqHJNnInRNF04lKqkuvLLO2KKXyPHd7P0laAMYVSFHuQ+EGJ4t3hFPc0KOnqmqanmq5xm4qG9T091e/L25Pn/kJcde3OsX058ft3mrzkLm1EQ8V72fV2lO3eIEdDDzLfWmh/5vnudvzaL6oCKwDAE19ZJam7jY5EXhnhLK21vTfNsLdf/o7NjCRuvnqY1f5oBa2t/uMWzzzl4R3BTdYNEUNM440vYqbOXRu8y4cPi7IDaeFd0+903junduzzR7pddyfd9azvGO2iFYAvAIxyBKIO02P/vBNDRpeP/Aqgck/2wy30q8+dpUPanFTKoSrY/jySBeV/lv+aGqn1efvvWrWdvkAAB9JSURBVMOq9KY6XeMqHYVW+5DuKZNyKqUGu6mwUkrf8Ni7I0opKXD7IukuTmvyuiZuz6M5hMvsS21aBwACOl30J3MXJw1fMuv1OjA1aBTylUiXX6NhsttETGd/pQ3Gehf0f5vaOdxGJv0Us7vQPNzdbsQ2vWPeTVmvKF3pLfdXblDVpM122jeYhaux5ZFvvju6Ds19b9qa2cS1d6cALE2nFqMm5s2MZeIo86e41SsnCzttP9Ac5csR//hi9K4/dpUPamF7O4F3122Fcs/Ze0OVG3GaBDayN1SZiUoXyW0xapMewoVss4WWocr8rnEf1fu+9xWbttNU/1aoNcdjAYDWKdwEeDs0AjGo08YPC1XeV0kmcNod0sL2drx313s5a9OY6wNaqry822kZqswCe6eGstKV/GBq2lq4FzJKS5UZAZt2Ta/Q/uWsnbK24I2YxCkATTqFmzbpRy6RaTNmfBRtvm9PycL2drx31zzfN5109/Ycuc0nVqZpcw1ap4vj9iYeK131lyfCocqMgE21p1PXMQ1j7vKmvNhp7wAsxNghZ2hL+zJknqqhNd1wJjEm3lyv196R13pMtPtoURRRJuc059LMsmzv/W307ow481NRFHo4eWCSUl1pnSrKvJIgcYbtJ4/XBMgUWWYlyCj79i8EAJi9sVPdoEbcX7ezLDyMOjDfQdJinqqAcEuV3uDe5hxvD2DgRXsaqL53Hgp31/a+1q71mCp37FT78gBYoLFbjoa2tG9CWqoGYrVn1HUtDUJZlpmNTGVZynK5F6awphTv71pWswHMaqHRu1BVlTuPQJ7n4UkNAlOie1+rPXPL1ka8sz+ER3e56rqWWemTJEnTVBZKcjJb9XaPeUsppSca7vRCAIDZGzvVDWo6+2uNWzcfCp+Mm8ZNR2mpajOcPFCYgE47Fd6Cd/JPr6zLjKCdCqxX21tdAGAau+VoaMlkTrvDWNzknxNRFMXO14OWJElVVU2Tf3aa5fIALcdF6Q67TlNuBo7CAWYT1cObujaJue/F7mknoDd+tZxZFABwStJdi7mOTkaazml/ZSZx+TswvB190BcNyMh3eS8C74J+s3inAECb12n3eAvb24W9uziYFarGLg4AzNLSTrt0/wEAAERAqAIAAIiAUAV4FEUhc2iNXRAAwGwsq7NTzzPUZFG1AQBAr5Y2pmpxk38u6t0FAACDofsPAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIjgp7ELMLQ0TQOP7na7wUoCTIH1ieAjAAAHW1yo4pwBmOQTkaYpHw0AOBLdfwAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAARTDFUKaWqqgqvUBSFtbAoCqWU9yEAAIC+Te42NVVV1XUdWEEpVde1UspdKH/XdV1VVTiWAQAAxDWhlqqqqoqiWK1W4XXcyFUURV3XeZ7vdrvdbpfnueSqHssKAADw1ITuopqmqf67LEurLUqvk2WZRCjdzSdPNHdEVnNzFXeNBbz4aADow9K+WybUUiXtTGVZNq0gMcvbBJVlmfXfcB8iAABAXBMKVWHSx9cUeN0hVgMUCQAAQJvcQHWvqqrW67W3Eavr2Cmzk7GTRTVgAgCAruYRqlarVZZl3vanro1SZCMAANCHGYQqGZBuTUAllwrKxFTuU7j0DwAADGwGoUqs12vzv3Vdm7NVkaIAAMC4ZjBQvSiK3VNJksisVBKq3Gv96rq2rgcEAADo1QxC1V66f1D+K39wsxoAADCk2XT/BSilyrJcrVb6yr6muUMBAAB6clJTncrIqkCcWtrUrkBLfDQA9GFp3y0L29uFvbtAS3w0APRhad8tpzCmCgAAYHSEKgAAgAgIVQAAABEQqgAAACIgVAHLdfD9xQEALkIVsFy73U7nKgIWAByJUAUsms5Vi7rsGQD6cAozqncS/jnOeQXLYX0WzP+aHwRrNT4jANBkcaGKUwIg5LMgU/MFGqvM1QYuIQDMC91/wHKZUYnMBABHIlQBy0WQAoCICFUAAAAREKoAAAAiIFQBAABEQKgCAACIYHFTKgAAAKag6wOhCgCAxWEKuj7Q/QcAABABoQoAACACQhUAAEAEhCoAAIAIGKiO2eMaFgDAFBCqMHtcwwIAmAK6/wAAACJYXEuV1VVkoakDAAAcZnGhitgEAAD6QPcfAABABIQqAACACAhVAAAAERCqAAAAIljcQHUsELODAgAGQKjC6WN2UADAAOj+AwAAiIBQBQAAEAHdf5PAoB8AAOaOUDUJDPoBAGDu6P4D0Gi73X748OHi4iJJkouLi/fv32+327ELBQATRagC4Lfdbt+9e3d7e7vZbG5vb8uyvL+/v7q6IlcBgNey+psm3r828eJN3N7ao3oDvJXz4cOH29vbT58+PX/+XJY8PDz88ssvZ2dnHz9+HLyMOEEMJx1d31+MS/viXdjeTvvdnXjxJo5QdQxv5VxcXGw2mzdv3pgLv379qpS6u7sbsHQ4cXw2R0Soimtx3X9p0NilA4amj3zrI3B/f393d3d5eWmtf3l5KQ8NWkoAmIPFhapd0NilA4bW9BE4Ozs7Pz//7bffrPW/f/8uDw1aSgCYg8WFKgAtXV9f39zcPDw86CUPDw9fvny5vr4esVQAMFnL6uyceOfuxIs3cYypim673V5dXb1+/frq6ury8vL79+9fvnz59u3bzc3Nq1evxi4dTgefzRExpiouWqoA+L18+fLm5ubs7EwpdX5+rpR68eIFiQoAmiwrQk48Mk+8eBNHS1WvqD30h6NrRLRUxUVLFQAAS8Q179ERqgAAWC6iVUSEKgAAFocs1QdCFQAAi0bAioVQBQDAspCiekKowrzx1QAAx+O7NIqfxi4AWuFe7k12u52+ZHdp1+4CwAHIT/0hVM0DoSFAclVC1gSAFsyvSk4rcRGqMGPW7y3zv3xNAAAGRqjCjJkNeDRWAQDGxUB1zJvZdk2iAgCMaHEtVeEBepyVZ4e3DAAwEYsLVZyDAQBAH+j+w4nbbrcfPny4uLhIkuTi4uL9+/fb7XbsQgEAThChCqdsu92+e/fu9vZ2s9nc3t6WZXl/f391dUWuAgBEt6wJKiY+Icfe4k28/OPyVs6HDx9ub28/ffr0/PlzWfLw8PDLL7+cnZ19/PjRfK75LCrZxbGH/nB0jajvyl/am7uwvZ32u0uoOoa3ci4uLjabzZs3b8yFX79+VUrd3d212QIElYP+cHSNiFAVF91/OFn39/d3d3eXl5fW8svLS3lolFIBAE4VoWp83IapJ2dnZ+fn57/99pu1/Pv37/LQGIUC0Fb61NjFAfYjVI1PzwaeELBiu76+vrm5eXh40EseHh6+fPlyfX09YqkAtLHb7aTnSP8BTNyyOjun3Lnb5i4rUy7/6LyVs91ur66uXr9+fXV1dXl5+f379y9fvnz79u3m5ubVq1dttgBB5aA/4aOLY69XjKmKi5aqkVkt27R1x/Xy5cubm5uzszOl1Pn5uVLqxYsX3kQFAMCRphghlVJFUSilrOVFUVRVpVfo9KiYbGRueUvgyZZ/Co6/dpLqDaByEJ0+qGipGhEtVXFNbm+rqlqtVmVZmqFKFiZJkmVZkiR1XSdPw4dSShaKLMskYFmm+e62/GZps8KSEap6ReWgD21+T3Ls9YpQFdeEuv+qqiqKQsKTRRbudruqqqqqKssySRLdHFUURV3XeZ7LYMY8z+u69oaqaVrUAQcAWpsWemBGJhQhzVFEVktVmqZ5npudemma6uYo9zNpPmq9xHT210VTyjFoqeoVlYO4AsNGrSONY69XtFTF9dPYBfgHqXfd02eSpimLmbqkW9D8r9kbCACYFHPMA41VOBkTClUBOj9J45M0WZmhyhrVbg2xAnAwaxI1znyIxWzD4LjCaZhHqNJ0I1ae5xKkuo6dOniqAj7zWCaOfPSEQwunZ2ahSo9VX6/XSZJ4Z17Yu4VeSgYAAJZtQlf/tSTTUGVZJrnKa0aX/gEAgNMwg1BVVZVSyspJVgMVKQrAKLjpLwBtBqFKRp1bscn8r3utX13X1vWAp42vdWAs3PQXgDaDUCXW67UOUnq2T/3fxGi7kj+a7lRzkvhaBwBgdPMYqC4TmZjzV2VZpmOTUqosy9VqpRtprLlDAQAA+janqU51S1VTYNI3VG7awsSndm0q3na7/ctf/vL58+e7u7vz8/Pr6+t//dd/ffnyZcunLwQzqmNEHDzH4IbKI2JG9bgWtrfTfne9xdtut+/evfvd7353dXV1eXn522+/3dzcyL9Wrpr43vWNUIURcfAcg1A1CmsAbk+VvLS3b2F7O+1311u8Dx8+3N7efvr06fnz57Lk4eHhl19+OTs7+/jx496nLwehCiPi4DkGoeqELe3tW9jeTvvd9Rbv4uJis9m8efPGXPj161el1N3d3d6nLwehCiPi4DkGoeqELe3tm83Vf8t0f39/d3d3eXlpLb+8vJSHRikVAABwEaom7ezs7Pz8/LfffrOWf//+XR4ao1AAAMCDUDV119fXNzc3Dw8PesnDw8OXL1+ur69HLNWk6ClPmfsUADCiZXV2Trxzt+nqv6urq9evX8vVf9+/f//y5cu3b99ubm5evXq19+nQGFOF/nDwHIMxVSdsaW8fLVVT9/Lly5ubm7OzM6XU+fm5UurFixduogIAAONaVoTc2zc0bm0c2ZSytB8EXdFShf5w8ByDb7bDDDPR1JGW9vbN4zY1ES3q3QUAnCo5nS0ttUwc3X9AkrRoxQQAIIxQBfwD0QqYiO12++HDh4uLiyRJLi4u3r9/v91uxy4UsAehCiBLAdMi9zy9vb3dbDa3t7dlWd7f319dXZGrMHHL6oudeN8zA9V7FaifWYz3xJTx6TuGW3vt73mKiR97Ey9edAvb22m/u4SqXjXVj9tMRTWiKz59x3Brr/09TzHxY2/ixYuO7j+cvq5TrtMbCIyIe55ivghVOH27p8yHyE/A1HDPU8wXoQqLZoUtb/ACMDDueYqZWlZn58Q7dxlTNS4qEAc7+ODhIonEV3vt73mKiX9xTbx40dFSNQldB/0AOA1WK+nYxZkK7nmKmVpWhJx7ZKalqldUIA525MGz8GOPb7ZjTLx+Jl686GipAgAAiIBQBQAAEMFPYxcAALBQegip/LGofiKcJEIVAGAcpCicmMWFqvC1dXzCAbS3tEG4AMIWF6pO9RuQiRiA4e12O52rCFgAGKh+UohWwMAkVyWn+4MNQHuLa6k6SWQpYHjW5878rw5YTJgOLAqh6tTQBwEMw+z1a2qsomcQWBS6/2aPZipgLGZaIjYBIFSdIGIWMAyCFAAToWreyE8AAEwEY6rmzfyhzLgNAABGREsVAABABIQqAACACOj+mwduOwoAwMQRquaBFAUAwMTR/QcAwMxw6fc0EaoAAJgZPY9/QsCaEkIVAADzw828J4gxVQAAzEmbm3ljFIsLVeFmUg5HAMDEtbmZN0axuFDFkQcAmDtu5j1NjKkCAGBmCFLTRKgCAACIYHHdf4CLCesBAMcjVAGkKACnxroqi2+5YRCqAGBMzNyIPphXCI5dlgVhTBUAjI9o5UrTVKpF/wFMHC1VADAaskIATSyYHVqqAGASCFjA3BGqAGAcpCjgxBCqAGAqiFnArBGqAGAE5Cfg9DBQHQBGYI7C5rp34DTQUgUAABABoQoAerHdbj98+HBxcZEkycXFxfv377fb7diFAtCjxYWqNGjs0gE4Edvt9t27d7e3t5vN5vb2tizL+/v7q6srchVwwpbVkX/aAxdOe++AKXM/fR8+fLi9vf306dPz589lycPDwy+//HJ2dvbx48e9TwdaCh88ox9aoxdgYAvb25N+d09774Apcz99FxcXm83mzZs35sKvX78qpe7u7vY+HWiJUDUpi+v+A4CIvPenu7+/v7u7u7y8tFa+vLyUh4YuJYBBEKpOAbcdBcaye0oWnp2dnZ+f//bbb9bK379/l4eGLiWAQRCqToH3ax3AiK6vr29ubh4eHvSSh4eHL1++XF9fj1gqAL1aVmfn0jp3AYxlu91eXV29fv366urq8vLy+/fvX758+fbt283NzatXr6yV+WrCwRhTNSm0VAFAfC9fvry5uTk7O1NKnZ+fK6VevHjhTVQATsayIuTSIjOAKdj7zcNXEw5GS9Wk0FIFAAAQAaEKAAAggimGKqVUVVXu8qIolFJKqaIouj4KAADQq8l1dlZVtVqtyrJUSpnLZfqlLMuSJKnrOkkScx2llCwUWZZ5Y9nSOncBTAFjqtAfxlRNyoRaqqqqKopitVq5D0l4KsuyqqqqquQd0msWRVHXdZ7nMktTnud1XXtDFQAAQE8mFCHNqcCtlqo0Ta3Gp6Io1uu1FF6eaO6Iu75ePp39BbAQtFShP7RUTcqEWqqknaksS/ehLMus3kB3Beu/Zm8gAABA334auwCtuG1O6/Xa/K8VuawhVgAAnBjdveN212AsE2qpaqmqKjmApE2r69ip9FB97AsAAIfhrq8TNI+WKk03QelBV+FuQRdHHgAA6MNsWqqkgUpf5RfOUlz6BwBYrO12++HDh4uLiyRJLi4u3r9/v91uxy7UIswjVMnkVVmW7XY779yepCgAAJIk2W637969u7293Ww2t7e3ZVne399fXV2RqwYwj1AlU1I1JSf3Wr+6rq3rAQEAWIK//OUvv/vd7z59+vTmzZvz8/Pf//73nz59ev369X/+53+OXbTTN4MxVTpLuW1UskSmDNU3t5GeQW5WAwBYoM+fP282m+fPn+slP//889XVlVLq48ePIxZsCWYQqjRrGoXkMTkppcqyXK1W+ho99y43AACcvPv7+7u7u8vLS2v55eWlPHR+fj5KwRbipKY6NVuqvJY2tSuAKWBGdQzp4uKiLMvf//735sKvX78qpe7u7gYuzNKO7YXt7cLeXQBTQKjCkN6/f39/f//p06eff/5Zljw8PPzpT3968eLFf/zHfwxcmKUd2wvb24W9uwCmIPDNY00szBcUjrfdbq+url6/fn11dXV5efn9+/cvX758+/bt5ubm1atXAxdmaafdeVz9BwAniUmxEd3Lly9vbm7Ozs6UUufn50qpFy9ejJKoFmhZEXJpkRnAFPDNg7GMfuyNXoCB0VIFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAEczpNjVRWLPCWBZ1kQIAAIhocaGK2AQAAPqwuFAFAIPRTePyBz/qgNNGqAKAvpCigEVhoDoAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIljcDZX1TeO9uPspAAA4zOJCFbEJAAD0ge4/AACACAhVAAAAERCqAAAAIiBUAQAARLC4geoAAJw8fam7/MFFWsMgVAEAcGpIUaOg+w8AACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAggsXNqK5n7vdiCloAAHCYxYUqYhMAAOgD3X8AAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIIKfxi7A0NI0DTy62+0GKwkAADgliwtVxCYAANAHuv8AAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIphiqFJKVVXV9dGiKJRSSqmiKHorGgAAgN/kblNTVVVd110fVUrp5XVdV1UViGUAAADRTailqqqqoihWq1XXR4uiqOs6z/Pdbrfb7fI8l1zVb3EBAAAM6XRuMJymqf67LEulVMtH5SFzR9I0zbLMzVVpOqH9BQDgtC3ttDuhlippZyrL8oBHsyyz/hvoQwQAAIhuQqHqGFazlvVfAACAvk1uoHpXXcdOmd2InSyqARMAAHQ1+1DVtVGKbAQAAPpwIt1/Fi79AwAAAzuRUEWKAgAA4zqFUOVe61fXtXU9IAAAQK9OIVTJfWn04Cr5g5vVAACAIc1+oHqSJEqpsixXq5W+ss+dOxQAAKBXJzXVqYysCsSppU3tCgDAiJZ22l3Y3i7s3QUAYERLO+2ewpgqAACA0RGqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAhO4TY1nehb2Xgtao4yAAAQ0eJCFbEJAAD0ge4/AACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABDBT2MXYGhpmgYe3e12g5UEAACcksWFKmITAADoA91/AAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEP41dgKGlaRp4dLfbDVYSAABwShYXqohNAACgD3T/AQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAR/DR2AYaWpmng0d1uN1hJAADAKVlcqCI2AQCAPtD9BwAAEAGhCm2Fe04RRu0diQo8BrV3DGoP7RGqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQARTDFVKqaqq3OVFUSillFJFUXR9FAAAoFeTm1G9qqq6rt3lSim9vK7rqqrM4BV+FAAAoG8TaqmqqqooitVq5T5UFEVd13me73a73W6X57kkpzaPAgAADCCdzr3wzFlry7JUSlkPmUVN0zTLMklO4Uetl5jO/s4OtXcMau9IVOAxqL1jUHvHWFrtTailStqZyrL0PpplmfVfs5cw/CgAAEDfJhSqwsyGq72PhlcGAACIbnID1V3e0VF6ZHrXsVPcGvMY1N4xqL0jUYHHoPaOQe2hpRmEqk5tVGGL6tkFAABDmk33nyXcQMWlfwAAYGCzCVWkKAAAMGXzCFXu1XwyMVXgUet6QAAAgF7NI1TJnWf08Cn5Q9+OJvwoAADAAGYwUD1JEqVUWZar1UpfgmFOZ+V9lFkVAADAkGY21amMnWoKTOFHZ0rfKFr+W1XVie1gf8yq082Z1F5L1oEnS2gAbq/pFu/UYRt8eA+w92QhNckR2Kt5dP9p4c/VSX7qrDH4q9WKUfkt6YpK01T+Xq1WfKG0ZB1mRVGs1+uRyjJLugLX6zWf2a68H16qMSx8spCarKqKObd6NbNQNUdyWBdFoY9v+UOOb1nitgHI7aXNj4RSamkfhihVp59SFEWe58v5Xo514MnCBSaquBVYPOq1zNMRpfbkUVkoYzz6L/iY+jtZFEUh98OtqirLsuUch8MjVPWrqio92Es3k0gjrf6CcH9AFEUhj65WK31hY1VV8+qrPVKsqlNKLareRMQDT/7bdFPOUxW3AuWJiwr0sWpPX8d9er0Qll5PFpXRFbioQ3EEO/RJTkX6v/J3lmVZlsmSLMvyPJe/8zyX5dZTyrJs+u8Ji1t1OhP0X/BJiFh7WZbJH8upvV28CpTt5HkuU8Dop5+2uLWnVzvtI7DXk4X1ZbiQ43AUtFQNrXo6mr6ua/MHRF3X1m8IJtzSjqk6uUQ0y7KT/73b5LDaK57OV7Jkh1WgNJRKl81ut7Nm1FuOg2svz/M0TZc28kFEPFlwHhkMoWpke4/1xX4L79Wy6vRwBGn0pj5F+9qr61qf1fSoYfDJPUb72pMwKsG092JN25GHnDVOCz0hVA1BDzZMnB/95hXXMpJaVuC4F8dXnYxUMLe2HFFqTzdrJ0my2+0W1Wp1fAUWRaFbWTj8tPbfe+bgoSU0t/R0sjBHsq/XaxJqj8boc1wQPZRHvg6kR9zsGt89HU4oS/TwC/l3yWOqolSduRq1p9fpdODtFjmmKkoF6v9y+B1Ze8PvyJD6PlnI8KxkMQP7xnLih+no9JDA8JepDMl0F/ZUqlmIW3Xe1U4YB96RolfgomqVD+8BBvjMLqEaRzeP29ScgHCniffRRfWzBESpusVWJgfekWJV4DJrlQ/vAfr7zC6tJkfBmKrecRwfjKo7BrV3JCrwGNTeAai0EzCze/8BAABMEy1VAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACCC/weBh0J9NSZGXgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "g['sigmaT_MS3_DRY'].Print()\n",
    "g['sigmaT_MS3_DRY'].Draw(\"AP*\")\n",
    "g['sigmaT_MS3_DRY'].GetHistogram().GetYaxis().SetLimits(100,150)\n",
    "g['sigmaT_MS3_DRY'].GetHistogram().GetYaxis().SetRangeUser(100,155)\n",
    "g['sigmaT_MS3_DRY'].GetHistogram().GetYaxis().SetTitle(\"#sigma_{t} (ps)\")\n",
    "g['sigmaT_OPT_DRY'].SetMarkerStyle(24)\n",
    "g['sigmaT_OPT_DRY'].SetMarkerSize(1.2)\n",
    "g['sigmaT_OPT_DRY'].Draw(\"PSAME\")\n",
    "g['sigmaT_OPT2_DRY'].SetMarkerStyle(22)\n",
    "g['sigmaT_OPT2_DRY'].SetMarkerSize(1.2)\n",
    "g['sigmaT_OPT2_DRY'].Draw(\"PSAME\")\n",
    "l.Draw()\n",
    "tL.DrawLatexNDC(0.11,0.92,'CMS MTD - POST MS3 Optimisation - Array Dry coupling')\n",
    "tL.DrawLatexNDC(0.15,0.83,'<MS3> = %.0f ps'%(R.TMath.Mean(6,g['sigmaT_MS3_DRY'].GetY())))\n",
    "tL.DrawLatexNDC(0.15,0.79,'<OPT> = %.0f ps'%(R.TMath.Mean(6,g['sigmaT_OPT_DRY'].GetY())))\n",
    "tL.DrawLatexNDC(0.15,0.75,'<OPT2> = %.0f ps'%(R.TMath.Mean(6,g['sigmaT_OPT2_DRY'].GetY())*1.5))#multiply by 1.5 because 4 of 6 elements are empty\n",
    "c1.Draw()\n",
    "for ext in ['png','pdf']:\n",
    "    c1.SaveAs('plots/LYSO_OPT_sigmaT_ARRAY_DRY_byProd.%s'%ext)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "id": "d31de83c-8d8d-4cd0-8303-fb9afdd98cf6",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x[0]=0, y[0]=0.0138488, ex[0]=0, ey[0]=0.01\n",
      "x[1]=1, y[1]=0.0455465, ex[1]=0, ey[1]=0.01\n",
      "x[2]=2, y[2]=0.0709276, ex[2]=0, ey[2]=0.01\n",
      "x[3]=3, y[3]=0.0865404, ex[3]=0, ey[3]=0.01\n",
      "x[4]=4, y[4]=0.0218292, ex[4]=0, ey[4]=0.01\n",
      "x[5]=5, y[5]=0.0186736, ex[5]=0, ey[5]=0.01\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file plots/LYSO_OPT_XT_ARRAY_DRY_byProd.png has been created\n",
      "Info in <TCanvas::Print>: pdf file plots/LYSO_OPT_XT_ARRAY_DRY_byProd.pdf has been created\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxwAAAI8CAIAAAD0vjrdAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElEQVR4nO3dsZLbSGI/YFL/Tawq1WgullyaN/CWH4CgnsSJR1N3odbBjgOSiTe5TX11Iwd+kSPBF7jSE0gKGHs00WTiP2gLhgAQ5JBNgmR/X6AagSDYaILEj92NRn+5XPYAANjNs64LAABwDoQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCoAgAiEKgCACH7qugDHKM/zPM/D31mWZVlWX6H83/YVGh9t335jkVpesf5oZckqYTuNK29Sqs2tKk/7qzyposor93q98Xi8SQE2Kc9TK3/zUtVfaItq3+Jweur2i40fz6aO6rW20/4tcVS2+/yesfrRFZYkWyH8nyU/GgwG9VoajUbldUajUfnR2WzW8mj5uZWHVm2/YjablVceDAZrH23ci8YXXVWkxl3bWntJ6q9S2alV+775yptXyNrCb/joqlLVtzCbzSrFaz8etqulrRVl2/FgKEobt3irrKrt41F+02N90PZki0/NGasfycW3aIK1QYXuv/+T53m/35/P5/WHJpNJv99f9cRK88NkMlm12qqHJpNJextGWaWEG7bBbGc4HO51+8WrlP+b53llSWE+n1feiCetHF1L5awqVa/Xq5RqOBxW3tPJZLL2J+8+dnw8Hvf7/X6/f1o/uLMsC8U+wLEaUflN3/zjf2w2OVYhHbr//k/5FDUajYqusSIJjcfjxu++xhxWV2xnNBoV28myLDz9SbmqrPFEUt5UeRdamk8Gg0HxrPJThsPhcl1T0+ZWvUq5y3LtG1FeudjaYDAoL6y8a1tUyCrlTrpNukVms1mxC8WuFbtQrBn2tCjbfD5v7w18Ui3tonys7rKdLMu2qO2tHfK1tlD5sG/4HdK58pGW53ko9tpj9eyNx+Ow+ylXAoFQ9b/Kp5/iLNj7PkglnMA2yT2rTmPF8nKq6H0/Q4fvps2/mMprNn4dl7eTZVn5RNuy2fJeh0aL8N9VaXIL5UE/4Y8iQ9RLWA5z5Tei/CVePLE+dClsuT7WYfMKabRJqCpKVTmWZrNZsQuVNUOZQxJaezw8tZYqQ0DCGXHrAVjlrRUjuspbe+r2GzfSuELLOptrebnGioryohX1Ruv6210pTDhzl7NyuZCVb5XGYue1kUDty+sqCbv49VL+6VUv9ng8LjJupXl11fKKTd6I9nUa97G+sPHY3vzbr/0T17Kp8kdm87eDY9Rt7+PxKMY3NA74qFRX0YNeH25SH6IUetnLv5srG599t6ps5eeGVywXslKSTcpfVt6XzR96qkptNO5aZc3GCqnvZsuutVTsk47/zSt/k72olKpew0W1txwSW9dSvf2mvs3KxusHeaVCCqPRqGX7jWOqGof0VXaqcTxcvTCVV6wXYJOXKw+OqawWcbhMuZZaPmXlEpbLUH5HKs9t3MGi5Ku+JRrrasMVKtVYLlLxd2V/i+dOp9NV+15WPwBavkZWrdO4C5WF5cK3DI3dZExVUeb6O1I5vOuDYjf5BuBoCVX/a9UR36geqiqfpfK3QPFQ5SO6+Qem8Su48lB7AGr8Qqnvy5Oe+CSN303LH88By9boGbSEsCfV6pP2a/PKLx6qHABP+nJcW7Yda6muXuZgbajaRHFc1U9FLT10xWu1XGHQWJiWtL12U8tDXbdROQfXi7GqniuhqrLXa+uz8cipJLZGq0q4/DE9rCrDbDZr/J7ZZHx3+x61r1P//dPbLFQ1Kgr5pFDVqF75a/eRU2GgetVTW1zDh6dosA0dN40bKX/MJpPJcDgMo2vH4/GqXqS1xdv8iUcoLw0AqlTjqq+kesWWv5XKtZpFGlHU+Ophyy2t9OVG/kqpWtr/sywrulxbvm23qKWy+skvFCnP88qZb+1nofFcHs4ExaZaRgsVVVFv6K18oHpNLWphneWPsW/Vm55/71TtNSW2xj2th8hYR1R5wF+v9D62bD80YNTLGZaH7aytz8aer1hj5hqNRqNQ7OJVysdDpR7qysvDHtW7C4ti18N66AHfruT1wLTqMqO16om2XvnFy21ynTJHLVo8O3FPqpDyuafSehH+XtXc3fI7uKUBvPL7sniJ4u/wzdWynZa927GlavajVatteBxu0uFYPKV4uVltSoL6lp+6X6tWLip8WWqVXNV01FKqel01RpNGW9RS8d/yU+rFbtxyS0tVvUjlo72yWv33fblU5XezfCyFY7uxWaKejRpXq6zT3ivUuCObVPiGH4TlxtVeLlulzKt68UKDUHvJ6ztYfzfrWtapvERLM2pR7Ol0Wt5sS602vrPFwpaybbiblYWrCl8srMSjDbv/6iXfsFRaqk6Rlqpdla9lWzvAcDweh3ofjUaVM+58Pn/SYPDyL7CuxjNmWTb80XbbaW8Db3n14o88z5ffTyqVio07q0K5Ra2lVTJoKVW9rrIsK68WazKLSvHKx1jnw2DLrVmhMS/7fglVUbYwzLloxtv93azs9SatRGsr6kkfhOItKHa/sRWn8Snty0ODUPFvv9+vN64U+1K+jKP3Yxf2k2w+qrrS5FOvh7rG61fCx2q5XPZWt9rWB57vot6AurlyyVf1M1QKr7HqpAlVVY0fm/y79ieGFTb5SITzRPhqKNbfsHm5OK/v+yKRffQthra9QvgpVo5H4Y9VZ5f2IoXTSajYrU8S7bIfLznsbXbGLUrV2AVQCKsVJV91bt6xlo7HeDyu/7QI6arcPxLCwY6TDjSenjtRlCTM8FS5VG2Xy2zzPO/3+8PhsKW6Vp3Xd7+8d1VaLRTHc3Fxa2V5RGs/I+26jTWdH6LsQqj6X8WnqPHLpf0HaCXlNH4k6pfLlh96UlFDCXcZLvCkF+q1jt0Z/GjtNkPCKFTqqnHMR3uRQq3W37XKFeZrC7ahcgtKZUlltcZShRorl6poiSnWWXt626KWjlaRNVe15JVnd6v8CHmSTVqktrb5B6Hy6vPvVq3wJEWlFd2pjW3A5RayrccJFZ6UVsvdeeGJezpE1447bNfttGGn8qOIZofqZzx25W+fSk92/dquyjCCSqNIZRTLqHY3mJaXXlu2sKT8cpUyxBpT1VIhW6jURosNK6plfE/9RZ/00NqVez9atr5B9aqrDBNpLMna4pVPFZWHyhVYr6XGMVX1oSGHGVNVH4FUqcnG8tT3bt9jquol39om5/hVZatvpHJ0tbwvq0ZrNdZJ3dqSDGqDxBvrqjyHQvijGF/V6AzGVJVL/tRSGVN1irRU/a9yE8JwOAyNB3mel3/MrepRamluKVTGshS/RfLSLNvb9VjFaiuez+dFZ0ToRAjLB4PBIZujyxVVXBoZ2vnai1TuM8r3Ob9z+by46i2rHEvlUq36ZV9cglcMG9qwwaNSS+ULKuuVUAzdy3+c3n3VC+3V5kPxKlOkbqF44mQyqfebH6xJr34xY6G9sfxJ8u9Xp27YELX5K2Yl5Zt6bbiFSg9gb13fX/EGFavlpas4g3KjdbFOy+y+G17tWJ76uF6e6OpDIDhVh89xx6yloso/vOq//+pVWvy3+FGyNjOtKlX991P5HDD78Qr2rVuq9nqE1GujRfu49foOtpe/8dfek/ausnK9KajxB+6GpWrZ2Yi11F6YVducrZunqnjidi1VLVFmtHrSo0rxKtupfySL8rS8XPuOxGqpammFbXx01cprW6rqKiWvfOrXlrxly0+tq8Y3a4uXbmyVbClYe26uFL79FSO2VG34ipwQLVU/WK747I1Go/Z+7nIz+Kp1xqV7NdSfvnzKnIrlXzP7+2Uz+nFmwoPJsqw+yKYoUv2NWLVyGFkSvX42rPxVx1Kvdu+axly79orIp9bSqs2W3+ID/0rO83xV+cffb9pTqZxwfUP4u2i9KLcltIyGKV8EUHjqR28XxZvS+I7v3kBVf3OLJZVqKb/WLpd0hNp7UsnLL7fJSzd+jiqf61VDx9oHVra8+qBpRvV9fDrqR/iGg1M5Wv1OzprHrzzqPPpnqWic39P2z8mT3oi9vmtbq1zBsKpguxR+7XPLc4pm329qtmrl9kf3ob38xaPtJ+8ih214nLS/Hadrw7evOCQOfwrY4qWL78yW/dpwnfYjPyT1wWCQf3eYD0L5tbLvgwT28ZuQfROq4PxVQlWnZeEohEMipIdDvu5sNnv79m0nL71WJVQd4BWLD2NllGT4w9n5FOn+A0hF/uPdkA6ZsEOgD4mqt1u349ko5tQIl5uEWdnCQzoBT9RPXRcAgMMpxlcNBoPdB3I96XXLL72POT9PznK5LFJU+VLNI2zGY0NCFZy/Yui6vj/CXBuHH3cYhmBnWXa0iSr7PrHOIWtmuVwWowaDQyZdojOmCgAggi7HVIW7lGRNd/N40sr9GkkfADiwzrr/stLU0uE2di1dyMWFKr1ebzKZTCYTFzEBAEelm5aq8Xg8n8/LU4233B44hKfZbBaCV+ivLN+1oFebeVZLFQBwYN2MqQotT+WX7vf7q653qD8U7scXnl7+GwCgK52NqapMwjEYDFbdX6L9hr5F2GrvQAQA2Ksux1RtuGY9KtVvvV5M9dEzwwcA0Im93ap5pXDzy/L9w5ffb9+94XN7pUFU4b/l4Vm91Tc/76iOASBR28eFE9TWUhXmL4h+kd3WG2y8zWTlMsDxeJznecud6peiFQAcRLkfKQVtY6paokl07R12eZ73+/3igsFyiqpHtLBEDyAAcEidjanaPPSEO4c/daSUWawAgEPq5uq/+rV+oRWqceUwJVVjogotWJVZqbRRAQCH1zZPVZggau1Emlu0CVUan8JgqaIk4dHRaBQGSIVQVY9coWDhucXIqjBtVXhu4x4ZUwUAh5HaaXdNqNpkE9vVV0g/xX/L483LkasIVS2vWynnqkTVS+/dBYAOpXba7aalqlC0VG29hWI7eZ6HOy63rJbauwsAHUrttLsmVJ3ZfYtTe3cBoEOpnXY7u00NAMA5EaoAACJom6eqcs9jAGAXSfWFJWjXzs4wPDxSYfYutc5dAI5KauegZ8+eJbXLbS1VaweqV+aXAgD2au1V88UM2PV1yvfGHQwGJ9QmcjJabrYcVhgMBvWHZrNZy6NHq31/AWCvvu2mmAd7NBqtWqc4v5cXTqfTxgwwnU4rT59Op8Xgn8FgUF/hSVI77bYNVF8ul+F+Mv1+v3zvl/F4XMxy7p4wAHBg5dmzy1bd0u3t27e9Xm80GoWQNJ1OQz57+/Zt+SlhzaI1az6fv337dtUU3DRYG7uKRqnRaHSiDVSFTfYXAPZkl1afb6WWql5TI9O3b9/KV5hVFtbXD1sbDAblhqXKmi2vtYnUTrvrp1TIsmy5XA4Gg8lkooEKALoVwlDRnlQ2n8/rV+6HNesjqCrbCWf2ylir0G/Y+FrUbTNP1UkPbeu36rp0ALCReg9gWFJuzVrrW2m4VUhOjU/XkrKh9aFqPB73+/35fF50/w2Hw9PNVe0Nd12XDgDWC9GnknVWXRgY2q6ePXvWOBKrWD8Mfq88PTzldE/6B7YmVGVZFip0NpuNx+PQFTgajcLo9bX3WgYANvTs2abdRyFUVUJSY99fr9crxkNPJpNnz56FdNXe+JTneRjzE7b5pNavlLW9f6GBKoxJL6fU8Xgc3qHJZCK9AsAuNs9SFeWhTu19f99+nCthMpm8fft2VdtVr9d7+/btZDJp6RCk0Zo3cjabNYbZYvS6wWsAsItv374VueqpjVXFOXrtpKBZls1ms2+l+RR6vV5xCVq9SOVpF1ZlLyrcpgYADqTlHBTi1LfS7J11k8lkMpmMRqMQd549ezYYDELfUeXvtZvqleavmk6nLafyDbe26rlJnXa3bHIsnFCiAoAj8axm1fL2tqvQX9TS9xfGUa3qdCrPqpDnuav8drRrqAIAnqo+SWbRFNQ4hWajogewpe8vjKPapP8uDLTaYl8oCFUA0KVnz56VE9XmTyyuAVx13V/ve9Kaz+eNrVDlJq6wZiV+abt6kp+6LgAAJG274UqFtdfoTafTt2/fvn37NsyWHuLXfD6vdBqORqMwZqv3vX2rWGfV/ZipEKoA4FSFJNRbd91fyFXz+bxyzX4x5j0Iq1Uaq9qHsVOW1tVwrv4DoEPdnoPyPC9CVeUef42rtayzodSu/ksrZAhVAHQotXNQaqHKQHUAgAiEKgCACIQqAIAIhCoAgAiEKgCACJKbp6rf77c8mtRFCgBARMmFKrEJANiH5EIVAJyH8tTnqybqbLx5X33NTe7xZ171tdKaDNPknwB0KNY5qLhJX0X9ljLPnq0cPF3coybP87dv37a/4mAwmM1mTy1napN/phUyhCoAOhTlHFTkpNFoVNz5uLi3TCVXhZUrt1uur1yJaJV7Lfe2vWWNUHXOhCoAOrT7OWg4HM7n88Z2o6L56tu3b8XCEKrKS1pWDkLD1XZNUxWphSpTKgDAaShamBrjTtGdNxwO126q0nZFFEIVAByFEJhahLal6XS6aoUQldZuhz0RqgCgS4vF4t///d//8Ic/DIfDP/zhD7e3t4vFonHNkJbaxzaFUVZrr+bb5HI/nkqoAoDOLBaLm5ub+/v7v/3tb//zP/8znU4fHh7evXu3KletFSJXpbEq/9FkMgmX++kEjMs8VQDQmb/85S//+I//+Oc///n58+e9Xu+f/umf/vznP//yyy//+Z//+R//8R/lNXdpW2qcMWEwGAhVcQlVANCZv/zlL3/7299Cogr+4R/+4V//9V+Hw2ElVO0y92Y9PG03RQLthCoAOLTynJz//M//3L5OfdaDFqFBK4ysKmiROgxjqgDg0L599/Lly7///e/ffvT3v//94uKi+G/xrBCVGudSL2wymJ09EaoAoDM3NzcfPnx4fHwsljw+Pv7Xf/3Xzc1NfeXQ4NQSqsIMVZVmKg5GqAKAzvzxj3/88uXLL7/88vHjx4eHh48fP/7bv/3bp0+f/vjHP9ZXzrIsBKbG6T0nk0nL1KAcgFAFAJ159erVX//614uLi+FweHl5ORwOX7x48de//vX169eN689ms8FgMJ/Pnz17NhwOJ5NJmCLh2bNna6cGZd/Suheee/8B0KH2c1Ce5xuOhSru3FfWeLe+Vff+ay+Ge/9tJ62QIVQB0KG456DiVoBHOz+CUHXO+v1++wpJ1QYAB5baWSa1UJXcPFVJvbsAwMEYqA4AEIFQBQAQgVAFABCBUAUAEIFQBQAQgVAFABBBclMqAMB5KE+qvmr+zzzP6wvrazautvZZVCQ3+WdS+wvAUYl1Dmq8TU2v15tOp5XoE25T02g0Go1Go973+9K0v2L5rjXhboPFZO6j0WhV3kpt8s+0QoZQBUCHopyDipw0Go0Gg0Gv15vP58Utayq5KqwcwlOhvnIlooX/lp9VtIRtftvBnlB13oQqADq0+zloOBzO5/PGEFPEnfLtk1fdULlx5aDlhspFm1Y5uhW5rRLdeumFKgPVAeA0FC1Mjc1CRawZDodrN1UPQJsIr17p75tOp73NRmWdvbRabrRUAdChHc9BoZmqPnCqrNI0taqlquWhlpaqUID6U1ZtSksVANCBlkHlQWgoar8KL4yyWttutF3D0mg0Cu1S9U2F102cKRUAoHtrE9WGsiybz+fz+bycvSoRaj6f14eib7jx+sIwymq7/sQzI1QBwLF49uxZY1ddb7dBS40zJoTZELbeZq80bn3VLFmpEaoAoGObNFPtklrq4Wn3GBTGV/VWXPeXpuRCVb/fb3k0qfF0AHSlJUXVH1rVdtWocYRT3NBTTMewanqqZCU3UH3ZquvSAZCEbyXtj5ZXCFGpcS71wiaD2XdRJKrpdCpRVSQXqgDgyK1qxwoNTi2hKsxQtb8L8cINagaDwbdv3wyiqhOqAKAzT7roL8uyEJgap/cs7se3vwakkOc0UK2S3JgqADgeTxov1ev1ZrNZGCH+7NmzMNh8MBgUUyT0vs9vvichtDUGQeOrekIVAJyW2WwWBjaF+aiK5d3GmnJJkpXWbVvcpgaADsU9BxW3AjzaaaJSu01NWiFDqAKgQ6mdg1ILVQaqAwBEIFQBAEQgVAEARCBUAQBEcNShajweZ1mWZdl4PI67MgBAXMd7NVyWZZXpN8JNIhuF2ySHeWaL+WTr15e6+g8ADia10+6RtlSNx+P5fD4ajcJ9jkej0Xw+XxWqQniazWZ5nud5Ht6/xin8AQD25EgjZGh5Kpet3++vaqyqPzQejyeTSX3XUovMANCh1E67R9pS1avdZDvc22jVmsc5kywAkI7jvfff5jmp3nxV3FcSAOAwjrGlqrGPb8OMled56DpcdVPJ/ra23x8AIAHH2FK1dV9eccFg46V/QVKduwDAwRxjS1WjlvkUet8bqIoLBg2xAgAO7BhbqoL2FFVZczgctk9kBQCwV0d6rWPoyKtMqTAajRpnS6/Pv7BKatd2AkCHUjvtHuneVhqfKhkrPBoyVvi71+uNRqPKRuoJLLV3FwA6lNpp90i7/7IsG41Gk8mkuOyufjVfpbOvPo2CmwACAAdz7BGyaKmKsrXUIjMAdCi1025ie5vYuwsAHUrttHsyUyoAABwzoQoAIAKhCgAgAqEKACACoQoAIAKhCgAgAqEKACACoQoAIAKhCgAgAqEKACCCI72h8v4Ud2hulNRs+gBARMmFKrEJANgH3X8AABEIVQAAEQhVAAARCFUAABEIVQAAEQhVAAARCFUAABEIVQAAEQhVAAARCFUAABEIVQAAEQhVAAARCFUAABEIVQAAEQhVAAAR/NR1AdhIv98v/3e5XHZVEgCgUXKhqpJOKo42rISC9fv9oy0hACQuuVAllAAA+2BMFQBABEIVAEAEQhUAQARCFQBABMkNVIc6M1YAsDuhCsxYAUAEuv8AACIQqgAAIhCqAAAiEKoAACIQqgAAIhCqTsBisbi9vb28vOz1epeXl7/++utisei6UADAD4SqY7dYLN69e3d/fz+dTu/v72ez2cPDw/X1tVwFAEclrYl5TnEiotvb2/v7+99///358+dhyePj4/v37y8uLn777bduy3ZmTvHwADhmqX2vJra3J/juXl5eTqfTn3/+ubzw48ePWZZ9/fq1q1KdpVM8PACOWWrfq7r/jtrDw8PXr1+vrq4qy6+ursJDnZQKAKhLLlT1W3VduqqLi4uXL19++fKlsvzz58/hoS4KBQA0SC5ULVt1XboGNzc3d3d3j4+PxZLHx8cPHz7c3Nx0WCoAoCKtzs5T7NxdLBbX19dv3ry5vr6+urr6/Pnzhw8fPn36dHd39/r1665Ld1ZO8fAAOGapfa8m11J1cl69enV3d3dxcZFl2cuXL7Mse/HihUQFAMcmrQh56pH51Mt/zMKIOtULEFFqpy0tVfB/jvBiBQBOhVAFshQAEQhV8AMBC4DtCFWkTooCIAqhCqrELAC2IFSRNPkJgFh+6roA0KXytb6pXfoLQFxaqgAAIhCqAAAiEKoAACIQqgAAIhCqAAAicPXfaSiu/HffXwA4TsmFqvZ5iY42rBxtwQCAILlQJZ0AAPtgTBUAQARCFQBABEIVAEAEyY2pOkuV0ffGjQHA4QlV5yCkKPcDBoAOCVVgGjAAIhCqQIoCIAID1QEAIhCqAAAiEKoAACIQqgAAIugyVI3H4yzLsiwbj8cbPiXLsjzPKwv7NZtvEAAgis6u/suybD6fh7/n83me5/W0VJHnefEUAICj0k1L1Xg8ns/no9FouVwul8vRaBRy1ar18zwfj8fD4bDxoV6vN5vNliVaqgCAA+tmDu76FIv9fn8wGKzKVeXbsMxmsyzLiv+Ox+PJZLLhXpz3nOPnvXd0xU2QgK2ldmLqbEzVYDCo/Lelay+0P81ms/pDRQ7bpAMReKrw6Sv/AUCjLsdURdxa+cd0S4sXAMCedNBS1Zh4ts5YoX2rMjyrZWv1SwU3tF3xAIBEdNBSFbeNqj7Eqv0iwXPtvxD7AKBbxzL559YddvWIFpak2QMoWgFAVzoLVfsOPXHbw46cLAUAnesmVNWv9QvTVj11O3me1+dPT7ONqiBgAUAnOpv8s1dqTAp/FNmoMSo1yrJsMBhMJpMiSBXTikYu8RGTogDgGHQzpUKWZaPRaDKZFIGgPgfVhg1OIYGVJ1sfjUaJz6ie2mRrAHAMOj77huS0+/inMPNnuD1zy2rnlzZWNVOd2W7SufP77AAHkNpXR2J7e9bv7nnvHd1ydAFbSO2r41imVAAAOGlCFQBABEIVsNJisbi9vb28vOz1epeXl7/++utisei6UABHSqgCmi0Wi3fv3t3f30+n0/v7+9ls9vDwcH19LVcBNEprBNl5j5g7773j8G5vb+/v73///ffnz5+HJY+Pj+/fv7+4uPjtt9+6LRtwElI7MSW2t2f97p733nF4l5eX0+n0559/Li/8+PFjlmVfv37tqlTACUntxKT7D2jw8PDw9evXq6uryvKrq6vwUCelAjhmQhXQ4OLi4uXLl1++fKks//z5c3ioi0IBHDWhCmh2c3Nzd3f3+PhYLHl8fPzw4cPNzU2HpQI4Wml1dq699/BJ10ZqXdfs22KxuL6+fvPmzfX19dXV1efPnz98+PDp06e7u7vXr193XTrgBKR2YkqupWrZquvSwRF59erV3d3dxcVFlmUvX77MsuzFixcSFcAqaUXI847M5713dMvRBWwhta+O5FqqAAD24aeuC0AExVix8EdSPwsA4EgIVedAigKAzun+AwCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIwG1qgDbuLAmwIaEKaCNFAWwouVBV/Oxu5PwBAGwnuVAlNgEA+2CgOgBABEIVAEAEQhUAQARCFWQdaEsAABHESURBVABABEIVAEAEQhUAQARCFQBABEIVAEAEQhUAQARCFQBABEIVAEAEQhUAQARCFQBABEIVAEAEQhUAQAQ/dV0A2Lt+v1/+73K57KokAJwxoYrzF1JUv98XpwDYn+RCVaXRosJJFwDYTnKhSmwCAPbBQHUAgAiEKgCACIQqAIAIkhtTBcBJMBkKJ0eoAuAYmQyFkyNUAUByNATug1AFAMnRELgPBqoDAEQgVAEARCBUAQBEIFRx5haLxe3t7eXlZa/Xu7y8/PXXXxeLRdeFAuAMCVWcs8Vi8e7du/v7++l0en9/P5vNHh4erq+v5SoAoktr2L/LHFJze3t7f3//+++/P3/+PCx5fHx8//79xcXFb7/91m3ZgE343t6rfVdvam9fYnub2LvL5eXldDr9+eefyws/fvyYZdnXr1+7KhWwOd/beyVUxaX7j7P18PDw9evXq6uryvKrq6vwUCelAuBcCVWcrYuLi5cvX3758qWy/PPnz+GhLgoFwNlKLlT1W3VdOiK7ubm5u7t7fHwsljw+Pn748OHm5qbDUgFwlpILVctWXZeOyP70pz99+fLl/fv3Hz9+fHh4+Pjx4y+//PLp06c//elPXRcNaGMyFE5RcqGKpLx69eru7u7i4iLLspcvX2ZZ9uLFi7u7u9evX3ddNGAlk6FwotIalp/aZQiUeffhVJgM5WBc/RdXYnub2LtLmXcfToXJUA4jjCQWqiLS/QfAETEZyoG5SCsioQqAI2IylMOQpfZBqALguJgM5cAErFjS6uxMrXOXMu8+nIrFYnF9ff3mzZvr6+urq6vPnz9/+PDh06dPLt2NpZ6i9vT1mNoXr5YqAI6LyVAOT2NVFGlFyNQiM2XefThFPrnRrcpP+6jn1N6+n7ouAABwOOWUk1ro2bfT6P4bj8dZlmVZNh6PN3xKlmV5nu+xTAAAJSfQUpVl2Xw+D3/P5/M8z9empTzPi6cAABzAsbdUjcfj+Xw+Go3CDY9Ho1HIVavWz/N8PB4Ph8MDlhEA4OgHqtcn0e/3+4PBYFWuKo+/m81mWZZVHj3y/WV/vPtwinxy98q9/+I69paqXq83GAwq/23p2gsNWrPZbP/l4mT0+/2Qtos/ACC60xhT1XUROG1J/U4CoCtH3VLV2Me3Y8bqb2uXFwUAzt5Rt1Tto41KowUAsA9H3VLVyOxTAMAROoFQJUUBAMfv2ENV/Vq/MG1VV+UBAGh07KEq3JemGFwV/ihuVpPneb/f3/zeNQAAe3LUA9V7vV6WZaPRaDKZFNff1eeg0j8IcH6Kr/36LNBwnE5mqtOQnHafT+FU9hc4CZX5VnzDcFrMqB5XYnub2LsLHIbvFk7OYX4PpPbROPbuPwAguqSyzsEc+0B1AICTIFQBAEQgVAEARCBUAQBEIFQBAEQgVAEARCBUAQBEIFQBAEQgVAEARJDcjOqVifkrzDALHJJbB8I5SS5U+c4Cjkf4Rkrt/mhwrnT/AQBEIFQBbKl9OAGQGqEKYEvL5bLIVQIWIFQBbK/IVQZFAckNVAeIotI0Vf6vgAVpEqoAtlG+cE9jFdDT/QewtfJUCBIVIFQBbEmQAsqEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIhCqALrm/DZwNk39C0ipndHMEdKU85RVwooQqSFp5WvCuy5IizVRwTpILVe1fYc4rQFdEWzh1yYUq31nAkdBMBWfGQHWAYyFmwUkTqgA6ID/B+Umu+w/gGJSHIhhNBedBSxUAQARaqgC2V/TihT80OEHKhCqA7UlRQEH3HwBABEIVAEAEQhUAQARCFQBABEIVAEAEQhWky6TeABEJVZCu5XJZmWYJgK0JVZC0IleZbwlgRyb/hERVmqbK/xWwALYgVEGiQnIKt/LVWAWwu+RCVfvAEScVkhISVfjbwQ+wo+RClTMHFHwcOud+zHBOkgtVAMdDioJz4uo/AIAIhCoAgAiEKgCACIQqAIAIDFQH2JfKHC6GpcN5E6oA9qU8w2rXZQH2TvcfAEAEQhUAQARCFQBABEIVAEAEQhUAQARCFQBABEIVAEAEQhUAQAQm/4SkFVN+hz/MUQmwteRCVeWuERXOKKTGMQ8QS3KhyikEANgHY6oAACIQqgAAIhCqAPZisVjc3t5eXl72er3Ly8tff/11sVh0XShgj4QqgPgWi8W7d+/u7++n0+n9/f1sNnt4eLi+vpar4Iz1kxq43e+ntb9AV25vb+/v73///ffnz5+HJY+Pj+/fv7+4uPjtt9+6LRscTGqn3cT2NrF3F+jK5eXldDr9+eefyws/fvyYZdnXr1+7KhUcWGqnXd1/AJE9PDx8/fr16uqqsvzq6io81EmpgH0Tqjh5/R91XRzoXVxcvHz58suXL5Xlnz9/Dg91UShg74QqTt5yuQzNy8Uf0Lmbm5u7u7vHx8diyePj44cPH25ubjosFbBXaXV2pta5mxRvLkdlsVhcX1+/efPm+vr66urq8+fPHz58+PTp093d3evXr7suHRxIat/MWqoA4nv16tXd3d3FxUWWZS9fvsyy7MWLFxIVnLcuI+R4PM7zvNfrZVk2Ho+3Xrk+jGY0GjVuMLXInBRvLkfLwUmyUjv4O7uhcpZl8/k8/D2fz/M8D5lp95UBAA6vm+6/8Xg8n89Ho1EYWTwajUJU2mLl8MdsNluWrG33AgCIq5t2udBhV37pfr8/GAwac1X7yuPxeDKZbLgXqbVDJsWby9FycJKs1A7+zgaqDwaDyn+LDr4nrVxustInCAB0pcsxVRFXLo9VX9XiBQCwPx20VDUmnlWxae3KocmqMuKqJYT1t/WUXQQAktNBS1XcNqrZbFZeJ8y80NKTmFTnLgBwMMcy+eeTOuzKK9dTV1iiBxAAOKTOQtXWKWpDT2oPAwDYUTehqn6tX5iJ6qkr53ne7/crs1JpowIADq+zyT97pcak8EeRjSpRqWXlLMsGg8FkMimCVDFT6J73gGPhGgIAjkQ3UypkWTYajSaTSXFGnM1mlXWKnJRl2Ww2Gw6HjSuHBDYcDoslq278x1laLpfF5HKpzTIHwFHp+CRU3CN5x5XDzJ9ZlrVvykn3XNWn3YdjUGlJdYiSmtROu4ntbWLv7tlr6fvzRgN0LrXTbmczqsPuyr1+GqsA6NaxzFMF2yn/DJKoAOiQUMVpE6QAOBJCFQBABEIVAEAEQhUAQASu/gOAc2OOtE4IVQBwbtxnohO6/wAAIkiupar9/rsSPQCwneRCldgEAOyD7j8AgAiEKgCACIQqAIAIhCoAgAiEKgCACIQqAIAIkptSgfNTzD0W/jBrBgCdEKo4eVIUAMdA9x8AQARCFQCclcVicXt7e3l52ev1Li8vf/3118Vi0XWhkiBUAcD5WCwW7969u7+/n06n9/f3s9ns4eHh+vparjqAflLjUfr9tPYXgNTc3t7e39///vvvz58/D0seHx/fv39/cXHx22+/HbgwqZ12E9vbxN5dAFJzeXk5nU5//vnn8sKPHz9mWfb169cDFya1067uPwC60f9R18U5Bw8PD1+/fr26uqosv7q6Cg91Uqp0CFUAdGO5XIZmjOIPdnRxcfHy5csvX75Uln/+/Dk81EWhEiJUAcD5uLm5ubu7e3x8LJY8Pj5++PDh5uamw1IlIq3OzrXNy0nVBsAxSG3Yzb4tFovr6+s3b95cX19fXV19/vz5w4cPnz59uru7e/369YELk9qbm1xL1bJV16UDgJ28evXq7u7u4uIiy7KXL19mWfbixYtOElWC0oqQqUVmgOPnm3l/Oq/bzgtwYO79BwCnpzKgJanscrSEKgA4PSFFpdYUdOSSG1MFALAPQhUAQARCFQBABEIVAEAEQhUAQARCFQBABEIVAEAE5qk6CuZwA1JjgiXOj5aqo1DcedAtCIFELJfL4vfk2rvdw0kQqgDoRpGr/JjkPOj+A+DQKk1T5f8KWFFUWgHV6mEIVQAcWvm+dc76+6A+O6H7D4AOlAeqSwCcB6EKgA4IUrswtP84Jdf9134g+pADcPxCt2m5F7XrEtHrJRiqHHkAnAHD0Y5QcqEKAE6aayePllAFAKfEtZNHy0B1ADgxrp08TkIVAJwYQeo4CVUAABEIVQAAEQhV3TOHGwCcAaGqe8XlGz0BCwBOllB1FFwWCwCnzjxVHTOHGwCcB6GqY+ZwA4DzoPuve+ZwA4AzIFR1T5ACgDOg+w+AblQufPYLk1MnVAHQDSmKM6P7DwAgguRaqtpn1/SzCQDYTnKhSmwCAPYhuVAFAGfAMP8jJFQBwOmRoo6QgeoAABEIVWyqfYw/7dTejlTgLtTeLtQemxOqAAAiEKoAACIQqgAAIhCqAAAiEKoAACIwT9VRMIcbAJy6swpV4/E4z/Ner5dl2Xg87rg0TyFFAcCp65/N6TzLsvl8Xvx3MBiEgFXW75/P/h6e2tuF2tuRCtyF2tuF2ttFarV3JmOqxuPxfD4fjUbL5XK5XI5Go/l8Xg9VAAB7ciYRsj4Uqd/v1xurUovMcam9Xai9HanAXai9Xai9XaRWe2fSUtXr9QaDQeW/5d5AAIC9Op9QlWVZ10UAANJ1Dlf/NY6dqoxbL7g15i7U3i7U3o5U4C7U3i7UHhs6h1C1eRtVUj27AMAhnU/3X4VL/wCAQzqfUCVFAQAdOpNQVb/WL0xb1VV5AIDUnEmoCjelKQZXhT9O6041AMBJO5NQlWVZmEW93+/3+/35fD6bzbouFACQkP93Ns054SbKWZb9y7/8y3//93+/efOm6xLFEd6gYnfyPD+bXdu3ctUVN9tWexuqHHhhidngNldUYDj2CupwEz68W1h7sgg16QjcqzNpqSpkWXZmR0xlAP5wODQkf0NFRfX7/fD3cDg8m18R+1Y5zMbj8WQy6agsJ6mowMlk4jP7VI0fXtXYrv1kEWoyz3Nzbu3VuYWqIxQO6+LHVrEkHN9hyXg8rpzs8zwvP6XX62VZltqHIUrVFU8Zj8ej0Sid7+VYB15YmGCiiluB4+/2WubjEaX2wqNh4Ww2Gw6H+y94l/Z3shiPx+FmuHmeDwaDdI7DwxOq9ivP8+FwGI7vopkkdKMUXxD1HxDj8Tg8OhwOi6sa8zxPavLSWFWXZVlS9RZEPPDCf1MbpBi3AsMTkwr0sWqvuKnrmXVB1O31ZJGXev2SOhQ7sGSfwqmo+G/4ezAYDAaDsGQwGIxGo/D3aDQKyytPmc1mq/57xuJWXZEJ9l/woxCx9gaDQfgjndpbxqvAsJ3RaBRmeCmeft7i1l6x2nkfgXs9WVS+DBM5DjuhperQwk+E4kfDfD4v/4CYz+eV3xDFDzV2qbosy8JXydn/3l1lu9ob/zhZScq2q8DQUBq6bJbLZeMNSVOwde2NRqNwTffhyno0Ip4snEcORqjq2NpjPdlv4bU2rLpiOEJo9Fafwea1V8xU0iuNGsYndxeb114IoyGY7r1Yx23HQ64yTos9EaoOoRhs2Kv96A8zQYS/w0jqsILjPti96sJIhfLW0hGl9opm7V6vt1wuk2q12r0Cx+Nx0cri8Cts/r1XHjyUQnPLnk4W5ZHsk8lEQt2jLvocE1IM5QlfB6FHvNw1vvxxOGFYUgy/CP+mPKYqStWVV1N7xTpPOvCWSY6pilKBxX8dfjvW3uF35JD2fbIIw7N6yQzs68qZH6adK4YEtn+ZhiGZ9YV7KtVJiFt1jaudMQfejqJXYFK16sO7hQN8ZlOoxs791OMg2jtNGh9Nqp+lRZSqS7YyHXg7ilWBadaqD+8W9veZTa0mO2FM1d45jrem6nah9nakAneh9rag0s5Af5nevIgAANFpqQIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCIQKgCAIhAqAIAiECoAgCI4P8D4qzqNdBe32kAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "g['XT_MS3_DRY'].Print()\n",
    "g['XT_MS3_DRY'].Draw(\"AP*\")\n",
    "g['XT_MS3_DRY'].GetHistogram().GetYaxis().SetLimits(0.000001,0.5)\n",
    "g['XT_MS3_DRY'].GetHistogram().GetYaxis().SetRangeUser(0.000001,0.25)\n",
    "g['XT_MS3_DRY'].GetHistogram().GetYaxis().SetTitle(\"XT\")\n",
    "g['XT_OPT_DRY'].SetMarkerStyle(24)\n",
    "g['XT_OPT_DRY'].SetMarkerSize(1.2)\n",
    "g['XT_OPT_DRY'].Draw(\"PSAME\")\n",
    "g['XT_OPT2_DRY'].SetMarkerStyle(22)\n",
    "g['XT_OPT2_DRY'].SetMarkerSize(1.2)\n",
    "g['XT_OPT2_DRY'].Draw(\"PSAME\")\n",
    "l.Draw()\n",
    "tL.DrawLatexNDC(0.11,0.92,'CMS MTD - POST MS3 Optimisation - Array Dry coupling')\n",
    "\n",
    "c1.Draw()\n",
    "for ext in ['png','pdf']:\n",
    "    c1.SaveAs('plots/LYSO_OPT_XT_ARRAY_DRY_byProd.%s'%ext)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "id": "a8100d21-d680-4df6-b6f8-038778a36647",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x[0]=0, y[0]=276.448, ex[0]=0, ey[0]=5.52896\n",
      "x[1]=1, y[1]=277.934, ex[1]=0, ey[1]=5.55869\n",
      "x[2]=2, y[2]=270.587, ex[2]=0, ey[2]=5.41174\n",
      "x[3]=3, y[3]=279.149, ex[3]=0, ey[3]=5.58298\n",
      "x[4]=4, y[4]=275.991, ex[4]=0, ey[4]=5.51983\n",
      "x[5]=5, y[5]=283.418, ex[5]=0, ey[5]=5.66836\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file plots/LYSO_OPT_SR_ARRAY_DRY_byProd.png has been created\n",
      "Info in <TCanvas::Print>: pdf file plots/LYSO_OPT_SR_ARRAY_DRY_byProd.pdf has been created\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxwAAAI8CAIAAAD0vjrdAAAABmJLR0QAAAAAAAD5Q7t/AAAgAElEQVR4nO3dwZLbNp7HcbArl0lVqu05O6n0G2xqH4CgT7sPsk7v5phc4j2QumzmMDluUmnvq0wTfIGUnyD2bvW93ae+mXv4jzFokIQoCRQh4fs5uNwkJYEQJf4EgGDR970CAADAYS7WLgAAAMA5IFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQwWdrFyBFxhhjjPxfa621Hm7g/hneYHRt+PlHixR4xeFab8kUeZ7RjeeUar6p8oRfZaeKcjdWSjVNM6cAc8qza+XPL9Xwhfao9j0Op12f3z55Ok+V1GvtJ/wtkZT9Pr9nbHh0yZJsKwT/0OOpsiyHtVTXtbtNXdfu2rZtA2vdx3qrpp7f07atu3FZllvXju7F6ItOFWl01/YWLsnwVbydmtr3+RvPr5CthZ+5dqpUw2do29YrXvh42K+W9mbLduDBYEsbt3hTpmo7He6bHuuDtpA9PjVnbHgk22/RDGsDHrr//sEYUxRF13XDVZvNpiiKqQd6zQ+bzWZqs6lVm80m3Ibh8ko4sw1mP1VVLfr89lXcP40x3hKr6zrvjdhp4+gClTNVKqWUV6qqqrz3dLPZbP3Ju8SON01TFEVRFKf1g1trLcU+wrEakfumz//4p2bOsQrkg+6/f3BPUXVd264xm4Saphn97hvNYUP2eeq6ts+jtZaH75SrXKMnEvep3F0INJ+UZWkf5T6kqqp+W1PTfFOv4nZZbn0j3I3ts5Vl6S703rU9KmSK20k3p1ukbVu7C3bX7C7YLWVPbdm6rgv3Bu5US4dwj9VDnkdrvUdt7+2Yr7UH78M+8ztkde6RZoyRYm89Vs9e0zSy+zlXAgSh6u/c0489C6pPg1TkBDYn90ydxuxyN1WoT2do+W6a/8Xkbjn6dew+j9baPdEGntbda2m0kD+n0uQe3EE/8h+bIYYldMOc+0a4X+L2gcOhS/LMw7EO8ytk1JxQZUvlHUtt29pd8LaUMksS2no87FpL3hAQOSPuPQDLfTY7ost9tl2ff/RJRjcIbDNf4OVGKyrKi3qGjdbDt9srjJy53azsFtL7VhktthmMBAovH/IStv314v70Gha7aRqbcb3m1anlnjlvRHib0X0cLhw9tud/+4U/cYGncj8y898OpGjd3sd02PENowM+vOqyPejD4SbDIUrSy+7+bvaevP1kqmzuY+UV3UJ6JZlTfpe7L/NX7cqrjdFd87YcrZDhbgZ2LVCxOx3/8yt/zl54pRrWsK32wCGxdy0N22+Gz+k9+fAg9yrEqus68PyjY6pGh/R5OzU6Hm5YGO8VhwWY83Lu4Bhvs4jDZdxaCnzK3BK6ZXDfEe+xoztoSz71LTFaVzM38KrRLZL9v7e/9rG3t7dT++4aHgCBr5GpbUZ3wVvoFj4wNHbOmCpb5uE74h3ew0Gxc74BkCxC1d9NHfGjhqHK+yy53wJ2lfcRnf+BGf0K9laFA9DoF8pwX3Z64E5Gv5v6p+eAPhg9RSCE7VSrO+3X/Mq3q7wDYKcvx61lO7CWhoZlFltD1Rz2uBqeigI9dPa1AlcYjBYmkLa3PlV/rOs2vHPwsBhT9eyFKm+vt9bn6JHjJbZRUyXsn6aHqTK0bTv6PTNnfHd4j8LbDH//qHmhapQt5E6hatSw8rfuI04FA9V9u7a4yofHNthKx83ok7gfs81mU1WVjK5tmmaqF2lr8eY/MEHGGQDkVePUV9KwYt1vJbdWdaQRRaOvLs8caKV3G/m9UgXa/7XWtss18G27Ry25hic/KZIxxjvzbf0sjJ7L5UxgnyowWshWxbCh1/tAqbEWNdmmfxr7pt5086lTVY0lttE9HYbIWEeUO+BPOe9j4PmlAWNYTlkuz7O1Pkd7vmKNmRtV17UU276Kezx49TDkLpc9GnYX2mIPw7r0gO9X8mFgmrrMaKthoh1Wvn25OdcpI2nR4tmJ26lC3HOP13oh/59q7g78Dg40gHu/L+1L2P/LN1fgeQJ7d2BLVfvU1GYzj8M5HY72Ifbl2sGUBMNn3nW/pja2Fd47rZJTTUeBUg3rajSajNqjluyf7kOGxR595kBL1bBI7tHubTb8fe+Wyn033WNJju3RZolhNhrdzNsm3Cs0uiNzKnzmB6GfXe1u2bwyT/XiSYNQuOTDHRy+m0OBbbyXCDSj2mLf3t66Txuo1dF31i4MlG3mbnoLpwpvF3rxaGb337DkM0tFS9UpoqXqUO61bFsHGDZNI/Ve17V3xu26bqfB4O4vsLXGM2qtq6f2e55wG3jg1e1/jDH9p5OKV7FxZ1VwW9QCrZIiUKphXWmt3c1iTWbhFc89xlYfBuu2Zkljnv50CZUtmwxzts14h7+b3l7PaSXaWlE7fRDsW2B3f7QVZ/Qh4eXSIGT/LYpi2Lhi98W9jEM97cLeyfxR1V6Tz7AehkavX5GPVd/3arrVdjjw/BDDBtT53JJP9TN4haex6qQRqnyjHxvzSfiBssGcj4ScJ+SrwW4/s3nZnteXvkhkib5Faduz5KeYG4/kP1Nnl3CR5HQiFbv3SSJMP73kUM0749pSjXYBWLKZLfnUufnAWkpH0zTDnxaSrtz+EQkHB046MHp6XoUticzw5F2qdshltsaYoiiqqgpU19R5/fDLe6fSqmWPZ3txq7c8oq2fkbB1Y83qhygOQaj6O/spGv1yCf8A9VLO6EdieLmsu2qnokoJDxkusNMLqeDYnfKprc8pCcPy6mp0zEe4SFKrw3fNu8J8a8FmcltQvCXeZqOlkhpzS2VbYuw2W09ve9RSsmzWnGrJc2d3836E7GROi9Te5n8QvFfvPpnaYCe20mx36mgbsNtCtvc4IWuntOp258kDFzpEt447DFt32rBT+VGEccfqZ0yd++3j9WQPr+3yhhF4jSLeKJZ6cDeYwEtvLZsscV/OK0OsMVWBCtmDVxsBMysqML5n+KI7rdq6sXqqD75Bw6rzhomMlmRr8dxThbfKrcBhLY2OqRoODTnOmKrhCCSvJkfLM9y7pcdUDUu+tznn+KmyDZ/EO7oC78vUaK3ROhnaWpJyMEh8tK7cORTkP3Z81agzGFPllnzXUjGm6hTRUvV3bhNCVVXSeGCMcX/MTfUoBZpbLG8si/0tYpxZtvfrsYrVVtx1ne2MkE4EWV6W5TGbo92KspdGSjtfuEhun5FZcn5n97w49ZZ5x5Jbqqlf9vYSPDtsaGaDh1dL7gWVw0qwQ/fM0+ndp15oUfOH4nlTpO7BPnCz2Qz7zY/WpDe8mNEKN5bvxHy6OnVmQ9T8V9QO96ZeM5/B6wFU2/r+7BtkNzPOVZzCbbS22wRm9515taM79fGwPNENh0DgVB0/x6UsUFHuD6/h779hldo/7Y+SrZlpqlTD30/uOaB9egX73i1Vix4hw9oICI9bH+5guPyjv/Z22jtv42FT0OgP3JmlCuxsxFoKF2bqOdtt81TZB+7XUhWIMvX0pEde8bznGX4kbXkCLxfekVgtVYFW2NG1Uxtvbaka8krufeq3ljzwzLvW1eibtcdLj7ZKBgoWzs1e4cOvGLGlauYr4oTQUvVEP/HZq+s63M/tNoNPbdM492oYPrzfZU5F99fMcr9s6qczEx6N1no4yMYWafhGTG0sI0ui18/Myp86ltTg3jWjuXbrFZG71tLU07pv8ZF/JRtjpsrffLppj1c5cn2D/N+2XrhtCYHRMO5FANauH71D2Ddl9B0/vIFq+ObaJV61uK91yCUdUns7ldx9uTkvPfo58j7XU0PHwgMrA69ejs2ovsSnY3iEzxycimQVq5w10+eOOo/+WbKN8ws9/znZ6Y1Y9F3bm3cFw1TBDin81se6c4rqTzc1m9o4vHYJ4fLbteGTt81hM4+T8Ntxuma+ffaQOP4pYI+Xtt+Zgf2auU34yJekXpal+eQ4HwT3tfSnQQJL/CbE0ghVwPnzQtWqZUES5JCQ9HDM123b9uXLl6u89FZeqDrCK9oPozdKUv7D2fkU0f0HALkwT++GdMyELYFeEpU6rNvxbNg5NeRyE5mVTVbRCXiiPlu7AACA47Hjq8qyPHwg106v6770EnN+npy+722Kci/VTLAZDzMRqoDzZ4eu0/cHmWvj+OMOZQi21jrZRKU/TaxzzJrp+96OGhTHTLqIjjFVAAAAESQ6pkqP3c5F7mqix+7+sXUtAADAolJsqZILSr1ZRrQzFbUadDmH1wIAACwtuZaq4V0IlFJN00jMkhlL67p2byccXgsAAHAEybVUFUVRlqXXUiXXR7hFlc0kOYXXAgAAHEFaLVWBMORN2iHBa+ZaAACApSUUqqRdaqp5KXyNq7eW68YBAMCRpTJPlTFms9mM3hdzNGbZkek79fHZadYAAMARpDbKaFGphKqqqmReuOGqndqotsrq3QUAYEW5tWUkEars7eXdORSMMXbqqeFDwg1UDFEHAABHlkSoEu6dj9SnG0WN3sR7iBQFAADWlcRA9aZp+qeUUjLvlISq4dV8MueC/H90Lbf4BgAAx7RPqJKOOembk//ELpXP9g/Kn/If+7rhtQAAAEewQ/efZKlhm5D61HPn3VgmIq1127ZVVdkhb+51gqNrmVUBAAAc09wZ1WUKA7lAb5icjDEyJ4JaMlqpT2OnpgJTeK1SqiiSm0EeAIBzldtpd9beFkUxPypJpklz5Hhu7y4AACvK7bSb2d5m9u4CALCi3E67Ea7+S7NRCgAA4Jh2DlVN0xRFIUFK/i8jxE8lWhVBa5cOAHDO+sysXd/Htlu7nDGmqir16WYvkkLatpWrAtOvvtzaIQEAScntHHRxcZHVLu82o7qMVZcKkqYpmbygaZqqqowxTGQAAMBytl7nbjuOhtsYY+y8SFP328Uhdr5NjZ2p3H1f7RV/vEMAACxks9nY2YvsbUU8L1++lP98/PjRLjTG2OWu29tb78QtEyRJ9irLsq5rzuzz7TamSmarkv+7t+rbGpwBAEAs3t1yrdHxzTZR1XV9e3v78ePH29tbyWQvX750HyJb2hN913UvX76UYT+YY+dQJf/Kf+QtsQOtCFUAABzHaH4aDVuyUIKUPY/bti73IZK9JHgJpVTXdadyLdrqdg5VdV13XSezq8sQK0lUU+2QAAAgIjnheneNE3J2Hi5UYw0f3vNIcvLGWt3e3k69Fob2mVJBrpO0ubVt277vuYExAABHM2yUssOt5j+JdAXK/yU5jT6clqqZIkz+Sa8fAADHZIffuAunxjdL29XFxcVo56Ddvq7rjx8/eg+Xh3Cin2nnUCWX+DFzJgAAcV1czD0pD4dDqYm+P6VU27byn81mc3FxIekq3Pgk1wBWVSXPyQifmXabUsGOSR992wAAwK4uLi7c6Q/mc4c6hfv+Pn786M6VYNPY1NQM7vwLJKr5dmupspN/mjFLlA8AgPP28eNH20a1a2OVPflundtIa922rTufglJKmqNGi+ROuzA1gwM8O3f/nXobFff+AwCkxuaq+U1WXg/gVN/fkFzIb4eoT82YYDdT09NiwbP/5J8nils/AgBWdzEwtTzcduV2543208k4qkBssk9Cp9Ph9un+Y/YEAAAO8fEp5bRRfRyYehLbAxjo+5PmqzlNTZvNZvRWNphvt1BljCnLcrPZ0HcGAEAU7kD1nUas2x7AQN+fJK2pPj63iUu29OIXbVc72fmGyur0h1UBAJCO/S79swKTdorb29uXL1++fPlSZkuXk3jXdV6nYV3X9obN3jZ2glCE7RaqtNaEVgAAEiFJSG277k9yldxlznu4m8ZkM6+x6vb2lsk/ZyrmjM4uikLm/jr1W1UXxaz9BQBgCeueg4wxNlR59/gb3SywzUwXFxdZnXbnhoymaaSWTzpdEaoAACvK7RxEqNqibVs7Jk5rfVpXAhKqAAAryu0clFuo2nnyz6qqjDF939d1bYwpiuLkohUAAEB0O4cqy01XSqmiKOwtGwEAAHIzqzvMGHMeI/+3TqaVVSslAODIcjvL0P03oqoq6eY7g/kUuE0NAABYwqxQZe8NdE7pCgAAIKLdroazEyvIn3Vda61PqGeQq/8AACvK7RyUW/ffniFjmK5O4gJAQhUAYEVxz0Hu1OdTE3WO9iwNt5zTAbVHGwqhajduujq84uTZlFKj0zQcslYQqgAAK4p1DrI36fMMbylzcTE5zsfeo8YY8/Lly/ArlmW5xzX+hKodGGOaprHtVQc+lczSbm/i6D2hXLg3tVZr7d7PqCzL0dBNqAIArCjKOcjmpLqu7WnRNnB4uUo29m63PNzYi2jevZbVvresyS1UbbkablTbtvIuCkmvezyPyyuMxOG6ruVPeTn7Kt5aede9P0eLtN/+AgAQxceDyQmxLMvhKpuB3IXDJYGNxe3t7dRL7Cq30+4Oe9u2rRd1bY6JUI7Bs8k7av/vvTHz13rLYxUYAIBdHZhRJO6MJiEhp2k3DwW2J1RFN2tKhaZpiqKoqsq2B8qDIw5Ob9t2+Gy2pdFtGBvlbVCWpdsbCABA+raeueQsbKPVkJ0CKW7BMNOsUCXvou3mW+JCP5ufjDFyS0GllH0h+U9RFE3TSMJTTy9V8Dp6T2iWBwBA5u7u7v7zP//zz3/+c1VVf/7zn1+/fn13dze6paSl8DlOWhm2Xs3HfJNL+GzORm3bHi2myHB19XR8nNZaGp/sSDq7dtfDYuudaqb0WQ21AwAcxd3d3fX19VdfffW3v/3t6urq/fv3b968+fbbb3/77bcXL17s8YRy5VbXde6J2ztX2vOpN6oHB5oVqo7Z8NP3vTRWyfstbVRyiNjZsJqm2Ww2stmuZSMbAQDS8euvv3711Vd//etfP//8c6XUP/3TP/31r3/94Ycffvnll//6r/9ytzykbWl0xoSyLAlVcc0KVS6ZRmFq1YGlETJLu+QqeS03UanB1KPLlQQAgEX9+uuvf/vb3yRRiT/96U//9m//VlWVF6oOaeAYhqf9pkhA2G6hys4mFZcENS8J2amn7JSew7W2pYoUBQA4Ie6cnP/8z/8c3sZeqTeHnBC967dokTqOWQPVLWkrmrqS8JBydF3nNYDZnDQam9ykNbzWr+u6rRcMAgCwFjvpwLNnz37//XdvJoLff//98vLSm/tAfYpKo3OpW3MGs2Mhu4UqNWN2gz3Y6VxtcpKJ2m2yLsvSdgXatbYkdtyV+2wncS9CAEDmrq+v37x58/j4aJc8Pj7+z//8z/X19XBjOS0GQpV7bxKsYKdZrWQq8z2nxNrGK5g3e+dwJqphwaypGd6XKzwAAFsNp8f8v//7v3/913/993//999///3+/v7333//j//4j3/5l3/53//939HpNGPNqB7A5J97221MlUxtUBRFXdfDpsUDGxv7vve6/FyyanR8lSyxD6fNEwBwKl68ePHbb7/98ssvVVU9PDxcXl5eX1//9ttvX3755ej2bdtWVdV13cXFhQw296YcCkwNiqXtfIPhwDxPuz7V8XFDZQDAisLnoPnzBG02m2EnoMzR7S2U0e4fdxzq/vLly9Fn21VuN1TeLWTIBFFlWY6OWEq/iYhQBQBYUdxzkJ1dKNn5EQhVITKRwelWEKEKALCi3M5BuYWqncdULVMMAACA07bPmCp3cvPTsvXGf1kFagDAkeV2lsmtpWq3UCVTn0/dHyb9iqP7DwCwotzOQbmFqp3v/aeYVQwAAGAgr5YbWqoAACvK7RyUW0vVzrepAQAAwNCsUFUUxfyR6U3TcJEgAADIzazuMGOMvUfjVGZqmkZmISvL0t5tJjV0/wEAVhT3HOROqj41/+foGXnqXnBhe7SY5Nb9t8OdDtu2DQ9R925ynKCd9hcAgLgOv0WxsPdO9tze3g5vaTylrmvZZs7tAt37K9/e3to8UJbl8EWzvaHyPi035hP5U3+y6/McHy1VAIAVRTkHye38lFJ1XUu46brO3rLm9vbWPSPLxl4IG27s3UlQ/nQfZVvC5t92UOXXUpVXyCBUAQBWdPg5qKoqGWkzDDE27rgNVFM3VB7dWARuqCyr1NPoZnPbsP0st1DF1X8AAJwG28I02ixkY40Mgw6b6kAMk1ev69ptDJPew2SHUx9TXi033KYGALCiA88y0kzldfB5vKapqZaqwKpAS5UUYPiQqaeiperMhYeYrV06AEC+7GCpKdJQFB7ELKOstrYb7dewVNf1cFS7PBV3W1H73aYGAADEtTVRzaS17rqu6zo3e3kRquu64VD0mU8+XCijrPbrTzwzhCoAAFJxcXExNQ/CIYOWJPd4yrI8MAnZcetTs2TlhlAFAMDK5jRTHZJahuHp8Bgk46vUxHV/edohVMmc6Vpr95Y1EpyNMZvNhjFJAADMEUhRw1XhOTw9oyOc4oYeOx3D1PRU2drtNjVW3/fSa+stjFy62JinCgCwotFz0PwUJY1D4Zah+Vf/TQlc/aecRBW+CNG+elan3VnD4iRR1XXdtq1UcVEU8r62n2RVawAALGeqHUuy1HBCc8veqHehgknHlNyyhkFUQ3O7/+q6tr1+bdtWVeUuAQAAe9jpoj+tdVmWXddVVTU6o3pgatAoJM/R5Tdln4HqEk6JqAAAHGinjjn1qV2j67qLiwsZbC4xyzZfzbk78t4ktI0GQcZXqfmhiggFAEAK2raVgU0yH5Vdvm6s8YZZ52nWwO2iKNq2dXPVcMlJYKA6AGBFcc9B9laAyU4TldtA9ezmqQrf/i+r9x4AcNK01mlmqWzNDVXDiVyHS07irSU2AQCAJczt/pvzXOnnFbr/AAAryu0cRPffCKafBwAACMur5YaWKgDAinI7B+XWUrXDnGNH0DSNDLsbnVZU7jw4tTb8WAAAgEXt0HIjYcVGlqZp7FRjh8+ubm8vKJPryzWibtnsy8ksZ95a70aEZVkOx9ErWqoAADii7E67/Qx2MrGyLGWJHWVV17XEoLqu5zzVFK8w8or2OUf/9Apj18qfcjvC4ascUkgAADBfbqfdHa7+c2f7lCX2sdJQNOepAi/hNXcVRWEbnIbP7zabeYXxHuu9yiGFBAAA8+V22t3hhsreNFTuTbCbpqmqSsY87VeO0fnZ7ZKu67x7bnu9jd5a20UIAABwHNsHqkt7jxtivMFV6lP6GR3GNJPNT8YYY4w0PnkvIUPRi6LQWnuv5QWyk5iGFAAAnJN9blNju+QCS/Ymw9WVMzmWPLkdpV7X9WazqapKGrd2TXIzJzIdyqoBEwAA7Gp7qLKtUIHOuIj6vpfGKklRbmOVjTVN0xRFUVVV3/e7JjmyEQAAWMLceapsA5KkHC/KyNpYnW7S01eWpeQqeVovxoUneT+kIxIAAGAPs0KVTGFQFEVRFF4DkjQa2W32MzrCfU5Es+GJFAUAANY1K1Rprfu+lympZHYoWW476Uav3dtJ13XeBX1uThpezeeO4hquXbSDEgAAYMQKc2ONkcLYGTu9+TzDs316ayVOMfknAADryu20m9CsXN51ed7sne5dcYZr7V1uxFTLWW6zkAEAsKLcTrtp7a3NSVOdieG5G7bO7JDbuwsAwIpyO+1mtreZvbsAAKwot9Pu3CkVAAAAEECoAgAAiIBQBQAAEME+9/47aeF7/2XV9QsAACLKLlQRmwAAwBLo/gMAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIILsZ1blNDQAAWEJ2oYrYBAAAlkD3HwAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAIsrtNDff+AwAAS8guVBGbAADAEuj+AwAAiIBQBQAAEAGhCgAAIIK0QlXTNFprrXXTNIHNRjeY+VgAAIAlFIkM3DbGVFWllCrLUinVdZ2aGFSute66rq5rNzzJQvtnWZbGmOFjiyKV/QUA4OzldtpNpaVKElXf98YYY0zbtkqpYZuTMcYNT6JpGolZfd/3fV/Xddd1o6EKAABgIamEKqVUXdf2/1prpdQwGFVVJU1Zrs1mo5wEJv+hExAAABxTKvNUtW0rQcrlLbFJaziBp5e0yrIcNmgBAAAsJ5WWKpufpPtPYpPb2iR9fFNds6PxCwAA4GhSaamyZHCVetobaIzZbDYy0Mqz69ip8G1qArIaagcAAHaVXKiyY9XdkVIylGq0/WnXRimyEQAAWEJyoUopJdNNSa5qmkZylTcBlTHGTkw1fAYu/QMAAEeWRKiShOQlIW/qKWm4srqu67rOHYm1dCEBAAACUhmo3nWdNwmCzUlN0/RPKaVkVioJVcNr/bquG868AAAAsJwkQpVko81m4wYpmc9zzsNt/6D7bMxTBQAAjimJ7j+lVN/3RVHYS/+UUmVZzgxGWuu2bauqslf2jc56BQAAsJy0bspjW6r2i0Ty8MBjc7sJEQAAK8rttJvZ3mb27gIAsKLcTrtJjKkCAAA4dYQqAACACFIZqH404dvUZNVKCQAAIsouVBGbAADAEuj+AwAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABBBdrep4d5/AABgCdmFKmITAABYAt1/AAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABFkN6M6t6kBAABLyC5UEZsAAMAS6P4DAACIgFAFAAAQAaEKAAAgAkIVAABABGmFqqZptNZa66Zp4q4FAABYVJHI1XDGmKqqlFJlWSqlujs038gAABqISURBVK5TT6/Uk6kQ3LVt22qtZa3WWhaKsiyNMcNXKYpU9hcAgLOX22k3lZYqSVR93xtjjDFt2yqlbJuThKe2bWWtvEPyENms67q6rvu+7/u+ruuu60ZDFQAAwEJSiZBFUdR17fbcFUVhG5zc/4umaTabjRReGrG8Zq3RxqrcIjMAACvK7bSbyuSfbl+eZZeUZTlc65JuQfdPtzcQAABgaclFSGlesr2BU5u5rVPDVi63Hct7VGr7CwDAucrttJtKS5VlR0rVdT26gR3SLuOudh07Fb73X0BWhwUAANhVcqHKjlXfbDbKGasu7FV+trsw3C04+vyRSgoAAPAPqVz955K5psqylFwljDFFUdir/MJZikv/AADAkSURqowxgVHq6lOXX1mWfd+Pzu1JigIAAOtKIlQppbqu89KSm5NkENVUchpe69d1nXc9IAAAwKKSGFMljVKbzUbuM6Oc+TyVk6WGbVSypGmaqqq01rKlfYYjlBwAAEAkdK2jd12enb3TXu43ZAvvbTM665XK79pOAABWlNtpN629tY1Su17T5z488Njc3l0AAFaU22k3s73N7N0FAGBFuZ12UxmoDgAAcNIIVQAAABEkcfXfMYVvU5NVKyUAAIgou1BFbAIAAEug+w8AACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABNndpoZ7/wEAgCVkF6qITQAAYAl0/wEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEkN2M6tymBgAALCG7UEVsAgAAS6D7DwAAIAJCFQAAQASEKgAAgAgIVQAAABGkFaqaptFaa62bpom7FgAAYFFFIlfDGWOqqlJKlWWplOq6Tj29Uk9rLQtFWZbGmJlrraJIZX8BADh7uZ12U2mpkkTV970xxhjTtq1SyrY5NU3TdV1d133f931f13XXdTY2hdcCAAAcQSoRsiiKuq7dnruiKGyDk8zY6RZ1/lrvVRLZXwAAzl5up91UJv9s21Zr7S10l0i3oPun198XWAsAALC0VLr/bH6S7j9pfHIbroaRa/ThczYGAACILpWWKksGVyml6rqW/4yOjrIj03cdOxW+919AVg2YAICT5p3sOIUdR3Khyo5V32w26tNECYHtd22U4sACAJw9OdnlNqRpdal0/7lkrqmyLCVXjQo3UHHpHwAAOLIkQpUxJjxKXZGiAABA2pIIVUqpruu8mdDdnDS8mk8mpgqs9a4HBAAAWFQSoUoapTabzXA+T/unctqu5D/u1KCBtQAAAEeQ0BA271IFb/ZOex8b4c1rFV7rvkQ6+wsAwKJWP+utXoAjS2tvbYqauqZPNthvrcrv3QUA5Gz1s97qBTiyzPY2s3cXAJCz1c96qxfgyJIYUwUAAHDqCFUAAJyVu7u7169fP3/+XCn1/PnzH3/88e7ubu1CZSG7UFUErV06AAAOcnd39+23397f39/e3t7f37dt+/Dw8OrVK3LVEeTV2Zlb5y4AIDevX7++v7//+eefP//8c1ny+Pj4/fffX15e/vTTT0cuTG6n3cz2NrN3FwCQm+fPn9/e3n7zzTfuwrdv32qtP3z4cOTC5Hbaza77DwCAc/Xw8PDhw4erqytv+dXVlaxapVT5IFQBAHAmLi8vnz179v79e2/5u3fvZNUahcoIoQoAlsKlMDi+6+vrm5ubx8dHu+Tx8fHNmzfX19crlioTeXV25ta5CyAFfPPgmO7u7l69evX111+/evXq6urq3bt3b968+eOPP25ubr788ssjFya3g5+WKgAAzseLFy9ubm4uLy+11s+ePdNaf/HFF6skqgzlFSFzi8wAUsA3D9ay+rG3egGOjJYqAACACAhVAAAAEXy2dgGOLXwBTlatlAAAIKLsQhWxCQAALCG7UAUAwBnwOl5oMkgBoQoAgNMjKSq3y+sSx0D108C8zMDJubu7e/369fPnz5VSz58///HHH+/u7tYuFIAFEapOQ9/38lvE/gdAyu7u7r799tv7+/vb29v7+/u2bR8eHl69ekWuwnHYX+D8FD+mvJoNT72Z9NTLD+Tj9evX9/f3P//88+effy5LHh8fv//++8vLy59++mndsuGcJH5eSLx40WW2tyf+7p56+YF8PH/+/Pb29ptvvnEXvn37Vmv94cOHtUqF85P4eSHx4kVH9x8ARPbw8PDhw4erqytv+dXVlaxapVQAlkaoAoDILi8vnz179v79e2/5u3fvZNUahQKwOEIVAMR3fX19c3Pz+Pholzw+Pr558+b6+nrFUgFYVF6dnafeuXvq5QfycXd39+rVq6+//vrVq1dXV1fv3r178+bNH3/8cXNz8+WXX65dOpyPxM8LiRcvuuxaqoqgtUs3jtlugGRNfYe8ePHi5ubm8vJSa/3s2TOt9RdffEGiQizJnrAyl1eEPMXILLPdfPXVV/J79/379zc3N/Lvixcv1i4dAKW2fbec4jcP0ifHlaSrZA+w3A7+zPb2BN9dZrsB0keowioST1Qqv4M/s709wXeX2W624q6iWB2hCscU6PtL7UjL7eDPbkzVaWG2mzm4hw+ArLhfeu4SvgBXl1aoappGa621bpom7toTxWw3AIAhtwWILJWOhNrlpD2zLEulVNd1Sqm2bbXWU2vdkmutZaEoy9IYM/oS6ezvTD/++OPDw8PPP//8pz/9SZY8Pj7+8MMPX3zxxV/+8pd1y5aUU3xzcTbo/sNaEj+6Ei9edKm0VEl4atvWGGOMkfegqqrRtW3bKqVsi1TTNF3X1XUtjZ91XXddNxqqkhWY2eG77757//79999///bt24eHh7dv3/7www9//PHHd999t1ZpAQDAUCoRsigKr3mpaZrNZiPFG17g4G4fXuu9SiL7O2qqeHd3d//93//966+/Pjw8XF5eXl9ff/fdd8x240n8zT1dXAcwBy1V++HoOlziR1fixYsulZaqsiy9sVBuJJJev/DDvT/d3sBT9+LFi59++kmGpX/48OEvf/kLiQpHw3UAh7BtzylPL7yi8NF1EpMzA65UQpUxxg6fEm4qkrxVFEXTNE3TyKfLTV3eY70/AWAJW8/0/VPHKdXZINDj5Hy2dgFGGGNkNJWMnVJKaa2l8Wmz2ciSuq7txjs9+d4/d/hIA/DIfNZ2lALfEkDmUmmpsrTWNlHZBie5uM8dir7ZbGTtro1S/b7i7iaA89Anf58QAEeTUKgyxhRFYcOTm5ZkoXu5X3jU1Gld+odDMNICa/HG+jD6B0AqoUq6/Mqy7Pt+dMT66KgpG55IUdmy7QSKgIXjckf8uEtosgKylUqoki6/0Wzk5SfhJq1hq1XXdVsvGMTZoP8Fa2FWawCuJAaq28A0vMOMLCnLUoaoy58y26eNTU3TVFWltXaT1gndrIbxrXvzmqbcP6lSHAGHGVbkNdJzNKYgidO5vdxvyBYvfCMa7xncQe6uZOOLFCzwwWCKvICttYe93d3d/fLLL7/++uuHDx+ePXsmE8++ePFi7XIlJ9nvlpSFL5zk2DsPuX00TmxvR8dXzVyr0n53yQT74YL25dzd3X377bdfffXVq1evrq6u3r9/f3NzI/9ybvNw+O1n6hcRx97ZyO2jkdnepvfuBsZWp1bUxAXeXNr59vP69ev7+/uff/75888/lyWPj4/ff//95eXlTz/9tG7ZUpPgd8upGP09ybF3NnL7aGS2t6m+u3RgHW7rm5vsu7+uQOJ8/vz57e3tN998427w9u1brbXcNAkWR9euwr8nOfbORm4fjVSu/ssZFxBhRe68AO7h9/Dw8OHDh6urK2/7q6srWXXUUuLsBCak4NjD6SJUrY8ghQRdXl4+e/bs/fv33vJ3797JqjUKhbMy9XuSYw+nK7tQVQStXTogIdfX1zc3N4+Pj3bJ4+Pjmzdvrq+vVywVzkbg9yTHHk5UXp2diXfuJl68xDGm6hBT17S/evXq66+/liuw3r179+bNmz/++OPm5ubLL79cpZzJ4ug6xLD2OPbORm4fjexaqgDM9OLFi5ubm8vLS631s2fPtNZffPEFZzUcAcceTlReETLxyJx48RJHS9UhqL29MWHH4cJHF8feScvt7UviNjUAcKKyOmEACKP7D2fu7u7u9evXz58/V0o9f/78xx9/vLu7W7tQAIAzRKjCybNXbg4v4ZSbXdzf39/e3t7f37dt+/Dw8OrVK3IVACC6vDo7E+/cTbx4p4ibXYTNv3MiByeWw5iqM5bb25fZ3qb97iZevFPEzS62mnmLJA5OLIdQdcZye/vo/ktCoAMLe+NmF3Nw00kAiIVQlYT+qbWLcya42UWYdyOBqVsLkPixHI4unJnsQhW3qckKN7sICNzRdrgZiR9LCB9dRC6cnLw6O3Pr3AU3uwibP1AdAPaQ23dLdi1VyAo3uwjL6ssOAJaWV4TMLTLDxbsfQOUAWEJu3y20VAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAIPlu7AMDi3DuxKCZnAgAsg1CF80eKCiBxAkAs2YWq8A2kOKMgNxzzABBLdqGKUwgAAFgCA9UBAAAiIFQBAABEQKgCAACIIK1Q1TSN1lpr3TTNcK0xJrA2/FgAAIBFFekM3Jbr8sqyVEp1XaeUattWay1rm6bZbDaygax1S661loWiLEtjzOhLpLO/AACct9xOu6m0VEl4atvWGGOMkfegqipZa4zZbDZ1Xfd9b4xp29Y+RCnVNE3XdbK27/u6rruuGw1VAAAAC0klQhZF4TUvSdOUFE8aotyiSh+f/DuctHD4bHZ5IvsLAMDZy+20m8o8VWVZemOh3EjUdZ10C1rext5a20UIAABwHKl0/8kgdHeJl4pkBLrWuigKrbXXCuU91vsTAABgaamEKpcxRnr0ZOyU5KfNZiMD1WXIVFVVsnzXsVPFviLvJAAAOC+pdP9Z9jo+99I/Yftlm6YpiqKqqr7vd22UyqpzFwAAHE1CLVXSQGWv47NpSf7jjZqq6zr8VAsVEgAAYFQqLVXGmKqqpuaXCjxKIhcpCgAArCuVax2H0yK4hlMquEuGa5lSAQCA1eV22k2ipcqmn+EdZux8VFVV2Yv+7Gyfdht3rbRdcbMaAABwTElESOn7G13lDk6Xq/+E1xDlPcNwkLvILTIDALCi3E67J7a3blvUrmtVfu8uAAAryu20m9neZvbuAgCwotxOuwlNqQAAAHC6CFUAAAARJHH13zGFbziTVSslAACIKLtQRWwCAABLoPsPAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQATZ3aaGe/8BAIAlZBeqiE0AAGAJdP8BAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARJDdjOrcpgYAACwhu1BFbAIAAEug+w8AACACQhUAAEAEhCoAAIAICFUAAAARpBWqmqbRWmutm6YJbDa6wczHAgAALKFI52o4meygLEulVNd1Sqm2bbXW3mZa667r6rp2w5MstH+WZWmMGX2JdPYXAIDzlttpN5WWKglPbdsaY4wx8h5UVeVtZoxxw5NomkZiVt/3fd/Xdd113WioAgAAWEgqEbIoCq95qWmazWbjFU8281qqpInL3XL4bHZ5IvsLAMDZy+20m0pLVVmW3lioYSSS1qzRJijpNHT/HDZoAQAALCeVGdWHUclLRdLHNxV4vaFX3hArAACApaXSUuUyxkiPXtu2dslms7F/ehvv9OTFvg7eLQAAcM5SaamybCOTe+lfVVVlWQ6vBFSDNqqtsurcBQAAR5NQqDLGyOV+3nQJ8n9vAipjjJ2YavSpFi0qAACAJ5VQJYlqan4ppdRms3H/7Lqu6zqbqEhRAABgXalc6zicFiG8sduaJT2GTKkAAEBScjvtJtFSZdPP6M1ntj68aZqqqrTW8jzSfMXNagAAwDElEaosr49PzctGWuu2bauqstfojd7fBgAAYDln1S7ntlSNyq0dEgCAFeV22s1sbzN7dwEAWFFup90UJ/8EAAA4OYQqAACACNIaqH4E4RvOZNVKCQAAIsouVBGbAADAEuj+AwAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABBBdrep4d5/AABgCdmFKmITAABYAt1/AAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABFkN6M6t6kBAABLyC5UEZsAAMAS6P4DAACIgFAFAAAQAaEKAAAgAkIVAABABGmFqqZptNZa66Zp4q4FAABYVJHO1XAy2UFZlkqpruuUUm3baq3nrNVay0JRlqUxZvQl0tlfAADOW26n3VRaqiQetW1rjDHGyHtQVdWctU3TdF1X13Xf933f13Xddd1oqAIAAFhIKhGyKAqvealpms1mI8XbulY9nYBquL1dnsj+AgBw9nI77abSUlWWpTcWyo1EZVnanr6ph3t/ur2BAAAAS0tlRvVhq5KbioZrN5uN+6cXubwhVgAAAEtLpaXKZYyRHr22bbeu3XXsVLGvCDsGAADOVyotVZZtZHIv7gusDXcLDmXVuQsAAI4moZYqaYKy1/F5aSm8dvhUS5YUAADAl0pLlTGmqqqp+aXCaxUpCgAArC2Vax2H0yLMXyt9gkypAABAUnI77SbRUmXTz/AOM03ThNfKv1VVaa1lS+kZ5GY1AADgmJKIkNK7N7qq7/vw2tFnGB3krvKLzAAArCi30+5Z7a3bUjUqt3cXAIAV5XbazWxvM3t3AQBYUW6n3YSmVAAAADhdhCoAAIAIkrj675jCN5zJqpUSAABElF2oIjYBAIAl0P0HAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEEF2M6pzmxoAALCE7EIVsQkAACyB7j8AAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqzBWejB5h1N6BqMBDUHuHoPYwH6EKAAAgguxuU8O9/wAAwBKyC1XEJgAAsAS6/wAAACIgVAEAAERAqAIAAIggrVDVNI3WWmvdNE3ctQAAAIsq0hm4LdfllWWplOq6TinVtq3WWtZqrWWhKMvSGGP/DK91XyKd/T051N4hqL0DUYGHoPYOQe0dIrfaS6WlSsJT27bGGGOMvAdVVcnapmm6rqvruu/7vu/ruu66zsam8FoAAIAjSCVCFkXhNS81TbPZbKR40ojlFtXdPrzWe5VE9vcUUXuHoPYORAUegto7BLV3iNxqL5WWqrIsvbFQXiSSbkH3T6+/L7AWAABgaalM/jlsVfJSkR1cNcpb6w2xAgAAWFoqocpljJHRVG3bqrG8pZzYtOvYKW6NeQhq7xDU3oGowENQe4eg9jBTcqHKpiV76d9ObVRhWfXsAgCAY0plTJVSyhhTFIW9ji+clsINVFz6BwAAjiyVlirp8puaX0qRogAAQNpSudZxOC2CS/oEvUkT6rqWCwZH1wbyGQAAQHRJhCo7Mr2ua2+VxCavHctLUaNr3dnYAQAAlpZWqBrykpNd7mWm8FoAAIClJRGq5rNtUXusPVH2RtHypzHmzHZwOW7V2Z5iam8m78CTJdytfL6pW7xTh3Pw4d3D1pOF1CRH4KISuvpvjvDn6iw/dd7IsKqqGCs2k62ooijk/1VV8YUyk3eYyW2jVirLSbIVuNls+MzuavTDSzWGhU8WUpNylf2RC5aVEwtVp0gO66Zp7PEt/5HjW5YM2wCMMe5DlFJa69w+DFGqzj6kaZq6rvP5Xo514MnCDBNV3ApsPlm0zOmIUnt2TG3TNG3bTo0SORvLnSyappFhx8aY4U3hEBGhalnGmKqq5Pi2zSTSSGu/IIY/IJqmkbVVVdn77RhjTquv9kCxqk5rnVW9iYgHnvwptzfIR9wKlAdmFehj1Z69r+v59UJ4Fj1ZGKcrMKtDcQU9liSnIvun/L8sy7IsZUlZljLZad/3dV3Lcu8hbdtO/XnG4ladzQTLFzwJEWuvLEv5Tz6118erQHmeuq7l0mb78PMWt/bsZud9BC56svC+DDM5DldBS9Wxmaej6buuc39AdF3n/YawP9RwSNVpreWr5Ox/707Zr/bsGOEjlTJh+1WgNJRKl03f99ne6H3v2qvruiiK3EY+iIgnC84jR0OoWtnWYz3bb+GtZladHY4gjd7Up5hfe13X2bOaHTUMPrmHmF97EkYlmC5erLQdeMh547SwEELVMdjBhmrwo9+94lpGUssGHPfi8KozzhxmudVqlNqzzdpKqX7bTTnPzOEV2DSNbWXh8LPmf++5g4dyaG5Z6GThjmTfbDYk1AWt0eeYETuUR74OpEfc7Rrvnw4nlCV2+IX8m/OYqihV525G7dltdjrw+izHVEWpQPsnh9+BtXf8HTmmpU8WMjxLZTOwby1nfpiuzg4JDH+ZypDM4cKFSnUS4lbd6GZnjAPvQNErMKta5cO7hyN8ZnOoxtV9pnAU4U6T0bVZ9bMERKm6bCuTA+9AsSowz1rlw7uH5T6zudXkKhhTtTiO471RdYeg9g5EBR6C2tsDlXYGTuzefwAAAGmipQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARECoAgAAiIBQBQAAEAGhCgAAIAJCFQAAQASEKgAAgAgIVQAAABEQqgAAACIgVAEAAERAqAIAAIiAUAUAABABoQoAACACQhUAAEAEhCoAAIAICFUAAAAREKoAAAAiIFQBAABEQKgCAACIgFAFAAAQAaEKAAAgAkIVAABABIQqAACACAhVAAAAERCqAAAAIiBUAQAARPD/6nFR/ecMJPEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "g['SR_MS3_DRY'].Print()\n",
    "g['SR_MS3_DRY'].Draw(\"AP*\")\n",
    "g['SR_MS3_DRY'].GetHistogram().GetYaxis().SetLimits(200,400)\n",
    "g['SR_MS3_DRY'].GetHistogram().GetYaxis().SetRangeUser(200,400)\n",
    "g['SR_MS3_DRY'].GetHistogram().GetYaxis().SetTitle(\"SR (mV/ns)\")\n",
    "g['SR_OPT_DRY'].SetMarkerStyle(24)\n",
    "g['SR_OPT_DRY'].SetMarkerSize(1.2)\n",
    "g['SR_OPT_DRY'].Draw(\"PSAME\")\n",
    "l.Draw()\n",
    "tL.DrawLatexNDC(0.11,0.92,'CMS MTD - POST MS3 Optimisation - Array Dry coupling')\n",
    "\n",
    "c1.Draw()\n",
    "for ext in ['png','pdf']:\n",
    "    c1.SaveAs('plots/LYSO_OPT_SR_ARRAY_DRY_byProd.%s'%ext)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "id": "412cc606-aa23-4fec-9311-40ae1e2f0c21",
   "metadata": {},
   "outputs": [],
   "source": [
    "fOut=R.TFile('LYSO_ARRAY_DRY_OPT_OPT2.root','RECREATE')\n",
    "for hn,h in g.items():\n",
    "    h.Write(hn)\n",
    "fOut.Close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8cc8348",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "217aa715",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "210ac3e6",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f2e08def",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ae8c2ed7",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}