{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "c04d00cf", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Requirement already satisfied: coffea in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (0.7.22)\n", "Requirement already satisfied: fsspec in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (2023.1.0)\n", "Requirement already satisfied: mplhep>=0.1.18 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (0.3.28)\n", "Requirement already satisfied: uproot3-methods>=0.10.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (0.10.1)\n", "Requirement already satisfied: numba>=0.56.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (0.56.3)\n", "Requirement already satisfied: correctionlib>=2.0.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (2.4.0)\n", "Requirement already satisfied: cachetools in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (5.3.2)\n", "Requirement already satisfied: scipy>=1.1.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (1.7.3)\n", "Requirement already satisfied: pandas in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (1.3.5)\n", "Requirement already satisfied: tqdm>=4.27.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (4.65.0)\n", "Requirement already satisfied: lz4 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (4.3.2)\n", "Collecting awkward<2,>=1.10.3\n", " Using cached awkward-1.10.5-cp37-cp37m-macosx_10_9_x86_64.whl (12.0 MB)\n", "Requirement already satisfied: packaging in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (23.0)\n", "Requirement already satisfied: toml>=0.10.2 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (0.10.2)\n", "Requirement already satisfied: cloudpickle>=1.2.3 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (2.2.1)\n", "Requirement already satisfied: typing-extensions in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (4.4.0)\n", "Requirement already satisfied: matplotlib>=3 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (3.5.3)\n", "Requirement already satisfied: pyarrow>=1.0.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (12.0.1)\n", "Requirement already satisfied: uproot3>=3.14.1 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (3.14.4)\n", "Requirement already satisfied: numpy<1.24,>=1.18.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (1.21.6)\n", "Requirement already satisfied: hist>=2 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (2.7.2)\n", "Requirement already satisfied: uproot!=4.2.4,!=4.3.0,!=4.3.1,==4.*,>=4.1.6 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from coffea) (4.3.7)\n", "Requirement already satisfied: setuptools in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from uproot!=4.2.4,!=4.3.0,!=4.3.1,==4.*,>=4.1.6->coffea) (67.4.0)\n", "Requirement already satisfied: importlib-resources in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward<2,>=1.10.3->coffea) (5.12.0)\n", "Requirement already satisfied: pydantic<2,>=1.7.3 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from correctionlib>=2.0.0->coffea) (1.10.13)\n", "Requirement already satisfied: rich in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from correctionlib>=2.0.0->coffea) (13.6.0)\n", "Requirement already satisfied: boost-histogram<1.5,>=1.3.1 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from hist>=2->coffea) (1.4.0)\n", "Requirement already satisfied: histoprint>=2.2.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from hist>=2->coffea) (2.4.0)\n", "Requirement already satisfied: fonttools>=4.22.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (4.38.0)\n", "Requirement already satisfied: cycler>=0.10 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (0.11.0)\n", "Requirement already satisfied: python-dateutil>=2.7 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (2.8.2)\n", "Requirement already satisfied: pyparsing>=2.2.1 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (3.0.9)\n", "Requirement already satisfied: kiwisolver>=1.0.1 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (1.4.4)\n", "Requirement already satisfied: pillow>=6.2.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from matplotlib>=3->coffea) (9.4.0)\n", "Requirement already satisfied: mplhep-data in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from mplhep>=0.1.18->coffea) (0.0.3)\n", "Requirement already satisfied: uhi>=0.2.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from mplhep>=0.1.18->coffea) (0.4.0)\n", "Requirement already satisfied: importlib-metadata in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from numba>=0.56.0->coffea) (4.11.4)\n", "Requirement already satisfied: llvmlite<0.40,>=0.39.0dev0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from numba>=0.56.0->coffea) (0.39.1)\n", "Requirement already satisfied: awkward0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from uproot3>=3.14.1->coffea) (0.15.5)\n", "Requirement already satisfied: pytz>=2017.3 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from pandas->coffea) (2022.7.1)\n", "Requirement already satisfied: click>=7.0.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from histoprint>=2.2.0->hist>=2->coffea) (8.1.7)\n", "Requirement already satisfied: six>=1.5 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from python-dateutil>=2.7->matplotlib>=3->coffea) (1.16.0)\n", "Requirement already satisfied: zipp>=0.5 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from importlib-metadata->numba>=0.56.0->coffea) (3.15.0)\n", "Requirement already satisfied: pygments<3.0.0,>=2.13.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from rich->correctionlib>=2.0.0->coffea) (2.14.0)\n", "Requirement already satisfied: markdown-it-py>=2.2.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from rich->correctionlib>=2.0.0->coffea) (2.2.0)\n", "Requirement already satisfied: mdurl~=0.1 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from markdown-it-py>=2.2.0->rich->correctionlib>=2.0.0->coffea) (0.1.2)\n", "Installing collected packages: awkward\n", " Attempting uninstall: awkward\n", " Found existing installation: awkward 2.2.4\n", " Uninstalling awkward-2.2.4:\n", " Successfully uninstalled awkward-2.2.4\n", "Successfully installed awkward-1.10.5\n", "Note: you may need to restart the kernel to use updated packages.\n", "Requirement already satisfied: awkward in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (1.10.5)\n", "Collecting awkward\n", " Using cached awkward-2.2.4-py3-none-any.whl (622 kB)\n", "Requirement already satisfied: awkward-cpp==17 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward) (17)\n", "Requirement already satisfied: importlib-resources in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward) (5.12.0)\n", "Requirement already satisfied: numpy>=1.17.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward) (1.21.6)\n", "Requirement already satisfied: typing-extensions>=4.1.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward) (4.4.0)\n", "Requirement already satisfied: packaging in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from awkward) (23.0)\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "Requirement already satisfied: zipp>=3.1.0 in /Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages (from importlib-resources->awkward) (3.15.0)\n", "Installing collected packages: awkward\n", " Attempting uninstall: awkward\n", " Found existing installation: awkward 1.10.5\n", " Uninstalling awkward-1.10.5:\n", " Successfully uninstalled awkward-1.10.5\n", "\u001b[31mERROR: pip's dependency resolver does not currently take into account all the packages that are installed. This behaviour is the source of the following dependency conflicts.\n", "coffea 0.7.22 requires awkward<2,>=1.10.3, but you have awkward 2.2.4 which is incompatible.\u001b[0m\u001b[31m\n", "\u001b[0mSuccessfully installed awkward-2.2.4\n", "Note: you may need to restart the kernel to use updated packages.\n" ] } ], "source": [ "%pip install coffea\n", "%pip install --upgrade awkward" ] }, { "cell_type": "code", "execution_count": 2, "id": "fb0fc16f", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/billyli/miniforge_x86_new/envs/phys239/lib/python3.7/site-packages/coffea/util.py:154: FutureWarning: In coffea version v2023.3.0 (target date: 31 Mar 2023), this will be an error.\n", "(Set coffea.deprecations_as_errors = True to get a stack trace now.)\n", "ImportError: coffea.hist is deprecated\n", " warnings.warn(message, FutureWarning)\n" ] } ], "source": [ "import itertools\n", "import logging\n", "from pathlib import Path\n", "import numba as nb\n", "\n", "import awkward as ak\n", "import click\n", "import h5py as h5\n", "import numpy as np\n", "import vector\n", "\n", "from coffea.hist.plot import clopper_pearson_interval\n", "import matplotlib.pyplot as plt\n", "\n", "# from src.data.cms.convert_to_h5 import MIN_JETS, N_JETS, N_FJETS\n", "\n", "vector.register_awkward()\n", "\n", "logging.basicConfig(level=logging.INFO)" ] }, { "cell_type": "code", "execution_count": 3, "id": "cad6b150", "metadata": {}, "outputs": [], "source": [ "# read test target file\n", "test_file = \"//Users/billyli/UCSD/hhh/reports/bv2/hhh_test.h5\"\n", "test_h5 = h5.File(test_file)\n", "\n", "# read baseline prediction\n", "baseline_file = \"//Users/billyli/UCSD/hhh/reports/bv2/pred_baseline.h5\"\n", "b_h5 = h5.File(baseline_file)\n", "\n", "# read spanet prediction\n", "spanet_file = \"//Users/billyli/UCSD/hhh/reports/bv2/dp_on/pred_v53.h5\"\n", "s_h5 = h5.File(spanet_file)" ] }, { "cell_type": "code", "execution_count": 4, "id": "cfb9df25", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "s_h5['INPUTS'].keys()" ] }, { "cell_type": "code", "execution_count": 5, "id": "458d96f8", "metadata": {}, "outputs": [], "source": [ "def sel_pred_bH_by_dp(dps, aps, bb_ps, dp_cut, ap_cut=1/13):\n", " # parse predicted bb assignment by DP\n", " dp_filter = dps>dp_cut\n", " ap_filter = aps>ap_cut\n", " ak8_filter = bb_ps>9\n", " filter = dp_filter&ak8_filter\n", " \n", " print(dp_filter)\n", " \n", " bb_ps_passed = bb_ps.mask[filter]\n", " bb_ps_passed = ak.drop_none(bb_ps_passed)\n", " \n", " return bb_ps_passed" ] }, { "cell_type": "code", "execution_count": 6, "id": "97f4563b", "metadata": {}, "outputs": [], "source": [ "def sel_target_bH_by_mask(bb_ts, bh_pts, bh_masks):\n", " bb_ts_selected = bb_ts.mask[bh_masks]\n", " bb_ts_selected = ak.drop_none(bb_ts_selected)\n", " \n", " bh_selected_pts = bh_pts.mask[bh_masks]\n", " bh_selected_pts = ak.drop_none(bh_selected_pts)\n", " \n", " return bb_ts_selected, bh_selected_pts" ] }, { "cell_type": "code", "execution_count": 7, "id": "81791964", "metadata": {}, "outputs": [], "source": [ "# A pred look up table is in shape\n", "# [event,\n", "# pred_H, \n", "# [correct, pred_H_pt]]\n", "def gen_pred_bH_LUT(bb_ps_passed, bb_ts_selected, fj_pts):\n", " LUT = []\n", " # for each event\n", " for bb_t_event, bb_p_event, fj_pt_event in zip(bb_ts_selected, bb_ps_passed, fj_pts):\n", " # for each predicted bb assignment, check if any target H have a same bb assignment\n", " LUT_event = []\n", " for i, bb_p in enumerate(bb_p_event):\n", " correct = 0\n", " predH_pt = fj_pt_event[bb_p-10]\n", " for bb_t in bb_t_event:\n", " if bb_p == bb_t+10:\n", " correct = 1\n", " LUT_event.append([correct, predH_pt])\n", " LUT.append(LUT_event)\n", " return LUT" ] }, { "cell_type": "code", "execution_count": 8, "id": "a56f264f", "metadata": {}, "outputs": [], "source": [ "# A target look up table is in shape\n", "# [event,\n", "# target_H, \n", "# target_bb_assign,\n", "# [retrieved, targetH_pt]]\n", "def gen_target_bH_LUT(bb_ps_passed, bb_ts_selected, targetH_pts):\n", " LUT = []\n", " # for each event\n", " for bb_t_event, bb_p_event, targetH_pts_event in zip(bb_ts_selected, bb_ps_passed, targetH_pts):\n", " # for each target fatjet, check if the predictions have a p fatject same with the t fatjet\n", " LUT_event = []\n", " for i, bb_t in enumerate(bb_t_event):\n", " retrieved = 0\n", " targetH_pt = targetH_pts_event[i]\n", " for bb_p in bb_p_event:\n", " if bb_p == bb_t+10:\n", " retrieved = 1\n", " LUT_event.append([retrieved, targetH_pt])\n", " LUT.append(LUT_event)\n", " return LUT" ] }, { "cell_type": "code", "execution_count": 9, "id": "f57e2d53", "metadata": {}, "outputs": [], "source": [ "# generate pred/target LUT\n", "# each entry corresponds to [recoH correct or not, reco H pt]\n", "# or \n", "# [targetH retrieved or not, target H pt]\n", "def parse_boosted_w_target(testfile, predfile, dp_cut=0.8):\n", " # Collect H pt, mask, target and predicted jet and fjets for 3 Hs in each event\n", " # h pt\n", " bh1_pt = np.array(testfile['TARGETS']['bh1']['pt'])\n", " bh2_pt = np.array(testfile['TARGETS']['bh2']['pt'])\n", " bh3_pt = np.array(testfile['TARGETS']['bh3']['pt'])\n", "\n", " # mask\n", " bh1_mask = np.array(testfile['TARGETS']['bh1']['mask'])\n", " bh2_mask = np.array(testfile['TARGETS']['bh2']['mask'])\n", " bh3_mask = np.array(testfile['TARGETS']['bh3']['mask'])\n", "\n", " # target assignment\n", " bb_bh1_t = np.array(testfile[\"TARGETS\"][\"bh1\"]['bb'])\n", " bb_bh2_t = np.array(testfile[\"TARGETS\"][\"bh2\"]['bb'])\n", " bb_bh3_t = np.array(testfile[\"TARGETS\"][\"bh3\"]['bb'])\n", "\n", " try:\n", " # pred assignment\n", " bb_bh1_p = np.array(predfile[\"TARGETS\"][\"bh1\"]['bb'])\n", " bb_bh2_p = np.array(predfile[\"TARGETS\"][\"bh2\"]['bb'])\n", " bb_bh3_p = np.array(predfile[\"TARGETS\"][\"bh3\"]['bb'])\n", " \n", " # boosted Higgs detection probability\n", " dp_bh1 = np.array(predfile[\"TARGETS\"][\"bh1\"]['detection_probability'])\n", " dp_bh2 = np.array(predfile[\"TARGETS\"][\"bh2\"]['detection_probability'])\n", " dp_bh3 = np.array(predfile[\"TARGETS\"][\"bh3\"]['detection_probability'])\n", "\n", " # fatjet assignment probability\n", " ap_bh1 = np.array(predfile[\"TARGETS\"][\"bh1\"]['assignment_probability'])\n", " ap_bh2 = np.array(predfile[\"TARGETS\"][\"bh2\"]['assignment_probability'])\n", " ap_bh3 = np.array(predfile[\"TARGETS\"][\"bh3\"]['assignment_probability'])\n", " except:\n", " # pred assignment\n", " bb_bh1_p = np.array(predfile[\"TARGETS\"][\"bh1\"]['bb'])+10\n", " bb_bh2_p = np.array(predfile[\"TARGETS\"][\"bh2\"]['bb'])+10\n", " bb_bh3_p = np.array(predfile[\"TARGETS\"][\"bh3\"]['bb'])+10\n", " \n", " # boosted Higgs detection probability\n", " dp_bh1 = np.array(predfile[\"TARGETS\"][\"bh1\"]['mask']).astype('float')\n", " dp_bh2 = np.array(predfile[\"TARGETS\"][\"bh2\"]['mask']).astype('float')\n", " dp_bh3 = np.array(predfile[\"TARGETS\"][\"bh3\"]['mask']).astype('float')\n", "\n", " # fatjet assignment probability\n", " ap_bh1 = np.array(predfile[\"TARGETS\"][\"bh1\"]['mask']).astype('float')\n", " ap_bh2 = np.array(predfile[\"TARGETS\"][\"bh2\"]['mask']).astype('float')\n", " ap_bh3 = np.array(predfile[\"TARGETS\"][\"bh3\"]['mask']).astype('float')\n", " \n", " # collect fatjet pt\n", " fj_pt = np.array(testfile['INPUTS']['BoostedJets']['fj_pt'])\n", " \n", " # convert some arrays to ak array\n", " dps = np.concatenate((dp_bh1.reshape(-1, 1), dp_bh2.reshape(-1, 1), dp_bh3.reshape(-1, 1)), axis=1)\n", " dps = ak.Array(dps)\n", " aps = np.concatenate((ap_bh1.reshape(-1, 1), ap_bh2.reshape(-1, 1), ap_bh3.reshape(-1, 1)), axis=1)\n", " aps = ak.Array(aps)\n", " bb_ps = np.concatenate((bb_bh1_p.reshape(-1, 1), bb_bh2_p.reshape(-1, 1), bb_bh3_p.reshape(-1, 1)), axis=1)\n", " bb_ps = ak.Array(bb_ps)\n", " bb_ts = np.concatenate((bb_bh1_t.reshape(-1, 1), bb_bh2_t.reshape(-1, 1), bb_bh3_t.reshape(-1, 1)), axis=1)\n", " bb_ts = ak.Array(bb_ts)\n", " fj_pt = ak.Array(fj_pt)\n", " bh_masks = np.concatenate((bh1_mask.reshape(-1, 1), bh2_mask.reshape(-1, 1), bh3_mask.reshape(-1, 1)), axis=1)\n", " bh_masks = ak.Array(bh_masks)\n", " bh_pts = np.concatenate((bh1_pt.reshape(-1, 1), bh2_pt.reshape(-1, 1), bh3_pt.reshape(-1, 1)), axis=1)\n", " bh_pts = ak.Array(bh_pts)\n", " \n", " # select predictions and targets\n", " bb_ts_selected, targetH_selected_pts = sel_target_bH_by_mask(bb_ts, bh_pts, bh_masks)\n", " bb_ps_selected = sel_pred_bH_by_dp(dps, aps, bb_ps, dp_cut)\n", " \n", " # generate correct/retrieved LUT for pred/target respectively\n", " LUT_pred = gen_pred_bH_LUT(bb_ps_selected, bb_ts_selected, fj_pt)\n", " LUT_target = gen_target_bH_LUT(bb_ps_selected, bb_ts_selected, targetH_selected_pts)\n", " \n", " # reconstruct bH to remove overlapped ak4 jets\n", " fj_eta = np.array(testfile['INPUTS']['BoostedJets']['fj_eta'])\n", " fj_phi = np.array(testfile['INPUTS']['BoostedJets']['fj_phi'])\n", " fj_mass = np.array(testfile['INPUTS']['BoostedJets']['fj_mass'])\n", " \n", " fjs = ak.zip(\n", " {\n", " \"pt\": fj_pt,\n", " \"eta\": fj_eta,\n", " \"phi\": fj_phi,\n", " \"mass\": fj_mass,\n", " },\n", " with_name=\"Momentum4D\"\n", " )\n", " fj_reco = fjs[bb_ps_selected-10]\n", " \n", " return LUT_pred, LUT_target, fj_reco" ] }, { "cell_type": "code", "execution_count": 10, "id": "3cf6e0d3", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[False, False, False], [False, False, False], ..., [False, False, False]]\n" ] } ], "source": [ "LUT_pred, LUT_target, fj_reco = parse_boosted_w_target(test_h5, s_h5)" ] }, { "cell_type": "code", "execution_count": 11, "id": "9d1e162d", "metadata": {}, "outputs": [], "source": [ "testfile = s_h5\n", "j_pt = np.array(testfile['INPUTS']['Jets']['pt'])\n", "j_eta = np.array(testfile['INPUTS']['Jets']['eta'])\n", "j_phi = np.array(testfile['INPUTS']['Jets']['phi'])\n", "j_mass = np.array(testfile['INPUTS']['Jets']['mass'])\n", "js = ak.zip(\n", " {\n", " \"pt\": j_pt,\n", " \"eta\": j_eta,\n", " \"phi\": j_phi,\n", " \"mass\": j_mass,\n", " },\n", " with_name=\"Momentum4D\"\n", ")" ] }, { "cell_type": "code", "execution_count": 12, "id": "1a5e9829", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
[[0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " ...,\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]]\n",
       "--------------------------------\n",
       "type: 42829 * 10 * int64
" ], "text/plain": [ "" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "js[:, np.newaxis].deltaR(fj_reco)[2]\n", "ak.sum(js[:, np.newaxis].deltaR(fj_reco)[2]<0.4, axis=-2)>0\n", "ak.local_index(js)" ] }, { "cell_type": "code", "execution_count": 13, "id": "5b622454", "metadata": {}, "outputs": [], "source": [ "overlapped = ak.sum(js[:, np.newaxis].deltaR(fj_reco)<0.4, axis=-2)>0" ] }, { "cell_type": "code", "execution_count": 14, "id": "48625f78", "metadata": {}, "outputs": [], "source": [ "def get_unoverlapped_jet_index(fjs, js, dR_min=0.8):\n", " overlapped = ak.sum(js[:, np.newaxis].deltaR(fjs)0\n", " jet_index_passed = ak.local_index(js).mask[~overlapped]\n", " jet_index_passed = ak.drop_none(jet_index_passed)\n", " return jet_index_passed" ] }, { "cell_type": "code", "execution_count": 15, "id": "b26eae3e", "metadata": {}, "outputs": [], "source": [ "jet_index_passed = get_unoverlapped_jet_index(fj_reco, js)" ] }, { "cell_type": "code", "execution_count": 16, "id": "f326ed3a", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
[[0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [2, 4, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [1, 2, 3, 5, 6, 7, 8, 9],\n",
       " [2, 3, 5, 7, 8, 9],\n",
       " [1, 2, 3, 4, 6, 7, 8, 9],\n",
       " ...,\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n",
       " [0, 2, 4, 5, 6, 7, 8, 9],\n",
       " [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]]\n",
       "--------------------------------\n",
       "type: 42829 * var * int64
" ], "text/plain": [ "" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "jet_index_passed" ] }, { "cell_type": "code", "execution_count": 17, "id": "659fd710", "metadata": {}, "outputs": [], "source": [ "def sel_pred_h_by_dp(dps, aps, b1_ps, b2_ps, dp_cut=0.0, ap_cut=0):\n", " # parse predicted bb assignment by DP\n", " dp_filter = dps > dp_cut\n", " ap_filter = aps > ap_cut\n", " b1_ak4_filter = b1_ps<10\n", " b2_ak4_filter = b2_ps<10\n", " filter = dp_filter & ap_filter & b1_ak4_filter & b2_ak4_filter\n", " \n", " b1_ps_passed = b1_ps.mask[filter]\n", " b1_ps_passed = ak.drop_none(b1_ps_passed)\n", " \n", " b2_ps_passed = b2_ps.mask[filter]\n", " b2_ps_passed = ak.drop_none(b2_ps_passed)\n", " \n", " return b1_ps_passed, b2_ps_passed" ] }, { "cell_type": "code", "execution_count": 18, "id": "89ab96fe", "metadata": {}, "outputs": [], "source": [ "def sel_target_h_by_mask(b1_ts, b2_ts, h_pts, h_masks):\n", " b1_ts_selected = b1_ts.mask[h_masks]\n", " b1_ts_selected = ak.drop_none(b1_ts_selected)\n", " \n", " b2_ts_selected = b2_ts.mask[h_masks]\n", " b2_ts_selected = ak.drop_none(b2_ts_selected)\n", " \n", " h_selected_pts = h_pts.mask[h_masks]\n", " h_selected_pts = ak.drop_none(h_selected_pts)\n", " \n", " return b1_ts_selected, b2_ts_selected, h_selected_pts" ] }, { "cell_type": "code", "execution_count": 19, "id": "a43a82f2", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "set((1,2))==set((2,1))" ] }, { "cell_type": "code", "execution_count": 20, "id": "50b0e6d4", "metadata": {}, "outputs": [], "source": [ "# A pred look up table is in shape\n", "# [event,\n", "# pred_H, \n", "# [correct, pred_H_pt]]\n", "@nb.njit\n", "def gen_pred_h_LUT(b1_ps_passed, b2_ps_passed, b1_ts_selected, b2_ts_selected, js, goodJetIdx, builder):\n", " # for each event\n", " for b1_ps_e, b2_ps_e, b1_ts_e, b2_ts_e, jets_e, goodJetIdx_e in zip(b1_ps_passed, b2_ps_passed, b1_ts_selected, b2_ts_selected, js, goodJetIdx):\n", " # for each predicted bb assignment, check if any target H have a same bb assignment\n", " builder.begin_list()\n", " for b1_p, b2_p in zip(b1_ps_e, b2_ps_e):\n", " correct = 0\n", " predH_pt = (jets_e[b1_p]+jets_e[b2_p]).pt\n", " for b1_t, b2_t in zip(b1_ts_e, b2_ts_e):\n", " if set((b1_p, b2_p)) == set((b1_t, b2_t)) and (b1_p in goodJetIdx_e) and (b2_p in goodJetIdx_e):\n", " correct = 1\n", " builder.begin_list()\n", " builder.append(correct)\n", " builder.append(predH_pt)\n", " builder.end_list()\n", " \n", " builder.end_list()\n", " return builder" ] }, { "cell_type": "code", "execution_count": 21, "id": "d38b23eb", "metadata": {}, "outputs": [], "source": [ "# A target look up table is in shape\n", "# [event,\n", "# target_H, \n", "# target_bb_assign,\n", "# [retrieved, targetH_pt]]\n", "@nb.njit\n", "def gen_target_h_LUT(b1_ps_passed, b2_ps_passed, b1_ts_selected, b2_ts_selected, targetH_pts, builder):\n", " # for each event\n", " for b1_ps_e, b2_ps_e, b1_ts_e, b2_ts_e, tH_pts_e in zip(b1_ps_passed, b2_ps_passed, b1_ts_selected, b2_ts_selected, targetH_pts):\n", " # for each target fatjet, check if the predictions have a p fatject same with the t fatjet\n", " builder.begin_list()\n", " for b1_t, b2_t, tH_pt in zip(b1_ts_e, b2_ts_e, tH_pts_e):\n", " retrieved = 0\n", " for b1_p, b2_p in zip(b1_ps_e, b2_ps_e):\n", " if set((b1_p, b2_p)) == set((b1_t, b2_t)):\n", " retrieved = 1\n", " builder.begin_list()\n", " builder.append(retrieved)\n", " builder.append(tH_pt)\n", " builder.end_list()\n", " \n", " builder.end_list()\n", " return builder" ] }, { "cell_type": "code", "execution_count": 22, "id": "5701eb0c", "metadata": {}, "outputs": [], "source": [ "def parse_resolved_w_target(testfile, predfile, dp_cut=0.5, fjs_reco=None):\n", " # h pt\n", " h1_pt = np.array(testfile['TARGETS']['h1']['pt'])\n", " h2_pt = np.array(testfile['TARGETS']['h2']['pt'])\n", " h3_pt = np.array(testfile['TARGETS']['h3']['pt'])\n", " \n", " # mask\n", " h1_mask = np.array(testfile['TARGETS']['h1']['mask'])\n", " h2_mask = np.array(testfile['TARGETS']['h2']['mask'])\n", " h3_mask = np.array(testfile['TARGETS']['h3']['mask'])\n", " \n", " # target assignments\n", " b1_h1_t = np.array(testfile[\"TARGETS\"][\"h1\"]['b1']).astype('int')\n", " b1_h2_t = np.array(testfile[\"TARGETS\"][\"h2\"]['b1']).astype('int')\n", " b1_h3_t = np.array(testfile[\"TARGETS\"][\"h3\"]['b1']).astype('int')\n", "\n", " b2_h1_t = np.array(testfile[\"TARGETS\"][\"h1\"]['b2']).astype('int')\n", " b2_h2_t = np.array(testfile[\"TARGETS\"][\"h2\"]['b2']).astype('int')\n", " b2_h3_t = np.array(testfile[\"TARGETS\"][\"h3\"]['b2']).astype('int')\n", " \n", " # predict assignments\n", " b1_h1_p = np.array(predfile[\"TARGETS\"][\"h1\"]['b1']).astype('int')\n", " b1_h2_p = np.array(predfile[\"TARGETS\"][\"h2\"]['b1']).astype('int')\n", " b1_h3_p = np.array(predfile[\"TARGETS\"][\"h3\"]['b1']).astype('int')\n", "\n", " b2_h1_p = np.array(predfile[\"TARGETS\"][\"h1\"]['b2']).astype('int')\n", " b2_h2_p = np.array(predfile[\"TARGETS\"][\"h2\"]['b2']).astype('int')\n", " b2_h3_p = np.array(predfile[\"TARGETS\"][\"h3\"]['b2']).astype('int')\n", " \n", " # resolved Higgs detection probability\n", " dp_h1 = np.array(predfile[\"TARGETS\"][\"h1\"]['detection_probability'])\n", " dp_h2 = np.array(predfile[\"TARGETS\"][\"h2\"]['detection_probability'])\n", " dp_h3 = np.array(predfile[\"TARGETS\"][\"h3\"]['detection_probability'])\n", " \n", " # ak4 jets assignment probability\n", " ap_h1 = np.array(predfile[\"TARGETS\"][\"h1\"]['assignment_probability'])\n", " ap_h2 = np.array(predfile[\"TARGETS\"][\"h2\"]['assignment_probability'])\n", " ap_h3 = np.array(predfile[\"TARGETS\"][\"h3\"]['assignment_probability'])\n", " \n", " # reconstruct jet 4-momentum objects\n", " j_pt = np.array(testfile['INPUTS']['Jets']['pt'])\n", " j_eta = np.array(testfile['INPUTS']['Jets']['eta'])\n", " j_phi = np.array(testfile['INPUTS']['Jets']['phi'])\n", " j_mass = np.array(testfile['INPUTS']['Jets']['mass'])\n", " js = ak.zip(\n", " {\n", " \"pt\": j_pt,\n", " \"eta\": j_eta,\n", " \"phi\": j_phi,\n", " \"mass\": j_mass,\n", " },\n", " with_name=\"Momentum4D\"\n", " )\n", " \n", " # convert some numpy arrays to ak arrays\n", " dps = np.concatenate((dp_h1.reshape(-1, 1), dp_h2.reshape(-1, 1), dp_h3.reshape(-1, 1)), axis=1)\n", " dps = ak.Array(dps)\n", " aps = np.concatenate((ap_h1.reshape(-1, 1), ap_h2.reshape(-1, 1), ap_h3.reshape(-1, 1)), axis=1)\n", " aps = ak.Array(aps)\n", " \n", " b1_ps = np.concatenate((b1_h1_p.reshape(-1, 1), b1_h2_p.reshape(-1, 1), b1_h3_p.reshape(-1, 1)), axis=1)\n", " b1_ps = ak.Array(b1_ps)\n", " b1_ts = np.concatenate((b1_h1_t.reshape(-1, 1), b1_h2_t.reshape(-1, 1), b1_h3_t.reshape(-1, 1)), axis=1)\n", " b1_ts = ak.Array(b1_ts)\n", " b2_ps = np.concatenate((b2_h1_p.reshape(-1, 1), b2_h2_p.reshape(-1, 1), b2_h3_p.reshape(-1, 1)), axis=1)\n", " b2_ps = ak.Array(b2_ps)\n", " b2_ts = np.concatenate((b2_h1_t.reshape(-1, 1), b2_h2_t.reshape(-1, 1), b2_h3_t.reshape(-1, 1)), axis=1)\n", " b2_ts = ak.Array(b2_ts)\n", "\n", " h_masks = np.concatenate((h1_mask.reshape(-1, 1), h2_mask.reshape(-1, 1), h3_mask.reshape(-1, 1)), axis=1)\n", " h_masks = ak.Array(h_masks)\n", " \n", " h_pts = np.concatenate((h1_pt.reshape(-1, 1), h2_pt.reshape(-1, 1), h3_pt.reshape(-1, 1)), axis=1)\n", " h_pts = ak.Array(h_pts)\n", " \n", " # select predictions and targets\n", " b1_ts_selected, b2_ts_selected, targetH_selected_pts = sel_target_h_by_mask(b1_ts, b2_ts, h_pts, h_masks)\n", " b1_ps_selected, b2_ps_selected = sel_pred_h_by_dp(dps, aps, b1_ps, b2_ps, dp_cut=dp_cut)\n", " \n", " # find jets that are overlapped with reco boosted Higgs\n", " if fjs_reco is None:\n", " goodJetIdx = ak.local_index(js)\n", " else:\n", " goodJetIdx = get_unoverlapped_jet_index(fjs_reco, js)\n", " \n", " # generate look up tables\n", " LUT_pred = gen_pred_h_LUT(b1_ps_selected, b2_ps_selected, b1_ts_selected, b2_ts_selected, js, goodJetIdx, ak.ArrayBuilder()).snapshot()\n", " LUT_target = gen_target_h_LUT(b1_ps_selected, b2_ps_selected, b1_ts_selected, b2_ts_selected, targetH_selected_pts, ak.ArrayBuilder()).snapshot()\n", " \n", " \n", " return LUT_pred, LUT_target" ] }, { "cell_type": "code", "execution_count": null, "id": "3e798507", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": 23, "id": "21e49987", "metadata": {}, "outputs": [], "source": [ "# calculate efficiency\n", "# if bins=None, put all data in a single bin\n", "def calc_eff(LUT_pred, bins):\n", " predHs = [predH for event in LUT_pred for predH in event]\n", " predHs = np.array(predHs)\n", " \n", " predHs_inds = np.digitize(predHs[:,1], bins)\n", " \n", " correctTruth_per_bin = []\n", " for bin_i in range(1, len(bins)):\n", " correctTruth_per_bin.append(predHs[:,0][predHs_inds==bin_i])\n", " correctTruth_per_bin = ak.Array(correctTruth_per_bin)\n", " \n", " means = ak.mean(correctTruth_per_bin, axis=-1)\n", " \n", " errs = np.abs(\n", " clopper_pearson_interval(num=ak.sum(correctTruth_per_bin, axis=-1),\\\n", " denom=ak.num(correctTruth_per_bin, axis=-1)) - means\n", " )\n", " \n", " return means, errs" ] }, { "cell_type": "code", "execution_count": 24, "id": "a597865b", "metadata": {}, "outputs": [], "source": [ "# calculate purity\n", "def calc_pur(LUT_target, bins):\n", " targetHs = [targetH for event in LUT_target for targetH in event]\n", " targetHs = np.array(targetHs)\n", " \n", " targetHs_inds = np.digitize(targetHs[:,1], bins)\n", " \n", " correctTruth_per_bin = []\n", " for bin_i in range(1, len(bins)):\n", " correctTruth_per_bin.append(targetHs[:,0][targetHs_inds==bin_i])\n", " correctTruth_per_bin = ak.Array(correctTruth_per_bin)\n", " \n", " means = ak.mean(correctTruth_per_bin, axis=-1)\n", " \n", " errs = np.abs(\n", " clopper_pearson_interval(num=ak.sum(correctTruth_per_bin, axis=-1),\\\n", " denom=ak.num(correctTruth_per_bin, axis=-1)) - means\n", " )\n", " \n", " return means, errs" ] }, { "cell_type": "code", "execution_count": 25, "id": "f2944545", "metadata": {}, "outputs": [], "source": [ "bins = np.arange(0, 850, 50)\n", "bin_centers = [(bins[i]+bins[i+1])/2 for i in range(bins.size-1)]\n", "xerr=(bins[1]-bins[0])/2*np.ones(bins.shape[0]-1)" ] }, { "cell_type": "code", "execution_count": 26, "id": "816f5d16", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[Text(0.5, 0, 'Gen H pT (GeV)'),\n", " Text(0, 0.5, 'Matching purity'),\n", " Text(0.5, 1.0, 'Matching Purity vs. Gen H pT, DP cut at 0.5, AP cut at 1/(13*13)')]" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABDUAAAHWCAYAAABnmLJqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAA9hAAAPYQGoP6dpAACwoElEQVR4nOzde1xT9f8H8Ne4bSAXQRgXBcW7eENRSIu0vJZiVpaaiqBZWqiJXTRL1Eo0zTAzyfJWavrNS2EX1DDNW2JapuL9EqlcJRli4zLO7w9/OzG2wYCNbfB6Ph57yD4757P3zubOZ+/zuUgEQRBARERERERERGRlbMwdABERERERERFRTTCpQURERERERERWiUkNIiIiIiIiIrJKTGoQERERERERkVViUoOIiIiIiIiIrBKTGkRERERERERklZjUICIiIiIiIiKrxKQGEREREREREVklJjWIiIiIiIiIyCrVq6TG+vXrIZFI8Ntvv1W5bd++fdG3b1/TB1VLycnJCA4Ohkwmg0QiwZ07dwAAX375Jdq3bw97e3s0btwYQM1fk0Qiwbx584wWMxERWa+7d+9CLpdj06ZN5g7FJEaNGoVnn31WvF8f2w61FRUVhRYtWpg7DCIiquDxxx/HpEmTzB2GSSQmJiIgIABFRUXV3rdaSQ31iV8ikeDQoUNajwuCAH9/f0gkEgwdOrTawQDAwoUL8c0339RoX2vRokUL8ThWvA0ePFjc7vbt23j22Wfh6OiIlStX4ssvv0SjRo1w/vx5REVFoVWrVvjss8+wevVqM76ahmv//v0a752trS3kcjlGjBiBc+fOmTs8vaKiouDs7Kz3cYlEgpiYGKM8l77PecXb/v37q1XvvHnzNPZ3cnJCQEAAIiIisG7dOp1fhlFRURr7uLq6omvXrvjggw9q9OVZG5s3b0ZCQoLR6/3kk0+wfv36au2TlJSE7t27QyaTISAgAHFxcSgtLa1yv+vXr+t9P7ds2VLDV/Cfc+fOQSKRQCaTicncivr27avxvB4eHujZsyfWrl2LsrKyWsdQHTU59oao7jlx1apVeOaZZxAQEACJRIKoqKhKt//nn39gZ2eH//3vf2LZ8uXL4eLiglGjRollGRkZmDVrFh555BG4uLhU+v924cKFeOCBB+Dl5QWZTAZvb+9K2w7Xrl0TE/cV2w579uzBxIkT0alTJ9ja2ur9oT1r1iz06dMH7dq1g4uLCxo3bozQ0FBs2LABgiBobPvGG29g+/btOHXqVKXHxtJUbDvI5XKEh4dj586dJn/ue/fuYd68edX+rrZk2dnZmDVrFjp37gxnZ2fIZDK0bt0a0dHROj+ndaGqBFvfvn3RqVMnozxXxe9PfbfqXvSq2DaSSqXw9vZG3759sXDhQuTk5GjtU/73hfp7v23btoiJiUFWVpZRXq+hjhw5gnnz5uk979TUDz/8UO1jee7cOQwePBjOzs7w8PDAuHHjdB4/XfT91pg8eXINotekUqng5+cHiUSCH3/8Uec2utppQUFBeOutt6BQKGodQ3XU5NgborptOUPPZ+WFhITgpZdeEu8fPnwYe/bswRtvvKGx3XvvvYdhw4aJ51t9r3fnzp0YNGgQ/Pz8IJVK0axZM4wYMQJnzpzRG4O+OC9cuIAZM2agd+/e4jn8+vXrOredMWMGunfvDg8PDzg5OaFDhw6YN28e7t69q7FdVFQUiouL8emnn+qNRx+7au8BQCaTYfPmzXjooYc0yg8cOIAbN25AKpXWpFoA9xtDI0aMwPDhw2tchyH27Nlj0vqrEhwcjJkzZ2qV+/n5iX8fP34cBQUFeOedd9C/f3+xfP/+/SgrK8Py5cvRunVrsbymr+nff/+FnV2NPgoEYNq0aejZsydKSkrw559/IjExEfv378eZM2fg4+Nj7vDM6ssvv9S4/8UXX2Dv3r1a5R06dKhR/atWrYKzszOKiopw8+ZN7N69GxMmTEBCQgK+++47+Pv7a2wvlUrx+eefAwDu3LmD7du349VXX8Xx48eN8kPcUJs3b8aZM2fwyiuvGLXeTz75BJ6enlX+kFX78ccfMXz4cPTt2xcrVqzA6dOn8e677yI7OxurVq0yqI7Ro0fj8ccf1yjr1atXdUPXsnHjRvj4+OCff/7Btm3b8Pzzz+vcrlmzZoiPjwcA5OTk4IsvvsDEiRNx8eJFLFq0qNZxGKq6x95Q1T0nLl68GAUFBQgNDUVGRkaV2+/evRsSiQQDBw4EAJSUlGD58uWYMWMGbG1txe0uXLiAxYsXo02bNujcuTOOHj2qt84TJ04gODgYo0aNgouLC77++mvxeTZs2ICHHnoIpaWlOHDgAPr164djx46hqKgIUqkUZWVlOHDgAPr06QPg/v+VrVu3onv37hrnx4o++ugjuLq6YuLEiQgICEBJSQn27t2LqKgoXLhwAQsXLhS37datG3r06IEPPvgAX3zxRZXHqDxLajvcunULn376KZ566imsWrXKKD9W1D777DONxOC9e/cwf/58AKgXPVVSU1MxZMgQFBQUYNSoUZg8eTKkUimuXbuGb775BuvXr8eBAwfw8MMPmztUk5kzZ47G9+rx48fx0Ucf4c0339Q4J3fp0qVG9avbRiqVCjk5OThy5Aji4uKwbNky/O9//8Ojjz6qtc+CBQsQGBgIpVKJQ4cOYdWqVfjhhx9w5swZODk51SiO6jpy5Ajmz5+PqKgosSe0Mfzwww9YuXKlwT+ub9y4gYcffhhubm5YuHAh7t69i6VLl+L06dNITU2Fg4NDlXXo+q3Rtm3bmoSvYd++fcjIyECLFi2wadMmPPbYY3q3VbfT7t69iz179uC9997Dvn37cPjwYUgkklrHYojqHntDVbctZ+j5TC0jIwO///47FixYIJYtWbIE/fr10/j9BwBvvfUWfHx80K1bN+zevVtvnadPn4a7uzumT58OT09PZGZmYu3atQgNDcXRo0fRtWtXAMBPP/2ERx55RKMdANxvMwwaNAgAcPToUXz00UcICgpChw4d8Mcff+h93uPHjyM8PBzR0dGQyWT4/fffsWjRIvz000/45ZdfYGNzv5+FTCbD+PHjsWzZMkydOrV6nxGhGtatWycAEJ566inB09NTKCkp0Xh80qRJQkhIiNC8eXNhyJAh1ala1KhRI2H8+PE12lcd3/Hjx2u0f10x9Phs2LBB5+uZP3++AEDIyckxVYhkgJ9//lkAIHz99dca5atWrRIACIsXLzZTZJUbP3680KhRI72PAxBefvllkzz3yy+/LFTza0enuLg4vf8HNm7cKNjY2AhhYWEa5bpet0qlEnr06CEAEG7evFnruAw1ZMgQoXnz5kavt2PHjkKfPn0M3j4oKEjo2rWrxnf5nDlzBIlEIpw7d67Sfa9duyYAEJYsWVLTcPUqKysTWrRoIcTGxgpPPvmk0LdvX53b9enTR+jYsaNGWWFhodCsWTOhUaNGQnFxsdFj06e6x95Q1T0nXr9+XSgrKzN433HjxmnEvWPHDgGAcPnyZY3tFAqFcPv2bUEQBOHrr78WAAg///yzQTGpz80ABBcXF6GkpET4+++/hW7dugnPPfecMHz4cMHZ2Vnw9vYWnJ2dhVdffVXc9+bNm+L7WNn/G32vdejQoUKjRo2E0tJSjfKlS5cKjRo1EgoKCqy67ZCRkSE0atRIaNu2rVGe4+7duzrLc3JyBABCXFycUZ7HnPLy8gRfX1/Bx8dH5/dcWVmZsHnzZiE1NbXOY6vqs6jrO89Yqvv/Wh99bSNBEIQ//vhDkMvlQuPGjYVbt26J5fped2xsrABA2Lx5c61iqo4lS5YIAIRr164Ztd7qtn+mTJkiODo6Cn/99ZdYtnfvXgGA8Omnn1a5f21+i1UlMjJS6N69u7B8+XKhUaNGOr839LXTnnrqKQGAcOTIEZPEpoux2p4VVbctZ+j5TG3NmjWCo6OjcO/ePUEQBCErK0uws7MTPv/8c61t1Z/XmnxXZ2ZmCnZ2dsKLL74oCML978Bnn31W6NGjh/DHH38IzZs3F65evSoMHDhQGDBggFBQUCAIgiDcvn1bUCgUgiDU7P/N0qVLBQDC0aNHNcp/++03AYCQkpJicF2CIAg1mlNj9OjRuH37Nvbu3SuWFRcXY9u2bXjuued07rN06VL07t0bTZo0gaOjI0JCQrBt2zaNbSQSCQoLC7Fhwwaxu1L5q143b97ExIkTxS4zgYGBmDJlCoqLizXqKSoqQmxsLLy8vNCoUSM8+eSTWt21Ko6LVXeX+9///of33nsPzZo1g0wmQ79+/XD58mWt17Ny5Uq0bNkSjo6OCA0NxcGDB4061rZv374YP348AKBnz57isWjRogXi4uIAAF5eXhpdjHQ9v1KpxLx589C2bVvIZDL4+vriqaeewpUrV8RtdHVTunnzJiZMmABvb29IpVJ07NgRa9eu1dimusfs2LFjePzxx+Hu7o5GjRqhS5cuWL58OQBg3bp1kEgk+P3337X2W7hwIWxtbXHz5k2dx2rbtm2QSCQ4cOCA1mOffvopJBKJ2K0qMzMT0dHRaNasGaRSKXx9ffHEE0/o7S5VE+Hh4QCgcYwBw44pYNh7VlhYiJkzZ8Lf3x9SqRTt2rXD0qVLtbpZm5J6mMqmTZvQrl07yGQyhISE4JdffqmzGHQZM2YMnn/+eRw7dkzjO0oXGxsb8f9MVZ+BO3fuYMaMGWjRooXYZS8yMhK5ubkA/us+W7Ee9f8Tdbftvn374vvvv8dff/0lfs9V1QVx3bp1ePTRRyGXyyGVShEUFKTVk6JFixY4e/YsDhw4INZb2fdRWloa0tLS8MILL2j01HrppZcgCILW93NlCgsLtb6Ha+Pw4cO4fv06Ro0ahVGjRuGXX37BjRs3DNrXyckJDzzwAAoLC6vsplvVOUXdfbaiiu91dY89YJxzoi7Nmzc3+MpGWVkZkpOTMWTIELHsm2++QYsWLdCqVSuNbV1cXODh4WFQvZW5e/cu9u7di2bNmuG3335D3759kZSUBKVSifz8fHTv3h1LliwRt9+8eTP69OmDJk2a4Mcff0RGRka1jlOLFi1QWFiICRMmaLzPx44dQ2FhocZ3hDW2HXx8fNChQwdcu3ZNI56Kw0TUw8XKD5FSD0W8cuUKHn/8cbi4uGDMmDHiY+rvpevXr8PLywsAMH/+fI1hCdZ47k5MTERGRgYSEhLQvn17rcclEglGjx6Nnj17apSbol1kKuphKidOnEDv3r3h6OiIwMBAJCYm1lkMunTt2hUJCQm4c+cOPv744yq3V/fmUH++9VH3Xu7cuTNkMhm8vLwwePBgcRiPrs+/Wvk28Lx58/Daa68BAAIDA8XPemWfs4MHD4pD/qRSKfz9/TFjxgz8+++/4jZRUVFYuXKl+HzqW2W2b9+OoUOHIiAgQCzr378/2rZtqzFcsCrFxcUoLCw0ePuq/Pvvv9i5c6c4N9G///6Lb7/91uD9DX1Pq2oLG/pdV5Nj/+2332LIkCHiOaNVq1Z45513oFKpxG1q0pbz8/ODvb19pduU9/333+ORRx6Bo6OjeL+0tFSj975abeZAksvlcHJyEodcSSQSbN26FcuWLcPkyZORkZGBp556Ci+//DL27NkjDmH38PCAi4tLjZ9XHXPFoV4hISHw8PCo1ucKqOHwkxYtWqBXr1746quvxC5HP/74I/Lz8zFq1Ch89NFHWvssX74cw4YNw5gxY1BcXIwtW7bgmWeewXfffSc2qL788ks8//zzCA0NxQsvvAAAYsPq1q1bCA0NxZ07d/DCCy+gffv2uHnzJrZt24Z79+5pdMOaOnUq3N3dERcXh+vXryMhIQExMTHYunVrla9t0aJFsLGxwauvvor8/Hy8//77GDNmDI4dOyZus2rVKsTExCA8PBwzZszA9evXMXz4cLi7u6NZs2YGHcOSkhLxx1B5jRo1gqOjI+bMmYN27dph9erVYne8Vq1aYfjw4fjiiy+wc+dOsUuXvq6BKpUKQ4cORUpKCkaNGoXp06ejoKAAe/fuxZkzZ7QarWpZWVl44IEHxB+tXl5e+PHHHzFx4kQoFAqtblaGHLO9e/di6NCh8PX1xfTp0+Hj44Nz587hu+++w/Tp0zFixAi8/PLL2LRpE7p166ZR/6ZNm9C3b180bdpUZ7xDhgyBs7Mz/ve//4ldltW2bt2Kjh07imNQn376aZw9exZTp05FixYtkJ2djb179yI9Pd1ok6KpT37u7u5imaHH1JD3TBAEDBs2DD///DMmTpyI4OBg7N69G6+99hpu3ryJDz/80KA4dX3+quvAgQPYunUrpk2bBqlUik8++QSDBw9Gamqq0cb91sS4ceOwevVq7NmzBwMGDKh0W/UJskmTJnq3uXv3LsLDw3Hu3DlMmDAB3bt3R25uLpKSknDjxg14enoaHNucOXOQn5+PGzduiO9VZXOcAPe/czp27Ihhw4bBzs4Ou3btwksvvYSysjK8/PLLAICEhARMnToVzs7OmDNnDgDA29tbb53qHyE9evTQKPfz80OzZs10/kjRZf78+XjttdcgkUgQEhKC9957TxzKUFObNm1Cq1at0LNnT3Tq1AlOTk746quvxMZmVa5evQpbW9tKuw5X55xSleoee6D250RjOH78OHJycjSGDx05cgTdu3evdd2CIOD27dsoLS3FxYsXAdxPInbt2lVsO6h/tJaVlek93uWP08qVK3HlypVKj1NUVBSKiorg7e2NDRs2YM2aNXBwcMC2bdu03meZTIbDhw+L31XW0HaoqKSkBH///Xel31+VKS0txaBBg/DQQw9h6dKlOrv4e3l5YdWqVZgyZQqefPJJPPXUUwDuD0sIDAy0unP3rl274OjoKL4OQ5iiXVSZ/Px8nefokpISg2P+559/8Pjjj+PZZ5/F6NGj8b///Q9TpkyBg4MDJkyYYHA9xjZixAhMnDhRHIpQGUPOzwAwceJErF+/Ho899hief/55lJaW4uDBg/j111+1znGVeeqpp3Dx4kV89dVX+PDDD8Vzuzqpp8vXX3+Ne/fuYcqUKWjSpAlSU1OxYsUK3LhxA19//TUA4MUXX8StW7d0Dr/V5ebNm8jOztYZe2hoKH744QeDXs++ffvg5OQElUqF5s2bY8aMGZg+fbpB++qTlJSEu3fvYtSoUfDx8UHfvn2xadMmvRe1KzLkPa3p7xddqnvsgfsXLpydnREbGwtnZ2fs27cPc+fOhUKhEJPuNWnLVUdJSQl++uknjaGTR44cQZMmTdC8efNa13/nzh2UlJQgMzMTCQkJUCgU6Nevn8Y2NjY2Ggmg2g4XKi0txZ07d1BcXIwzZ87grbfegouLC0JDQ7W27d69Ow4fPly9J6hOt47y3cM+/vhjwcXFRewS88wzzwiPPPKIIAi6uzypt1MrLi4WOnXqJDz66KMa5fq6kEZGRgo2NjY6u+Spu9qq4+vfv79YJgiCMGPGDMHW1la4c+eOWNanTx+NLrfq7nIdOnQQioqKxPLly5cLAITTp08LgiAIRUVFQpMmTYSePXtqdNlev369AMCg7sfNmzcXu+JWvMXHx4vb6euOp69LV8XXtHbtWgGAsGzZMr3HTBAErW5KEydOFHx9fYXc3FyNfUaNGiW4ubmJ76Whx6y0tFQIDAwUmjdvLvzzzz964xg9erTg5+cnqFQqsezkyZMCAGHdunVar6G80aNHC3K5XKOLcUZGhmBjYyMsWLBAEARB+Oeff4zaXV79+teuXSvk5OQIt27dEpKTk4XWrVsLEolEo+uqocfUkPfsm2++EQAI7777rsbjI0aMECQSiVa38YrGjx+v9/Onvhky/ES97W+//SaW/fXXX4JMJhOefPJJnfvUxfATQfjvvS4fh3r4SU5OjpCTkyNcvnxZWLhwoSCRSIQuXbpU+nxz584VAAg7duzQeqzi90/Frnfqz0n5Lr3V7bJY8ftTEARh0KBBQsuWLTXKqjMEQt1VMD09Xeuxnj17Cg888ECl+//111/CwIEDhVWrVglJSUlCQkKCEBAQINjY2AjfffedQTHoUlxcLDRp0kSYM2eOWPbcc88JXbt21dq2T58+Qvv27cX39Ny5c8K0adMEAEJERESlz2PIOUX9OatI13td3eEntT0nGqKqfd9++22Nz2FJSYkgkUiEmTNnVlqvId3UMzIytL5XFi5cKLYdLl26JPTo0UPw9/cXQkJChD59+gje3t5Co0aNhNdee02sp/xxGjJkiBAQEFDpcYqPj9d4Th8fH73vc5s2bYTHHnvMqtoOAwcOFD/vp06dEkaNGiUAEKZOnaoRT8X3Rj1crPy5VH0umDVrltZzjR8/XuOzUVmXZms7d7u7uwvBwcFa5QqFQjy2OTk5Gl3qjd0u0qf8UC19N0OGn/Tp00cAIHzwwQdiWVFRkRAcHCzI5XKdQ/PqYviJWteuXQV3d3fxvvp1//TTT0JOTo7w999/C1u2bBGaNGkiODo6Cjdu3NBb1759+wQAwrRp07QeU/9/1vX5V6v4ua5uN3pd5+f4+HhBIpFoDB2pTvvn+PHjAgDhiy++0HrstddeEwAISqWy0joiIiKExYsXC998842wZs0aITw8XAAgvP766wbFoM/QoUOFBx98ULy/evVqwc7OTsjOztbYTn3+vHDhgpCTkyNcu3ZN+PTTTwWpVCp4e3sLhYWFep/DkLZwdb7rqtv21PWevvjii4KTk5PGca/NUOKq9k1JSdH6HD700ENCSEhIpfUaOvykXbt24neKs7Oz8NZbb4nf4WVlZcJzzz2nNfxkwIABGsNPyjPk/83Ro0c1vsvatWun9/vmhRdeEBwdHSt9DRXVeElXdZej7777DgUFBfjuu+8qzdKpu84A97PH+fn5CA8Px8mTJ6t8rrKyMnzzzTeIiIjQmbWsmDl64YUXNMrCw8OhUqnw119/Vflc0dHRGleM1EMJrl69CgD47bffcPv2bUyaNEmjy/aYMWM0rsxXJSwsDHv37tW6jR492uA6qrJ9+3Z4enpi6tSpWo/py7YJgoDt27cjIiICgiAgNzdXvA0aNAj5+fla71lVx+z333/HtWvX8Morr2hdOS0fR2RkJG7duoWff/5ZLNu0aRMcHR3x9NNPV/paR44ciezsbI1uaNu2bUNZWRlGjhwJ4P5n0MHBAfv378c///xTaX3VMWHCBHh5ecHPzw+DBw9Gfn4+vvzyS7HranWOqSHv2Q8//ABbW1tMmzZN4/GZM2dCEAS9M1GXJ5PJdH7+qhquUVGvXr0QEhIi3g8ICMATTzyB3bt3a3TTq2vqbHlBQYFGeWFhIby8vODl5YXWrVvjzTffRK9evapcPWD79u3o2rUrnnzySa3H6mKiq/Lfn+oreH369MHVq1eRn59fozrVXWN1Tewsk8k0us7qEhAQgN27d2Py5MmIiIjA9OnT8fvvv8PLy0vnJMiG+vHHH3H79m2N78LRo0fj1KlTOHv2rNb258+fF9/TDh06YMWKFRgyZIjOoV1q1T2nmEJtzonG8sMPP2gMPcnLy4MgCNU6l+nj4eGBvXv3YteuXeL/m3///VdsO/z22294++23kZubixEjRgC4/7l7+OGHNWIqf5xKSkpQVlZW6XEaPXo09u7di82bN2P06NHIycnBI488ovN99vDw0Lgabg1thz179oif965du+Lrr7/GuHHjsHjxYoPrqGjKlCk13hewvnO3QqHQeUV13Lhx4rH18vISVxcwRbuoKitXrtR5fq7OpJ12dnZ48cUXxfsODg548cUXkZ2djRMnThhcjyk4OztrnZ+B+8MrvLy84O/vj1GjRsHZ2Rk7d+7U29sHuH9+lkgk4tDs8ur6u7ywsBC5ubno3bs3BEEwuMdjRVWdn8tvo09SUhJef/11PPHEE5gwYQIOHDiAQYMGYdmyZQYP56zo9u3b2L17t8b5+emnnxaHXOnSrl07eHl5ITAwEC+++CJat26N77//vtKJX2vy+8WYyr+nBQUFyM3NRXh4OO7du4fz58+b/PmB++fnoKAgjZ5ot2/fNsr5Gbg/rDk5ORmffPIJOnTogH///Vdst6uHcP7666/ixKGBgYHYs2eP2HulJoKCgrB371588803eP3119GoUSOt1U/U3N3d8e+//+LevXsG11/jJS+8vLzQv39/bN68Gffu3YNKpRIbJrp89913ePfdd/HHH39oLJ9oyIczJycHCoXC4O7s5cefAf8NAzDkZFjVvurGTcVZZ+3s7KrVBdLT01PnmChjunLlCtq1a1etlU1ycnJw584drF69Wu9SsdnZ2Rr3qzpm6q5mVb1/AwYMgK+vLzZt2oR+/fqhrKwMX331FZ544okqx2wNHjwYbm5u2Lp1q9h9auvWrQgODhZnepZKpVi8eDFmzpwJb29vPPDAAxg6dCgiIyNrtUrJ3LlzER4ejrt372Lnzp3YsmWLOIsvUL1jash79tdff8HPz0/rmKhnKzekAW5ra2uUz1+bNm20ytq2bYt79+4hJyfHbKu/qL8kKx4jmUyGXbt2AYA4tt6Qbt9XrlypsnFuSocPH0ZcXByOHj2q9QWfn58PNze3atepPmnrWs5WqVRqnNQN5eHhgejoaCxatAg3btyoUZf6jRs3IjAwEFKpVByD3qpVKzg5OWHTpk0aXTGB+8MhP/vsM3EZwDZt2kAul1f6HNU9p5hCbc6JxpCZmYmTJ09qzKquJhhhbh4HBwfxOyY3Nxc7d+7EO++8g9DQUPTv3x9bt27F8OHDoVKp8PjjjyM5ORnA/e6u5Yci6DpOq1at0nucmjdvLnbNffTRR/HVV1/h+PHj+Pfff7U+04IgaNRjDW2HsLAwvPvuu+ISiR06dKjVCg12dnY1HvqiZm3nbhcXF50N6QULFojLmZcftmiKdlFVQkNDdSbi3N3dDR466ufnh0aNGmmUqY/p9evX8cADDxhUjyncvXtX52dj5cqVaNu2Lezs7ODt7Y127dpptKd0uXLlCvz8/Iwy509NpKenY+7cuUhKStJ6f2t60aGq83P5bQwlkUgwY8YM7N69G/v378fYsWOrHdfWrVtRUlKCbt26acwRExYWhk2bNonDYcvbvn07XF1dYW9vj2bNmhk0dKQmv1+M6ezZs3jrrbewb98+reVna/qeVtf333+PiIgIrXJjnJ8BzVXqRo0aJf6GWLp0KQDoHbo9ePDgGj+nq6ur2C544oknsHnzZjzxxBM4efKkmDxRU7/O6rSJavVpee655zBp0iRkZmbiscce03tiPXjwIIYNG4aHH34Yn3zyCXx9fWFvb49169Zh8+bNtQlBp4rLz6gZ8kGozb71gXoJt7Fjx4oTlVZU8UqBsY6Zra0tnnvuOXz22Wf45JNPcPjwYdy6dcugL16pVIrhw4dj586d+OSTT5CVlYXDhw9r/QB65ZVXEBERgW+++Qa7d+/G22+/jfj4eOzbt09rPLChOnfuLP4nHT58OO7du4dJkybhoYcegr+/f42OKdWOenK5ij8gjJXM0UXfF29te6xcuXIF/fr1Q/v27bFs2TL4+/vDwcEBP/zwAz788EONZRerw9fXF8D9JcMqLn2bkZGhc4yjIdR15eXlVfvHkkKhwK5du6BUKnUmzDZv3oz33ntP41g3atTI6t7Tuj4n6vLjjz9CJpPhkUceEcs8PDwgkUiM2pOtvCZNmohjr8u3Hbp06YL9+/dr/biveJwSExPx119/iRdUDKVQKPDLL7+Iy9Cp/fPPPxqfM2toO1R1QaS6n1mpVFrlj8aqWNu5u3379jh16hRKSko0Ju3Tdx42Z7uoPiopKcHFixd1JpX1JXNqy1Tf5SqVCgMGDEBeXh7eeOMNtG/fHo0aNcLNmzcRFRVllPNzRRkZGfDw8NDZi6Mq5c/PNbFp0yYAwIMPPqjz8atXr6Jly5YaZQ8//HC15h0zlKne0zt37qBPnz5wdXXFggUL0KpVK8hkMpw8eRJvvPFGjd/T6rh27RrOnz+vNSF8kyZNTHJ+dnd3x6OPPopNmzaJSY3yjLmYQnlPPfUUxo0bhy1btmglNf755x84OTlVK3lXq6TGk08+iRdffBG//vprpRNpbd++HTKZDLt379b4T7hu3TqtbXV9SL28vODq6ir+UDEn9RWgy5cvazQGS0tLcf36dYv6cdqqVSscO3ZM68RdGS8vL7i4uEClUhnth4I6K3vmzJkq64yMjMQHH3yAXbt24ccff4SXl5dWQ1SfkSNHYsOGDUhJScG5c+cgCILYfbViPDNnzsTMmTNx6dIlBAcH44MPPsDGjRur/+J0WLRoEXbu3In33nsPiYmJ1TqmhrxnzZs3x08//YSCggKNKx3qLnHGmEDIUJcuXdIqu3jxIpycnCqdWMvU1JNBGfrZqUqrVq2q/P5RX4mrOIuzrp4z1ck879q1C0VFRUhKStK4+le+q3dN6g0ODgZwv1t8+QTGrVu3cOPGDXFiyupSd6+uyfu/Y8cOKJVKrFq1SqsRdOHCBbz11ls4fPgwHnrooRrFpmboOaX8e1o+aV/b97S250RjqDirOnD/qn2rVq2qnJW+poqLi5Gfn1/jtsOOHTu0roir6Ws7ODk54d69e1pX10pLS/H3339j2LBhtXtRBqqrtkN1voeqo6rPoTWdu4cOHYpff/0VO3fuxLPPPlvl9qZoF9WFW7duobCwUKO3hnrSXmNNjF4T27Ztw7///mvU8/Pu3buRl5ent7eGqc7Pp0+fxsWLF7FhwwZERkaK5bqG8lan3qZNm8LLy0tcvaW81NRU8fxdXbU5P1+7dg1HjhxBTEyM1sS+ZWVlGDduHDZv3oy33nqrRrGVZ0hb2FTv6f79+3H79m3s2LEDDz/8sFiu67xoyvOzm5ubVlunffv22L59u0me899//62zXihqRUVFKCsr0/m8165dE3uPGKpW6XlnZ2esWrUK8+bN09lFRs3W1hYSiUQje3b9+nV88803Wts2atRI6wNqY2OD4cOHY9euXTr/g9dl5rtHjx5o0qQJPvvsM5SWlorlmzZtMtnVrZp6+umnkZubq3PZLH3HzNbWFk8//TS2b9+us8Ff1RKJunTv3h2BgYHiMl6VxdGlSxd06dIFn3/+ObZv345Ro0YZ3P2sf//+8PDwwNatW7F161aEhoYiMDBQfPzevXtitz21Vq1awcXFRaOLX0ZGBs6fP1+tWcYr1vn0009j/fr1yMzMrNYxNeQ9e/zxx6FSqbS2+fDDDyGRSMQVierC0aNHNcYS//333/j2228xcOBAvVeqTG3z5s34/PPP0atXL62ZnGvq6aefxqlTp3TOvaF+X9TJu/JL2qpUKp3dlRs1amTwyUN9HMv/X8nPz9f5A1jX96c+HTt2RPv27bF69WqN72Z11/7ywwnz8/Nx/vx5jZh1fRfcvHkTa9euRZcuXcQrTdWxceNGtGzZEpMnT8aIESM0bq+++iqcnZ3FK0W1Yeg5Rdd7ql46tKLqHPvanhNrq6SkBHv37tWYu0KtV69eOo+JoQoLC/WOgS0oKECPHj1q1XYoLS3VeZwcHR11th3UP3IqNozT0tKgVCrRu3dvA19Z7dRV26F58+awtbXVWlr7k08+qVW96vHv+j6L1nTunjJlCry9vTFjxgzxR355FdslpmgX1YXS0lJ8+umn4v3i4mJ8+umn8PLy0pgLqy6dOnUKr7zyCtzd3XUOVaiJp59+GoIgYP78+VqPqd9LV1dXeHp6GvT/Qp0EMuR7V9f5WRAELF++vFb1Avdf13fffYe///5bLEtJScHFixfxzDPPiGUlJSU4f/68Rq+OvLw8rR4LJSUlWLRoERwcHDQSq4ZSn3tff/11rfPzs88+iz59+hjl/AwY1hauznddbd/T4uJivfWaIhHwww8/YODAgVrfob169cI///xj8Nw8uui6MHD9+nWkpKSYpJcU8N9KKxV9/vnnALRX4QOAkydPVvv8XOvBSvq64pU3ZMgQLFu2DIMHD8Zzzz2H7OxsrFy5Eq1bt8aff/6psW1ISAh++uknLFu2DH5+fggMDERYWBgWLlyIPXv2oE+fPnjhhRfQoUMHZGRk4Ouvv8ahQ4dqNaa0OhwcHDBv3jxMnToVjz76KJ599llcv34d69evR6tWrQzO2t28eVPn1QVnZ2cMHz7cKLFGRkbiiy++QGxsLFJTUxEeHo7CwkL89NNPeOmll/DEE0/o3G/RokX4+eefERYWhkmTJiEoKAh5eXk4efIkfvrpp2p3W7OxscGqVasQERGB4OBgREdHw9fXF+fPn8fZs2exe/durbhfffVVAKjWmD97e3s89dRT2LJlCwoLC7W6UF28eBH9+vXDs88+i6CgINjZ2WHnzp3IysrCqFGjxO1mz56NDRs24Nq1azW+mvHaa6/hf//7HxISErBo0SKDj6kh71lERAQeeeQRzJkzB9evX0fXrl2xZ88efPvtt3jllVeMuuRjVTp16oRBgwZpLOkKQGfjwhB9+/bFgQMHDE5Ubtu2Dc7OziguLsbNmzexe/duHD58WJxAz1hee+01bNu2Dc888wwmTJiAkJAQ5OXlISkpCYmJiejatSs6duyIBx54ALNnzxavGG3ZskXjB4xaSEgItm7ditjYWPTs2RPOzs56f9wNHDgQDg4OiIiIwIsvvoi7d+/is88+g1wu1+qaGhISglWrVuHdd99F69atIZfLxTXhdVmyZAmGDRuGgQMHYtSoUThz5gw+/vhjPP/88xoZ8p07dyI6Ohrr1q1DVFQUgPsNG/XQGD8/P1y/fh2ffvopCgsLtRp069ev19q/IvVEgxUnwFWTSqUYNGgQvv76a3z00UfVWutdF0POKQMHDkRAQAAmTpyI1157Dba2tli7di28vLyQnp6uUV91jr0xzon67Nq1C6dOnQJwvxH7559/4t133wUADBs2DF26dMGhQ4egUCh0JjWeeOIJfPnll7h48aI4/l5NXY96wtYvv/wShw4dAgDx6tylS5fQv39/jBw5Eu3bt4eNjQ2++uorABCX9Aaq13YICgpCly5d8OuvvyIvLw9yuRz37t3Drl27xP83jRo1wq5duzBgwAC0aNECDg4OOH78OG7cuAEnJyeMHTtW431etWoVHB0dMWDAAJNd+SrPWG2Hqri5ueGZZ57BihUrIJFI0KpVK3z33Xd6e7gYytHREUFBQdi6dSvatm0LDw8PdOrUSWMIgbWcuz08PLBz505ERESga9euGDVqFHr27Al7e3v8/fff4rmjfM84Y7eL6oKfnx8WL16M69evo23btti6dSv++OMPrF69ukbfn/PmzcP8+fPx888/o2/fvlVuf/DgQSiVSqhUKty+fRuHDx9GUlIS3NzcsHPnTqPNufXII49g3Lhx+Oijj3Dp0iUMHjwYZWVlOHjwIB555BFxnpTnn38eixYtwvPPP48ePXrgl19+0ZnUUid85syZg1GjRsHe3h4RERFa85MA96+ct2rVCq+++ipu3rwJV1dXbN++XWeiUl3vtGnTMGjQINja2mp8fit688038fXXX+ORRx7B9OnTcffuXSxZsgSdO3dGdHS0uN3NmzfRoUMHjB8/HuvXrwdwf5LQd999FyNGjEBgYCDy8vKwefNmnDlzBgsXLtQ49tevX0dgYKDG/rps2rQJwcHBWsNV1YYNG4apU6fi5MmTtV4W3JC2cHW+66pz7Hv37g13d3eMHz8e06ZNg0QiwZdffqmzXVqdthwA/Pnnn0hKSgJwv9defn6+eF7t2rUrIiIi8O+//+Lnn39GYmKi1v5DhgyBnZ0dfvrpJ63etF9++SX++usv8aLCL7/8ItY9btw4sbdg586d0a9fPwQHB8Pd3R2XLl3CmjVrxKRXdeTn52PFihUAIC6/+vHHH6Nx48Zo3Lix+H9v//79mDZtGkaMGIE2bdqguLgYBw8exI4dO9CjRw+t88WJEyeQl5en93eqXtVZKkXfEqMV6VrSdc2aNUKbNm0EqVQqtG/fXli3bp3O5fLOnz8vPPzww4Kjo6MAQGM5ur/++kuIjIwUvLy8BKlUKrRs2VJ4+eWXxWWz9MWna9kffcuyVVyCSt8yUB999JHQvHlzQSqVCqGhocLhw4eFkJAQYfDgwZUeG/XxgZ6lusov71PbJV0F4f6yRHPmzBECAwMFe3t7wcfHRxgxYoRw5coVcRvoWPonKytLePnllwV/f39xv379+gmrV6+u8TE7dOiQMGDAAMHFxUVo1KiR0KVLF2HFihVaxycjI0OwtbUV2rZtq+8Q6rV3714BgCCRSIS///5b47Hc3Fzh5ZdfFtq3by80atRIcHNzE8LCwoT//e9/Gtupl7mrajmvqpYt69u3r+Dq6iouB2jIMRUEw96zgoICYcaMGYKfn59gb28vtGnTRliyZInGcoT6qJc21QfVWNL15ZdfFjZu3Cj+3+7WrVuly8FVtaxWSEiI4OPjU+Vzq/8PqG8ymUxo1qyZMHToUGHt2rU6lzqr6nVX5fbt20JMTIzQtGlTwcHBQWjWrJkwfvx4jSX+rly5IvTv319csuzNN98UP5Plj8vdu3eF5557TmjcuLHW/3tdkpKShC5duggymUxo0aKFsHjxYnHJs/Kf08zMTGHIkCGCi4uLAAOXidy5c6cQHBwsSKVSoVmzZsJbb72ltdyf+ruo/P/pzZs3Cw8//LDg5eUl2NnZCZ6ensKTTz4pnDhxQus5VqxYIQAQkpOT9cbxwQcfCACElJQUvduol7/89ttvBUG4/51nyPKG+lR1ThEEQThx4oQQFhYmODg4CAEBAcKyZct0Lula3WNvjHOiLpUt2ax+/1599VUhKChI5/5FRUWCp6en8M4772g9pq/e8jHn5OQIL7zwgvg96+DgIHh7ewsAhL1791Yau662Q3R0tN7nLH8s1qxZI3h4eAgSiUQAINjZ2QkPPvigsG7dOuH69eta77O3t7cwevRoQRCsq+1Q8fjokpOTIzz99NOCk5OT4O7uLrz44ovCmTNntOKp7Dux4pKugiAIR44cEUJCQgQHBwed7QZrOXeXj/e1114TgoKCBEdHR/GzERkZKfzyyy9a25uiXVRRVe1sQ7/z1Nv99ttvQq9evQSZTCY0b95c+Pjjj/XuU9WSrjNnzhQkEolw7ty5Sp9bfQzUN3t7e8HLy0t4+OGHhffee09r6U9BMPz3hT6lpaXCkiVLhPbt2wsODg6Cl5eX8Nhjj2mcj+7duydMnDhRcHNzE1xcXIRnn31WyM7O1vlZfuedd4SmTZsKNjY2VX6m0tLShP79+wvOzs6Cp6enMGnSJOHUqVNa73dpaakwdepUwcvLS/yeqsqZM2eEgQMHCk5OTkLjxo2FMWPGCJmZmRrbqD9b5b8Pf/vtNyEiIkJsrzg7OwsPPfSQ1v8ZQRCE06dPC4DupZ3VTpw4IQAQ3n77bb3bXL9+XQAgzJgxQxAE/b9VDGVIW9jQ77rqHvvDhw8LDzzwgODo6Cj4+fkJr7/+urB79+5at+UqW7JZ/f599913gkQiEbKysnTWMWzYMKFfv35a5eplnHXdysccFxcn9OjRQ3B3dxfs7OwEPz8/YdSoUcKff/5Zaey6qD97Vf2evXz5shAZGSm0bNlScHR0FGQymdCxY0chLi5OY+lstTfeeEMICAgw6DdNedVKapB+KpVK8PDwEJ5//nlzh2L1cnJyBDs7O3GNerJMhiY/DKVQKAQ7O7tKG11knZ555hmhZ8+e5g6D/l+HDh2E1157Te/jCxYsEAIDA4XS0tI6jKru/P7774JEIhF+//13c4dSr9oOPHdbjtomfHXp2bOnMGLECKPWSea3cuVKoVGjRlrJEjKPKVOmVNpe+uWXXwQbGxvh4sWLdRhV3VEqlYKPj4+QkJBQ7X1rN+V1A6VUKrW6IX3xxRfIy8szqEseVW79+vVQqVQYN26cuUOhOvTLL7+gadOmmDRpkrlDISMSBAH79+8Xu0GSeRUXF2PkyJEa3ZcrmjFjBu7evYstW7bUYWR1Z9GiRRgxYkSNJ9urqfreduC5u/5SKBQ4deqUziWgybqph316e3ubOxTC/UncKxvCHR4ejoEDB+L999+vw6jqzrp162Bvb4/JkydXe1+JUPEMS1Xav38/ZsyYgWeeeQZNmjTByZMnsWbNGnTo0AEnTpyAg4ODuUO0Svv27UNaWhrefvttPPLII9ixY4e5Q6JKSCQSvPzyyzonciIiIk31te3Ac7fl6du3L3Jzcy1i1UAiorpQ64lCG6IWLVrA398fH330kTghYGRkpDirMNXMggULcOTIETz44IPixDNERET1QX1tO/DcTURE5saeGkRERERERERklTinBhERERERERFZJSY1iIiIiIiIiMgqcU4NIykrK8OtW7fg4uICiURi7nCIiIggCAIKCgrg5+cHGxtexzAGnu+JiMiS8FzPpIbR3Lp1C/7+/uYOg4iISMvff/+NZs2amTuMeoHneyIiskQN+VzPpIaRuLi4ALj/YXJ1dTVzNERERIBCoYC/v794jqLa4/meiIgsCc/1TGoYjboLqqurKxs5RERkUThMwnh4viciIkvUkM/1DXPQDRERERERERFZPSY1iIiIiIiIiMgqMalBRERERERERFaJSQ0iIiIiIiIiskpMahAREZHVWrlyJVq0aAGZTIawsDCkpqbq3bZv376QSCRatyFDhtRhxERERGRMTGoQERGRVdq6dStiY2MRFxeHkydPomvXrhg0aBCys7N1br9jxw5kZGSItzNnzsDW1hbPPPNMHUdORERExsKkBhEREVmlZcuWYdKkSYiOjkZQUBASExPh5OSEtWvX6tzew8MDPj4+4m3v3r1wcnJiUoOIiMiKMalBREREVqe4uBgnTpxA//79xTIbGxv0798fR48eNaiONWvWYNSoUWjUqJHebYqKiqBQKDRuREREZDmY1CAiIiKrk5ubC5VKBW9vb41yb29vZGZmVrl/amoqzpw5g+eff77S7eLj4+Hm5ibe/P39axU3ERERGReTGkRERNTgrFmzBp07d0ZoaGil282ePRv5+fni7e+//66jCImIiMgQduYOgIiIiKi6PD09YWtri6ysLI3yrKws+Pj4VLpvYWEhtmzZggULFlT5PFKpFFKptFaxEhERkemwpwYRERFZHQcHB4SEhCAlJUUsKysrQ0pKCnr16lXpvl9//TWKioowduxYU4dJREREJsaeGkRERGSVYmNjMX78ePTo0QOhoaFISEhAYWEhoqOjAQCRkZFo2rQp4uPjNfZbs2YNhg8fjiZNmpgjbCIiIjIiJjWIyGpkK5TILiiq9n5yFynkrjITRERE5jRy5Ejk5ORg7ty5yMzMRHBwMJKTk8XJQ9PT02Fjo9kp9cKFCzh06BD27NljjpCJqAo81xNRdUkEQRDMHUR9oFAo4Obmhvz8fLi6upo7HKJ66cO9F7E85VK195verw1mDGhrgoiILBvPTcbHY0pkWjzXE1UPz0vsqUFEVmRMWAAGBGku36gsUWFE4lEAwLbJvSCzt9XaT+7CSf6IiIgsSbZCiU3H0jEmLECjh4WxzvX66iei+odJDSIyGlN3GZW7yrS2u1dcKv4d5OcKJwd+rREREVm67IIiLE+5hAFB3hrndmOd6/XVT0T1D1v/RGQ0m46ls8soERERGZWq7L/R8qnX8hDexgu2NhIzRkREloRJDSKqNlN3GU27lY/5u9IQFxGEID+3OoufiIiILEvymQzEJZ0V70etOw5fNxniIoIwuJOvGSMjIkvBpAYRVZupu4xezLqLY9fycDHrrmmSGuySSkREZPGSz2RgysaTqLiqQWa+ElM2nsSqsd2Z2CAi2FS9CRFRzVXsMlr+vjXUT0RERKajLFHhXnGp1q1AWYK4pLNaCQ0AYtm8pDQUKEt07q8sUdXlyyAiM2JPDSIyGVN3GWWXVCIiIuumHqJaXQKATIUSneftMW5ARGR12FODiExC3WU0S6G5Goq6y2jymQyLrp+IiIiIiCwfe2oQUY2pu4xWpCoTKu0yKsH9LqMPtvbUOXt5calK/NcU9bNLKhERkWXYNrkXgvxctcpTr+Uhat3xKvdfH90ToYEeWuVptxQ17gVCRNaFSQ0iqjFTdxl9fftpvL79tMnqJyIiIvOS2dvqnEQ8vI0XfN1kyMxX6ryIIQHg4ybTu7yrrtXWiKh+4vATIiIiIiKyKLY2EsRFBAG4n8AoT30/LiJIZ0KDiBoW9tQgohozVZfR707dwuvbT+P9pztjaFc/o9fPLqlERESWb3AnX6wa2x1xSWc15tDy4aTgRFQOkxpEVGOm6jLqYGcr/ssuqURERA3X4E6+eLC1pzikdH10T73ndyJqmDj8hIiMztRdRtkllYiIqOEofz4PDfTg+Z2INDCpQUQmoe4yKneVapT7uMmwamz3WncZNXX9REREZDpyFymm92sDuYu06o0tsH4ishwcfkJEJmPqLqPskkpERGSd5K4yzBjQVqs8W6FEdkGRRln5pdjTbil0DiOVu0ghd5VVWT8R1T/1OqmxcuVKLFmyBJmZmejatStWrFiB0NBQnduuX78e0dHRGmVSqRRKpbIuQiWqt0zdZZRdUomIiOqPTcfSsTzlkt7H9U30Pb1fGyYxiBqoepvU2Lp1K2JjY5GYmIiwsDAkJCRg0KBBuHDhAuRyuc59XF1dceHCBfG+RMIfR0S66OvSaayrK229nREW6IG23s5Gjvy/52OXVCIiIsszJiwAA4K8q70fz+lEDZdEEARdiwdYvbCwMPTs2RMff/wxAKCsrAz+/v6YOnUqZs2apbX9+vXr8corr+DOnTs1ej6FQgE3Nzfk5+fD1VV7iUuihuDDvRcrvbqiT22urtwrLkXQ3N0AgLQFg3SulkLUUPHcZHw8pkREZEl4XqqnPTWKi4tx4sQJzJ49WyyzsbFB//79cfSo7i5rAHD37l00b94cZWVl6N69OxYuXIiOHTvq3LaoqAhFRf9dkVYoFMZ7AURWildXiIiIiIioLtXLpEZubi5UKhW8vTV/XHl7e+P8+fM692nXrh3Wrl2LLl26ID8/H0uXLkXv3r1x9uxZNGvWTGv7+Ph4zJ8/3yTxE1kruatMYxiJsRlreAsREREREdUP9TKpURO9evVCr169xPu9e/dGhw4d8Omnn+Kdd97R2n727NmIjY0V7ysUCvj7+9dJrEQNFScPIyIiIiKi8uplUsPT0xO2trbIysrSKM/KyoKPj49Bddjb26Nbt264fPmyzselUimkUnaZJ6pLHN5CRERERETl1cukhoODA0JCQpCSkoLhw4cDuD9RaEpKCmJiYgyqQ6VS4fTp03j88cdNGCkRVYeph7cQEREREZF1qZdJDQCIjY3F+PHj0aNHD4SGhiIhIQGFhYWIjo4GAERGRqJp06aIj48HACxYsAAPPPAAWrdujTt37mDJkiX466+/8Pzzz5vzZRARERERkZHomp/LEJyfi8hy1dukxsiRI5GTk4O5c+ciMzMTwcHBSE5OFicPTU9Ph42Njbj9P//8g0mTJiEzMxPu7u4ICQnBkSNHEBQUZK6XQERERERERlTV/Fz6cH4uIsslEQRBMHcQ9QHXByYiIkvDc5Px8ZgSWTd9K6mpJxvfNrkXV1Ijq8LzUj3uqUFERERERFServm57hWXin8H+bnCyYE/kYisiU3VmxARERERERERWR4mNYiIiIiIqF7JVijx4d6LyFYorbJ+IjIc+1YREVGDwBnviYgajuyCIixPuYQBQd4m+Q43df1EZDgmNYiIqEHgjPdERKSLquy/dRNSr+UhvI0XbG0kZoyIiKqDSQ0iImoQxoQFYECQt0aZoTPeExFR/ZR8JgNxSWfF+1HrjsPXTYa4iCAM7uRrxsiIyFBMahARUYPAGe+JiAyXrVBi07F0jAkL0PjuNNZQPn3116XkMxmYsvEkhArlmflKTNl4EqvGdmdig8gKsPVGREREREQa9M0ZYayhfHU1J4WyRKWRwFZTlQmISzqrldAAAAGABMC8pDQ82NpT51AUZYnK6LESUc0wqUFERERERAbRNZSvsKgUI1f/CgCYFxGEbgHuWokAcw3lUw8xrC4BQKZCic7z9hg3ICIyOiY1iIiIiIjIIBWH8lWck2LerjTOSUFEdYpJDSIiIiIiqjZrmJNi2+ReCPJz1SpPvZaHqHXHq9x/fXRPhAZ6aJWn3VLUuBcIERkXkxpERERERKSTtc9JIbO31TkJdHgbL/i6yZCZr9T5GiQAfNxkepd31bVaFhGZB5MaRERERESkU32dk8LWRoK4iCBM2XgSEkAjsaFOYcRFBOlMaBCRZbExdwBERERERER1bXAnX6wa2x1yV81JTH3cZBYxdIaIDMOeGkREREREpFN9n5NicCdfPNjaU+xRsj66p94hJ0RkmZjUICIiIiIinRrCnBTl4wsN9GBCg8jKcPgJERERERFVi3pOCuC/OSjULGFOCrmLFNP7tYHcRVr1xhZYPxEZjj01iBqQbIUS2QVF1d5P7iLVWJOeiIiISD0nRVzSWWQp/mtf+LjJEBcRZNY5KeSuMswY0NZq6yciwzGpQVQPZSuU2HQsHWPCAjSSEZuOpWN5yqVq1ze9XxuNE7e++omIiKh+q3iBpJm7Ez4a1Q0jV/8KAJgXEYRuAe6wtZHgzM18cTteICEiU2FSg6geyi4owvKUSxgQ5K3RgBgTFoABQd4a2ypLVOJEXdsm99I5xrVi10p99RMREVH9VtUFknm70nSWV7xAQkRkLExqEDUgcleZVhKiQFki/n23qFS8ukJEZA1WrlyJJUuWIDMzE127dsWKFSsQGhqqd/s7d+5gzpw52LFjB/Ly8tC8eXMkJCTg8ccfr8OoiSyfvjkjdF0gMbQ+Q+onIqouJjWIGrDkMxmISzor3o9adxy+FjAOlojIEFu3bkVsbCwSExMRFhaGhIQEDBo0CBcuXIBcLtfavri4GAMGDIBcLse2bdvQtGlT/PXXX2jcuHHdB09k4fTNGaHrAokx6zc1XfOLKUtU4t9ptxR6e62ydyqRZWJSg6iBSj6TgSkbT2otw5aZr8SUjSexamx3JjaIyKItW7YMkyZNQnR0NAAgMTER33//PdauXYtZs2Zpbb927Vrk5eXhyJEjsLe3BwC0aNGiLkMmIjOraviMekhuRRw+Q2S5mNQgqseUJSrcKy7VKleVCYhLOqtzXXkB95dim5eUhgdbe+ocilL+igYRkTkUFxfjxIkTmD17tlhmY2OD/v374+hR3T9KkpKS0KtXL7z88sv49ttv4eXlheeeew5vvPEGbG21r8wCQFFREYqK/ruqq1AojPtCiKhOGWv4DBFZDiY1iOoxfVcbqiIAyFQo0XneHuMGRERkJLm5uVCpVPD21vxx4u3tjfPnz+vc5+rVq9i3bx/GjBmDH374AZcvX8ZLL72EkpISxMXF6dwnPj4e8+fPN3r8RGQexho+o4+u4S2G4PAWoppjUoOIiIgahLKyMsjlcqxevRq2trYICQnBzZs3sWTJEr1JjdmzZyM2Nla8r1Ao4O/vX1chE5GVqWp4iz4c3kJUc0xqENVj2yb3QpCfq1Z56rU8RK07XuX+66N7IjTQQ6s87Zaixr1AiIiMwdPTE7a2tsjKytIoz8rKgo+Pj859fH19YW9vrzHUpEOHDsjMzERxcTEcHBy09pFKpZBK2e2ciAyja3iLskQltpu2Te6ldyJSIqoZJjWI6jGZvS2cHLT/m4e38YKvmwyZ+Uqd82pIAPi4yRDexkvnnBq6Tsb1AbuMElkPBwcHhISEICUlBcOHDwdwvydGSkoKYmJidO7z4IMPYvPmzSgrK4ONjQ0A4OLFi/D19dWZ0CAiqi5dw1vKz28W5Oeqs21GRDXH/1FEDZCtjQRxEUGYsvEkJIBGYkOdwoiLCNKZ0KjP2GWUyLrExsZi/Pjx6NGjB0JDQ5GQkIDCwkJxNZTIyEg0bdoU8fHxAIApU6bg448/xvTp0zF16lRcunQJCxcuxLRp08z5MoiIiKgWmNQgaqAGd/LFqrHdEZd0FlmK/3on+LjJEBcR1CCXc2WXUSLrMnLkSOTk5GDu3LnIzMxEcHAwkpOTxclD09PTxR4ZAODv74/du3djxowZ6NKlC5o2bYrp06fjjTfeMNdLICIiolpiUoOoAak4vKKZuxM+GtUNI1f/CgCYFxGEbgHusLWR4MzNfHG7hjK8gl1GiaxPTEyM3uEm+/fv1yrr1asXfv31VxNHRURERHWFrXOiBqSq4RXzdqXpLOfwCiIiIiIiskRMahDVQ3IXKab3a6M1LELX8ApD6zOkfiIiIiIiorrEpAZRPSR3lensWaFreIUx6yciIiIiIqpLNlVvQkRERERERERkeZjUICIiIiIiIiKrxKQGETU42QolPtx7EdkKpVXWT0RERERE9zGpQUQNTnZBEZanXNJY3taa6iciIiLz4gUSIsvBpAYRUSVUZYL4d+q1PI37RERE1DDxAgmR5WBSg4hIj+QzGei/7IB4P2rdcTy0eB+Sz2SYMSoiIiKyJrxAQmRaTGoQEemQfCYDUzaeRJZC8wpJZr4SUzaeZGKDiIiIqsQLJESmx6QGETVYyhIV7hWXat0KlCWISzoLXddR1GXzktJQoCzRub+yRFWXL4OIiIgsEC+QENUNO3MHQERkLiMSj9ZoPwFApkKJzvP2GDcgIiIisirqCyQVqcqESi+QSHD/AsmDrT1hayPRWS8RGYZJDSIiIiIiohrgBRIi82NSg4garG2TeyHIz1WrPPVaHqLWHa9y//XRPREa6KFVnnZLUeNGDhERERERGY5JDSJqsGT2tnBy0P4aDG/jBV83GTLzlTq7jUoA+LjJEN7GS2eXUZm9rfGDJSIiIovDCyRE5sekBhFRBbY2EsRFBGHKxpOQABqJDXUKIy4iSGdCg4iIiBqOur5Akq1QIrugSOdjlZG7SCF3lVV7PyJrwKQGEZEOgzv5YtXY7ohLOqsxa7mPmwxxEUEY3MnXjNERERGRJTPVBZJNx9KxPOVSteOZ3q8NZgxoW+39iKwBkxpERHoM7uSLB1t7ipN4rY/uqfeKChEREVF5prhAMiYsAAOCvDXKlCUqcajKtsm9dPbykLtIq/1cRNaCSQ0iokqUT2CEBnowoUFEREQGM/YFErmrTGsYSfklZYP8XHUOhyGqz2zMHQARUV2Tu0gxvV8bk121MHX9REREZD14gYTItJjGI6IGR+4qM+m4UlPXT0RERObFCyREloNJDSIiIiIiomrgBRIiy8HhJ0RERERERERklZjUICIiIiIiIiKrxKQGEREREREREVklJjWIiIiIiIiIyCpxolAiov+XrVAiu6BIo0xZohL/TrulgMzeVms/uYtUa814IiIiIiIyPSY1iIj+36Zj6Viecknv4yMSj+osn96vDWcoJyIiIiIyAyY1iIj+35iwAAwI8q72flxDnoiIiIjIPJjUICL6f3JXGYeREBGRVchWKLHpWDrGhAWY5Nxl6vqJiIyFSQ0iIiIiIiuTXVCE5SmXMCDIWyPpoGt+KENUnB9KX/1ERJaGSQ0ionrCWA1ZIiKyXlXND6UP54cyDk46TlT3mNQgIqon2JAlIqKK80MduZyLT3+5ituFxWKZp7MDXghvid6tPcUyzg9lHJx0nKjuMalBRFRP6JroVFmiEhtQ2yb30nt1qD7hOHMiasjKzw+VfCYD8T+eh1Bhm9t3ixH/43msGtsdgzv51n2Q9RgnHSeqe0xqEBHVE7omOr1XXCr+HeTnCieH+v+1b+px4BxnTkSWRFmi0viuV1OVCYhLOquV0AAAAYAEwLykNDzY2hO2NhKd9VL1cdJxorpX/1u3REREeqjK/mvup17LQ3gbL52NeyIiS6VvOENVBACZCiU6z9tj3ICIiOqYjbkDMKWVK1eiRYsWkMlkCAsLQ2pqqkH7bdmyBRKJBMOHDzdtgEREZDbJZzLQf9kB8X7UuuN4aPE+JJ/JMGNURERERFQd9banxtatWxEbG4vExESEhYUhISEBgwYNwoULFyCXy/Xud/36dbz66qsIDw+vw2iJiKguJZ/JwJSNJ7W6ZWfmKzFl40mOMyciq7Ftci8E+blqladey0PUuuNV7r8+uidCAz20ytNuKWrcC4Rqj/NDERmu3vbUWLZsGSZNmoTo6GgEBQUhMTERTk5OWLt2rd59VCoVxowZg/nz56Nly5Z1GC0RERmbepx5xVuBsqTScebA/XHmBcoSnftznDkRWRKZvS2cHOy0buFtvODrJoO+AXUSAL5uMoS38dK5v66JpanuqOdvqslS7ZZQP1Fdqpc9NYqLi3HixAnMnj1bLLOxsUH//v1x9Kj+jPOCBQsgl8sxceJEHDx4sNLnKCoqQlHRf18CCoWi9oETEZHRcJw5ETVktjYSxEUEYcrGk5AAGolcdaIjLiKI8wjVA5wfihq6etlTIzc3FyqVCt7emsspeXt7IzMzU+c+hw4dwpo1a/DZZ58Z9Bzx8fFwc3MTb/7+/rWOm4iIiIioNrIVSpy5mY8zN/PRzN0Jsx9rD49GDhrbNHF2wOzH2qOZu5O4bbZCaaaIqTY4PxRRPe2pUV0FBQUYN24cPvvsM3h6ehq0z+zZsxEbGyveVygUTGwQEVkQjjMnooZo07F0LE+5VOk2uXeLsfDH8xpl0/u1wYwBbU0ZGhkZ54ciuq9eJjU8PT1ha2uLrKwsjfKsrCz4+PhobX/lyhVcv34dERERYllZWRkAwM7ODhcuXECrVq009pFKpZBKpSaInoiIjEE9zrwi9TjzzHylznk1JAB8/n+cua7uuxxnTkSWbExYAAYEeVe9YQVyF7ZrLZF6fqiKVGVCpfNDSXB/fqgHW3vqPJdxfiiqT+plUsPBwQEhISFISUkRl2UtKytDSkoKYmJitLZv3749Tp8+rVH21ltvoaCgAMuXL2cPDCIyimyFskYTcsldpJyZ3Ig4zpyI6gO5ixTT+7XRSkbIXWVGOWfoq5/qFueHIqpavUxqAEBsbCzGjx+PHj16IDQ0FAkJCSgsLER0dDQAIDIyEk2bNkV8fDxkMhk6deqksX/jxo0BQKuciKimDOkSrAu7BBvf4E6+WDW2O+KSziJL8V+iycdNhriIIHbXJSKLJ3eVmfTcYOr6yTLxAgxZo3qb1Bg5ciRycnIwd+5cZGZmIjg4GMnJyeLkoenp6bCxqZfzpBKRhdLVJVhZohKvwmyb3Evn0AZeJTONwZ188WBrT/Eq1vronpwxnoiILEpdzw/FCzCVY9LHMtXbpAYAxMTE6BxuAgD79++vdN/169cbPyAiatB0dQkuP042yM9V5xwQZDrlExihgR5MaBARkUWp6/mheAGmckz6WCa2nomIqF4x9ThwjjO3LCtXrsSSJUuQmZmJrl27YsWKFQgNDdW57fr168VhqGpSqRRKJZeyJCLrYqr5oXgBpnJM+limhvuJJCKieonjzBuOrVu3IjY2FomJiQgLC0NCQgIGDRqECxcuQC6X69zH1dUVFy5cEO9LJOydQ0TWifND1T0mfSwTJ5UgIiIiq7Rs2TJMmjQJ0dHRCAoKQmJiIpycnLB27Vq9+0gkEvj4+Ig39VxbRETWaHAnX/wU20e8vz66Jw698SgTGtSgMKlBREREVqe4uBgnTpxA//79xTIbGxv0798fR4/qXwLx7t27aN68Ofz9/fHEE0/g7NmzlT5PUVERFAqFxo2IyJJwfihq6CwmqdGnTx988cUX+Pfff80dChEREVm43NxcqFQqrZ4W3t7eyMzM1LlPu3btsHbtWnz77bfYuHEjysrK0Lt3b9y4cUPv88THx8PNzU28+fv7G/V1EBHpwvmhiAxnMUmNbt264dVXX4WPjw8mTZqEX3/91dwhERERUT3Sq1cvREZGIjg4GH369MGOHTvg5eWFTz/9VO8+s2fPRn5+vnj7+++/6zBiImqo1PM3mWoZUFPXT1SXLCapkZCQgFu3bmHdunXIzs7Gww8/jKCgICxduhRZWVnmDo+IiIgsiKenJ2xtbbXaCFlZWfDx8TGoDnt7e3Tr1g2XL1/Wu41UKoWrq6vGjYiIiCyHxSQ1AMDOzg5PPfUUvv32W9y4cQPPPfcc3n77bfj7+2P48OHYt2+fuUMkIiIiC+Dg4ICQkBCkpKSIZWVlZUhJSUGvXr0MqkOlUuH06dPw9eWEekRERNbKopIaaqmpqYiLi8MHH3wAuVyO2bNnw9PTE0OHDsWrr75q7vCIiIjIAsTGxuKzzz7Dhg0bcO7cOUyZMgWFhYWIjo4GAERGRmL27Nni9gsWLMCePXtw9epVnDx5EmPHjsVff/2F559/3lwvgYiIiGrJYhbRzc7Oxpdffol169bh0qVLiIiIwFdffYVBgwaJa8hHRUVh8ODBWLp0qZmjJSIyn2yFEpuOpWNMWIBJxsKaun4iYxk5ciRycnIwd+5cZGZmIjg4GMnJyeLkoenp6bCx+e/6zT///INJkyYhMzMT7u7uCAkJwZEjRxAUFGSul0BERES1ZDFJjWbNmqFVq1aYMGECoqKi4OXlpbVNly5d0LNnTzNER0RkOFMnBc5nFmB5yiWENHc3TVKjoAjLUy5hQJA3kxpk8WJiYhATE6Pzsf3792vc//DDD/Hhhx/WQVRERERUVywmqZGSkoLw8PBKt3F1dcXPP/9cRxEREdWMqZMCeYXFGv9WRlUmiH+nXstDeBsvrl9PRERERPWGxcyp0axZM1y6dEmr/NKlS7h+/XrdB0REVAcqJh3K36+t5DMZ6L/sgHg/at1xPLR4H5LPZBjtOYiIiIiIzMlikhpRUVE4cuSIVvmxY8cQFRVV9wEREZmYKZMOyWcyMGXjSWQpijTKM/OVmLLxJBMbRERERFQvWExS4/fff8eDDz6oVf7AAw/gjz/+qPuAiIhMyBhJh+JSFe4Vl2rdCpQliEs6C119PtRl85LSUKAs0bm/skRV+xdIRERERFQHLGZODYlEgoKCAq3y/Px8qFRsYBOR9VGW3E86VKQqEypNOkhwP+nwYGtPnfNfFJfe/058fftpvL79dLXjEgBkKpToPG9PtfclIiIiIrIkFpPUePjhhxEfH4+vvvoKtra2AACVSoX4+Hg89NBDZo6OiKj6RiQerdF+TDoQERERERnGYpIaixcvxsMPP4x27dqJq6AcPHgQCoUC+/btM3N0RMZl6iU/TV0/WYb3n+6MoV39tMpTr+Uhat3xKvdfH90ToYEeWuVptxQ1TsgQVaZPnz6YOHEinnnmGTg6Opo7HCIiIqoHLCapERQUhD///BMff/wxTp06BUdHR0RGRiImJgYeHtqNbiJrpm/Jz2yFEtkFRZXsqZvcRapZj4mXFCXDbJvcC0F+rlrltU06fHfqFl7ffhoOdrZwctD+Gg9v4wVfNxky85U6h7hIAPi4yfQu7yqzt60yNqKa6NatG1599VVMnToVzz77LCZOnIgHHnjA3GERERGRFbOYpAYA+Pn5YeHCheYOg8hsNh1Lx/IU7aWNqzK9XxvMGNDWBBFRbcjsTZN0cLCrPOlgayNBXEQQpmw8CQmg8Rzq2uIignTWTWRKCQkJWLp0KZKSkrBhwwY8/PDDaN26NSZMmIBx48bB29vb3CESERGRlbGopMadO3eQmpqK7OxslJWVaTwWGRlppqiI6s6YsAAMCPqvUX/kci4+/eUqbhcWi2Wezg54Ibwlerf2FMvkLtI6jZNqpy6SDoM7+WLV2O6ISzqrscKKj5sMcRFBGNzJt8Z1E9WGnZ0dnnrqKTz11FPIzs7G6tWr8fbbb+PNN9/E448/jmnTpuHRRx81d5hERERkJSwmqbFr1y6MGTMGd+/ehaurKySS/xrzEomESQ1qEOSuMnG4SPKZDMT/eF7rSv7tu8WI//E8Vo3tzh+mVqwukg6DO/niwdae4oSj66N76u39QVTXUlNTsW7dOmzZsgVyuRxRUVG4efMmhg4dipdeeglLly41d4hERERkBSwmqTFz5kxMmDABCxcuhJOTk7nDIaoTplryU1nCZZCtQV0kHcrXFRrowYQGmVV2dja+/PJLrFu3DpcuXUJERAS++uorDBo0SLyYERUVhcGDBzOpQURERAaxmKTGzZs3MW3aNCY0qEHhkp9Uk6SDRyN7jX+NTe4ixfR+bTisiYyuWbNmaNWqFSZMmICoqCh4eXlpbdOlSxf07NnTDNERERGRNbKYpMagQYPw22+/oWXLluYOhYioVkydFGjv44rp/dqgvY/2yirGIHeVceJZMomUlBRx2XZ9XF1d8fPPP9dRRERERGTtLCapMWTIELz22mtIS0tD586dYW+veQVy2LBhZoqMyHRMteRn2i1FjXuBUO2ZOinApANZq7i4OOzYsQONGzfWKFcoFBg+fDj27dtnnsCIiKxEtkKJ7IIijbLyw47Tbil0Ls0ud5GK87YR1TcWk9SYNGkSAGDBggVaj0kkEqhUnCOA6h9TLfmp62RGRGRuBw4cQHFxsVa5UqnEwYMHzRAREZF12XQsHctTLul9XN9Fren92vCCCNVbFpPUqLiEK1FDVhdLfhIR1ZU///wTACAIAtLS0pCZmSk+plKpkJycjKZNm5orPCIiqzEmLAADgryrvR/nyaL6zGKSGuUplUrIZOweRQ1P+S6FzdydMPux9vj0l6u4Xfjflc0mzg54Ibwlmrk74czNfADsUkhEli04OBgSiQQSiQSPPvqo1uOOjo5YsWKFGSIjIrIuclcZ23xEFVhMUkOlUmHhwoVITExEVlYWLl68iJYtW+Ltt99GixYtMHHiRHOHSGRyVXUpBIDcu8VY+ON5jTJ2KSQiS3bt2jUIgoCWLVsiNTVVY9UTBwcHyOVy2Npy2BwRUV3JViix6Vg6xoQFmCRJYur6icqzmKTGe++9hw0bNuD9998X59cAgE6dOiEhIYFJDWoQ2KWQiOqj5s2bA+BQUyIiS5FdUITlKZcwIMjbNEkNE9dPVJ7FJDW++OILrF69Gv369cPkyZPF8q5du+L8+fOV7ElkffQt+WmsLoWmXlKUiMhQSUlJeOyxx2Bvb4+kpKRKt+VKZ0RElkdV9t/sbqnX8vROVE9kLhaT1Lh58yZat26tVV5WVoaSkhIzRERkOlzyk4gaiuHDhyMzMxNyuRzDhw/Xux1XOiMisjzJZzIQl3RWvB+17jh83WSIiwjC4E6+ZoyM6D825g5ALSgoSOdybtu2bUO3bt3MEBERERHVVllZGeRyufi3vhsTGkREliX5TAambDyJLEWRRnlmvhJTNp5E8pkMM0VGpMliemrMnTsX48ePx82bN1FWVoYdO3bgwoUL+OKLL/Ddd9+ZOzwiIiKqhZKSEgwePBiJiYlo06aNucMhIiIAyhIV7hWXapWrygTEJZ2FoGMfAYAEwLykNDzY2lPnUBRlCRPVVHcsJqnxxBNPYNeuXViwYAEaNWqEuXPnonv37ti1axcGDBhg7vCIiGqt/JK9auVP+mm3FJDZa68AwSV7qT6wt7fHn3/+ae4wiIionBGJR2u0nwAgU6FE53l7jBsQUQ1YTFIDAMLDw7F3715zh0FEZBJVLdmrr2HBJXupvhg7dizWrFmDRYsWmTsUIiIiqicsKqlBRFSfccleauhKS0uxdu1a/PTTTwgJCUGjRo00Hl+2bJmZIiMiapi2Te6FID9XrfLUa3mIWne8yv3XR/dEaKCHVnnaLUWNe4EQVZdZkxoeHh64ePEiPD094e7uDolE/9JAeXl5dRgZEZHxGWvJXiJrdebMGXTv3h0AcPHiRY3HKmsDEBGRacjsbeHkoP2TMLyNF3zdZMjMV+qcV0MCwMdNpnd5V13DaYlMxaxJjQ8//BAuLi7i32zQEBER1V8///yzuUMgIiID2NpIEBcRhCkbT0ICaCQ21L/Y4iKCdCY0iOqaWZMa48ePF/+OiooyXyBEREREREQkGtzJF6vGdkdc0lmNZV193GSIiwjC4E6+ZoyO6D8WM6fGDz/8AFtbWwwaNEijfM+ePVCpVHjsscfMFBkREREZwyOPPFJpr8x9+/bVYTRERFSVwZ188WBrT3GVk/XRPfUOOSEyFxtzB6A2a9YsqFTa6xmXlZVh1qxZZoiIiIiIjCk4OBhdu3YVb0FBQSguLsbJkyfRuXNnc4dHREQ6lE9ghAZ6MKFBFsdiempcunQJQUFBWuXt27fH5cuXzRARERERGdOHH36os3zevHm4e/duHUdDRNRwyV2kmN6vjclWWDN1/UTlWUxPDTc3N1y9elWr/PLly1pLvhEREVH9MXbsWKxdu9bcYRARNRhyVxlmDGhrslXZTF0/UXkW01PjiSeewCuvvIKdO3eiVatWAO4nNGbOnIlhw4aZOToiIsuXrVAiu6BIo0xZ8t+wvrRbCp1LrMldpGx0kFkdPXoUMhk/g0RE9Z2utooh2FahylhMUuP999/H4MGD0b59ezRr1gwAcOPGDYSHh2Pp0qVmjo6IyPJtOpaO5SmX9D4+IvGozvLp/dpgxoC2pgqLSPTUU09p3BcEARkZGfjtt9/w9ttvmykqIiKqK1W1VfRhW4UqYzFJDTc3Nxw5cgR79+7FqVOn4OjoiC5duuDhhx82d2hERFZhTFgABgR5V3s/jneluuLm5qZx38bGBu3atcOCBQswcOBAM0VFRER1RVdbRVmiEi+8bJvcS2+vUiJ9LCapAQASiQQDBw5kw4aIqAbkrjJ2zSSLtm7dOnOHQEREZqSrrXKvuFT8O8jPFU4OFvUTlayAWT8xH330EV544QXIZDJ89NFHlW47bdq0OoqKiIiITOm3337DuXPnAABBQUEICQkxc0RERERkrcya1Pjwww8xZswYyGQyvcu8Afd7cDCpQUREZN1u3LiB0aNH4/Dhw2jcuDEA4M6dO+jduze2bNkizqlFREREZCizJjX++OMPcXzttWvXzBkKERERmdjzzz+PkpISnDt3Du3atQMAXLhwAdHR0Xj++eeRnJxs5giJiIjI2tiY88k9PDyQnZ0NAHj00Udx584dc4ZDREREJnTgwAGsWrVKTGgAQLt27bBixQr88ssvZoyMiIiIrJVZkxrOzs64ffs2AGD//v0oKSkxZzhERERkQv7+/jrP9SqVCn5+fmaIiIiIiKydWYef9O/fH4888gg6dOgAAHjyySfh4OCgc9t9+/bVZWhERERkZEuWLMHUqVOxcuVK9OjRA8D9SUOnT5+OpUuXmjk6IiIiskZmTWps3LgRGzZswJUrV3DgwAF07NgRTk5O5gyJiIjqqWyFEtkFRRplyhKV+HfaLQVk9rZa+8ldpFwq10iioqJw7949hIWFwc7ufhOktLQUdnZ2mDBhAiZMmCBum5eXZ64wiYiIyIqYNalRUlKCyZMnA7h/pWbx4sXibOhERETGtOlYOpanXNL7+IjEozrLp/drgxkD2poqrAYlISHB3CEQERFRPWPWpIa7uzsyMjIgl8shkUjMGQoREdVzY8ICMCDIu9r7yV2kJoimYRo/fry5QyAiIqJ6xmImCj1w4AAnCiUiIpORu8rQqalbtW8cemLZVq5ciRYtWkAmkyEsLAypqakG7bdlyxZIJBIMHz7ctAESERGRSVnMRKGCIHCiUCIiIjLY1q1bERsbi8TERISFhSEhIQGDBg3ChQsXIJfL9e53/fp1vPrqqwgPD6/DaImIiMgUOFEoERERWaVly5Zh0qRJiI6OBgAkJibi+++/x9q1azFr1iyd+6hUKowZMwbz58/HwYMHcefOnTqMmIiIiIzNrEkNR0dHThRKRERE1VZcXIwTJ05g9uzZYpmNjQ369++Po0d1T/oKAAsWLIBcLsfEiRNx8ODBKp+nqKgIRUX/rZqjUChqFzgREREZlVnn1Cjv559/RuPGjVFcXIwLFy6gtLTU3CERERGRhcrNzYVKpYK3t+bkr97e3sjMzNS5z6FDh7BmzRp89tlnBj9PfHw83NzcxJu/v3+t4iYiIiLjMmtPjfL+/fdfxMTEYMOGDQCAixcvomXLlpg6dSqaNm2qtxspERERWYcnn3xS52pnEokEMpkMrVu3xnPPPYd27doZ/bkLCgowbtw4fPbZZ/D09DR4v9mzZyM2Nla8r1AomNggIiKyIBbTU2PWrFk4deoU9u/fD5nsv5nm+/fvj61bt5oxMiIiIjIGNzc37Nu3DydPnoREIoFEIsHvv/+Offv2obS0FFu3bkXXrl1x+PDhKuvy9PSEra0tsrKyNMqzsrLg4+Ojtf2VK1dw/fp1REREwM7ODnZ2dvjiiy+QlJQEOzs7XLlyRefzSKVSuLq6atyIiIjIclhMUuObb77Bxx9/jIceekjjKk7Hjh31NjSIiIjIevj4+OC5557D1atXsX37dmzfvh1XrlzB2LFj0apVK5w7dw7jx4/HG2+8UWVdDg4OCAkJQUpKilhWVlaGlJQU9OrVS2v79u3b4/Tp0/jjjz/E27Bhw/DII4/gjz/+YO8LIiIiK2Uxw09ycnJ0Lr9WWFios6sqERERWZc1a9bg8OHDsLH575qKjY0Npk6dit69e2PhwoWIiYkxeKnV2NhYjB8/Hj169EBoaCgSEhJQWFgoroYSGRmJpk2bIj4+HjKZDJ06ddLYXz05ecVyIiIish4W01OjR48e+P7778X76kTG559/rvOKiyFWrlyJFi1aQCaTISwsDKmpqXq33bFjB3r06IHGjRujUaNGCA4Oxpdfflmj5yUiIiJtpaWlOH/+vFb5+fPnoVKpAAAymczgixkjR47E0qVLMXfuXAQHB+OPP/5AcnKyOHloeno6MjIyjPcCiIiIyOJYTE+NhQsX4rHHHkNaWhpKS0uxfPlypKWl4ciRIzhw4EC169u6dStiY2ORmJiIsLAwJCQkYNCgQbhw4YLOHiEeHh6YM2cO2rdvDwcHB3z33XeIjo6GXC7HoEGDjPESiYiIGrRx48Zh4sSJePPNN9GzZ08AwPHjx7Fw4UJERkYCAA4cOICOHTsaXGdMTAxiYmJ0PrZ///5K912/fr3Bz0NERESWyWKSGg899BD++OMPLFq0CJ07d8aePXvQvXt3HD16FJ07d652fcuWLcOkSZPELqiJiYn4/vvvsXbtWp0rqfTt21fj/vTp07FhwwYcOnSISQ0iIiIj+PDDD+Ht7Y33339fnODT29sbM2bMEOfRGDhwIAYPHmzOMImIiMiKWExSAwBatWpVrbXj9SkuLsaJEycwe/ZssczGxgb9+/fH0aNHq9xfEATs27cPFy5cwOLFi3VuU1RUhKKiIvG+QqGoddxERET1ma2tLebMmYM5c+aI582Kq4kEBASYIzQiIiKyUhaV1DCW3NxcqFQqcUytmre3t86xvGr5+flo2rQpioqKYGtri08++QQDBgzQuW18fDzmz59v1LiJiIgaCi6NSkRERMZgMROFWgIXFxf88ccfOH78ON577z3ExsbqHY87e/Zs5Ofni7e///67boMlIiKyMllZWRg3bhz8/PxgZ2cHW1tbjRsRERFRddXLnhqenp6wtbUVx+uqZWVlwcfHR+9+NjY2aN26NQAgODgY586dQ3x8vNZ8GwAglUohlUqNGjcREVF9FhUVhfT0dLz99tvw9fXlku1ERBYmW6FEdkGRRpmyRCX+nXZLAZm9dhJa7iKF3FVm8viIdKmXSQ0HBweEhIQgJSUFw4cPBwCUlZUhJSVF7wzpupSVlWnMm0FEREQ1d+jQIRw8eBDBwcHmDoWIiHTYdCwdy1Mu6X18RKLu+Qmn92uDGQPamiosokrVy6QGAMTGxmL8+PHo0aMHQkNDkZCQgMLCQnE1lMjISDRt2hTx8fEA7s+R0aNHD7Rq1QpFRUX44Ycf8OWXX2LVqlXmfBlERET1hr+/PwRBMHcYRESkx5iwAAwI8q56wwrkLuzBTuZjMUmNJ598Umc3VIlEAplMhtatW+O5555Du3btDKpv5MiRyMnJwdy5c5GZmYng4GAkJyeLk4emp6fDxua/KUUKCwvx0ksv4caNG3B0dET79u2xceNGjBw50jgvkIiIqIFLSEjArFmz8Omnn6JFixbmDoeIiCqQu8o4jISsjkSwkEsmUVFR+Oabb9C4cWOEhIQAAE6ePIk7d+5g4MCBOHXqFK5fv46UlBQ8+OCDZo5Wm0KhgJubG/Lz8zmjOxERWQRLOze5u7vj3r17KC0thZOTE+zt7TUez8vLM1NkhrO0Y0pEZO3uFZciaO5uAEDagkFwcrCY6+4GMXf8PC9ZUE8NHx8fPPfcc/j444/FHhRlZWWYPn06XFxcsGXLFkyePBlvvPEGDh06ZOZoiYiIqLoSEhLMHQIRERHVMxaT1FizZg0OHz6sMSTExsYGU6dORe/evbFw4ULExMQgPDzcjFESERFRTY0fP97cIRAREVE9YzFJjdLSUpw/fx5t22rOmnv+/HmoVPeXEZLJZFz+jYiIyIooFAqxO6xCoah024babZaIiIhqzmKSGuPGjcPEiRPx5ptvomfPngCA48ePY+HChYiMjAQAHDhwAB07djRnmERERFQN7u7uyMjIgFwuR+PGjXVenBAEARKJRLyIQURERGQoi0lqfPjhh/D29sb777+PrKwsAIC3tzdmzJiBN954AwAwcOBADB482JxhEhERUTXs27cPHh4eAICff/7ZzNEQERFRfWMxSQ1bW1vMmTMHc+bMEbunVuyGGhAQYI7QiIiIqIb69Omj828iIiJjy1YokV1QVO395C5SLmVrxSwmqVEex9QSETU81t4Qsfb468qdO3eQmpqK7OxslJWVaTymHm5KRERUE5uOpWN5yqVq7ze9XxvMGNC26g3JIllMUiMrKwuvvvoqUlJSkJ2dDUEQNB7nOFsiovrN2hsi1h5/Xdi1axfGjBmDu3fvwtXVVWN+DYlEwqQGERHVypiwAAwI8tYoU5aoMCLxKABg2+RekNnbau0nd5HWSXxkGhaT1IiKikJ6ejrefvtt+Pr6cpUTIqIGxtobItYef12YOXMmJkyYgIULF8LJycnc4RARUT0jd5Vp9X68V1wq/h3k5wonB4v5CUxGYjHv6KFDh3Dw4EEEBwebOxQiIjIDa2+IWHv8deHmzZuYNm0aExpERERkNDbmDkDN399fa8gJERER1R+DBg3Cb7/9Zu4wiIiIqB6xmEtGCQkJmDVrFj799FO0aNHC3OEQERGRESQlJYl/DxkyBK+99hrS0tLQuXNn2Nvba2w7bNiwug6PiIiIrJzFJDVGjhyJe/fuoVWrVnByctJq6OTl5ZkpMiIiIqqp4cOHa5UtWLBAq0wikXBScCIiIqo2i0lqJCQkmDsEIiIiMrKKy7YSERERGZPFJDXGjx9v7hCIRNkKJTYdS8eYsACtif+soX4iIiIiIqKGwKwThSoUCo2/K7sR1aXsgiIsT7mE7IIiq6yfiMgSTZs2DR999JFW+ccff4xXXnml7gMiIiIiq2fWpIa7uzuys7MBAI0bN4a7u7vWTV1OZKlUZQKOXrmNb/+4iaNXbkNVxlV8iIh02b59Ox588EGt8t69e2Pbtm1miIiIiIisnVmHn+zbtw8eHh4AgJ9//tmcoRAZJFuh1OhdceRyLlYfvIrcu8VimaezA14Ib4nerT3FMrmLlMNMiKjBu337Ntzc3LTKXV1dkZuba4aIiIiIyNqZNanRp08fnX8TWapNx9KxPOVSpdvk3i3Gwh/Pa5RN79cGMwa0NWVoREQWr3Xr1khOTkZMTIxG+Y8//oiWLVuaKSoiIiKyZhYzUSgA3LlzB6mpqcjOztaaLT0yMtJMUVFDpixR4V5xqXj/qe5NEd7GE6oyAS9tOonbhcV6923SyAGfjOkOWxsJvFykGvUoS7hsIRE1PLGxsYiJiUFOTg4effRRAEBKSgo++OADroJGRERENWIxSY1du3ZhzJgxuHv3LlxdXSGRSMTHJBIJkxpkFiMSj9Z439uFxRi5+lcjRkNEZN0mTJiAoqIivPfee3jnnXcAAC1atMCqVat4niciIqIasZikxsyZMzFhwgQsXLgQTk5O5g6HiIiITGDKlCmYMmUKcnJy4OjoCGdnZ3OHRERERFbMYpIaN2/exLRp05jQIIuybXIvBPm5apWnXstD1LrjVe6/PronQgM9tMrTbilq1QuEiMgaPfroo9ixYwcaN24MLy8vsVyhUGD48OHYt2+fGaMjIiJjyVYoselYOsaEBZhksnxT10/WxaxLupY3aNAg/Pbbb+YOg0iDzN4WTg52WrfwNl7wdZNBomc/CQBfNxnC23jp3F9mb1uXL4OIyCLs378fxcXacxEplUocPHjQDBERmU62QokP915EtkJplfUT1UZ2QRGWp1zSWDXQmuon62LWnhpJSUni30OGDMFrr72GtLQ0dO7cGfb29hrbDhs2rK7DI9JSfknX6N4ttFY5URP+//FzGQoAXNKViBq2P//8U/w7LS0NmZmZ4n2VSoXk5GQ0bdrUHKERmYz6R9eAIG/TXKk2cf1EdUVVJoh/p17LQ3gbL9ja6Lt0SKTNrEmN4cOHa5UtWLBAq0wikUCl4moRZH6GLOmqVj7hwSVdiaghCw4OhkQigUQiEVc9Kc/R0RErVqwwQ2REda/8BRI1VZmAs7fy8c+9Erg72aOjn5vWjzpeIKH6KPlMBuKSzor3o9Ydh6+bDHERQRjcydeMkZE1MWtSo+KyrUSWbkxYAAYEeWuUGdoQISJqqK5duwZBENCyZUukpqZqzKfh4OAAuVwOW1sOy6OGoToXSMrjBRKqb5LPZGDKxpMQKpRn5isxZeNJrBrbnYkNMojFTBRKZA3krjKdV0m6+jeu+2CIiKxE8+bNAfBiBhGgeYHkyOVcvUNZAeDNx9qjd2tPALxAQtZJWaLCveJSrXJVmYC4pLNaCQ3g/jBuCYB5SWl4sLWnzqEoyhL24qf/WExSY9q0aWjdujWmTZumUf7xxx/j8uXLSEhIME9g1CDJXaSY3q+NyRoQpq6fiMiSpaWlIT09XWvSUM6fRfVRxR91zjI7OMvsoCoTMHHDNb37SQCsPXwdo8MCxB915evhjzqyBjVd7U8AkKlQovO8PcYNiOoli0lqbN++XWPiULXevXtj0aJFTGpQnZK7ykzaxdPU9RMRWaKrV6/iySefxOnTpyGRSCAI96/RSST3f7Bx/iyqj/ijjojItCwmqXH79m24ublplbu6uiI3N9cMEREREZExTZ8+HYGBgUhJSUFgYCBSU1Nx+/ZtzJw5E0uXLjV3eEREZGTbJvdCkJ+rVnnqtTxErTte5f7ro3siNNBDqzztlqLGCUOqfywmqdG6dWskJycjJiZGo/zHH39Ey5YtzRQVERERGcvRo0exb98+eHp6wsbGBjY2NnjooYcQHx+PadOm4ffffzd3iERGxx911JDJ7G3h5KD9kzO8jRd83WTIzFfqnFdDAsDHTaZ3eVeZPSeXpv9YTFIjNjYWMTExyMnJEZd7S0lJwQcffMChJ0RERPWASqWCi4sLAMDT0xO3bt1Cu3bt0Lx5c1y4cMHM0RGZRsUfdeolXd2dHODp7IDcu8V69/V0doC7kwOu5hRqLenKH3VkzWxtJIiLCMKUjSchATQSG+oURlxEkM6EBlFFFpPUmDBhAoqKivDee+/hnXfeAQC0aNECq1atQmRkpJmjIyIiotrq1KkTTp06hcDAQISFheH999+Hg4MDVq9ezV6Z1GBUZ0nX3LvFeGLlYQBc0pXqn8GdfLFqbHfEJZ1FlqJILPdxkyEuIojLuZLBLCapAQBTpkzBlClTkJOTA0dHRzg7O5s7JCIiIjKSt956C4WFhQCABQsWYOjQoQgPD0eTJk2wdetWM0dHVDfKL+kK3F/WdfXBqxo9NjydHfBCeEtxOVeAS7pS/TS4ky8ebO0pToi7Prqn3iEnRPpYTFLj0UcfxY4dO9C4cWN4eXmJ5QqFAsOHD8e+ffvMGB0RERHV1qBBg8S/W7dujfPnzyMvLw/u7u7iCihE9Z3cVaYxjKRTUzdMDG+J1Gt5yC5QQu4iQ2igB3/UUYNR/rNuKZ/9bIUSm46lY0xYgMb/V2upv6GxmKTG/v37tdarBwClUomDBw+aISIiIiIyNQ8P7QkQiRoaWxsJerVqYu4wiIxG7iLF9H5tTNbDyNT1ZxcUYXnKJQwI8jZNUsPE9Tc0Zk9q/Pnnn+LfaWlpyMzMFO+rVCokJyejadOm5giNiIiIjGDChAkGbbd27VoTR0JUd6z9Rx1RbchdZSadA8bU9VeHquy/aU5Tr+Vx+IwZmD2pERwcDIlEAolEIq56Up6joyNWrFhhhsiIiIjIGNavX4/mzZujW7duEARdi/cR1T8N6UcdUUOVfCYDcUlnxftR647DlxOd1jmzJzWuXbsGQRDQsmVLpKamasyn4eDgALlcDltbLllFRERkraZMmYKvvvoK165dQ3R0NMaOHWu0YScrV67EkiVLkJmZia5du2LFihUIDQ3Vue2OHTuwcOFCXL58GSUlJWjTpg1mzpyJcePGGSUWIiJqOJLPZGDKxpOomKrPzFdiysaTWDW2OxMbdcTG3AE0b94cLVq0QFlZGXr06IHmzZuLN19fXyY0iIiIrNzKlSuRkZGB119/Hbt27YK/vz+effZZ7N69u1Y9N7Zu3YrY2FjExcXh5MmT6Nq1KwYNGoTs7Gyd23t4eGDOnDk4evQo/vzzT0RHRyM6Ohq7d++ucQxERFR/KUtUuFdcqnUrUJYgLumsVkIDgFg2LykNBcoSnfsrS1R1+TLqPbP31KgoLS0N6enpWpOGDhs2zEwRERERUW1JpVKMHj0ao0ePxl9//YX169fjpZdeQmlpKc6ePVujZdyXLVuGSZMmITo6GgCQmJiI77//HmvXrsWsWbO0tu/bt6/G/enTp2PDhg04dOiQxsos5RUVFaGoqEi8r1Aoqh0nERFZpxGJR2u0nwAgU6EUl6ol07KYpMbVq1fx5JNP4vTp05BIJOKVG/USbyoVs1lERET1gY2NjXiur+n5vbi4GCdOnMDs2bM16u3fvz+OHq26ESoIAvbt24cLFy5g8eLFereLj4/H/PnzaxQjERERmZ7FJDWmT5+OwMBApKSkIDAwEKmpqbh9+zZmzpyJpUuXmjs8IiIyEq793jAVFRVhx44dWLt2LQ4dOoShQ4fi448/xuDBg2FjU/3RsLm5uVCpVPD29tYo9/b2xvnz5/Xul5+fj6ZNm6KoqAi2trb45JNPMGDAAL3bz549G7GxseJ9hUIBf3//asdLRETWZ9vkXgjyc9UqT72Wh6h1x6vcf310T4QGas8hlXZLUeNeIKTNYpIaR48exb59++Dp6QkbGxvY2NjgoYceQnx8PKZNm4bff//d3CESETVo2QolsguKqt6wArmLVCO5wLXfG56XXnoJW7Zsgb+/PyZMmICvvvoKnp6eZonFxcUFf/zxB+7evYuUlBTExsaiZcuWWkNT1KRSKaRSLplJRNQQyext4eSg/ZM5vI0XfN1kyMxX6pxXQwLAx02md3lXmT3njTQmi0lqqFQquLi4AAA8PT1x69YttGvXDs2bN8eFCxfMHB0REW06lo7lKZeqvd/0fm1qvOwg136vHxITExEQEICWLVviwIEDOHDggM7tduzYYXCdnp6esLW1RVZWlkZ5VlYWfHx89O5nY2OD1q1bA7i/rPy5c+cQHx+vN6lBRERUka2NBHERQZiy8SQkgEZiQ91KiYsIYpuljlhMUqNTp044deoUAgMDERYWhvfffx8ODg5YvXo1WrZsae7wiIgavDFhARgQpNnVX1miErtPbpvcS+eVB7lLza5yc+33+iMyMlKcI8tYHBwcEBISgpSUFAwfPhwAUFZWhpSUFMTExBhcT1lZmcZEoERERIYY3MkXq8Z2R1zSWWQp/juP+LCtUucsJqnx1ltvobCwEACwYMECDB06FOHh4WjSpAm2bt1q5uiIiEjuKtMaznGvuFT8O8jPVWcXzZrg2u/1y/r1601Sb2xsLMaPH48ePXogNDQUCQkJKCwsFFdDiYyMRNOmTREfHw/g/qSfPXr0QKtWrVBUVIQffvgBX375JVatWmWS+IiIqH4b3MkXD7b2FFc5WR/dk71KzcBikhrll1Jr3bo1zp8/j7y8PLi7uxv96g4REZmfeu33ilRlQqVrv0twf+33B1t76mw0cO33hmPkyJHIycnB3LlzkZmZieDgYCQnJ4uTh6anp2tMQlpYWIiXXnoJN27cgKOjI9q3b4+NGzdi5MiR5noJRERk5cq3RUIDPZjQMAOLSWro4uGhPVMsERHVD1z7nYwhJiZG73CT/fv3a9x/99138e6779ZBVEREZM3kLlJM79emxkNozV1/Q2P2pMaECRMM2m7t2rUmjoSIiIiIiIgaOrmrrMaTnFtC/Q2N2ZMa69evR/PmzdGtWzcIgq7OxkREVB9x7XciIiIiqi2zJzWmTJmCr776CteuXUN0dDTGjh3LYSdERA0A134nIiIiotqyqXoT01q5ciUyMjLw+uuvY9euXfD398ezzz6L3bt3s+cGEVEDpF77HfhvrXc1rv1OREREROWZPakBAFKpFKNHj8bevXuRlpaGjh074qWXXkKLFi1w9+5dc4dHRER1TL32u9xVcwItHzcZl3MlIiIiIpHZh59UZGNjA4lEAkEQoFJxWT4iooaKa78TERERUVUsoqdGUVERvvrqKwwYMABt27bF6dOn8fHHHyM9PR3Ozs7mDo+IiMyEa78TERERUWXM3lPjpZdewpYtW+Dv748JEybgq6++gqenp7nDIiIiE+Ha70RERERkLGZPaiQmJiIgIAAtW7bEgQMHcODAAZ3b7dixo44jIyIiU+Da70RERERkLGZPakRGRkIiYXdiIiIiIiIiIqoesyc11q9fb+4QiIiIiIiIiMgKWcREoURERERERERE1cWkBhERERERERFZJSY1iIiIiIiIiMgqMalBRERERERERFaJSQ2yStkKJT7cexHZCqVV1k9ERERERES1x6QGWaXsgiIsT7mE7IIiq6yfiIiIiIiIas/sS7qa0sqVK7FkyRJkZmaia9euWLFiBUJDQ3Vu+9lnn+GLL77AmTNnAAAhISFYuHCh3u3JMmUrlFqJCFWZgLO38vHPvRK4O9mjo58bbG0kGtvIXaSQu8rqMlQiIiIiIiKqpXqb1Ni6dStiY2ORmJiIsLAwJCQkYNCgQbhw4QLkcrnW9vv378fo0aPRu3dvyGQyLF68GAMHDsTZs2fRtGlTM7wCqolNx9KxPOVStfeb3q8NZgxoa4KIiIiIiIiIyFTqbVJj2bJlmDRpEqKjowEAiYmJ+P7777F27VrMmjVLa/tNmzZp3P/888+xfft2pKSkIDIysk5iptobExaAAUHeAIAjl3Ox8Mfzerd987H26N3aE8D9nhpERERERERkXeplUqO4uBgnTpzA7NmzxTIbGxv0798fR48eNaiOe/fuoaSkBB4eHjofLyoqQlHRf8McFApF7YKmGlGWqHCvuFS87yyzg7PMDqoyARM3XNO7nwTA2sPXMTosQByKUr4eZYnKZDETERERERGRcdTLpEZubi5UKhW8vb01yr29vXH+vP4r9+W98cYb8PPzQ//+/XU+Hh8fj/nz59c6VqqdEYmGJakqEgBkKpToPG+PcQMiIiIiIiKiOsPVT3RYtGgRtmzZgp07d0Im0z155OzZs5Gfny/e/v777zqOkoiIiIiIiKhhq5c9NTw9PWFra4usrCyN8qysLPj4+FS679KlS7Fo0SL89NNP6NKli97tpFIppFLOw2Bu2yb3QpCfq1Z56rU8RK07XuX+66N7IjRQe4hR2i1FjXuBEBERERERUd2olz01HBwcEBISgpSUFLGsrKwMKSkp6NWrl9793n//fbzzzjtITk5Gjx496iJUqiWZvS2cHOzE211lKa7mFMLdyQGezg6V7uvp7AB3JwdczSnEXWWpRj0ye9s6egVERERERERUU/WypwYAxMbGYvz48ejRowdCQ0ORkJCAwsJCcTWUyMhING3aFPHx8QCAxYsXY+7cudi8eTNatGiBzMxMAICzszOcnZ3N9jqoeqqzpGvu3WI8sfIwAC7pSlRetkKJTcfSMSYsAHJX3UPwLLl+IiIiImo46m1SY+TIkcjJycHcuXORmZmJ4OBgJCcni5OHpqenw8bmv44qq1atQnFxMUaMGKFRT1xcHObNm1eXoVMtlF/SFbi/rOvqg1eRe7dYLPN0dsAL4S3F5VwBLulKVF52QRGWp1zCgCBv0yQ1TFw/ERERETUc9TapAQAxMTGIiYnR+dj+/fs17l+/ft30AZHJyV1lGj+SOjV1w8Twlki9lofsAiXkLjKEBnqIy7gSUe2oygTx79RreQhv48X/X0RERERUZ+p1UoMIAGxtJOjVqom5wyCqd5LPZCAu6ax4P2rdcfi6yRAXEYTBnXzNGBkRERERNRT1cqJQqv/kLlJM79fGZMNGTF0/kbVLPpOBKRtPIktRpFGema/ElI0nkXwmw0yREREREVFDwp4aZJXkrjKTTuxp6vqJrIGyRIV7xaVa5aoyAXFJZyHo2EcAIAEwLykND7b21DkURVmiMnqsRERERNQwMalBREQ6jUg8WqP9BACZCiU6z9tj3IAsBFeHISIiMi+ei6k8JjWIiIiqgavDEBER1Uy2QonsAs2hq+V7cKbdUkBmb6u1n9xFqnFO5LmYymNSg4iIdNo2uReC/Fy1ylOv5SFq3fEq918f3ROhgR5a5Wm3FDXuBWJtuDoMERHRfzYdS8fylEt6H9fXPpjer02Nh4bzXFz/MalBREQ6yext4eSgfZoIb+MFXzcZMvOVOufVkADwcZPpbTTougJTH3F1GCIiIk1jwgIwIMi72vvVdPJ+nosbBiY1iIioWmxtJIiLCMKUjSchATQSG+oURlxEUIO+CqJeHaZi0ke9Osyqsd3ZmCIiogZH7iqrs+EcPBc3HFzSlYiIqm1wJ1+sGtsdclfNKyc+brIG00hQrw5T8VagLKl0dRjg/uowBcoSnftzdRgiIiLD8FxMAHtqEBFRDQ3u5IsHW3uKq5ysj+5Zq3Gqxpo8rK5wdRgiIiLz4rmYACY1iIioFsonMEIDPWo15MQck4cRERERkXVjUoOIiDTIXaSY3q9NjSflqmn9dT15WG1xdRgiIiLz4rmYACY1iIioArmrzKQ9H/TVX5eThxkDV4chIiIyL56LCeBEoUREREalXh0G+G81GDWuDkNERGR6PBc3LExqEBERGRlXhyEiIjIvnosbDg4/ISIiMgFjrw5DRMaVrVBi07F0jAkLMMnQN1PXT0RV47m4YWBPDSIiIhMx5uowRGRc2QVFWJ5ySWspaWupn4gMw3Nx/ceeGkRERNVgrtVhiMh8VGUCUq/lIbtACbmLjD+MiMyM52Iqj0kNIiKiajDX6jBEVDeyFUqN3hVHLudi9cGryL1bLJZ5OjvghfCW6N3aUyyTu0g5zISojvBcTOUxqUFERERWa+XKlViyZAkyMzPRtWtXrFixAqGhoTq3/eyzz/DFF1/gzJkzAICQkBAsXLhQ7/bUMG06lo7lKZcq3Sb3bjEW/nheo2x6vzb8EUREZAZMahAREZFV2rp1K2JjY5GYmIiwsDAkJCRg0KBBuHDhAuRyudb2+/fvx+jRo9G7d2/IZDIsXrwYAwcOxNmzZ9G0aVMzvAKyBMoSFe4Vl4r3n+reFOFtPKEqE/DSppO4XVisd98mjRzwyZjusLWRwMtFqlGPskRl0riJiOg+JjWIiIjIKi1btgyTJk1CdHQ0ACAxMRHff/891q5di1mzZmltv2nTJo37n3/+ObZv346UlBRERkbqfI6ioiIUFf03FEGhUBjxFZAlGJF4tMb73i4sxsjVvxoxGiIiqi6ufkJERERWp7i4GCdOnED//v3FMhsbG/Tv3x9Hjxr2I/XevXsoKSmBh4eH3m3i4+Ph5uYm3vz9/WsdOxERERkPe2oQERGR1cnNzYVKpYK3t7dGube3N86fP69nL01vvPEG/Pz8NBIjFc2ePRuxsbHifYVCwcRGPbNtci8E+bmK93MKipBTUIQzN/Mxb1dalfvPiwhCp6Zu8HKRwqvcSglptxS16gVCRESGYVKDiIiIGpxFixZhy5Yt2L9/P2Qy/StWSKVSSKVc0q8+k9nbwsnhvybxjpNXq5wotDx14qPiRKEye1vjBUlERHoxqUFERERWx9PTE7a2tsjKytIoz8rKgo+PT6X7Ll26FIsWLcJPP/2ELl26mDJMskJjwgIwIOh+D6Ajl3O1Vjkp783H2ovLuspdmPwiIjIHJjWIiIjI6jg4OCAkJAQpKSkYPnw4AKCsrAwpKSmIiYnRu9/777+P9957D7t370aPHj3qKFqyJnJXGeSu93vvdGrqhoAmTpi/Kw0Z+UpxG183GeIigjC4k6+5wiQiov/HpAYRERFZpdjYWIwfPx49evRAaGgoEhISUFhYKK6GEhkZiaZNmyI+Ph4AsHjxYsydOxebN29GixYtkJmZCQBwdnaGs7Oz2V4HWbbBnXwxIMgHqdfykF2ghNxFhtBAD9jaSMwdGhERgUkNIiIislIjR45ETk4O5s6di8zMTAQHByM5OVmcPDQ9PR02Nv8t9LZq1SoUFxdjxIgRGvXExcVh3rx5dRk6WRlbGwl6tWpi7jCIiEgHJjWIiIjIasXExOgdbrJ//36N+9evXzd9QGQ15C5STO/XxmRzYZi6fiIiuo9JDSIiIiJqcOSuMo3VSqytfiIius+m6k2IiIiIiIiIiCwPkxpEREREREREZJWY1CAiIiIiIiIiq8SkBhERERERERFZJU4USkREBslWKJFdUKRRpixRiX+n3VJAZm+rtZ/cRQq5q8zk8RERERFRw8OkBhERGWTTsXQsT7mk9/ERiUd1lk/v14YrABARERGRSTCpQUREBhkTFoABQd7V3k/uIjVBNERERERETGqQiWQrlNh0LB1jwgJM0u3c1PUTkTa5q4z/34iIiIjIonCiUDKJ7IIiLE+5pDX+3lrqJyIiIiIiIsvHnhpkdqoyAanX8pBdoITcRYbQQA/Y2kjMHRYRERERERFZOCY1yKySz2Rg/q40ZOQrxTJfNxniIoIwuJOvGSMjIiIiIiIiS8ekBtWp8ktCHrmci4U/ntfaJiNfickbT+LNx9qjd2tPAFwSkoiIiIiIiLQxqUEmpSxR4V5xqXh//ZHr+GT/FYP2LZ/weKlvK8Q82lqjXiIiIiIiImrYmNQgkxqReNQo9Xyy/4rByRAiIiIiIiJqGLj6CRERERERERFZJfbUIJPaNrkXgvxcxfs5BUXIKSjCmZv5mLcrrcr950UEoVNTN3i5SOHlIhXL024pjNYLhIiIiIiIiKwTkxpkUjJ7Wzg5/Pcx23HyKpanXDJ4f3XiY3q/NpgxoK1GvURERERERNSwMalBdWpMWAAGBHkD0L/6iVrF1U+IiIiIiIiIymNSg+qU3FUmLs3aqakbApo4Yf6uNGTkK8VtfN1kiIsIwuBOvuYKk4iIiIiIiKwAkxpkVoM7+WJAkA9Sr+Uhu0AJuYsMoYEesLWRmDs0IiIiIiIisnBMapDZ2dpI0KtVE3OHQURERERERFaGS7qSSchdpJjer43J5sIwdf1ERERERERk+dhTg0xC7irTWK3E2uonIiIiIiIiy8eeGkRERERERERkldhTg4iIiIiIiKgK2QolsguKNMqUJSrx77RbCsjsbbX2k7tIxRUgyfiY1CAiIiIiIiKqwqZj6Viecknv4yMSj+osn96vDYfOmxCTGkRERERERERVGBMWgAFB3tXej4sbmBaTGkRERERkcbIVSmw6lo4xYQEm6bZt6vqJqP6Ru8r4fWGBOFEoEREREVmc7IIiLE+5pDV+3VrqJyKiusGeGkRERERk1VRlAlKv5SG7QAm5iwyhgR6wtZGYOywiIqoDTGoQERERkdVKPpOB+bvSkJH/f+3de1xU1fo/8M9wmQGEAZU7iaCYZipeIfSYnUTNPBanjkc9WKhEWVKo3VCPopVSPy21jmE30dJSM7WOpmkopnlLFBVRvItHRTNTbiLCPN8//LFjCyYCw8yGz/v1mlfO2muvefau8dk9e89aRUqbj6sDEga2xSPtfCwYGRER1QUWNYiIiIhIM8ovqbjt2CVMX3u4Qp/zV4swatEeTOjfBt2D3AFwSUUiovqKRQ0iIiIislpFN0pRWFyivF+w7RQ+TD1epX3LFzxeeKglYh8OUo1LRETax6IGEREREVmtf8zbXivjfJh6vMrFECIi0g6ufkJEREREREREmsQnNYiIiIjIai0fFYa2vkbl/X82HqvWExe3/vwk81xurT0FQkREllNvixpz587FjBkzkJOTg+DgYHzwwQcICQmptO/BgwcxefJkpKWl4fTp05g1axbGjBlTtwETEZGmlZ+8sEz53+xnnsuFg71thf04eSHRn3Owt4WT/o9L1uHdA/Boex+UmgTRC3/Bpfzi2+7r7qzHZ1HdYGujg6eLQTVOZd9HIiLSnnpZ1Fi6dCnGjRuHefPmITQ0FLNnz0a/fv2QlZUFT0/PCv0LCwvRokULDBo0CGPHjrVAxEREpHWLd2ZjTsrR226/3R3huN6tMLbPveYKi6je8TQ6KIXAtyLa4flFewAAUq6P7v//862Idghu5lan8RGR5fAGQ8NUL4sa7733HmJiYjBixAgAwLx587BmzRrMnz8f8fHxFfp369YN3bp1A4BKtxMREd1JZKg/+rT1uuv9PF0MZoiGqGF4pJ0PkoZ1xtT/ZuL81SKl3dvVAQkD2+KRdj4WjI6I6hpvMDRM9a6oUVxcjLS0NIwfP15ps7GxQXh4OLZvr73fTV6/fh3Xr/9RBczNza21sYmISHvK3z0morrzSDsf9GnrjV0nL+NiXhE8XRwQEtgEtja6O+9MRPUKbzA0TPWuqHHp0iWUlpbCy0v9H7OXlxcOHz58m73uXmJiIqZOnVpr4xEREVHDcjG3CIt3ZiMy1N8sBTFzj29NbG10CGvZ1NJhEJGF8QZDw8QlXatp/PjxuHr1qvI6c+aMpUMiIiIiDbmYdx1zUo5W+P23VsY3N08XA+J6tzLbHVRzj09ERHWj3hU13N3dYWtriwsXLqjaL1y4AG9v71r7HIPBAKPRqHppxcXcIszacAQXc4vu3NmKP4OIiKihKDUJth//Dd+mn8X247+h1CR33snMtJ7rPY0OGNvnXt7VJSLSuHpX1NDr9ejSpQtSUlKUNpPJhJSUFISFhVkwMutRF3dutH53iIiIyFqsyziPv7yzEUM/2YG4JekY+skO/OWdjViXcd6icfFJEyIisgb1bk4NABg3bhyioqLQtWtXhISEYPbs2SgoKFBWQ3n66afh5+eHxMREADcnF83MzFT+fPbsWaSnp8PZ2RlBQUEWOw5rUWoSTr5FRERUR8ovSbjt2CVMX1txTrDzV4swatEeTOjfBt2D3AFY95KEvJYgIiJzqZdFjcGDB+PXX3/F5MmTkZOTg44dO2LdunXK5KHZ2dmwsfnjIZVz586hU6dOyvuZM2di5syZ6NWrF1JTU+s6fKuyLuN8hWXSfLhMGhERUa0pulGKwuIS5f2CbafwYerxKu1bvuDxwkMtEfvwHzdjim6U1l6QNcBrCSIiMqd6WdQAgNjYWMTGxla67dZCRUBAAEQs/9tUa7Mu4zyeX7QHt56ZnKtFeH7RHiQN68yLESIisqi5c+dixowZyMnJQXBwMD744AOEhIRU2vfgwYOYPHky0tLScPr0acyaNQtjxoyp24Ar8Y95tbPk/Iepx6tcDKkrvJYgIiJzq7dFDbqzW+8MAcCvedfxa951lJoEE1dmVLgIAaC0TVyZgcZOetja6ODhYoBHudnDreXuEBER1V9Lly7FuHHjMG/ePISGhmL27Nno168fsrKy4OnpWaF/YWEhWrRogUGDBmHs2LEWiLh+uvV6gtcSRERUl1jUsELmXlf+ckExgJrfGfqtoBiDP95RGyERERHdtffeew8xMTHKnFnz5s3DmjVrMH/+fMTHx1fo361bN3Tr1g0AKt1uKctHhaGt7x+rqJUVBTLOXsWU/2becf8pA9uinZ9rhaJA5rncWnsK5M/U5DN4LUFERDXFooYVKpvtu09bL7MWNYiIiLSquLgYaWlpGD9+vNJmY2OD8PBwbN9ee/8jf/36dVy//sfqG7m5ubU2dhkHe1s46f+4JFux5wTmpByt8v5lhY+43q0wts+9SnvZ0xOXC8yzegivJ4iIyBqwqKFxNZlN/P892R5/C/ZVtWnt7hARETVMly5dQmlpqTIJeBkvLy8cPlxxtZDqSkxMxNSpU2ttvKqIDPVHn7Y3j+t2q5+UuXX1k/IuF9xQ/bO2lRU1br2e4LUEERHVJRY1NKyms4nr7dR3hoDauzvkYG9b5TGIiIis1fjx4zFu3DjlfW5uLpo1a2bWz/Q0OihParbzc4V/Uyezrx5Sk5skt15P8FqCiIjqEosaGlIX69bX1t0hIiIic3J3d4etrS0uXLigar9w4QK8vb1r7XMMBgMMBsvmuEfa+aBPW+9qFx3upLaXXOW1BBER1SUWNayYudatLy65/Wzilrg7REREdLf0ej26dOmClJQUREREAABMJhNSUlJuu6S7ltna6BDWsmmtj2uOJVd5LUFERHWJRQ0rZg3r1pv77hAREVF1jRs3DlFRUejatStCQkIwe/ZsFBQUKKuhPP300/Dz80NiYiKAm5OLZmZmKn8+e/Ys0tPT4ezsjKCgoNt+jrl4uhgQ17uV2Z5QaNJIDwBopLc1y5Krf3aTpDxeSxARkTmxqEF3ZK67Q0RERDUxePBg/Prrr5g8eTJycnLQsWNHrFu3Tpk8NDs7GzY2Nkr/c+fOoVOnTsr7mTNnYubMmejVqxdSU1PrOnx4Gh1Uc0jUtrKiRswXadUeo7aWXOW1BBERmQuLGlbMXOvW/3zsEmI+T0OTRvZmiRsw/90nIiIiAIiNjb3tz01uLVQEBARApLLnEqgmzHU9wWsJIiKqChY1rJi51q33cXUEADRpZL6LBHPffSIiIqKqufUmyX82HqvWz1JvnaOrbMlVc11P8FqCiIiqgkUNDeFs4kRERHS3br1JMrx7AB5t74NSkyB64S+4lF98233dnfX4LKobbG108HQxqMbhkqtERGQNWNTQEM4mTkRERDVV/nrirYh2eH7RHgBQTRhaNoXnWxHtENzMrU7jIyIiuhssamgYZxMnIiKimniknQ+ShnWucJPEmzdJiIhII1jU0DjOJk5EREQ1wZskRESkZSxqWCFzz/bN2cSJiIioPN4kISIirWJRwwqZe7ZvziZORERU//EmCRERNQQsahARERHVQ7xJQkREDYGNpQMgIiIiIiIiIqoOFjWIiIiIiIiISJNY1CAiIiIiIiIiTWJRg4iIiIiIiIg0iUUNIiIiIiIiItIkFjWIiIiIiIiISJNY1CAiIiIiIiIiTWJRg4iIiIiIiIg0iUUNIiIiIiIiItIkO0sHUF+ICAAgNzfXwpEQERHdVJaTynIU1RzzPRERWRPmehY1ak1eXh4AoFmzZhaOhIiISC0vLw+urq6WDqNeYL4nIiJr1JBzvU4ackmnFplMJpw7dw4uLi7Q6XQ1Gis3NxfNmjXDmTNnYDQaaynCuqX1Y2D8lsX4LYvxW1Ztxi8iyMvLg6+vL2xs+IvT2sB8/wfGb1mM37IYv2Ux/j8w1/NJjVpjY2ODe+65p1bHNBqNmvySlqf1Y2D8lsX4LYvxW1Ztxd9Q79qYC/N9RYzfshi/ZTF+y2L8NzX0XN8wSzlEREREREREpHksahARERERERGRJrGoYYUMBgMSEhJgMBgsHUq1af0YGL9lMX7LYvyWpfX4qeq0/u+a8VsW47csxm9ZjJ/K40ShRERERERERKRJfFKDiIiIiIiIiDSJRQ0iIiIiIiIi0iQWNYiIiIiIiIhIk1jUICIiIiIiIiJNYlHDCs2dOxcBAQFwcHBAaGgodu3aZemQAAA//fQTBg4cCF9fX+h0OqxatUq1XUQwefJk+Pj4wNHREeHh4Th69Kiqz+XLlxEZGQmj0Qg3NzdER0cjPz/f7LEnJiaiW7ducHFxgaenJyIiIpCVlaXqU1RUhNGjR6Np06ZwdnbGk08+iQsXLqj6ZGdnY8CAAXBycoKnpydeffVVlJSUmD1+AEhKSkKHDh1gNBphNBoRFhaGtWvXaib+8t5++23odDqMGTNGabP2+KdMmQKdTqd6tWnTRjPxA8DZs2cxbNgwNG3aFI6Ojmjfvj12796tbLfm73BAQECF86/T6TB69GgA1n/+S0tLMWnSJAQGBsLR0REtW7bEm2++ifJzdVvz+afax1xvHlrP98z1zPU1xVzPXN8gCVmVJUuWiF6vl/nz58vBgwclJiZG3Nzc5MKFC5YOTb7//nuZOHGirFixQgDIypUrVdvffvttcXV1lVWrVsm+ffvksccek8DAQLl27ZrS55FHHpHg4GDZsWOHbNmyRYKCgmTo0KFmj71fv36SnJwsGRkZkp6eLo8++qj4+/tLfn6+0mfUqFHSrFkzSUlJkd27d8sDDzwg3bt3V7aXlJRIu3btJDw8XPbu3Svff/+9uLu7y/jx480ev4jId999J2vWrJEjR45IVlaWTJgwQezt7SUjI0MT8ZfZtWuXBAQESIcOHSQuLk5pt/b4ExIS5P7775fz588rr19//VUz8V++fFmaN28uw4cPl507d8qJEyfkhx9+kGPHjil9rPk7fPHiRdW537BhgwCQTZs2iYj1n/9p06ZJ06ZNZfXq1XLy5En5+uuvxdnZWebMmaP0sebzT7WLud58tJ7vmeuZ62uCuZ65vqFiUcPKhISEyOjRo5X3paWl4uvrK4mJiRaMqqJbL3RMJpN4e3vLjBkzlLYrV66IwWCQr776SkREMjMzBYD88ssvSp+1a9eKTqeTs2fP1lnsIjf/0gQgmzdvVmK1t7eXr7/+Wulz6NAhASDbt28XkZsXejY2NpKTk6P0SUpKEqPRKNevX6/T+Ms0btxYPv30U83En5eXJ61atZINGzZIr169lAsdLcSfkJAgwcHBlW7TQvyvv/66/OUvf7ntdq19h+Pi4qRly5ZiMpk0cf4HDBggI0eOVLU98cQTEhkZKSLaO/9UM8z1dac+5Hvmeub6qmKuZ65vqPjzEytSXFyMtLQ0hIeHK202NjYIDw/H9u3bLRjZnZ08eRI5OTmq2F1dXREaGqrEvn37dri5uaFr165Kn/DwcNjY2GDnzp11Gu/Vq1cBAE2aNAEApKWl4caNG6r427RpA39/f1X87du3h5eXl9KnX79+yM3NxcGDB+sw+puPty1ZsgQFBQUICwvTTPyjR4/GgAEDVHEC2jn/R48eha+vL1q0aIHIyEhkZ2drJv7vvvsOXbt2xaBBg+Dp6YlOnTrhk08+UbZr6TtcXFyMRYsWYeTIkdDpdJo4/927d0dKSgqOHDkCANi3bx+2bt2K/v37A9DW+aeaYa6vW1rO98z1lomfud46vsPM9dbxd6hW2Fk6APrDpUuXUFpaqvoiAoCXlxcOHz5soaiqJicnBwAqjb1sW05ODjw9PVXb7ezs0KRJE6VPXTCZTBgzZgx69OiBdu3aKbHp9Xq4ubmp+t4af2XHV7atLhw4cABhYWEoKiqCs7MzVq5cibZt2yI9Pd3q41+yZAn27NmDX375pcI2LZz/0NBQLFiwAK1bt8b58+cxdepU9OzZExkZGZqI/8SJE0hKSsK4ceMwYcIE/PLLL3jppZeg1+sRFRWlqe/wqlWrcOXKFQwfPlyJy9rPf3x8PHJzc9GmTRvY2tqitLQU06ZNQ2RkpCoGLZx/qhnm+rqj1XzPXM9cX13M9cz1DRWLGtTgjB49GhkZGdi6daulQ7lrrVu3Rnp6Oq5evYrly5cjKioKmzdvtnRYd3TmzBnExcVhw4YNcHBwsHQ41VJWZQeADh06IDQ0FM2bN8eyZcvg6OhowciqxmQyoWvXrpg+fToAoFOnTsjIyMC8efMQFRVl4ejuzmeffYb+/fvD19fX0qFU2bJly7B48WJ8+eWXuP/++5Geno4xY8bA19dXc+efSCu0mu+Z6y2Hud56MNfT3eDPT6yIu7s7bG1tK8zie+HCBXh7e1soqqopi+/PYvf29sbFixdV20tKSnD58uU6O77Y2FisXr0amzZtwj333KO0e3t7o7i4GFeuXFH1vzX+yo6vbFtd0Ov1CAoKQpcuXZCYmIjg4GDMmTPH6uNPS0vDxYsX0blzZ9jZ2cHOzg6bN2/G+++/Dzs7O3h5eVl1/JVxc3PDvffei2PHjln9+QcAHx8ftG3bVtV23333KY/VauU7fPr0afz444945plnlDYtnP9XX30V8fHxGDJkCNq3b4+nnnoKY8eORWJioioGaz//VHPM9XVDy/meuZ65vrqY65nrGyoWNayIXq9Hly5dkJKSorSZTCakpKQgLCzMgpHdWWBgILy9vVWx5+bmYufOnUrsYWFhuHLlCtLS0pQ+GzduhMlkQmhoqFnjExHExsZi5cqV2LhxIwIDA1Xbu3TpAnt7e1X8WVlZyM7OVsV/4MAB1V80GzZsgNForJBA6orJZML169etPv7evXvjwIEDSE9PV15du3ZFZGSk8mdrjr8y+fn5OH78OHx8fKz+/ANAjx49KixreOTIETRv3hyA9X+HyyQnJ8PT0xMDBgxQ2rRw/gsLC2Fjo065tra2MJlMALRz/qnmmOvNqz7me+Z65vqqYq5nrm+wLD1TKaktWbJEDAaDLFiwQDIzM+XZZ58VNzc31Sy+lpKXlyd79+6VvXv3CgB57733ZO/evXL69GkRublEkZubm3z77beyf/9+efzxxytdoqhTp06yc+dO2bp1q7Rq1apOlih6/vnnxdXVVVJTU1VLRRUWFip9Ro0aJf7+/rJx40bZvXu3hIWFSVhYmLK9bJmovn37Snp6uqxbt048PDzqbJmo+Ph42bx5s5w8eVL2798v8fHxotPpZP369ZqI/1blZ0QXsf74X375ZUlNTZWTJ0/Kzz//LOHh4eLu7i4XL17URPy7du0SOzs7mTZtmhw9elQWL14sTk5OsmjRIqWPNX+HRW6uEOHv7y+vv/56hW3Wfv6joqLEz89PWeZtxYoV4u7uLq+99prSx9rPP9Ue5nrz0Xq+Z65nrq8J5nrm+oaKRQ0r9MEHH4i/v7/o9XoJCQmRHTt2WDokERHZtGmTAKjwioqKEpGbyxRNmjRJvLy8xGAwSO/evSUrK0s1xm+//SZDhw4VZ2dnMRqNMmLECMnLyzN77JXFDUCSk5OVPteuXZMXXnhBGjduLE5OTvL3v/9dzp8/rxrn1KlT0r9/f3F0dBR3d3d5+eWX5caNG2aPX0Rk5MiR0rx5c9Hr9eLh4SG9e/dWLnK0EP+tbr3Qsfb4Bw8eLD4+PqLX68XPz08GDx6sWvfd2uMXEfnvf/8r7dq1E4PBIG3atJGPP/5Ytd2av8MiIj/88IMAqBCTiPWf/9zcXImLixN/f39xcHCQFi1ayMSJE1VLzFn7+afaxVxvHlrP98z1zPU1xVzPXN8Q6URE6uqpECIiIiIiIiKi2sI5NYiIiIiIiIhIk1jUICIiIiIiIiJNYlGDiIiIiIiIiDSJRQ0iIiIiIiIi0iQWNYiIiIiIiIhIk1jUICIiIiIiIiJNYlGDiIiIiIiIiDSJRQ0iIiIiIiIi0iQWNYioQZg0aRKeffZZi31+cXExAgICsHv3bovFQEREVN9kZWXB29sbeXl5FoshPj4eL774osU+n6ihY1GDSKOGDx8OnU4HnU4He3t7BAYG4rXXXkNRUZHF44qIiKjQnpqaCp1OhytXrlRr3AULFijHe7vXqVOnKt03JycHc+bMwcSJEyu0x8XFISgoCA4ODvDy8kKPHj2QlJSEwsLCKsX17rvvonHjxpWe98LCQhiNRrz//vvQ6/V45ZVX8Prrr9/1sRMREQG1k7eq68/yeEBAAGbPnl3tse+U36dMmXLbfcePH48XX3wRLi4uSpuI4JNPPkFYWBiMRiOcnZ1x//33Iy4uDseOHatSTBcuXIC9vT2WLFlS6fbo6Gh07twZAPDKK69g4cKFOHHiRNUPmohqDYsaRBr2yCOP4Pz58zhx4gRmzZqFjz76CAkJCZYOyywGDx6M8+fPK6+wsDDExMSo2po1a1bpvp9++im6d++O5s2bK20nTpxAp06dsH79ekyfPh179+7F9u3b8dprr2H16tX48ccfqxTXU089hYKCAqxYsaLCtuXLl6O4uBjDhg0DAERGRmLr1q04ePBgNc4AERE1ZLWVt6xR+Vw+e/ZsGI1GVdsrr7xS6X7Z2dlYvXo1hg8frrSJCP71r3/hpZdewqOPPor169cjMzMTn332GRwcHPDWW29VKSYvLy8MGDAA8+fPr7CtoKAAy5YtQ3R0NADA3d0d/fr1Q1JS0t0fPBHVGIsaRBpmMBjg7e2NZs2aISIiAuHh4diwYYOy3WQyITExEYGBgXB0dERwcDCWL1+uGuPgwYP429/+BqPRCBcXF/Ts2RPHjx9X9n/jjTdwzz33wGAwoGPHjli3bp3Zjken0yEpKQn9+/eHo6MjWrRoocTr6OgIb29v5aXX6+Hk5KRqs7W1rXTcJUuWYODAgaq2F154AXZ2dti9ezf++c9/4r777kOLFi3w+OOPY82aNar+V65cwTPPPAMPDw8YjUY8/PDD2LdvHwDA09MTAwcOrPSiZ/78+YiIiECTJk0AAI0bN0aPHj1ue9eHiIjodmorbwHAlClT0LFjR3zxxRcICAiAq6srhgwZYrafcJQ9xTl16lQlplGjRqG4uBgAVLnc1dUVOp1O1ebs7FzpuMuWLUNwcDD8/PyUtqVLl2LJkiVYunQpJk2ahAceeAD+/v544IEH8M477yA5OVk1xqeffor77rsPDg4OaNOmDT788ENlW3R0NFJSUpCdna3a5+uvv0ZJSQkiIyOVtoEDBzK/E1kIixpE9URGRga2bdsGvV6vtCUmJuLzzz/HvHnzcPDgQYwdOxbDhg3D5s2bAQBnz57Fgw8+CIPBgI0bNyItLQ0jR45ESUkJAGDOnDl49913MXPmTOzfvx/9+vXDY489hqNHj5rtOCZNmoQnn3wS+/btQ2RkJIYMGYJDhw5Ve7zLly8jMzMTXbt2Vdp+++03rF+/HqNHj0ajRo0q3U+n0yl/HjRoEC5evIi1a9ciLS0NnTt3Ru/evXH58mUANy96Nm7ciNOnTyv7nDhxAj/99JNyF6dMSEgItmzZUu3jISKihqe28xYAHD9+HKtWrcLq1auxevVqbN68GW+//bbZjiElJQWHDh1CamoqvvrqK6xYsQJTp06t0ZhbtmxR5XcA+Oqrr9C6dWs89thjle5T/jwtXrwYkydPxrRp03Do0CFMnz4dkyZNwsKFCwEAjz76KLy8vLBgwQLVGMnJyXjiiSfg5uamtIWEhOB///vfbX8KS0RmJESkSVFRUWJrayuNGjUSg8EgAMTGxkaWL18uIiJFRUXi5OQk27ZtU+0XHR0tQ4cOFRGR8ePHS2BgoBQXF1f6Gb6+vjJt2jRVW7du3eSFF16oUlzlXw4ODgJAfv/999vuC0BGjRqlagsNDZXnn3++Qt9evXpJXFzcbccqs3fvXgEg2dnZStuOHTsEgKxYsULVt2nTpkq8r732moiIbNmyRYxGoxQVFan6tmzZUj766CMRESkpKRE/Pz9JSEhQtk+aNEn8/f2ltLRUtd+cOXMkICDgjnETERGVqe28lZCQIE5OTpKbm6tsf/XVVyU0NPS2MWzatEkAVMjvjRo1Ep1OJ7NmzbrtvlFRUdKkSRMpKChQ2pKSksTZ2blCnkxOThZXV9c/PR9lgoOD5Y033lC1tWnTRh577DFVW1xcnBKrn5+f0t6yZUv58ssvVX3ffPNNCQsLU97Hx8dLYGCgmEwmERE5duyY6HQ6+fHHH1X7Xb16VQBIampqlWInotpjZ7lyChHV1F//+lckJSWhoKAAs2bNgp2dHZ588kkAwLFjx1BYWIg+ffqo9ikuLkanTp0AAOnp6ejZsyfs7e0rjJ2bm4tz586hR48eqvYePXqoHmH9s7jK27lzpzK3xJ8JCwur8D49Pf2O+93OtWvXAAAODg537Ltr1y6YTCZERkbi+vXrAIB9+/YhPz8fTZs2rTBu2c90bG1tERUVhQULFiAhIQEigoULF2LEiBGwsVE/EOfo6Gj2ydyIiKhhqG7eAm5O7ll+ck0fHx9cvHjxjp+5ZcsW1X4A8NBDD91xv+DgYDg5OSnvw8LCkJ+fjzNnzqjmvLob165dq1J+nzhxImJjY7FixQpMnz4dwM15MY4fP47o6GjExMQofUtKSuDq6qq8HzlyJN5++21s2rQJDz/8MJKTkxEQEICHH35Y9RmOjo4AwBxPZAEsahBpWKNGjRAUFATg5vwNwcHB+OyzzxAdHY38/HwAwJo1a1S/NQVuzsUB/JGAzRlXmf/9739m+aw7cXd3BwD8/vvv8PDwAAAEBQVBp9MhKytL1bdFixYA1OclPz8fPj4+SE1NrTB2+cdOR44cicTERGzcuBEmkwlnzpzBiBEjKuxz+fJlJQ4iIqKqMEfeuvWGhk6ng8lkumMsgYGBqnEAwM7OMv9L4e7ujt9//13V1qpVqwrnycPDAx4eHvD09FTayq6TPvnkE4SGhqr6l5+jq1WrVujZsyeSk5Px0EMP4fPPP0dMTIzqZywAlJ/2MMcT1T3OqUFUT9jY2GDChAn497//jWvXrqFt27YwGAzIzs5GUFCQ6lW2SkiHDh2wZcsW3Lhxo8J4RqMRvr6++Pnnn1XtP//8M9q2bWu249ixY0eF9/fdd1+1x2vZsiWMRiMyMzOVtqZNm6JPnz74z3/+g4KCgj/dv3PnzsjJyYGdnV2F81hWMCn7nF69emH+/PlITk5GeHh4pXeeMjIylCdliIiIqsIceauu7du3T3l6EriZ352dnW+7cllVdOrUSZXfAWDo0KHIysrCt99++6f7enl5wdfXFydOnKhwngIDA1V9o6Oj8c033+Cbb77B2bNnVautlMnIyIC9vT3uv//+ah8PEVUPixpE9cigQYNga2uLuXPnwsXFBa+88grGjh2LhQsX4vjx49izZw8++OADZQKs2NhY5ObmYsiQIdi9ezeOHj2KL774QrnD8eqrr+Kdd97B0qVLkZWVhfj4eKSnpyMuLs5sx/D1119j/vz5OHLkCBISErBr1y7ExsZWezwbGxuEh4dj69atqvYPP/wQJSUl6Nq1K5YuXYpDhw4hKysLixYtwuHDh5W7NOHh4QgLC0NERATWr1+PU6dOYdu2bZg4cSJ2796tGjM6OhorVqzAypUrK0wQWmbLli3o27dvtY+HiIgaJnPkrbpUXFyM6OhoZGZm4vvvv0dCQgJiY2Mr/EzzbvTr1w/bt29HaWmp0jZkyBD84x//wJAhQ/DGG29g586dOHXqFDZv3oylS5eqnsKYOnUqEhMT8f777+PIkSM4cOAAkpOT8d5776k+Z9CgQbC3t8dzzz2Hvn37VlqI2bJlC3r27Gm2p2CJ6E9YelIPIqqeqKgoefzxxyu0JyYmioeHh+Tn54vJZJLZs2dL69atxd7eXjw8PKRfv36yefNmpf++ffukb9++4uTkJC4uLtKzZ085fvy4iIiUlpbKlClTxM/PT+zt7SU4OFjWrl1brbjKJhi700Shc+fOlT59+ojBYJCAgABZunRppX2rOlGoiMj3338vfn5+FSYjO3funMTGxkpgYKDY29uLs7OzhISEyIwZM1STmeXm5sqLL74ovr6+Ym9vL82aNZPIyEjV5KMiIoWFheLq6ipNmjSpMEGbiMi2bdvEzc1NCgsLqxQ3ERFRebWVtxISEiQ4OFg19qxZs6R58+a3/ew/y+PNmze/40Shjz/+uEyePFmaNm0qzs7OEhMTU2muvJuJQm/cuCG+vr6ybt06VXtpaanMmzdPQkNDpVGjRqLX66VFixYSExMjmZmZqr6LFy+Wjh07il6vl8aNG8uDDz5YYUJWEZFnn31WAMiyZcsqjaV169by1VdfVSluIqpdOhERy5ZViIhu0ul0WLlyJSIiImp1XBFBaGgoxo4di6FDh9bq2Hdj8ODBCA4OxoQJEywWAxERUV0bPnw4rly5glWrVtX62HPnzsV3332HH374odbHrqq1a9fi5Zdfxv79+y02vwhRQ8afnxBRvafT6fDxxx+jpKTEYjEUFxejffv2GDt2rMViICIiqm+ee+45PPjgg8jLy7NYDAUFBUhOTmZBg8hC+KQGEVkNcz2pQURERJZjzic1iIhY1CAiIiIiIiIiTeLPT4iIiIiIiIhIk1jUICIiIiIiIiJNYlGDiIiIiIiIiDSJRQ0iIiIiIiIi0iQWNYiIiIiIiIhIk1jUICIiIiIiIiJNYlGDiIiIiIiIiDSJRQ0iIiIiIiIi0qT/A1VFeMkCjxXPAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "dp_cut = 0.5\n", "LUT_pred_spanet, LUT_target_spanet = parse_resolved_w_target(test_h5, s_h5, dp_cut=dp_cut, fjs_reco=None)\n", "eff_s, efferr_s = calc_eff(LUT_pred_spanet, bins)\n", "pur_s, purerr_s = calc_pur(LUT_target_spanet, bins)\n", "\n", "fig, ax = plt.subplots(1, 2, figsize=(12, 5))\n", "ax[0].errorbar(x=bin_centers, y=eff_s, xerr=xerr, yerr=efferr_s, fmt='o', capsize=5, label='Trust if a fatjet is assigned')\n", "ax[0].set(xlabel=r\"Reco H pT (GeV)\", ylabel=r\"Matching efficiency\", title=f\"Matching Efficiency vs. Reco H pT, DP cut at {dp_cut}, AP cut at 1/(13*13)\")\n", "\n", "ax[1].errorbar(x=bin_centers, y=pur_s, xerr=xerr, yerr=purerr_s, fmt='o', capsize=5, label='Trust if a fatjet is assigned')\n", "ax[1].set(xlabel=r\"Gen H pT (GeV)\", ylabel=r\"Matching purity\", title=f\"Matching Purity vs. Gen H pT, DP cut at {dp_cut}, AP cut at 1/(13*13)\")" ] }, { "cell_type": "markdown", "id": "ef22b104", "metadata": {}, "source": [ "### Gen-level Overlap Removal Study\n", "Below is some code intended to provide some info about how often a resolved Higgs and a boosted Higgs overlap\n", "\n", "* N1: number of resolved H without overlapping removal\n", "* N2: number of resolved H with overlapping removal\n", "* N3: number of Higgs that are bothe resolved and boosted\n", "\n", "N1-N2-N3 = Number of Higgs that were removed - Number of Higgs that are both resolved and boosted" ] }, { "cell_type": "markdown", "id": "f3d94ea3", "metadata": {}, "source": [ "#### Getting N1" ] }, { "cell_type": "code", "execution_count": 27, "id": "8585ac1c", "metadata": {}, "outputs": [], "source": [ "def get_unmasked_arr(arr, mask):\n", " return ak.drop_none(arr.mask[mask])" ] }, { "cell_type": "code", "execution_count": 28, "id": "107d34b0", "metadata": {}, "outputs": [], "source": [ "testfile = test_h5\n", "# N1\n", "h1_pt = ak.Array(testfile['TARGETS']['h1']['pt'])\n", "h2_pt = ak.Array(testfile['TARGETS']['h2']['pt'])\n", "h3_pt = ak.Array(testfile['TARGETS']['h3']['pt'])\n", "\n", "# mask\n", "h1_mask = ak.Array(testfile['TARGETS']['h1']['mask'])\n", "h2_mask = ak.Array(testfile['TARGETS']['h2']['mask'])\n", "h3_mask = ak.Array(testfile['TARGETS']['h3']['mask'])\n" ] }, { "cell_type": "code", "execution_count": 29, "id": "d9d53d14", "metadata": {}, "outputs": [], "source": [ "N1 = ak.count(get_unmasked_arr(h1_pt, h1_mask))+ak.count(get_unmasked_arr(h2_pt, h2_mask))+ak.count(get_unmasked_arr(h3_pt, h3_mask))" ] }, { "cell_type": "markdown", "id": "3db03a68", "metadata": {}, "source": [ "#### Getting N2" ] }, { "cell_type": "code", "execution_count": 30, "id": "25a07b92", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[False, False, False], [False, False, False], ..., [False, False, False]]\n" ] } ], "source": [ "# assume perfect reconstruction\n", "_, _, fjs_reco = parse_boosted_w_target(test_h5, test_h5)" ] }, { "cell_type": "code", "execution_count": 31, "id": "29893773", "metadata": {}, "outputs": [], "source": [ "# target assignments\n", "b1_h1_t = np.array(testfile[\"TARGETS\"][\"h1\"]['b1']).astype('int')\n", "b1_h2_t = np.array(testfile[\"TARGETS\"][\"h2\"]['b1']).astype('int')\n", "b1_h3_t = np.array(testfile[\"TARGETS\"][\"h3\"]['b1']).astype('int')\n", "\n", "b2_h1_t = np.array(testfile[\"TARGETS\"][\"h1\"]['b2']).astype('int')\n", "b2_h2_t = np.array(testfile[\"TARGETS\"][\"h2\"]['b2']).astype('int')\n", "b2_h3_t = np.array(testfile[\"TARGETS\"][\"h3\"]['b2']).astype('int')\n", "\n", "\n", "b1_ts = np.concatenate((b1_h1_t.reshape(-1, 1), b1_h2_t.reshape(-1, 1), b1_h3_t.reshape(-1, 1)), axis=1)\n", "b1_ts = ak.Array(b1_ts)\n", "\n", "b2_ts = np.concatenate((b2_h1_t.reshape(-1, 1), b2_h2_t.reshape(-1, 1), b2_h3_t.reshape(-1, 1)), axis=1)\n", "b2_ts = ak.Array(b2_ts)\n", "\n", "# reconstruct jet 4-momentum objects\n", "j_pt = np.array(testfile['INPUTS']['Jets']['pt'])\n", "j_eta = np.array(testfile['INPUTS']['Jets']['eta'])\n", "j_phi = np.array(testfile['INPUTS']['Jets']['phi'])\n", "j_mass = np.array(testfile['INPUTS']['Jets']['mass'])\n", "js = ak.zip(\n", " {\n", " \"pt\": j_pt,\n", " \"eta\": j_eta,\n", " \"phi\": j_phi,\n", " \"mass\": j_mass,\n", " },\n", " with_name=\"Momentum4D\"\n", ")\n", " \n", "goodJetIdx = get_unoverlapped_jet_index(fjs_reco, js, dR_min=0.8)\n" ] }, { "cell_type": "code", "execution_count": 32, "id": "c76afb0d", "metadata": {}, "outputs": [], "source": [ "N2 = 0\n", "for tb1_e, tb2_e, goodJetIdx_e in zip(b1_ts, b2_ts, goodJetIdx):\n", " for b1Hx, b2Hx in zip(tb1_e, tb2_e):\n", " if (b1Hx in goodJetIdx_e) & (b2Hx in goodJetIdx_e):\n", " N2 += 1" ] }, { "cell_type": "code", "execution_count": 33, "id": "15cc3c00", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "81810 73866\n" ] } ], "source": [ "print(N1, N2)" ] }, { "cell_type": "markdown", "id": "78f8ee7b", "metadata": {}, "source": [ "#### Getting N3" ] }, { "cell_type": "code", "execution_count": 34, "id": "aade4838", "metadata": {}, "outputs": [], "source": [ "# mask\n", "# mask\n", "h1_mask = np.array(testfile['TARGETS']['h1']['mask'])\n", "h2_mask = np.array(testfile['TARGETS']['h2']['mask'])\n", "h3_mask = np.array(testfile['TARGETS']['h3']['mask'])\n", "h_masks = np.concatenate((h1_mask.reshape(-1, 1), h2_mask.reshape(-1, 1), h3_mask.reshape(-1, 1)), axis=1)\n", "\n", "bh1_mask = np.array(testfile['TARGETS']['bh1']['mask'])\n", "bh2_mask = np.array(testfile['TARGETS']['bh2']['mask'])\n", "bh3_mask = np.array(testfile['TARGETS']['bh3']['mask'])\n", "\n", "bh_masks = np.concatenate((bh1_mask.reshape(-1, 1), bh2_mask.reshape(-1, 1), bh3_mask.reshape(-1, 1)), axis=1)" ] }, { "cell_type": "code", "execution_count": 35, "id": "cce02929", "metadata": {}, "outputs": [], "source": [ "N3=np.sum(h_masks&bh_masks)" ] }, { "cell_type": "markdown", "id": "449fd705", "metadata": {}, "source": [ "#### Summary" ] }, { "cell_type": "code", "execution_count": 36, "id": "c1542006", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.017112822393350446" ] }, "execution_count": 36, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# this percents of resolved Higgs should not be removed\n", "(N1-N2-N3)/N1" ] }, { "cell_type": "code", "execution_count": 37, "id": "62add9d4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "7944" ] }, "execution_count": 37, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# this number of resolved Higgs removed for overlapping with a boosted Higgs\n", "N1-N2" ] }, { "cell_type": "code", "execution_count": 38, "id": "b3cb315b", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "6544" ] }, "execution_count": 38, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# this number of resolved Higgs can also be reconstructed as a boosted Higgs\n", "N3" ] }, { "cell_type": "code", "execution_count": 39, "id": "a38cdea2", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.07999022124434665" ] }, "execution_count": 39, "metadata": {}, "output_type": "execute_result" } ], "source": [ "N3/N1" ] }, { "cell_type": "markdown", "id": "aac2b02f", "metadata": {}, "source": [ "How good is the correlation between removal and bi-category?" ] }, { "cell_type": "code", "execution_count": 40, "id": "58036d92", "metadata": {}, "outputs": [], "source": [ "N_overlap_bicat = 0\n", "N_overlap_resolve = 0\n", "N_unoverlap_bicat = 0\n", "N_unoverlap_resolve = 0\n", "\n", "bAndr = h_masks&bh_masks\n", "for tb1_e, tb2_e, goodJetIdx_e, h_masks_e, bh_masks_e, in zip(b1_ts, b2_ts, goodJetIdx, h_masks, bh_masks):\n", " for b1Hx, b2Hx, resolve, boost in zip(tb1_e, tb2_e, h_masks_e, bh_masks_e):\n", " \n", " if (b1Hx in goodJetIdx_e) & (b2Hx in goodJetIdx_e):\n", " overlap = False\n", " else:\n", " overlap=True\n", " \n", " if resolve & boost:\n", " bicat = True\n", " else:\n", " bicat = False\n", " \n", " if overlap & bicat:\n", " N_overlap_bicat += 1\n", " elif ~overlap & bicat:\n", " N_unoverlap_bicat += 1\n", " elif overlap & resolve:\n", " N_overlap_resolve += 1\n", " elif ~overlap & resolve:\n", " N_unoverlap_resolve += 1" ] }, { "cell_type": "code", "execution_count": 41, "id": "350bb51e", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " Resolved+Boosted | Resolved only\n", "Removed: 6544 1400\n", "Not Removed: 0 73866\n" ] } ], "source": [ "print(\" \", \"Resolved+Boosted |\", \"Resolved only\")\n", "print(\"Removed: \", N_overlap_bicat, \" \",N_overlap_resolve)\n", "print(\"Not Removed: \", N_unoverlap_bicat, \" \", N_unoverlap_resolve)" ] }, { "cell_type": "markdown", "id": "eaeb4aab", "metadata": {}, "source": [ "### Reco-level Overlap Removal Study" ] } ], "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.9.10" } }, "nbformat": 4, "nbformat_minor": 5 }