{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "dfad2597",
   "metadata": {},
   "source": [
    "# df104_HiggsToTwoPhotons\n",
    "The Higgs to two photons analysis from the ATLAS Open Data 2020 release, with RDataFrame.\n",
    "\n",
    "This tutorial is the Higgs to two photons analysis from the ATLAS Open Data release in 2020\n",
    "(http://opendata.atlas.cern/release/2020/documentation/). The data was taken with the ATLAS detector\n",
    "during 2016 at a center-of-mass energy of 13 TeV. Although the Higgs to two photons decay is very rare,\n",
    "the contribution of the Higgs can be seen as a narrow peak around 125 GeV because of the excellent\n",
    "reconstruction and identification efficiency of photons at the ATLAS experiment.\n",
    "\n",
    "The analysis is translated to a RDataFrame workflow processing 1.7 GB of simulated events and data.\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "**Author:** Stefan Wunsch (KIT, CERN)  \n",
    "<i><small>This notebook tutorial was automatically generated with <a href= \"https://github.com/root-project/root/blob/master/documentation/doxygen/converttonotebook.py\">ROOTBOOK-izer</a> from the macro found in the ROOT repository  on Tuesday, May 19, 2026 at 08:10 PM.</small></i>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "52a0e1e2",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:33.012580Z",
     "iopub.status.busy": "2026-05-19T20:10:33.012466Z",
     "iopub.status.idle": "2026-05-19T20:10:34.001253Z",
     "shell.execute_reply": "2026-05-19T20:10:33.998978Z"
    }
   },
   "outputs": [],
   "source": [
    "import ROOT\n",
    "import os"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2662ed1d",
   "metadata": {},
   "source": [
    "Enable multi-threading"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "4478ec1b",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:34.013448Z",
     "iopub.status.busy": "2026-05-19T20:10:34.013300Z",
     "iopub.status.idle": "2026-05-19T20:10:34.139803Z",
     "shell.execute_reply": "2026-05-19T20:10:34.139406Z"
    }
   },
   "outputs": [],
   "source": [
    "ROOT.ROOT.EnableImplicitMT()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "257363a9",
   "metadata": {},
   "source": [
    "Create a ROOT dataframe for each dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d7a5b3b2",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:34.154622Z",
     "iopub.status.busy": "2026-05-19T20:10:34.154461Z",
     "iopub.status.idle": "2026-05-19T20:10:36.714702Z",
     "shell.execute_reply": "2026-05-19T20:10:36.708679Z"
    }
   },
   "outputs": [],
   "source": [
    "path = \"root://eospublic.cern.ch//eos/opendata/atlas/OutreachDatasets/2020-01-22\"\n",
    "df = {}\n",
    "df[\"data\"] = ROOT.RDataFrame(\"mini\", (os.path.join(path, \"GamGam/Data/data_{}.GamGam.root\".format(x)) for x in (\"A\", \"B\", \"C\", \"D\")))\n",
    "df[\"ggH\"] = ROOT.RDataFrame(\"mini\", os.path.join(path, \"GamGam/MC/mc_343981.ggH125_gamgam.GamGam.root\"))\n",
    "df[\"VBF\"] = ROOT.RDataFrame(\"mini\", os.path.join(path, \"GamGam/MC/mc_345041.VBFH125_gamgam.GamGam.root\"))\n",
    "processes = list(df.keys())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e44bab0",
   "metadata": {},
   "source": [
    "Apply scale factors and MC weight for simulated events and a weight of 1 for the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "5d292aa0",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:36.726223Z",
     "iopub.status.busy": "2026-05-19T20:10:36.726038Z",
     "iopub.status.idle": "2026-05-19T20:10:36.958279Z",
     "shell.execute_reply": "2026-05-19T20:10:36.947569Z"
    }
   },
   "outputs": [],
   "source": [
    "for p in [\"ggH\", \"VBF\"]:\n",
    "    df[p] = df[p].Define(\"weight\",\n",
    "            \"scaleFactor_PHOTON * scaleFactor_PhotonTRIGGER * scaleFactor_PILEUP * mcWeight\");\n",
    "df[\"data\"] = df[\"data\"].Define(\"weight\", \"1.0\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1c634fd1",
   "metadata": {},
   "source": [
    "Select the events for the analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "8a2dce74",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:36.979354Z",
     "iopub.status.busy": "2026-05-19T20:10:36.979177Z",
     "iopub.status.idle": "2026-05-19T20:10:37.351800Z",
     "shell.execute_reply": "2026-05-19T20:10:37.336535Z"
    }
   },
   "outputs": [],
   "source": [
    "for p in processes:\n",
    "    # Apply preselection cut on photon trigger\n",
    "    df[p] = df[p].Filter(\"trigP\")\n",
    "\n",
    "    # Find two good muons with tight ID, pt > 25 GeV and not in the transition region between barrel and encap\n",
    "    df[p] = df[p].Define(\"goodphotons\", \"photon_isTightID && (photon_pt > 25000) && (abs(photon_eta) < 2.37) && ((abs(photon_eta) < 1.37) || (abs(photon_eta) > 1.52))\")\\\n",
    "                 .Filter(\"Sum(goodphotons) == 2\")\n",
    "\n",
    "    # Take only isolated photons\n",
    "    df[p] = df[p].Filter(\"Sum(photon_ptcone30[goodphotons] / photon_pt[goodphotons] < 0.065) == 2\")\\\n",
    "                 .Filter(\"Sum(photon_etcone20[goodphotons] / photon_pt[goodphotons] < 0.065) == 2\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7d18b67a",
   "metadata": {},
   "source": [
    "Compile a function to compute the invariant mass of the diphoton system"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "34f8eb84",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:37.378328Z",
     "iopub.status.busy": "2026-05-19T20:10:37.378098Z",
     "iopub.status.idle": "2026-05-19T20:10:37.562754Z",
     "shell.execute_reply": "2026-05-19T20:10:37.551718Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ROOT.gInterpreter.Declare(\n",
    "\"\"\"\n",
    "using namespace ROOT;\n",
    "float ComputeInvariantMass(RVecF pt, RVecF eta, RVecF phi, RVecF e) {\n",
    "    ROOT::Math::PtEtaPhiEVector p1(pt[0], eta[0], phi[0], e[0]);\n",
    "    ROOT::Math::PtEtaPhiEVector p2(pt[1], eta[1], phi[1], e[1]);\n",
    "    return (p1 + p2).mass() / 1000.0;\n",
    "}\n",
    "\"\"\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0aee2a6f",
   "metadata": {},
   "source": [
    "Define a new column with the invariant mass and perform final event selection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "7abbe325",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:37.565165Z",
     "iopub.status.busy": "2026-05-19T20:10:37.565033Z",
     "iopub.status.idle": "2026-05-19T20:10:37.929862Z",
     "shell.execute_reply": "2026-05-19T20:10:37.909676Z"
    }
   },
   "outputs": [],
   "source": [
    "hists = {}\n",
    "for p in processes:\n",
    "    # Make four vectors and compute invariant mass\n",
    "    df[p] = df[p].Define(\"m_yy\", \"ComputeInvariantMass(photon_pt[goodphotons], photon_eta[goodphotons], photon_phi[goodphotons], photon_E[goodphotons])\")\n",
    "\n",
    "    # Make additional kinematic cuts and select mass window\n",
    "    df[p] = df[p].Filter(\"photon_pt[goodphotons][0] / 1000.0 / m_yy > 0.35\")\\\n",
    "                 .Filter(\"photon_pt[goodphotons][1] / 1000.0 / m_yy > 0.25\")\\\n",
    "                 .Filter(\"m_yy > 105 && m_yy < 160\")\n",
    "\n",
    "    # Book histogram of the invariant mass with this selection\n",
    "    hists[p] = df[p].Histo1D(\n",
    "            ROOT.RDF.TH1DModel(p, \"Diphoton invariant mass; m_{#gamma#gamma} [GeV];Events\", 30, 105, 160),\n",
    "            \"m_yy\", \"weight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0837308",
   "metadata": {},
   "source": [
    "Run the event loop"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "594e0d39",
   "metadata": {},
   "source": [
    "RunGraphs allows to run the event loops of the separate RDataFrame graphs\n",
    "concurrently. This results in an improved usage of the available resources\n",
    "if each separate RDataFrame can not utilize all available resources, e.g.,\n",
    "because not enough data is available."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "304d2be9",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:37.940430Z",
     "iopub.status.busy": "2026-05-19T20:10:37.940272Z",
     "iopub.status.idle": "2026-05-19T20:10:50.042499Z",
     "shell.execute_reply": "2026-05-19T20:10:50.042005Z"
    }
   },
   "outputs": [],
   "source": [
    "ROOT.RDF.RunGraphs([hists[s] for s in [\"ggH\", \"VBF\", \"data\"]])\n",
    "\n",
    "ggh = hists[\"ggH\"].GetValue()\n",
    "vbf = hists[\"VBF\"].GetValue()\n",
    "data = hists[\"data\"].GetValue()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "142c1328",
   "metadata": {},
   "source": [
    "Create the plot"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "491b6b96",
   "metadata": {},
   "source": [
    "Set styles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "3a03b1a5",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.054838Z",
     "iopub.status.busy": "2026-05-19T20:10:50.054691Z",
     "iopub.status.idle": "2026-05-19T20:10:50.160694Z",
     "shell.execute_reply": "2026-05-19T20:10:50.160238Z"
    }
   },
   "outputs": [],
   "source": [
    "ROOT.gROOT.SetStyle(\"ATLAS\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e509683",
   "metadata": {},
   "source": [
    "Create canvas with pads for main plot and data/MC ratio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "059d62b9",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.170370Z",
     "iopub.status.busy": "2026-05-19T20:10:50.170215Z",
     "iopub.status.idle": "2026-05-19T20:10:50.334128Z",
     "shell.execute_reply": "2026-05-19T20:10:50.333398Z"
    }
   },
   "outputs": [],
   "source": [
    "c = ROOT.TCanvas(\"c\", \"\", 700, 750)\n",
    "\n",
    "upper_pad = ROOT.TPad(\"upper_pad\", \"\", 0, 0.35, 1, 1)\n",
    "lower_pad = ROOT.TPad(\"lower_pad\", \"\", 0, 0, 1, 0.35)\n",
    "for p in [upper_pad, lower_pad]:\n",
    "    p.SetLeftMargin(0.14)\n",
    "    p.SetRightMargin(0.05)\n",
    "    p.SetTickx(False)\n",
    "    p.SetTicky(False)\n",
    "upper_pad.SetBottomMargin(0)\n",
    "lower_pad.SetTopMargin(0)\n",
    "lower_pad.SetBottomMargin(0.3)\n",
    "\n",
    "upper_pad.Draw()\n",
    "lower_pad.Draw()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f6eb6600",
   "metadata": {},
   "source": [
    "Fit signal + background model to data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "b90d80a0",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.343243Z",
     "iopub.status.busy": "2026-05-19T20:10:50.343087Z",
     "iopub.status.idle": "2026-05-19T20:10:50.525867Z",
     "shell.execute_reply": "2026-05-19T20:10:50.525371Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "****************************************\n",
      "Minimizer is Minuit2 / Migrad\n",
      "Chi2                      =      19.9699\n",
      "NDf                       =           26\n",
      "Edm                       =  3.10817e-10\n",
      "NCalls                    =          161\n",
      "p0                        =        94325   +/-   72.0526     \n",
      "p1                        =     -1777.22   +/-   0.778157    \n",
      "p2                        =      11.5606   +/-   0.00536058  \n",
      "p3                        =   -0.0256281   +/-   2.66824e-05 \n",
      "p4                        =        119.1                      \t (fixed)\n",
      "p5                        =          125                      \t (fixed)\n",
      "p6                        =         2.39                      \t (fixed)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<cppyy.gbl.TFitResultPtr object at 0x5564afb53680>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n",
      "Error in <TFitResultPtr>: TFitResult is empty - use the fit option S\n"
     ]
    }
   ],
   "source": [
    "fit = ROOT.TF1(\"fit\", \"([0]+[1]*x+[2]*x^2+[3]*x^3)+[4]*exp(-0.5*((x-[5])/[6])^2)\", 105, 160)\n",
    "fit.FixParameter(5, 125.0)\n",
    "fit.FixParameter(4, 119.1)\n",
    "fit.FixParameter(6, 2.39)\n",
    "fit.SetLineColor(2)\n",
    "fit.SetLineStyle(ROOT.kSolid)\n",
    "fit.SetLineWidth(2)\n",
    "data.Fit(\"fit\", \"0\", \"\", 105, 160)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8aa18962",
   "metadata": {},
   "source": [
    "Draw data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "590f375e",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.542576Z",
     "iopub.status.busy": "2026-05-19T20:10:50.542424Z",
     "iopub.status.idle": "2026-05-19T20:10:50.677950Z",
     "shell.execute_reply": "2026-05-19T20:10:50.677549Z"
    }
   },
   "outputs": [],
   "source": [
    "upper_pad.cd()\n",
    "data.SetMarkerStyle(20)\n",
    "data.SetMarkerSize(1.2)\n",
    "data.SetLineWidth(2)\n",
    "data.SetLineColor(\"black\")\n",
    "data.SetMinimum(1e-3)\n",
    "data.SetMaximum(8e3)\n",
    "data.GetYaxis().SetLabelSize(0.045)\n",
    "data.GetYaxis().SetTitleSize(0.05)\n",
    "data.SetStats(0)\n",
    "data.SetTitle(\"\")\n",
    "data.Draw(\"E\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18e11ef9",
   "metadata": {},
   "source": [
    "Draw fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "99e722d1",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.682553Z",
     "iopub.status.busy": "2026-05-19T20:10:50.682421Z",
     "iopub.status.idle": "2026-05-19T20:10:50.797207Z",
     "shell.execute_reply": "2026-05-19T20:10:50.796716Z"
    }
   },
   "outputs": [],
   "source": [
    "fit.Draw(\"SAME\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4030c166",
   "metadata": {},
   "source": [
    "Draw background"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "6cc9b596",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.799581Z",
     "iopub.status.busy": "2026-05-19T20:10:50.799387Z",
     "iopub.status.idle": "2026-05-19T20:10:50.918647Z",
     "shell.execute_reply": "2026-05-19T20:10:50.918173Z"
    }
   },
   "outputs": [],
   "source": [
    "bkg = ROOT.TF1(\"bkg\", \"([0]+[1]*x+[2]*x^2+[3]*x^3)\", 105, 160)\n",
    "for i in range(4):\n",
    "    bkg.SetParameter(i, fit.GetParameter(i))\n",
    "bkg.SetLineColor(4)\n",
    "bkg.SetLineStyle(ROOT.kDashed)\n",
    "bkg.SetLineWidth(2)\n",
    "bkg.Draw(\"SAME\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d25ddd1d",
   "metadata": {},
   "source": [
    "Scale simulated events with luminosity * cross-section / sum of weights\n",
    "and merge to single Higgs signal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "fca64e84",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:50.943523Z",
     "iopub.status.busy": "2026-05-19T20:10:50.943370Z",
     "iopub.status.idle": "2026-05-19T20:10:51.108813Z",
     "shell.execute_reply": "2026-05-19T20:10:51.099650Z"
    }
   },
   "outputs": [],
   "source": [
    "lumi = 10064.0\n",
    "ggh.Scale(lumi * 0.102 / 55922617.6297)\n",
    "vbf.Scale(lumi * 0.008518764 / 3441426.13711)\n",
    "higgs = ggh.Clone()\n",
    "higgs.Add(vbf)\n",
    "higgs.Draw(\"HIST SAME\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2a276931",
   "metadata": {},
   "source": [
    "Draw ratio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "ea94a438",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:51.114245Z",
     "iopub.status.busy": "2026-05-19T20:10:51.114096Z",
     "iopub.status.idle": "2026-05-19T20:10:51.291440Z",
     "shell.execute_reply": "2026-05-19T20:10:51.290853Z"
    }
   },
   "outputs": [],
   "source": [
    "lower_pad.cd()\n",
    "\n",
    "ratiobkg = ROOT.TH1I(\"zero\", \"\", 100, 105, 160)\n",
    "ratiobkg.SetLineColor(4)\n",
    "ratiobkg.SetLineStyle(ROOT.kDashed)\n",
    "ratiobkg.SetLineWidth(2)\n",
    "ratiobkg.SetMinimum(-125)\n",
    "ratiobkg.SetMaximum(250)\n",
    "ratiobkg.GetXaxis().SetLabelSize(0.08)\n",
    "ratiobkg.GetXaxis().SetTitleSize(0.12)\n",
    "ratiobkg.GetXaxis().SetTitleOffset(1.0)\n",
    "ratiobkg.GetYaxis().SetLabelSize(0.08)\n",
    "ratiobkg.GetYaxis().SetTitleSize(0.09)\n",
    "ratiobkg.GetYaxis().SetTitle(\"Data - Bkg.\")\n",
    "ratiobkg.GetYaxis().CenterTitle()\n",
    "ratiobkg.GetYaxis().SetTitleOffset(0.7)\n",
    "ratiobkg.GetYaxis().SetNdivisions(503, False)\n",
    "ratiobkg.GetYaxis().ChangeLabel(-1, -1, 0)\n",
    "ratiobkg.GetXaxis().SetTitle(\"m_{#gamma#gamma} [GeV]\")\n",
    "ratiobkg.Draw(\"AXIS\")\n",
    "\n",
    "ratiosig = ROOT.TH1F(\"ratiosig\", \"ratiosig\", 5500, 105, 160)\n",
    "ratiosig.Eval(fit)\n",
    "ratiosig.SetLineColor(2)\n",
    "ratiosig.SetLineStyle(ROOT.kSolid)\n",
    "ratiosig.SetLineWidth(2)\n",
    "ratiosig.Add(bkg, -1)\n",
    "ratiosig.Draw(\"SAME\")\n",
    "\n",
    "ratiodata = data.Clone()\n",
    "ratiodata.Add(bkg, -1)\n",
    "for i in range(1, data.GetNbinsX()):\n",
    "    ratiodata.SetBinError(i, data.GetBinError(i))\n",
    "ratiodata.Draw(\"E SAME\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "33d5883f",
   "metadata": {},
   "source": [
    "Add legend"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "809f1b11",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:51.293137Z",
     "iopub.status.busy": "2026-05-19T20:10:51.293013Z",
     "iopub.status.idle": "2026-05-19T20:10:51.421351Z",
     "shell.execute_reply": "2026-05-19T20:10:51.420878Z"
    }
   },
   "outputs": [],
   "source": [
    "upper_pad.cd()\n",
    "legend = ROOT.TLegend(0.55, 0.55, 0.89, 0.85)\n",
    "legend.SetTextFont(42)\n",
    "legend.SetFillStyle(0)\n",
    "legend.SetBorderSize(0)\n",
    "legend.SetTextSize(0.05)\n",
    "legend.SetTextAlign(32)\n",
    "legend.AddEntry(data, \"Data\" ,\"lep\")\n",
    "legend.AddEntry(bkg, \"Background\", \"l\")\n",
    "legend.AddEntry(fit, \"Signal + Bkg.\", \"l\")\n",
    "legend.AddEntry(higgs, \"Signal\", \"l\")\n",
    "legend.Draw()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1482f03e",
   "metadata": {},
   "source": [
    "Add ATLAS label"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "01c84bd8",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:51.423225Z",
     "iopub.status.busy": "2026-05-19T20:10:51.423107Z",
     "iopub.status.idle": "2026-05-19T20:10:51.538917Z",
     "shell.execute_reply": "2026-05-19T20:10:51.537353Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<cppyy.gbl.TLatex object at 0x5564b0632000>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "text = ROOT.TLatex()\n",
    "text.SetNDC()\n",
    "text.SetTextFont(72)\n",
    "text.SetTextSize(0.05)\n",
    "text.DrawLatex(0.18, 0.84, \"ATLAS\")\n",
    "text.SetTextFont(42)\n",
    "text.DrawLatex(0.18 + 0.13, 0.84, \"Open Data\")\n",
    "text.SetTextSize(0.04)\n",
    "text.DrawLatex(0.18, 0.78, \"#sqrt{s} = 13 TeV, 10 fb^{-1}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c90531a",
   "metadata": {},
   "source": [
    "Save the plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "a5ad27d6",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:51.540534Z",
     "iopub.status.busy": "2026-05-19T20:10:51.540399Z",
     "iopub.status.idle": "2026-05-19T20:10:51.713129Z",
     "shell.execute_reply": "2026-05-19T20:10:51.712324Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file df104_HiggsToTwoPhotons.png has been created\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Saved figure to df104_HiggsToTwoPhotons.png\n"
     ]
    }
   ],
   "source": [
    "c.SaveAs(\"df104_HiggsToTwoPhotons.png\")\n",
    "print(\"Saved figure to df104_HiggsToTwoPhotons.png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0b45b013",
   "metadata": {},
   "source": [
    "Draw all canvases "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "c5daf506",
   "metadata": {
    "collapsed": false,
    "execution": {
     "iopub.execute_input": "2026-05-19T20:10:51.714571Z",
     "iopub.status.busy": "2026-05-19T20:10:51.714438Z",
     "iopub.status.idle": "2026-05-19T20:10:51.889270Z",
     "shell.execute_reply": "2026-05-19T20:10:51.888421Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "\n",
       "<div id=\"root_plot_1779221451878\" style=\"width: 700px; height: 750px; position: relative\">\n",
       "</div>\n",
       "\n",
       "</div>\n",
       "<script>\n",
       "   function process_root_plot_1779221451878() {\n",
       "      function execCode(Core) {\n",
       "         Core.settings.HandleKeys = false;\n",
       "         \n",
       "Core.unzipJSON(73911,'WkwIyFsAtyABeAHsvXuTHLe14PlVGDX3D8sGc3DwRtXuRlCvkWf1YFiyLQ2H90ZJLFI9Irvb3U2ZskPffeN3gMxENdl06+nru7TVLACJTLwPzvv8ffMfV9+dH073zw6b7eazd/an3+4v/3z48tPT/fnl12dXG7N5/MfTk788P/z+3c3Wms3jt0+uLlvqky//z+GrK8o3VPvk/Ork7LRn/t+T00ebrTebx8uXtn9/VVs3NeB9cD4Es3n84cnp4Z2zp2cXm6307KdX3z09rNk/nzy6+rpl3z95+rRXprNk58rW6uuHx1cf7S+enJxutnaSZDaP/3Dy5OuhzEZGeXZ1dfbsWsXPzs7XEq32+eMTOuLM5vEXa/JeS9Lc55dX+ytaqpU6R7l7LcfL71/snx2u952ya4Nf6h0Paqk6F8+f5PVxdvjk22cXjw4Xn578rc/gUPjR2aNDW9vPpf1+0X8/d21+v+i/n1+d3fvy8v7Ji8PTb/obV2evym62qTLJXxy9kJ1O2frKccFme7cV/PHoLVb0j+s7Y3Zu5k9HL7SP/Gl95bhgaWYeytXZ5300OpSXs6yktRJ8ykWSrVIjW3R9/4tv2jz198nSzTm72d7VD/gSbZAYxYWkX/j86dlfP373nT7pY+bzP57rA9b0iyH956X0gyV178vLow/d+/Ly6Fv3vrxcX7v35eX65h9fPNMTwQx/tyZfPNu/aAP643dL8rOvD1f7zdbrwL4+6al7l+eHr67+sL86OWuj+Pj5sy8PFy392clX3/QPkfyuffPDsycv2vMPz558t6T+1lL394/u709OOT1m8/idi7PLy6/3J/2DS/b+WYdG47Zmc7X8uqM/Ont08vjk8Gizfbx/enkwm8f/4+Lk0Yvj7Hdr9t6Xl++cnV0M9d97dHK1/5LzfnXxnA+8f/Li8Oho3POn71+cPDu5Ovn2cPkS2Pvw5BKoOkPcnt1fXGy2Dx6azdn5FYnvzebxey8OX11utqfPnz41m8cfNxD9FQDzs5MruqFg9+Pnz+7vnx6urmYgyWx9fHhx9XLpu7//9P6H977YbDf/NifN5vG7Z8+/fHp4+/njx/Nq/eFwtT85Zar6QD+/PPnb4Y+X8/MvjrP69A+H/dPNVsqUUkpZYd3J33qpo0t/Pjl9dPbXz87OOURj/osx34HVWuGDA9C574G/zsf8na83/Sy/s7+6emmO711dtZuM8X3+9uHqr4fDaQfTRzmdyPcvzp59dna+2crExvn80f4KKKiZL+YMV9m9lpHvzeabj86+PXxyvv/L82U7fPOHA0M+Lnz8wcmTrz9kCP1W0q25v/rq63lyv/n067O/vvft4fTq06v91fPLZQd+c+/51Rl7YKn50eH0+dv7i5Znh9z7ij22vPH4D4f9o09On343v/H4zydXX589vxq347xFP9hf9g02l4y1Hly7rH82lABweSNK8OfDl3rkT06f3IQXsDPeebq/vOzngXoNERkLztmqG7t1MZr+t5OtNdbYndNSUn5r9andBX3mYtzF5R27S70m5bnXJV22pRonYkrYVdIlGCdpJzQYTP/biWylOiPFGcluJ26bk9H/duK3Yp3pfzsJW3HZ9L+dxK1Eb/rfTtJWcjH9byd562ww/W8nZetcNf1vJ3XrgjP9b+esVpbK+3bn5DjrtlK9kSpGUtk5v5WSjBSyfudCazdFI9HuXNxK8TwxEuPO0aveyWB3Lm8lRiPBGfF258pWfNQBlrRzdSs5GmcTPdl5uxXvjDad/M7LVnymB0aS3Xm3FW+NxGL4sPetnVLau0yVa+1Ev/NxK5KNvh/izqet2EAPtMs+b8VWoy94t/NlKy7ovDpbd76usx7jLtitZGt0uFV2QbYssfbZx11wW+fE6PyQZapK66S3uxC2OsWsmgu7ELc6p3RD6i6kLZ+pwdS0C7k1w0SFvAtlqy2KdnMX6lZCNjoIcbtotRO6xXZRthKS0T27i07noWf8tpY5HbY5zOm4jTKn09blOZ23fk4WNnb/TN1aE9jAcZfs1ppUWpqjU6Wl3ZZB6bHaJU6Q+NqecIjYHfp6JFNyyyTakNZI0pPkeyvauh5S2SWa14zLu0z7mrFhl+kAGSl2l+mBZmLa5fkMs76ZHugTW3aZHpApcZe1AzGaJLs8n2Sfd3luXvwu123robG7Yretg6RlG3UYpN226DSQ5gg36GJ3hV2p00CmHV/as7vCMUlzhqPbemV3hYOr00CGY6uzbeyudsDlCxlpEEkCmQa8pAK+qtcnkhKZ0DLBkWlQjE1ld7WBsapfzlorC+VF0yGTrpp23tidWLsVEROiAgWxojs0OBPcTqxrQKGyvLIT67clG9Fdw+Og45CYjCcbt84V47OJYSc2bSVakzmgdScWkJG0riQ+xfF0xnvjHFmmxCRmVXaisDVyhtxOBDAmxotx5JzOKsAo1Z0AW0nz5+JOAK6sJhPPdyW2Hrf9LEBXK6ZWI9buRPIWAMd4pOwE6Mq+S9UIEL1upVQalp04u5UAzApGEvAdONbmrZB1W33NZlPqTpzfOq/t2J0AWkvkWjAccwG2ijcKyhyVkx51J6X1Geg6XAHOyrYKR7LKzrEkgMT2t3MW2ORN/9s51kQYj/7tnI0N3gLWd45FKV5BgN05m7eczSwm8+Wydb62bkneOVu3vAWwjjvHkrALeHnnRLZZTHt558QpXO5v7xyLQnu8vHMSeLO9vHMSebO9vHOsBxBdMm/vnOTW27nRsgxOW63a37lZrjvAUev0zumS0NvWLkvSxirF75wuiQ7W+bBzLiyj1ZZZlD5aOu2A4uusujyPt71clgHzsp/PMLDNhp2fjzGbiXuLk9xv8ZZvh1mfW7fznOf5ueYHnEDz7VRTP8rOc7B7dc2WBYPQLLCtoRtRdsGuiJHdBdtwAl4n55Y3ya14CrnAdcKNsguW2WjQi/EFPcwK+3U8gcNMnts4+l1g61i3jC+IbXmeW7cLIEt0kRtU8462+vCCeM3ZoP0X7QdVGY3oJcfLmptBvd0FUVDfxiWAel6hnIuGF+wuOO6ZNiTXl6ThjLsAOtSHxJUWOLR9SDpk9kcfUnsOJHE65PYc7KMPiSlwZR5Sq17nMWltb8nSgFb2otm2arvg50tPX/V667X13gWvF+/yot68y3tMhy7aLrIl+x7hK7HKstQtvy49PYp1Xfz2PLTxtDHuIluytaNdjjXR5bZqu1jznNNv6djbAHex6tDJWbdLliXQW3GXbMPSeSdZxqyLtkuWEeui7ZJt2DoYA7uwnTDtYgJ8tRPW80WH3Jdtl7hTeN6WbZcUhOkObE2yC9sO1CElrhY9Ye1z4vVzDJEOcrPoJuyV23zMw5KRfEgdawfq8WY7HX1w3CrsxjYdTjcCq7ZLbMl5fLwHxOrjY5ckQFY7YT3fL/6+/ElhVhufTrVi6XrCer5db22DJFfn4dE2WHobneba+ew1/UpEJd8mBYwkgZ1zxjTdpqOldWuQzDMgbN3e5Q4I+7TvMruudzHKLrPpepaX4zJefbguP8PLHb3RLuxyR3BoSSu30bWHZYaBCh0KMLDPa5RdsQ0S6KLsSgeCCiSKXdacJ707DfaUvht7f3dl3o0KaQp7se89bWOZbHDMBpFbE30TUpmKXKUdrFBTIWH/5AoIedLBURtu4R6dAbJ1u8I92ra9Vi5LTttYYBFIrm5AmifTb4c2O7vCFTrDdT4LkGS3NWBcgJE9y8sdE+5d6ruxz8Ku9N04t9TOKjm6BM3YGwLZ7rCrdYrN2MfGsSiQjH2KOBbFt0lq+2ZXVgBJj73eF/Sp1e1Qu1/Uxc/0Qau8wmxtSEmU+dUw0iglKMTmO/rdoBC7Hb5dCUop9X0SZoDNOMPaHXoX1u60zyzdYWBh6U6rPHeHD8W5N5rpM9Su6F2BXOzrogOJyzbS70Z61KFzaVRjn9Bdia1HfcJ3JdIjfQorJGqP1iw96vCxKA3Z74OS2oXWlm1XGh3ZFm1XlJDsgLM0QrL1dldSu97bVtiVNFyytK/k5PrVNj/zqJWiXF5tJOX8LHdWyNyh3G+/uUu5bex5UnLb2A1il9w29jzQ3JChZRYylM1MYZfcNvY8ERCa7SE7qrSN3a6/UtqyzVVLW7V2aZSia9ZuiVIUB+rdgeLUvdUWF5qzw0ydo9KO/txGO/pLG3r0exO13z36WofU8wx0QN3b78iBXlOlA+meaUBxXs0Oo5f17DB6/sxyAemG7jBaP1Q7iO6ZFUIzb3UA0ezgCpDud7F1u9rBdGuldjjdR1JnON2XunY43SazdpxhzjWMoa907YC6P+zcNd1TdApyg4lvf7sKsGQft79ddbqC7Q6tLmyjh6wWz6u6zQOkcd5Vl7bRKgUeIPrzNlQIdQl8tGxDMRCeIe2qq9uQTLRGQtlVb7cB8Gkkyq562YZgtDNhV73bBm8UPOVd9X4bnIH8TnASwjZYA8ctuV31ceuriRVe2K76tPXFpMYirD5vfTaQyJkWytZHAwGcaaHCNknJSE67GuzWe5OykVx3NcjWO6h5aOMa3NYryScl7mrwW1cNtHDJuxrC1hVlElTYGhG2VY5wC3cVTloyORuptJDhEORqpNICxKop7GTZ1VC3zpnS8LqqfDSjpF/eVYBjNSXCk9pVQGMxJRknflcj3ANTinESdzWGrSRT4V6VXY1KEMIIc7KrwMVoqus5jryyJxzfBL+D5dJySjbX1HIJvMrU3HNsGVOh82HVwDOAF+Gk7mry5ATGhWa1M7prNRvbU9+fppZtjIuacssyyLKrqTAssYySLAxZA2tCsxnWp4FV6iTvahbNwg7TLCwfA0uoZXWGRFzPwiU2cC9ggVcgYzGivAaycBcav1CzuWVpl6cQ80bgcjDbQEaydCPuarGahVOgWWlZusFTeNfK9GhZ37J0I+wqDDmeMnyyMGSNOIZPNrUs3SCbW5Zu+F0FQFKZbpCtmvV0w+9qbb2C9yRuVyGiCtzHnm1zBb9JZFeBknCn6QbZtoKebpBtK+jpBjw43VDi6QZZ3VHKXdWsbikJ7O+6q5D0wQhcLAsvDaYdeTZ8IQ8zCoAx5+GqGwmwI/S53woc6GSczdSHgwlHes7Hli9zHraMkVDhiFM/az7SHc1zqwB2jLOR58pOlUh/yIvd1n4Va1Y0q2IMnrotzPTGKxErXrP0RZ8GzdIVz9O4haVDTzSbyCY6AscRBh78xDlbNEs39Gndws1jVmAxOqvvKneQrNAQ4Kwyp87RScQflSlTVFfZrVKZMnBd1/OMCWwXHiTQiW4rE6/n6aiiu03MUekLLCN6DvzSzrTtrzCS5iHAKmIbI4X2PRhLz9M+RBiiJZVKkIfD2HiKhf54eLE9T388xGHji8KftF4lAsonbZzcDOAV2LOF2fAFsIzgBgguFvBO++TpDwCe9hu8FwuIp33ytA+Qp31uANoPnsui8TGZjxC4Slqe/oTIRSPA80z7IXENSWHOaT9kLqmW1/YLV1iTVjEfoXLBtTzzES3Xn4rfEv2JwuWobNFE+9FxkbY87XecWGVcrEcM28j7XOi0D16cjRQIB9qPacsdqXltL28T/eG5tle27B3NM/5Yt4i44PsqqzxZOIgtT3tJtvCseR5oL6nAsOVpL/ktc4OELtBeCltEi5qnvRS3tKV52ktJGd2aZ3wpb5VrDx6h7ZUtwj/yyo5PdVuR/pGnvWyRKrU888u1YOkw8I4GuRggcLSAFhVnZooUIlIQtlwjrYA2uR44NKWqmFAsFwRcVgqctgpgoVkKtFnEAjRLgTaL3I4zalXYJpaLwiFapICV5apg62oBzXJZqMSRGuwtrgs2I5JP+PyWC4PdSA2hWa4Mtt9aAE8eqSU1mH+uDTaYFtAxLo7Q+4G0wHJ1sMW0Bh3j8gA70wImaObBaQE9VSZcnw+hp8oP6fOhkF4vEeYUiQVd12ukz6nCci4S3XnUoKdcJYAlLaCnXCYNziBqpACEoK8tBSoEUkhQjABmVQzE0dYCpC6zJEgLVEKD3GH+BqIU7hQVRBQDPEUapLuZPLIVFT0wtGIAqCoQYiJKRlAqKhBSAXo2AFgkQlLZb9kAYJEJSWVY2QBgEQsJchzyKiWCjUOH5nzjVGmeDqsYgiGB4pGHcqE/SSU2iIaQWlIfAIxsCJGw5ukfwiHdzAUhoyAdQsjCeAC4CIic7u1iCv2DlwLER5xL+3D3+lYH4KqICAhfrQHgIiVSmTN52oe7B8SvYjLzBT9FD4IYAC5SIgeEJ6/tI2FkTzsDwBX4e0B48swPTBU9Jt4AcIULBtFW9aYJ0tzW6akJBoAresGw4YMB4AoXjB4ikH7yces8+x+SgHzaIsPUPO1zweiRSgaAK1wwgfYSglkR33QXWFMArgS7dXrAsgHgCiJ5IDh52kcmr+etGACwQEfocSsGACxQEkGpAxNpH1pCT181AGCBmlAERUlMEegJ5TNYE5kPKArNOwNAFmgKzXsTaB+qQvPRAJAlQqOCWiQDQJYIl1QRMANAlgj3FkykmkB/IjJbEEdR0apE+kfeGwC0ILLXfDAAaIn0j+fJAKAl0j/yGWG2SKR/5Kvx9CfRP6gVa4DXkugfeWcA15LoH/lgPP1J9I98NABrSfSPfDae+Uht/thDCCkl0T9rHKgx/Un0j7wYILck+kfeG6f9oX/kEZsjaqV/5JMBbAssGc1DDZGnfzwHEyfP+lbDngGIS27ri66L0J8c+3POIPnU82Ak5KEXeR+sizz7jzzCXRSJ3n7+GF2yf1NtvM37T8/2V95tzOapao/FaDbfbrYPqgumOki/ZKoDthRTXTXVW1MR8Xpnqvem+mCqj6b6ZKoHRhVTfTU1WFODmBqcqcGbGoKpAeIRHRFgWzE1VFOjNTWKqdGZGr1B17QiA4/J1AhMLKbGamqypiYxNTlTkzc1oWgSTU3J1AQsLaamamq2pmYxNTtTszc1B1MzpGgyNWdTM3C4mlqsqUVMLc7U4k0twdQCkZpMBZYq/K6mVmsqMujqTK3eVBRcajS1QsJmUytwX+Xq8GGAxFY5IRYJhHJGLbIYQLAFbbAAX4vCjgXsWmCtBcBa8GoLaLXAU6tUMZDUAj4tWLIFcFrwBAuItGDCFuBogYgW9MACCy0A0AL1LMisBd5ZgJwFd7WAN6uKSwAyC0ZgAWEWuGVBPy0Qy7LlLNimZbNZoJJlm9mmvMMb7CoLY8YCc2zkDeXKgVRa7ngL+8WCU1pueAs2abnaLdwWy51u4bRYKBLLjW5hrliucgsqabnDLQSC5fK2sFAsOKTl6rZg6pY724KDWxQILBikBaW2XNUWBNJyR1tQRwvuaUEYLVig1atdkQauZQuyZLmQLXij5Sq2lTdAGy03sQVhtKw5XAuQBv5ByYE1b1pNrDmMCO57/oHiZ80VUYTLIIohNn0I1lzxQ2HNm2IAa676UrAGuFX5RxUpaIM1h+7neuQf2mDNFSFsGlmsuaKDwppDrHN38Q9vsOYQ4dxHMDJogzVX5A+CmouFf3iDNVeFLcXzII+5FviHN1hzZQYKa676carxIqy5KrhBxgKo+Yc3WHPVsFMqtGmGseaKvqmGnOJtShQqwtZUSFhzRdeUOlM8TekuyeXh999/b34pbU5MJW7U5mzGCK8xIOkam/c++/Dep9Sb1Zq14M7y+jUbCLKDncN1E4fVAEQwuXjZAuSj/cU3h4vBoqQV9E+qrnIvaSYSqgf82eHF1b3TJ2hdo4FKtj22zRBEnz89eXK6gS3SKgxN8Pj9M3TZg6oY71+cvKwZfu/q6h7l6G8/Ovn25PLk7PRys41CizwZPvjh/svDbOhCe5ofWtD8J48fXx7U+gRA2ysdd/vkq28+PJw+wXzGThYlZ12E+U2Z0NfVomvvXT2d9c11uOTn5tFM/uK//Aj/17/MCH/cCr69vxgsh97eX8x7QjWwscfCtuLp/U/biXj3Yv/XZnLR8p+cX63mHS3TLTxapht5fHJ+9W5Tt28WZOjALx/o8KA/el9NrZZMr/aS9jwVtPDRyRU2Y3P+s7Ozp6o8T0GzTXnn7PTq7PnF5aad+ntXvS/XoOW9qyvOr0Ko14ABxxG7NRxIekwYezMqYpXIYQ0xae6900fvXVycdUsteq1ZrU5T7z8//WqACGQHsEi2rx+nGluXXpkvke2V56fDYn94eHI4fTQa1KgdnpYO8JUPrYVz27OZHZ9oECEB8lpFpqcDErN5/AF2EIfLaxC8l356vv8KOwCFL4vd2zCGxehtnITZju64N0vVufgn28d9cHLJhhz7Q9G1q2quNTfb5mauOJcyE/3dvghMGEUfnZyePHv+7H8dLs5WIw8evGyJ2Oxd7l8cHh8u/seHa+1WPkxbKxgXl9bG0tVmqpW+e3j8wWYbLUu2lPx5s83HJZ9vkAAMVb7oBff349a7vz/aWTS+FK0ta9ErTDDv7x+9PPb7+0evMOu8v3/0KsvO+/tH7P7P1ynqJV8clXAzMqLWPXIMR3PnzXry8w427s8FX2wQTm4ef/rVxeFw+v7+K4U/vANUG5aALCdh2LgUXV8TyoYTNGfX80PJvIna1mp1Lp4BRDdpCmpGQ6FCjVp6XxTIQKq1rv15tnai5gdk1C7206uLk/N3D1+dPNs/vVxMihQoz6jOgjDMQwHl0hrXxqdl1weohcMIl/w6RC06HqMWrUBTs8OANM+I5p4woAaAzn5/enq4+AMDVHRoRiQvN9sHGBLduaP/iOM/kuFO4YeCdCfckUTOWQrCHYEzMOS0aCx39hUVKCq2Pyv2TrB3wp2gNVvrVPj5/3kIRDnsHx0uuLDVakonasm9f3L1/rxxYt84aovFQi5PdD/pvH21f6ovswP+59nJKYUzLvDO/nzMfnbybEEocylVSlDI8ftn+ycHPrSA+Xf2p4+eHv789cnlN4eLP+xPn3QT5Vb+9tmLXtbWr5VqTwYrzT+dnD09OZ1LuwFjq/rOycVXT6/bRPdHmJfS6eEa/Byc+70X55+PmM9c+MVY+MWras6FRzWp+NH+xbsnT9SMnm34ycXV12fv7J8dLvYdAr1MsN3fP/rZLPBYsRtptvv7RzdRbBxtZmg+6j07H86eHebv2v34MiV2DWiDZLwKYl+zxedc///LEL/GKVRfXJZQLEytTbPLxxQdg1Vs8pOffI02SU4FhtxGTfVLcDLZ6GK1kf9T+8hC/m7NZSqu2mqtxBzZHIMl/+sec7tPLlopUmDJZbX6Pfp6SHXy1toYs2MA9Go1xn/N02Yn31+TUmoqcF9/qj8AbdFsRo8A14s227tz0WrVj1eAVyzC4CTgNU8Voa85xupClBCV9XndZ8ArFmr2GYALgb7Sc9Fme1fyVIMvMVYfQ0i6sqPbgNXuf/L6cHAd8Do3AnbSb73Wk0D74qudCfT3uz8BmH2rR4E+itmrQOIod78CxeqV8ONdC7TrXl0L6G3/xrXALVwLXCN0b++aoHPunp+fHy7+47zdGYpP/ECnBD+fEXu/6H8VG3ZF6WcHOdcnUWnNm+5QEPCf9w7FUUs7Zeqq5ehK4CDgsKUfrhH9Zggt32mul3GOn4hwvMccNDdAr52vD+Tdm2YrKjO277VH+6s9Fcdtduup5Ob5h+jIj2IMf/zV4SnUkaeJG/i6nal7kx+lPsAXylAcRnjn2X/8/b892T97tm//fn/nwf84/OkhVX4N/rBeG/+QQTz3d2EP99dewx+m/1+ewNpWFzKfq7uZtoc/b05mdOd+3urgFOXx+ycXl53h9+F+TuGJys2E3rPDuyeX50/3gwcX6I6FcmHWlCm+elT56OzRh/sve/41TOvbLd9315dPfXr8atz8wEX7j1arXZJXT9fVam/darFY6rZWfKavFKlffZ1uYL3fbp3+dn2d/gudp/88S3RNXjAIExqi997p1cUJjJecfFUZ0+XzZ3+9lnXX8i/AwqPkIFOtNcasu5cXX0ACOWS8kpObYk0FauSj/Qv4p8sN2PmpK0L98dnFs5lRx67uEoHmiunxp8+f0QeFQEfe7WCu4+zuZRne61Coa9jW+8LOe9WlINHDR+z3wuOTJnrglG62m988sA9/90Ae/vbF7x64h7998e/udw88v/6t3z0ID397eHH+m7t2ir/9zW9e3H0QH771v//7g/TwrX93bykIXAl44DDZVxPw127M67JUJutH3Zc3QfuPz3FsxKJ9/OjkWWO1fnyOz7BGGHx33kWvH58/VtaJ3h0fv/t+Z7q+8/XJ5V+e7y/gGNWpplqyBHEOLYXNuvB3VcNg2Blzwf39hYpXdElHBaCEZpAqAOWm/pMdFHXyPgfJxSdjp5yLxBRLDLX66AxoWPTJxlxqKd5mHKlMKRXnc3FOpIo/3I3q5ev+/uIjXK8dqx0dtXq+QQCD6pFInaSZNUy+Lu/jpO3Hvv/pHg9SN72N3oa2HKMNE/oi3opHqzH6mqbsgq/JZZ+cia7KFFy0IQYbcjJRSpmkhFiDK7lkE21JU805cD6jBBNqyUyejcXHKsmEUu0Uky022+xTMCHXCAtBnIs1V2dCtjLVbIOT7GJF807sFEMSJz6ivhkiPA8ryZesFYL3k0ebw6WIXkvwIU8h1GJTEp+TCU65J+KdtRllmyBFJtRVnM0JjY0g1k6xSPbRJnRtgnUy5WJdcjWkZHwNYUpFaiipuhKML6lOLoWa0ZaJaOjVOOHwLxQpGRW+7PwUJKeUoo1JjE/RTbV6712yGSW8WMLkSvKpROeqNz5KnnwtZGMp3vgQ3RTFS8hwGKLxvtSpSsjWcwKM957NaG2uPlicsmQ/uehKyjkkR4Gzk1C9Bp/QDJQsk3ifk0hVVUJxacIXp0/FowTjbYmTrSX47ND58jakKWcp1jtVZ/LW1in7UGyShDakqxyaEr0EVqUaV71MJaeaoosOO7xS8hTpfYo5YA5VfJmcs9ZJiQWjn1zCZIOUJMln7FSyS5OVaLPytrJxKcWp2OS99UVdFCUbJ+9LyjG5oJ57QprISwUORFT/7OQyl0r1qJ+54POUQ0g21BxQ3vSlTJW5QsELsxIf/JRdjLlmj6WWtzJFySnanFDXdS7JlGuMIZXA5eaco+fe22BjUiOrkqcSXGEXNzuj6KZSsgvVFZTTnEiZasrRV/EZ/0HMeQ25JHR/ilrV+CkHccVzXtS+xU2u2hCsc+h21SpTjMEDCtHNqkmmaG10rqAdVr2bxFofY7JNa9f6ybuMclz0Te06ThJD9rmkoJYVIU/cus6qTlJxdoq2Fgneq6ZSrn6qJVlnU/XoouWUp5psjOhNoZiVg5uCzcWhiIYOWGZDW+dCqtahRpaqm5IUF6IrGVWxlMrkvXWeTYS6VQphCuKSS1lwnskWm4pnP7qiynqx5imWxOZ0KFTFHKbkUgDuJvTNYpSp+hBrjqIa2NHVKZTosqs5oVQXbdZjEqPNqBeGEieHYpWNEcEZG3nyPgWbnPUosoUgU7LibMnJU4FrCMgWrfXoG/paphBKjZJKQLPO50QN3UuqguhjnGK03vmgvfI+TADb5ABnmIaKm7K1KRZ6FbFRk8kxf8XXwJq5VKdUoys1e4cinwtlstkCLkJGbc+5NHlXY3SRecceIkwh5FhyiEXV7wqXS8nFhpxVsyzZCeZwKN551cPzeUqViy1HMC1U6ZEnhhxjUP3N6ibrvJ6ArKqTqUw5JsZpi2oyhjhJqrGmJEkVCx1DSTHF4AIatDVN1Tvrkq8eJUlgiPOW+6ZEYFkNeUJfskpJnldUtTBZvRbB7PZPT/42ioYo+/3p1eHJBb4qQWEgD58/3b/MxOnlNyFojQn/SyBo7yANu79Hj+HqcAtspAbv4mRtZhumWsxdyTlPwBwY59waIlypMVnhkgzB4IPWxeSKLamGHJsm4kuYxb2nT9eOfIqH0dUjqD5Abvr3zb+dq1fWDT//12cIjU+fmJPTK/NZn0Wt+wlOh/8fJrSR8ptzZJ6Xh6/OcBRtUVH8sR8Cj54/JD/lQwhn5w+5n/IhP3zI/5QPgXvOPQo/5UNx+FD8KR9CBD73KH2vvJn5DG1+8+Ac2uS8ESfn7uFvP/tof/X1dvvpX37z4q3fPTj3D397fvbX37ww/i2yI7Wy1lS65bwRLufp4VtvKeEyUAWDX2Hol/3F/f0FPsmVTvvTAbWL4dzD2Jl3a3el+87Zs/Ozy5PmrpwyxtEc75JT+NH84a6KCLN/XG3k7ZNTVCXeu7j4BHe9QBLyn3x7uHj89OyvKJQh4bm4gB11E2oNxxCnsx2/9lVMwBQgZIkGjNkEcR5EMoAJJnA70LMiBooCpMoalzG+iNiWezyoYavh1MBYtb8LV0DmQstcF3qVJdSkkQVK4N5yehWhTI9dSdNjbgrn1v8CusOf3vvovd25qin9x+X+24OKm2/BJr6ZWoan/ktA4jekMsr4b0hls5Dab0jlN6TyG1L5Daksb0jlN6TyG1K5hQya0f83pPJMs78hlX8AF+ANqfzid/9SpDLMg6OYYT9RT+gXJQi//EY11m8hO70mHEVoOAhHZ1npYJ8zqBP9SsJR+jSwQV4nHIUjobJREoNolOwiBJ9Fn4uUfC4YZKEPNGaKsTB6Zknl9TKkj2PZL0Um4djwSKIYr0kU8csxShRLxE/iIFHEcn6UKBYf4ihRTLFKHSWKNuGEYJUoKg88mUWkWHIpCXZNFynmUoqPg0jRxSwRBxGzSNHZmqTUuogUbXI1Z1pRkWIuzvta+UaTKeYkviJqc12mGKrUZBPsHpUphowkJvnsfZMpeinB1epwYaAyxZx8CDZlmD4+Fj/lgqAgJR8qMsU0+ZA80gy4RD5EO6m0qvqK1wEkaAjaivcp+lyN985PNlW8JyAM8Ui8Ysouh1yZMO+sTKGEnKKUnOFaBYdMrFZM4SMywxKnkKNPVSU9xltEMDV7pCMBT401h6n45KK1FTGOq06mkBCz11jw24GI1SeHFMwGfLoVKVOoCHJjRd7kcipTiRal+IKsyGVnJ4+oyZZAnCqXsp9s9aWK8xFhXXLIHaKvoUjCXV4sMpVSnRQf8fbgos9TZdFKQjyJEDFO1kuozruMmDFEP0kqOScELR6PJW7ykmzILkX82vnsELpkJJy4s3De+yn6glQ74TzRIS1OUkMt6hAAgS9WAKySOOSOzrnJ25TVMgDBpKh0PWVfimQNSRTdFHKOpUiIGbeCUqZcfI4lZo9OAqtQ2Amp6IZyFqmircjBcsFNo0WqiJpDsLnii0PFiqHWkqIr+ItArhhqteLR06eGd5OlgqtIsmZyKdVofXXqQ6pJFhOCHTyXyCxZRFMiIZRaZIu45IBDOcsWpSK6VjdzXbZoa0jINwfZog8ac2mRLZbq8ccwyhZRPziWLSK3OpItqpuYI9miuh18I1t8I1v81WSLPwyBeUl4+HOJB7n+34gDf6yA8pjG+QkyM5DUAQP8BQRhPzd2/8HvP/3sDig+Xb+FmOc11gCDnOfJkw/43ozWv3ty/vXZ1dnpnZPTb/cXJ/vTqzvP9peqej0IcZKdsfhXazwi6xiQelDll00WX6/z2MQluJdpgkPybwwE2jq9MRC4rnj+6xoI3Mqcg02vp+pHrRbb/Y2FwBynHA9Ns0uk2fhkMRp/pUemWy0Rk/xfYYlubSGAvWy3Dm02Ar7kyQaX0QLLONjcPNYn2AJNTlK2qkJqHc79+rMXm20oEOloyYYoZTU6wHYgoeprw5R9SWg3vkrs/RLf5FbGA7fT/pDJorFaq0c9OMvhLmJnF5yElNCjTU7L/OSKRSu0wkeIKJGbBBVYUEu1Sl1TJpP3LpRS0DINIR7uBuMnb2OpLqDyVUQoU3tBm33M3gX+56qW+WJDzLn6nELyWk+cUutZ1fWKp5qHGrdQeWjYalHy1SmRyiQTmHKysVpX0Ai3MUSnRa647GPIAbPwqi/aCgkcfA42VNzv0TUXXVaaOaVcUDFXjUovwcNncYEhyIResnUWFkiyOqwAwa6agjWhftmmJOQsyTvrcyhSW1kK0Xs6n1FPpsxPaIn7EAXV5Hq4m4xDmzQmlGvhrSTK0K21tZQsvrAYlMXJOx9DybgT1BnOxk0lo0WZQglV+F4xMgWPN0C0KVEIpaxONdhSxKEw7lixYtDLDJJqDk5X+HA3G5kkinMlwMCRrGUORXNU1WEmFHZPVh2mn8FmZVGM+qcoReEJq0oM3geXo/q3ZNvl5JNXdfOMRjZ7LPlcxAdJ0Udh8wjuEVOVKCHEikmG88W6XEosqBIb9PmxyIjeR8GoyMhUJXsfBJMCX4tJU3G6QUItvmgszimj4M8Oj3BUkptQ/C24OChFfVz6qVoPAypXdU6NCrqzoWIFgc/IqNYKFUhjAwqbkqaYJUuoUUqEVThFjwI6DBaYTEbQ4LWe4+CqJwTw5AIsjRpczrTGHJScghPnnVSr0WInizY3PA1HD4l0ht6nyzm4UBwKo+3o5WIZIJpmIbt29FQvH8uC4mppZ9u6It6lGiRXB6dHoYcPmUlzNeUC0wkN/FR9RGM7utyOqNTkbbWsTwp6bP2E0rUkL8k5FxP1KiwxmJch+pStvhvQKI+YbyjjSt9Nk5Oc0VILtsYEdPsFNNVuSbE0x2ZUfpUp189o8A1uptbe3qHujq+Pgsp19wYSkpcpArLR2s40rI5BYpyqi/yXbCpclpiEZ4l5KjWUKABs/fTH777DDaqONr44yn3u+rOCE5gvlhz2d6NVOV38/Smu+bh4P/16/+jsr4O7mHfOLvC/tH90Qph7Ki8Q6ssLWoBc6mLl+5glcXu/xuOemv/xvDneVBPs636njh1vju4BromxfoDlXltyTBi/u2ndFZG+sevzuFq/52yfqTk74o6r1R4TOwim5uyrBVN8aqBhyb6egL3B++knX/6fw1dXqsp6cXiMhR1atYrhbrabd7tl/rKYTw/nm5cYCf/MSftB0rwfPmksQpvXley/PmVefV3MU/b2/qtvnlycPT9VwLFO3H+uaXs9v4SnN/NLbrHXbjNtLo877dOTJ6f7p3d+d+ftb55MHL3/rDP303hNP9Pcef/y3F2ftAW/U9X4438U5VMo9+nhfH+xb/rydurr9gSrVTxAAbXfOXv6/Bn20DxsmaOXfnbGIgO5BT/xw/3V4cVNQFowLVrvHKrNLMXVJ/WNQPwHOn7OnBd1BiXczeqxsMyAabgkXw/cPzx5dnLVLNYbrKH/n1ycPDk5XVg2EXtPTD5DtLHIv9rUf3J+OL0z3yk/2/Q3v9tsWGUKvJn+m3b+f7v8y8XV3y+/v/N/3xF/57PDn8wdsXcef/nvf78r3+sRec2JYD/r89u7QmdJ1hOheOsPQnducyJ8bf8jxIK3/vvvH76EnrzxsPjso/0M0wEqo49FINdjdb8NEPtiTd5rSYD+55fNdfdUlUY4yt1rOcUYZg/NA7a7uGcewOBSbxbWNF+HS9W5GCSET7JnRt+TlI2kyVLxyOuV0j2v9O2HN627rsQpS3A1FI9RONVvdrZIrJp/5F9RL4HVBeJ1/4rj41v4V4wFS8+UUkix+KL3yvpxwUfg+rS7Z1QfgupfMWTJFn2aBDcDAvbIeSOD+eP6tTG72abKZPzp6AVUlI58K7roJlr9U/8KM9qLfjnfigxk/TreE/G5MJWa4V06rLs3iyPFb169yIOfRTeJ1Jok1iDeC1+/wc8iT17hZLF7S/yzUtDswg86La3c/tf6WWwjWd+89+Xl0cuvcrJ4l9gnm1f6WGz7842Lxe4q8P2TF4dH9y7PD19d/QHcdnD/fQGWdfIt/ol+AJdgwaUb4vzi8NVPeb8zQp6e/fWNi8Vri6DA/Sak/mfkuPU759jFYj9g7R7giMJOa0fr+n1zdNX87Gj4vc9/r1FzbkEFfSC/v2m6jnws/o1AAwMVRHpAxmZqpd+9QPvh2iV7Mz+gz+WPYj/NKhREkvq5nSz+c30ssln/kdM+7owjga8wC1oy4EuaH1iG3cNi8551k9OtX9913w1xgV7jui8QXWzhC7+kRAGheOfuwhAaHWPeRSH8l4icdJtVs4r5HK2bBfO47bqx0dUxpuIVzS0mZb+as8VrMPcHsMdZzeZb8yagw5gGMvKubvGj0DADEbk87bt9ybfdrtkP919+jJ+/OfOn/dPn3V/ch/sv+RiOiY+5XQ+/hz32r+lLEhTvB8ONf01nn7fWE9FN1TREluTsHlZ9RHZXsbO/SCotpi/tAh9UPBpb89d0D3mMPv4C/kd+fzpEMxXrfv7AgD9I0fP9m6DDEUqibOfLkyMzrrFsQE9mfGRgDfyq6EmMvwB+wiSNl9q/XDRA+t8urdhCSP0rIyPAjE6YvYSJ/BdaqX8+9qFemd5czzei9f909HBeoF/nep4tU4cLujkn+q91Q1+LON4ulPPNNtbmexi9RlfEoSy1OZfNNmVrNt/K5toTFHwsr+C9enypW/OqguT8odeVvfrZq754/M2hxkOzOfebbQ7VbL711/oz1DtKH3+t9WKsS8mYH9Nzn6/XufmN47ePc9e/Mj4d063V45Ixd1P6uFfHudu9M9Zqltrfhmsb4jRsCN91HjfbQozGb+P1ChFpg9l8m4hvP2W0XVVp+DbpsQNj/bH8hvRDs/k2b7Z58ug6V93Xp3mzTWbzbbllV4aX+VzdbMdObE5r+5xYrNRxDO2czV3zcU6Pb/xMdeiLCNLNtcnNKSUMTjiZ+iigGdw1sW9Kjz0e62RbannFu2P5mL7p3bHObdI1eggk7f9Yfyy/Ie2kZtcUw4/Svb7Omm+zNlfdnFKisxb6kxAED7CbU0r0Sezz6YgUFlDvtTenUV8tXcX8Fmkv2ZeuFzvWH8s9WrGpa88P9Y/KxzpjGvV+29+9KV1t6lP+6mQQzKtVG3qsMRbjjTr6VuXG8oR/7Kbbe2O63UI6xzelIx4GYtMvvjHtcY7e9u9ROkl0zeoAx+hrusTqc9t3Yxp17ShqWXC7tMczQa8fq+Aem/2S8ELfTQ/GdLbF2tz2y1G6QS199ygd9e1WPqbLfGrWVLExpdDW5CgdQ6iuzfGN6VzYFdrOmFYF9A7hqsvZd/hSo5PQ4cVRusF7vkN4jpp1Oo6SeMvvFYZkBmi311AS19WWISnO23YWRUIOTVcfxfzSVlEEZ87tCyuME+LkSeuDbiKdHcFUJLQ+DEksA3Jr2Kcs7VQfJfH5rXMkQfCBrP1dIeGYjLzY6kaPd2Sty13YEiVV0e0qScQ3iCDJxxzaS6mE2Pq9prL32BIwt8R+cO05SvJt1aV4sbE1AOhqMESqDc1KwUp1WJToF2rG+IckLsH7+XI2sqFaaY44BqeCuGCdzoyuetvXM0ylgnOStGPO4UK/Jb1UrzEtcN/dD5vzCV8Y+hL2KGpOYR1KBzpuF3US9XlU6KTJxPK3viTWpvV7hc8ury3gfrxtBAycUlt9V7AWaF9YwbSrufTT7vUb2jFc4kpbu/lu8lKLayvjXQAhxnjE1VjaeniPuUMr9dUTPqEbNXWAgGFIX3uvPdHZ9NgZtTn2CS8YrTRHXMfrF4rXWAJ8bIHavmb10G4nG6wrzWQjiKUFLVzBN9YZ/dgE76QfkOAzR0fr4n8m6nSqlVXb3iFJcFFXIaR+wEJOHFJ9qbAPda5C9dJ3X7RAcB2BbqK206Pk2HdU1FXSr2KV0jZM1B7oYAmdgCWMGuEsoDom3M1rZ2MR2y9kQlb0VU5WITivDdCb6KW5DSG5yi2jFXx2bfFSqEHascFtvGunImGDlbSxAWqnyhzpeLONAFm+lSV2c5mMPVu7L7MnwIMeoIwVTpvlnIAXOncZi6O2uTAg6jOql4JoBey2+jQWp5VprAUl0DEMIFyBiNXJK9gstWnCRUmfpqqWLboZq3PzXq9E6GjQp8Yw7+GaqvNt/YDABH6w6nwltjg0bYPxLbEEd2ipkL2a9LSGWH2xBTsrNZyyYjOTRYSMrLtbgZ2OVHSJ29MZjIvMp0CcJUAC7zoX8VlNaobm4hJgV8tKaFZa3voG28Q7blMeEtSkdRSPRGrupAusAESCTiD1ggeIaioCDDWVAxE2eFprA38Sna1VuxIxU9KnMQEbqLfc/7LsKNGwDvq9lHLrcapY2/ACRoStiRw4R1qWre4UcMB29IRYPzqfBS1cbb7kvkmkWtz08Gb1CbtfUqlvTak1JDXE1OgHurPdsoUJMBG0MY0jpKdhBfBES7EKGpyTgFqXnQgFELUNh+GlTi2Wh4Q3wSwOwzpG4XyJLmoqSA2KgLlA6AD9SqguZcaBWWS7m9yC3DgCpjAgh1GWbiRHuBiF/y7PV41DBUtbxRavtVBCcDoXrgA4WVpXncb6INXn01tuZ46VXzAmL3oDUCb4AOLDen/pBHgWQNsnEoYodoKJKQ6HLGFhxOvSYjgrOjCP+yo1k/MRn0fsRh9zRe/ATj7pZtVUCYRcsZPPvh80AFvDxXyZURJfaqlq8eorgUSYn2DBijTBFtCUYC3IMIIURzwYOxEDo+2t4G2QqmX4ulJrWZ0xBY4EFGlHJGgLdC8kSam9kUofWjsXrG0o4qye5VCw8KWohqKTh6cs3dqxASyMJEWxHlLO9RssemsxGrNT1GuPriu6qic0RqLqMDtRoTxbJqbccRjCiHg9SHhfa2sRq8SoF3GsNTWKIOFuSr+ix1J7lRZQQoyfqLetAmlF6hL0j8LVFadPiQPJvKfMjmIceJhq4041FoL+2An8p114WQKdoAyvXXrDZO9Lu2BywKEaI8cYOjm+RwQ77UDOihxRtOxpvfR1WoqNfR3BpdpWaSG9WRW9GLXzJVTaoAzIoF8uuHhTPEGPiHa51BnoNWBCpyrRtliNqqCR3mnsLR1F24JalrBIptlaNPYOqY4TEDoricUaVpFQT8Ir3iqTbfcOIcWrd80w1VVsWKd2rWGe3K4P4meFWvkFW8ckXYLPraRtKZnU3Jbfhi3IhM86oUmHjzwedSxUiDTugpYkkG9KGjotk6vF86RBATV+18tMJh+KL3RCl03rFIAWJbWmqLbYkjX+0RS8gk6Z8EhHBCU0ZIXYNBNW2tZh9A2RRhNREVkSgZA1JBQ7IVGCI+DPlHBcwKPkfM5aAkjmrdS2gEzEiNKZbFvQCLGlCsPLIXuclE0ZlJmXgE4E05n6EZap43QylYjfAhK5JJywTVWJFBLAcE3gC44P1nlYtR18Ny+1I0yVE09CTyyJCrw36MFmguxg8czMk0ilVZYabKXEORsilTtxgO84D+mBqzprcdpHjCmMD9yEBwCPdX1fEzcFNp6bgoYBI5FyazvU5PUr0UWiKzkiW1XtTSzoMRNeT/ebw6TfJT6TGhJLyKSMKwiHewcwJDflBCWlxsqEGCPhddHcVJICTzeVtkEcUYgK/vWmmhLe/DwB1diLHhcIwWkiVWbAT9LW3E94bbQ8khygnf3kBKSUhBIj+HLwFR+NGIgXKyQiaBiJSpwt43HEWHwhofDI+CnqBeLxP6Ax4SZwQe2N4sw0kLJG6PKTBp/jbSULeUQcOKLRTYW48ZVExume8VNtcb38VJMXH01oLJFEooHHMNlKRDtDXCodXiDAXyGUX4tDxVv9kg+Tdy5S16fqC0+CKDoVpkAwLh7hy8BmEgrwTZgSDjIriVgFY/jOYAwT0fwcVYpYcJwwca2rwXwVPCJgLI7zhmJip7ZJwLERE0FEXeIRzhFTNnEiphhPiE4YqEsUKXUz0XHtyK2Ko+8WuS9qggBcfKbdIpH7hCshTsSLS1QpEGBUKQ1/VQts1gr/FcmWYNJkQy7RmTQ1fgkJwBAJ13CKNLlabaaOb5h9moLLNVcShPPDQ0GMFVQkTclFHCWmqVMEacoJmwmcGGTAWZpKhddsErNUnMmT9TEYgmjFLJ4EDjsp6RQorHUdLVEH9dBnPHASHoZofRV3oZmDB92Zp6SeMEye8BoSKGl4UsaxAoh3JkgYc0RMMh9KMgX0nHxDSE0hUqPLlQS3GTXwEmIKJ55ZLVOnB8uUG7QpTDgoUpkK97QzROzz8KPqfFPVCReY0ZnKqqeUTGVs+OKsU4D6oiTUjBOHOkWcKVCSsk448gG8f5jKrVxwVdD52wJBnzTEGbQJABtMOgCDLYHWgOXgcxDJRKok1h2JhE8EEgRO1MpFkkZGm6oQmVG4gZPXhMJeSvReJAFjQBN4S9A6EXYJj0CONJF9ifpWCSEUHtVIAEchginOGgQYDoEjAOqEf0ugsXq+dABdOu8mHTGh5ziqAqz0xJ4jHiE9xxes5ExCMnErxU+ueickgq3coH6KLkei4E0pqIdPj5SJ6fJTKVxvCm74cJjYAZ4E0USJvgeizfwFXSJ+E04jSORQmIgwQXcJodAcIf8m8eBYEieihGrC4/RH4hQV5EucUhQuxohP1Paowt3CI4eNRD6VNAluOkh4J/QgTSEWxpimWMEEkqKH+lJRekPyZPGbQUJiINBSnhx+KkgEX5hvnIgKKEGeskvErstgeaAoBdjlPQlJweMddfKE7SMBdat1krDGhUNOjwsgkMmrE24xNIGTDOI+Tj7jsEcqSAsDrpNGoiXRbiBCFepiEwtUUXHc5ErAH67FTa/wmzyho5ydcsWzsrNTTRK8ceB4kbsdVMw6HM1OSjeSSK4V5PlJJdKgcW4Spx5miVoJ/ua47yHiHHe4lEiitE64CVc7+M2diOvqSDR2luNCTDjQnZLSrI4bDbqPoJuQIiQapHfcQMyoC4SIxpdvmNoBd1wdQcS4COsAhILQocAmx45xIDwR17iE2ooTlJ03BCAtOCYiQipb1qUpZFsjiYRLZRIlZnzHZDYVMwl45WQ5otAS58vlKSnocHlSPyzGFdzcgOgAGXPURPCF+S9UBhkqxBBlUiogCP/HdXJcqCRCRMbp6tyNOpVsCTILE6zqr9LbFIRa1CVz8/xLSXWSvPEEgg2gIQ2TriRi8kQlE+L2EoqMWNzEq/ZuarxKDwJHcFXvoMn4Dk6QIOE8vpbBjTUib3DEp4WGxGU08ICt7IkhGomgG+BwSCERfM5aknCKQ4nyy4znWqcgTj7izsZH6MqqJcRkzsZz17K/fFKKwZPQLeNTQ4tINGajx5wuUpd7DoyTaLAWbhQ+n4OUajwxjZ3PJGBcEKFtSklj6ZapEo3S+DpJ9FYT7Zb0dVLcn0Rh85gApyuDz9hJGS8kkg/gM3YqJeDNXBcgEwmOOpnIxFNSV0RBpkIQWxPUBVEuJBRtIZEII0qCKt4Ejg1M7eC5iWIgAeWhj2qA4x8C1JfPJEIKraTRIgEnPMUGE4CmthJjeYq4riKhTITALBMklgC7Ks8IuDXSUHap3z+E6xQX8bg++SJgcBkYLER1Zsv5akKZBI4mCY0VTCJ7oraGOuHdqpJwxYVMIsJpJlE84DTCyeS0RSYV7CbaKQHBI0QxmyhyiYL+RNxW4zY8spchSSL3oq+JhAoHSCSJoZAoFQNHPHtnwvRFP+lFpmF4I77ow4RsopLA37mQwKVbJlFwQ24i2KdIINHwkQjeGAiOnaCzgybATBIlBDf1JBqKH/MEE7GSCIqqRJAtrs9YQHViJqEcPRIJSo0E/piSiXVyDtdAkatA+CmwAQ3xXRUthT/lYqUgVWQjCc6AC94k6O1iK4mINSwJZWiqE6ycqIy7cWjW5KacrUsmgRtksFrfSHASDTAm6AyujxTgHoKphgkvWoFEge1mEgGIYUEnpiuB+UbCyNJD7mjYvylNPmLVmdKUYOuSaNRTIoYvoC0Rt7h6T0K5KiZxuUpyJPCc5ElkroTEDZokkfCJkOOE1E01FRI1E0A1EzQaB+rZwtBqicZNI3h2rtWTCIUwvplA4SE4ozGwcTqfQbJK1ZJURRIl+P4TQ3RuRRWyn5pAKXso0Woy51P0twn6cY+vTDCcqdvMdyPoFk1itgyikHHc70Ik4QE3JBqhTPBhlW8QLVayDSQQ5hQSReD35wJzn/4VeJNBE8qLMrlOVorNJAgZHkgkAWvORCAHTBfcMRIavsBHlJpIQLwlQ3x17XORKViAVZEpM77ChcKiFRADghoXN8VqYyHBbetN8ZNohO8CXQy7hKAGyoMnSHjiDJcA7zY4EsoaJEFAaG+I5awYeuFqT+JI4HYxmAJBpm8lIDf0R5pSyVZIEGO9Goyu8XNIAvFXIoH3NGcIT+GiCySgcD0JopxrohZsmYrGhS+VBCGVA4nibRFTEbpA/VSiFMA+qHAuQfCrnWoMPhqN+OwDv21OK/urejHEXNadW93kg+NtwFeGbAF8gVJVbgKYGBXWAvhEVS6BExLE8o6mgncL5E9QMFHD3EDkzJViqno6DJVE5BoiwZrVBKmSCglfQMBrAgpFLVG9BFMzOBzEWZ6Cc6IlMKoLJRUxp6l6DZRMgtjUWtLkn7VAepVkakUqopGiIdOh5Kpy0QoJJc+MWDs1bN/aCQN+fmHVaKJ62JuWEwtLkJgVoO1WJuUMWq7YlI1YBx6kvz6DdltoImgWq67YeOSRLfIybjY1FjWYDZ8HsRHqMr80DF5DOALibnsiF9gwpeR4mXsC8sBypLUBZhp03sYpBugOCwQkIoKNSuUasVAphESw7FntDYccdqBNUxFIKZvVG6AmHNuHEjaAJmCia6Kww8QW+IB0rOCDlMqFrtJ5yG/tWJmKBmy3FVUC3qqcdpqAotbhVA4uX65Tzrk9wmulEYFaZlxKJEMYQhvD1SGRShBNNOILShhxJglnIWUggJWXTGzrQNx4kSkrgSZE6/BQu3AzmSiIW+5cEh7/hiSih0yEuCUYC4kC51I5yxCyeCol0RYTQpYbkwRBTDSRA/xopV+h5CBgHYxuCNjEMkC3WoYN3Yo8jQRuVTWBn09NwCTRRNU4FtCwlWjv0K44ciSBqoImULrQRIR9xKNMMHISwGxNVE/IeqhYjYoOFYvaFwkH/Qc1yz1EQlmRJGKbmTQl2C2U6IpJmpSlQqICKvEqYREnkhDwRBLojWhC0TBKCGfTEjAnKYnKVxCIHe1NnjQeMb+1f68k9hZ+K7wOs0yN400CuEBCjyK/uuwFISQ7AiIZ3gaJuaoKFinhktAE/OOWyBxSKfBh2ncVhackw+jRhB4CjZmiAy8qutFHpW06OE4AE6W/dTnBtTiLJJJusQI9QxPgr7rFODG6nBBuiDJIcLNoAmasJsCANFEyXa3EqemVERe2Et28PFoSyyOd04pjUPpcJ4vuhSaWJpZGAxCKOkvHWleXzo/DWQY4D7n0YIhMS5sN1Rc6mrp1MufpbTBvXIIIyUvJskzLwgWnIGFY3GW5lw0w7whRuD7umb6Lln0FGgjzZtx7WWrb3cv2XDbssoVVo4q3NAw3Caehs8ej0A/Hclw4QETgHk7UcsSWQwfXQKv4ua5LMNvGs2sVOEvkfLNFhhOvyDYlECL8LmDCNyiroITvKXAhAcODywfO2QyAFpC0AKkFbHl0BKi8gDaVPnTwx3f8ChAXEBmigoUBjC6AFRnKCHcHSLzA5gVaO3RyANIWpvEA4weov9wDy82w3BX99lCRKqPghtErYrhzlPPPo3Yv2dpvKr27WL/xNgsQN5S4hhjoHchdutyKwz253JwuuHapook73reAxH4DL3eygygf7u3lJte7nTt5uO1VoXLACAYcYcEalLikTscslCY9Rj5UMt8QFMVCZoxlQGFmnGZBcjrasyBCiho1pGlBltqBXPEpECz2jGJcPdFxMMj063jajLgpJrfgdortESmsDvhfxwgVR0zHWOOCR6pC5YBrIrjo2GdwiIJXDDWBs4J9D1jsgtd2TLdjvgMm7F3HjTu2vODP4NM2jxi24tzuCAlvHOKGqLsyoO6KzINbK3rvKomO8C8kgBIF2a1UwkI2DITEQlp0UctMfowEiUtKogD/lWhRMgY5zEDYdFJnIX6UHIoyEEgLyQQRlXKGiFrIqk4ALqTXQIx18mwh2JSEk9iJOikkOpnnO6thJAU7cbiQi0pAFqUkXQ2pDETmSnUuZGjv10KqDsRrI2cXAreRvE5ILERwJ4vLTCgPpHMnphfyeiC4VRd3pckHKn2h2zslD20PW7RR+/AVFvpfOQJVYA0sPILGNJiZCMpVqAObYWE8uJUVsTAnbGNXLAyMgaXRmRx2ZnsMjJDGGmnMksY1QaKufJSFodJZLDPPZWDCdLbMwqgZWDeqGNG4Ozk1fo8yhwYOUOcJLVyigW/UOUkLb2ngNnX+08KRGnhUTd+n8bG88ro6Z2vgdXXul/LD4EgphywOLLOFiaZstSiw1VQLeOG8Lay4lTXXWXUL825g5y0MPvRw/MAEHNiCnVG4sA6VmeiLiYiyOnuxMxyVBQkjdGBKdjblwrgcWJkLc7OzO7n40TBpHFFJA4904ZoOfNSFs9p5rXBfg+QjfqyTBL9y4dmuTNyFq9v4vAvnd+QFq5YjJRZTiIGDPPCUpUJarXzngRO98KYbeGz866ys7Qij+IjHbQPwpLHB4asPjPHOKl+Y52Flpy8M9s5y9zDhU+fKw+tWPj1M/YVzv/Dy4e4H2OADv1+CSgAWmQBSglxVXOBL1GijXMLRDbKFQdrg4SAvAolFRDEKLRYxRlOgXkUdg/AjuGBVLtIFJIvIZBCiLGIVFITrIHpBGEOMUxXOFBXKoHGg4hobEDh6bgoV6QxCnkXs0+RAi2AIlkoXFS3Coy5O4s6N6UjitMig2t5uciq0d5BcWXqjsqxSSARY3Ku4axGADSKxRUiGApqKzbogbRGtqbAtEPd0Fb8tAjkl/RaR3SjEa1zwUdDXRX/QXh5NJhUPEkEUlg3MpCZCVMFhlykuUkaVOyJKRBLpxZFYZJNdWon8ErWlJtEkdukg45yFnl0KushFlf98TXa6SFO7fHWRuA4y2EUqu8hpF8ltl+Uu0l2V99raBcDokQ0iYdWwRH5cS7XSBMmQIqNoeRE2L+LnLpBGTUJF1Kpx16TXTd1glGsvku6GBa3S8EE+vkjMZxH6IlNvUvZF7t4k8TLL5jWK6iqt91lgX0BtdIm+VVbKIvQf1ABmvYBFUWBRHSCabURToKsXLAoHqoIAsTkoJSxqCq7Clxg0GbpuQ5q1HQb9h1khQsVhTWUC/ZNBiaKrVcRZ0UJVLyCuBmWMRT2DoLyJR4sGh5oLUNKUPND5CMSTVS0QVBtVL2TWE2n6IYvmyKJLsmiXEDJItVUKiu9oouQUoQtVSYWfRWkFSx0taarvTdWFn676osowUHerdgzaSYmChG0FiUWhBlmBatYsSjeSMly8RTFHVXWU5l2VdxZ1nkXBZ1H5WZSAvNMpHxSFMLeDtTArE6l6EaQ3+G1XOFpUkJLUAl2GxoPya1fFpUWVaVFuElVDWBSgFpWoVUdqUZpa1KgWxSqsCENXvkpVRnWsduQHja2ubbxodQ16Xovm16ILlnxGUU/1xdAQXfTHFn2yRcNs0TnrOmjNiLNpqTH5qreGQhOKbNjKNNU2ND8GZbdF/W1RiGuGcaozxx0+aNEtenWLpt2ie6e6eKqaR+zvQVlvUd/rFmaDht+i87doATbTrlFTcNEdXLQJF/3CReNw0UHsWomqp4hm5aC5uOgyFuzH0JZc9B2bRv2gErnoSEa1/hn1KBfNykXXctG+7Ip5g4Kmvn2kxNmMKFdFzzBV8RhqDsqgBbMQ3xRG0Y8ZVEhnndKk1lKj2mnv6aCaSrAidB5W9dVZn3VRcHVeRR+qBBtpSDwqlk1PFnAwas5inKZqtslrEPVV33bRwF10chct3VVvV4hSPur2Ltq+i/5vjJhRq4ow4xuVhhc14kWxeFU1bmqQqo7MSRz0k2eF5UWFeVVqntWcF8VnjQZHB7tytKpLYyg1KFAvKtWrkvWidj0rYjcx7qirPWtvp5ysQ2l71fDuKt+LDnhs8li36okvmuMhZYdOGVdbTZpo2uaL9vmqj75oqM866y5jV3Sk1+5ULcituu+rNvyiHz9rzC869KtW/aJnj6Uv6v5dF3/Uzl/09dUYiEdq2Txq+a96/x6djME0YDEWWM0HmqwHhUNl5q4mB81aebRKmE01VssFwq6Nxg0xtbDrGEBgUDuaRMxGEovZBOxtJB6DacVibMHmavYYWA2h6DRbaKRmQOJRMMf+wXuNfK+GHhVxwGL5sZqCLMYhi7lIMzLuNyJTtBiZLGYnix1KLCpnG0xVFuOV2ZplNm9ZDF7sYgKzGMXMM7IazjSe42Jcs9rbqA2OKjajTrnY5RAcUO2MBvud1aZnMfNZbMhK3yxojs3Wx4PV0GpJRBBCNQcaLI5c7m4BSrOB5CvKSFHzJnQw1IJpsGpaLZ0W46fF8Dv3A2ZHu6nVlmq1r1psrlKNGOge22at9lpYmahR3GDXtdp6oZeiNk+DTdhqJ4aRcLMdg6ZqVml1tjEb7M5WWzR12NBt1lBqIiBa04sghfY2sxJDYLFINRtATS1mcaup3GxBPljUNZcDze6uuz8IBA/MfHgwz1NTVjXZW834FhvvwdxvNQEMrvW4Wanpq6v14GpROM/YYHi4bKTFPnE1WRzMGMtsmjuYO2Zsm7A+86tZ5GoqqQYt+jSqtT71OM9q3zyYXq7mmARPVNvK1WpzteTEO0Gz7hwsPhcr0GXTequeI0YD0sGmdDHhHm1PF3vUxdZ7MFtdTVkH89bZ4nU1gl0NYxcz8tGA1nW/F4Oh7WB8K6kZGA9GuoPh7mLMiyWgGqg7iaoOiZKV1G4IvBoHN8N7nlp1v6AmxrNh8Whs7LulvjTdVaKxLHbKRRdNy2iVZRlMnFer59lYdbUJTwVfM9QnRqz6rBjNqmfj/cH8Os72sKOZNi4xGATQt5tzDybe0WGzwNPVFHwxD8cnQDMZH8zI5/07mpsvJuiLVbpLapV0bLy+GrQvRu5SFJVAa2s1hsddDUePILDNnHo0pF+M61eDe1XY16lafKuMxvozBLGoxTWr/iNb/8EDwOAXoOPZOANgW/L50bFAqYBBAM6REwK1fWOVjxwW4DeiuzEYnBvoMjH7oyOEwT3C6DRhcKWAf5Lm0SVjn93cw4zOGFYPDYvN++jMIaH4p1ejTcsms+vVcOQkQtFxFuPIn8TgZWK5FUaPFIOfitF7xeDTYvB0Mfi/GLxirK4yBv8Zg1cNsXhKoWODB46AuVdzSDJ66xh8eHRj+IZgNMcnITbeBHZIg1+Q1dUTPg9mdyCIN5ubHkW7mneWYDMIgH5hcU2ygv/RicmCW1ifuwsHTo+SJrg7WZ1X+ai2WFoag+pfXHOusjrAGh2xDO5ZnHPdS9PgUEs5aM0NxujyZXUE4yoTont4dBpTii3NZcroYAYj/+YoSZepeSUaXdQMjmvUQ42eGBejtGUaPd+o9Z+es9FJjp8v9tGfjuq1ti4i01KfB9Y5cOT2AYFWaRVQoO2FTXyozn9WNz+2GQZoqOTZF4Bd4fux+6DVqVBJuXlzOHZAtLolmoG7cgsXD0aLwzkYVs2JCBGGQRK1C0mZC5ocvCXNXnSOvCqBBrWaodTuuEYWkG4bFaUjV0aTYnFWfFZrcobrY3M5YWXxBmCx+pI2c8rNbWBLnAY417cW3HBx/se3BkdUCyg/dlo1A3MrorryvKVwuHVx9IU1eMjCP4KeLFG+aevC6mNLI4h3Z0LHLrlW+H3swksZx9rkkZuv0f3XAMUtpuh9J9ni8WKmO9gWmR0y2CO3Y6tjss67wy/amAwzymqPXZ5BDjcYYle3CtfSo3u11WePVUKynafj9Oimbbl+j128KWrefOUcuYSLerno9jp2IbdCbosFeXdrdex+Tv1XtDm+yXXdTW7vBmB97D4Pi6MGlUafTccu+dTpSVufo/TqGmrwHXVDcoXLR8kVMI9AGp7jDJpvTq+AevSeNQLl26XVrZQehhEyH6dX518jnB5B8nE6skn0lrG3SnfHkeodckyvbs/sUfoGb5VH5SswtrdJ3+QZc3AE98PTK2DGhcvqFfSm9AqcbU+r583U/ZV293ibU0rUv2ae/ZWuDvfsUXr13zcC1n+U1lZLb7WD080pJbgvlnpLH7Cj/9YxXZpfGnWY/EPTg1/Zo3fH8jE9fJ9ROXvdr62W8AQfsKD4UTTs/eaUksKT6/6ax2qoYsw+m8fyIa3t+uu+fSnRr4fN9gGqLbBEcLU7poev61fiZjs+3pxSwk5wuAh2Uzs8fGVMj6/80PLbvHuLOtr7vNmOzW9OKdE5wKuw4KckW3+99+MrN6XHd29K3+bd29S56fs3lQ/f1Fmom+1YdXOqJQWH3HazxQ0CLrnxizxWG9P2KDc+GdPHtY5zx/Xwkn1ccnNufDKmr3/j+NmYG9PNP/fLJa//Wnvr5ndf/b35rfHpmJ6f8zuWj+m5zvWy6/m53vGX/lHp/Hz4mvppl83WB+KLfUtyeMqecTwMQDCSw8OfP5DMe3d+UCiZd28VSubR/mpPRQ3BRTyo4+h2xFQgXsw/K3yM//mD293550a30winPzxMVX/t6ukt4tt5IrL8V40o8963h9Ory18rrky4TVAxjazF0qzR5ucTdYtghMtatVXbv9i0UDO/WlC7OWrJv2bUt1/hPP3nWaJbB5ZBBKlx539s7LefFlzm75t/219cbLYteEsKG7N5ejgFvTKEbhGz+XazfaCuhgN+I4L6xMC2OeAIFG9LxuMg1mONq3ZlMDaNq9kal701LuJzBvGgumGEI2jUjEQKxlBqnIQnrqbopZplEa9DKOOKwwJWbdTQzJZm14N1KUbexA2HQqxZQxZ/cn51cnbab8b3n59+RfaHxTC+FjXyfbnpZhZ8xizBtB6fXI03828e2Ie/eyAPf/vidw/cw9+++Hf3uweeX//W7x6Eh789vDj/zV07xd/+5jcv7j6ID9/63//9QXr41r+7t67d6T8oJJwQn5MIth0DAEj9qHi1N91HH5/vL5SwJaLbswb7Pj4HCFqFq9+d95iVH58/Prm63LTb7eN3399sHZP1ztcnl395vr84bLZojqZasgQIb93/L8VGWoIbzhv8/v7ivYuLswtd0lfvWrC7zfYBrioQrnk081Fbsrj7kZhwp1Frl0LZ6JONueDV2GKZ5CaE9T4XpyJJf7gbdWPd3198dHJKMO0bWz3fEF6YtjFDFYMnchTcl/e5LH7s+5/uvz3c/DaevtohjTbgCMV5KyhSGpzB4AfW14RRizPRVUElEqkLvD4TBZcPmI5zpHHjYtXoKYdYU8G/WaglM3k2wv6UZEKp6hq32GyzT0AE/CkoqyRWHI6EbFEDsgE5QKzRhIRKZkjiYENjHxAduiT4wsIrokHwg16jrw4HiDj6Ceh0Vxj6guYOji9QgcShXwYABMEHh9jg8Jwfkwn4RsArGb5iMFvHgnvKxbqE0VQyvoYwpSI1lFQdHppKUl8WVSVAeJVkJGjShSIlSzI+Oz8FyQkhMDYimNPg6sR7lJAcHpxKmFxJPhW0JvDghHltRSXXIYgzyGynKF5Cxv9GBFDiYSRk6zkBxnucWllr0ZqzeKLCi0p0GJwF9cvk8GtGdawsMDIQ3LR4n5NgB4KCPh68fAgeSX/AB1bBUrYEj3ZPMN4GNYgv1jtUNo23GIb6UCyybKvw2iE2RTHB4xDAVY9FEN4/XHSqmo+FHL1PMSO6cMUXWBEWr90FdfmMXVuQkgQ1VjyD4WJMIi7sUZcxLiX8S6Lr4tUDiku4iEDWHZMLODeNuGBS2Tcx0vCDVnFXg6Dy/6vt7FojKaIw/F8GLzZr91B1qk7VqSCCrizuhbiQ6E0IMia9ScN89HYmGJX8d3nOGDMoghj2KsxXTU2lu86pqvd9TvO6BJIxZOSMgAENZycJnx9jFXETxA7NFGYVra2CfeCGgA1TNNQidFQK4nhHMnLUSGSi5ymFHBQCO1oH+GTGVYx8RyJuIKtwxC34E24qdwp/wtAljHnL1Upk31I6zJzLmqNY0nJ4wmlxIecggFxjw/CoOTEVYmZueK9CUBE3vbcky8jOthaKWHSxgecQbM5NIbggyEBvAF2wgK+KljGcVT2wUKMJFmH0PylFJ/i5ZbME4UwaX34tddmggcSYExrummWZQzWJoTr9onJBO2e5BUdmFLdymsAZrkTqUmyZUkD47NzzWIBwIUMoNTpduAiMaK5HMRwbURtq/cLF6XYBrWjfQdhB9eENGrHSaqsawRBFhQpnKlUoVIHmPsBHEgMcjkQ8w8XB1AoWgHwiM2UkJxVz2NrFjLEoRFQmB1IjYYiZTf1ooIvJeS3WlAIonnlgaTTudHh0wCJV0alikfVepYRqVviBHCSSqwDsDUWNXqGbd+IN//DkOCCUPUjGxFpFI9RFyY6HYLrIDo4UweTWVFEBuq4e1gosgZqd5xqjEVysIkKnfkyMMLCAk1uSg70+1WVpBDYUYbjPscw3zVUV30kM2G8RwzNJuJ8bO6MWfmcwJ3lmxS6vDSoPMAcmVEkFYLSA5mgYjJMLdFpiVmlVEBUTb0yZyxrQtAAGzAqKaQzcXSzBwyLFD1fr8bfhenH6YbW+G0iodvPm3XY/3Myr9aHg8dvdvLlfr4h/P+1/nYbtasPux/nT893iww/b8eP98O6bwwe+9twjghX5NAnam/W4vXm/mlebYT/8h2yk5SR4PiqXIYSvHh/MkjnHohlRAw5MYeaJBMmcu55jci1CFZ+G/cYBIP/ILL5ar587cjbs94zlfr5nKP2Fu8Xpxe+Lz6bVSI7Pny/Oz/bzuL3pxu2+expFf+/38/Uwf8mAjvMdFcensOgWd8PVbnu9OA2P3f9viDz6qaH4kobkqCF5SUPpqKH0koZYMT39tPyShvSoIX1JQ+WooULR+Od7aPHqYmJtMh0WJ5Ncvj7/brW/PT09+/jq4eTziyldvp52v7x66NIJD49XK8/v9HXLdFi4TOXy5MQXLkergvvNz8N8uB3Zk1zN71czK4ILevPjcLXfzUf3PZVyn65WtjhYK+w20+5uPCzofNvjr/riPPL54xPUGPcKpv+2CC44MbB3NUMaH7ueNOPP1Njw1uJ3rSR60iwX6OG92dKKVXLIlNmbhhMbDU6YoOR0zxwZlpEG5Fpi1wOrdT0CPEZCWZ/JoiylokT6yPRQi0kmvwuWoyoAAbRwaKzE3FamnkxYqchj4eLxzVCfVNVrkVlZVuZpikNl0sQ+EmYqOR3RTLse/C79jkEFswPo7ypW1KKkCARQAStTRaxJQLzd9aksc0bCkBGxJO16/OKlSMrmGkX/kDaIxK2QlxWnPZpwtkzJJpif1KiIKUd/XXPXy3O0twaXuhdAnOTJLBGod9KTRUAkpK5Jo/IUVk5g+bGQf6p0fWJp6Ll6Ni2kaeDhMhlBbGqM9VEUjLgS/r7X8Ph46bfV2dU8TlzWvhf/7Xhzux5vbvdvdtvtcLV/Dmpvx4fh+rDj6HHu8Q/RV4Xs').then(json => {\n",
       "   const obj = Core.parse(json);\n",
       "   Core.draw('root_plot_1779221451878', obj, '');\n",
       "});\n",
       "\n",
       "      }\n",
       "      const servers = ['/static/', 'https://root.cern/js/7.11.0/', 'https://jsroot.gsi.de/7.11.0/'],\n",
       "            path = 'build/jsroot';\n",
       "      if (typeof JSROOT !== 'undefined')\n",
       "         execCode(JSROOT);\n",
       "      else if (typeof requirejs !== 'undefined') {\n",
       "         servers.forEach((s,i) => { servers[i] = s + path; });\n",
       "         requirejs.config({ paths: { 'jsroot' : servers } })(['jsroot'],  execCode);\n",
       "      } else {\n",
       "         const config = document.getElementById('jupyter-config-data');\n",
       "         if (config)\n",
       "            servers[0] = (JSON.parse(config.innerHTML || '{}')?.baseUrl || '/') + 'static/';\n",
       "         else\n",
       "            servers.shift();\n",
       "         function loadJsroot() {\n",
       "            return !servers.length ? 0 : import(servers.shift() + path + '.js').catch(loadJsroot).then(() => execCode(JSROOT));\n",
       "         }\n",
       "         loadJsroot();\n",
       "      }\n",
       "   }\n",
       "   process_root_plot_1779221451878();\n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from ROOT import gROOT \n",
    "gROOT.GetListOfCanvases().Draw()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
