retrieve_bulk.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "2021-08-13 17:29:46.477 INFO    rdkit: Enabling RDKit 2021.03.4 jupyter extensions\n"
     ]
    }
   ],
   "source": [
    "import streamlit as st\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import re\n",
    "from PIL import Image\n",
    "import webbrowser\n",
    "import json\n",
    "import pickle\n",
    "import sys \n",
    "import joblib\n",
    "\n",
    "sys.path.append('./CC/')\n",
    "\n",
    "import chemaxon\n",
    "from chemaxon import *\n",
    "from compound import Compound\n",
    "from compound_cacher import CompoundCacher\n",
    "from rdkit.Chem import rdChemReactions as Reactions\n",
    "from rdkit.Chem import Draw\n",
    "from rdkit import Chem"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def load_smiles():\n",
    "    db = pd.read_csv('./data/cache_compounds_20160818.csv',\n",
    "                     index_col='compound_id')\n",
    "    db_smiles = db['smiles_pH7'].to_dict()\n",
    "    return db_smiles\n",
    "\n",
    "def load_molsig_rad1():\n",
    "    molecular_signature_r1 = json.load(open('./data/decompose_vector_ac.json'))\n",
    "    return molecular_signature_r1\n",
    "\n",
    "\n",
    "def load_molsig_rad2():\n",
    "    molecular_signature_r2 = json.load(\n",
    "        open('./data/decompose_vector_ac_r2_py3_indent_modified_manual.json'))\n",
    "    return molecular_signature_r2\n",
    "\n",
    "\n",
    "def load_model():\n",
    "    filename = './model/M12_model_BR.pkl'\n",
    "    loaded_model = joblib.load(open(filename, 'rb'))\n",
    "    return loaded_model\n",
    "\n",
    "\n",
    "def load_compound_cache():\n",
    "    ccache = CompoundCacher()\n",
    "    return ccache\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def count_substructures(radius, molecule):\n",
    "    \"\"\"Helper function for get the information of molecular signature of a\n",
    "    metabolite. The relaxed signature requires the number of each substructure\n",
    "    to construct a matrix for each molecule.\n",
    "    Parameters\n",
    "    ----------\n",
    "    radius : int\n",
    "        the radius is bond-distance that defines how many neighbor atoms should\n",
    "        be considered in a reaction center.\n",
    "    molecule : Molecule\n",
    "        a molecule object create by RDkit (e.g. Chem.MolFromInchi(inchi_code)\n",
    "        or Chem.MolToSmiles(smiles_code))\n",
    "    Returns\n",
    "    -------\n",
    "    dict\n",
    "        dictionary of molecular signature for a molecule,\n",
    "        {smiles: molecular_signature}\n",
    "    \"\"\"\n",
    "    m = molecule\n",
    "    smi_count = dict()\n",
    "    atomList = [atom for atom in m.GetAtoms()]\n",
    "\n",
    "    for i in range(len(atomList)):\n",
    "        env = Chem.FindAtomEnvironmentOfRadiusN(m, radius, i)\n",
    "        atoms = set()\n",
    "        for bidx in env:\n",
    "            atoms.add(m.GetBondWithIdx(bidx).GetBeginAtomIdx())\n",
    "            atoms.add(m.GetBondWithIdx(bidx).GetEndAtomIdx())\n",
    "\n",
    "        # only one atom is in this environment, such as O in H2O\n",
    "        if len(atoms) == 0:\n",
    "            atoms = {i}\n",
    "\n",
    "        smi = Chem.MolFragmentToSmiles(m, atomsToUse=list(atoms),\n",
    "                                       bondsToUse=env, canonical=True)\n",
    "\n",
    "        if smi in smi_count:\n",
    "            smi_count[smi] = smi_count[smi] + 1\n",
    "        else:\n",
    "            smi_count[smi] = 1\n",
    "    return smi_count\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def decompse_novel_mets_rad1(novel_smiles, radius=1):\n",
    "    decompose_vector = dict()\n",
    "\n",
    "    for cid, smiles_pH7 in novel_smiles.items():\n",
    "        mol = Chem.MolFromSmiles(smiles_pH7)\n",
    "        mol = Chem.RemoveHs(mol)\n",
    "        # Chem.RemoveStereochemistry(mol)\n",
    "        smi_count = count_substructures(radius, mol)\n",
    "        decompose_vector[cid] = smi_count\n",
    "    return decompose_vector\n",
    "\n",
    "\n",
    "def decompse_novel_mets_rad2(novel_smiles, radius=2):\n",
    "    decompose_vector = dict()\n",
    "\n",
    "    for cid, smiles_pH7 in novel_smiles.items():\n",
    "        mol = Chem.MolFromSmiles(smiles_pH7)\n",
    "        mol = Chem.RemoveHs(mol)\n",
    "        # Chem.RemoveStereochemistry(mol)\n",
    "        smi_count = count_substructures(radius, mol)\n",
    "        decompose_vector[cid] = smi_count\n",
    "    return decompose_vector\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_reaction_formula_side(s):\n",
    "    \"\"\"\n",
    "        Parses the side formula, e.g. '2 C00001 + C00002 + 3 C00003'\n",
    "        Ignores stoichiometry.\n",
    "\n",
    "        Returns:\n",
    "            The set of CIDs.\n",
    "    \"\"\"\n",
    "    if s.strip() == \"null\":\n",
    "        return {}\n",
    "\n",
    "    compound_bag = {}\n",
    "    for member in re.split('\\s+\\+\\s+', s):\n",
    "        tokens = member.split(None, 1)\n",
    "        if len(tokens) == 0:\n",
    "            continue\n",
    "        if len(tokens) == 1:\n",
    "            amount = 1\n",
    "            key = member\n",
    "        else:\n",
    "            amount = float(tokens[0])\n",
    "            key = tokens[1]\n",
    "\n",
    "        compound_bag[key] = compound_bag.get(key, 0) + amount\n",
    "\n",
    "    return compound_bag\n",
    "\n",
    "\n",
    "def parse_formula(formula, arrow='<=>', rid=None):\n",
    "    \"\"\"\n",
    "        Parses a two-sided formula such as: 2 C00001 => C00002 + C00003\n",
    "\n",
    "        Return:\n",
    "            The set of substrates, products and the direction of the reaction\n",
    "    \"\"\"\n",
    "    tokens = formula.split(arrow)\n",
    "    if len(tokens) < 2:\n",
    "        print(('Reaction does not contain the arrow sign (%s): %s'\n",
    "               % (arrow, formula)))\n",
    "    if len(tokens) > 2:\n",
    "        print(('Reaction contains more than one arrow sign (%s): %s'\n",
    "               % (arrow, formula)))\n",
    "\n",
    "    left = tokens[0].strip()\n",
    "    right = tokens[1].strip()\n",
    "\n",
    "    sparse_reaction = {}\n",
    "    for cid, count in parse_reaction_formula_side(left).items():\n",
    "        sparse_reaction[cid] = sparse_reaction.get(cid, 0) - count\n",
    "\n",
    "    for cid, count in parse_reaction_formula_side(right).items():\n",
    "        sparse_reaction[cid] = sparse_reaction.get(cid, 0) + count\n",
    "\n",
    "    return sparse_reaction\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def draw_rxn_figure(rxn_dict, db_smiles, novel_smiles):\n",
    "    # db_smiles = load_smiles()\n",
    "\n",
    "    left = ''\n",
    "    right = ''\n",
    "\n",
    "    for met, stoic in rxn_dict.items():\n",
    "        if met == \"C00080\" or met == \"C00282\":\n",
    "            continue  # hydogen is not considered\n",
    "        if stoic > 0:\n",
    "            if met in db_smiles:\n",
    "                right = right + db_smiles[met] + '.'\n",
    "            else:\n",
    "                right = right + novel_smiles[met] + '.'\n",
    "        else:\n",
    "            if met in db_smiles:\n",
    "                left = left + db_smiles[met] + '.'\n",
    "            else:\n",
    "                left = left + novel_smiles[met] + '.'\n",
    "    smarts = left[:-1] + '>>' + right[:-1]\n",
    "    # print smarts\n",
    "    smarts = str(smarts)\n",
    "    rxn = Reactions.ReactionFromSmarts(smarts, useSmiles=True)\n",
    "    return Draw.ReactionToImage(rxn)  # , subImgSize=(400, 400))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_rule(rxn_dict, molsig1, molsig2, novel_decomposed1, novel_decomposed2):\n",
    "    if novel_decomposed1 != None:\n",
    "        for cid in novel_decomposed1:\n",
    "            molsig1[cid] = novel_decomposed1[cid]\n",
    "    if novel_decomposed2 != None:\n",
    "        for cid in novel_decomposed2:\n",
    "            molsig2[cid] = novel_decomposed2[cid]\n",
    "\n",
    "    molsigna_df1 = pd.DataFrame.from_dict(molsig1).fillna(0)\n",
    "    all_mets1 = molsigna_df1.columns.tolist()\n",
    "    all_mets1.append(\"C00080\")\n",
    "    all_mets1.append(\"C00282\")\n",
    "\n",
    "    molsigna_df2 = pd.DataFrame.from_dict(molsig2).fillna(0)\n",
    "    all_mets2 = molsigna_df2.columns.tolist()\n",
    "    all_mets2.append(\"C00080\")\n",
    "    all_mets2.append(\"C00282\")\n",
    "\n",
    "    moieties_r1 = open('./data/group_names_r1.txt')\n",
    "    moieties_r2 = open('./data/group_names_r2_py3_modified_manual.txt')\n",
    "    moie_r1 = moieties_r1.read().splitlines()\n",
    "    moie_r2 = moieties_r2.read().splitlines()\n",
    "\n",
    "    molsigna_df1 = molsigna_df1.reindex(moie_r1)\n",
    "    molsigna_df2 = molsigna_df2.reindex(moie_r2)\n",
    "\n",
    "    rule_df1 = pd.DataFrame(index=molsigna_df1.index)\n",
    "    rule_df2 = pd.DataFrame(index=molsigna_df2.index)\n",
    "    # for rid, value in reaction_dict.items():\n",
    "    #     # skip the reactions with missing metabolites\n",
    "    #     mets = value.keys()\n",
    "    #     flag = False\n",
    "    #     for met in mets:\n",
    "    #         if met not in all_mets:\n",
    "    #             flag = True\n",
    "    #             break\n",
    "    #     if flag: continue\n",
    "\n",
    "    rule_df1['change'] = 0\n",
    "    for met, stoic in rxn_dict.items():\n",
    "        if met == \"C00080\" or met == \"C00282\":\n",
    "            continue  # hydogen is zero\n",
    "        rule_df1['change'] += molsigna_df1[met] * stoic\n",
    "\n",
    "    rule_df2['change'] = 0\n",
    "    for met, stoic in rxn_dict.items():\n",
    "        if met == \"C00080\" or met == \"C00282\":\n",
    "            continue  # hydogen is zero\n",
    "        rule_df2['change'] += molsigna_df2[met] * stoic\n",
    "\n",
    "    rule_vec1 = rule_df1.to_numpy().T\n",
    "    rule_vec2 = rule_df2.to_numpy().T\n",
    "\n",
    "    m1, n1 = rule_vec1.shape\n",
    "    m2, n2 = rule_vec2.shape\n",
    "\n",
    "    zeros1 = np.zeros((m1, 44))\n",
    "    zeros2 = np.zeros((m2, 44))\n",
    "    X1 = np.concatenate((rule_vec1, zeros1), 1)\n",
    "    X2 = np.concatenate((rule_vec2, zeros2), 1)\n",
    "\n",
    "    rule_comb = np.concatenate((X1, X2), 1)\n",
    "\n",
    "    # rule_df_final = {}\n",
    "    # rule_df_final['rad1'] = rule_df1\n",
    "    # rule_df_final['rad2'] = rule_df2\n",
    "    return rule_comb, rule_df1, rule_df2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_ddG0(rxn_dict, pH, I, novel_mets):\n",
    "    ccache = CompoundCacher()\n",
    "    # ddG0 = get_transform_ddG0(rxn_dict, ccache, pH, I, T)\n",
    "    T = 298.15\n",
    "    ddG0_forward = 0\n",
    "    for compound_id, coeff in rxn_dict.items():\n",
    "        if novel_mets != None and compound_id in novel_mets:\n",
    "            comp = novel_mets[compound_id]\n",
    "        else:\n",
    "            comp = ccache.get_compound(compound_id)\n",
    "        ddG0_forward += coeff * comp.transform_pH7(pH, I, T)\n",
    "\n",
    "    return ddG0_forward"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_dG0(rxn_dict, rid, pH, I, loaded_model, molsig_r1, molsig_r2, novel_decomposed_r1, novel_decomposed_r2, novel_mets):\n",
    "\n",
    "    # rule_df = get_rxn_rule(rid)\n",
    "    rule_comb, rule_df1, rule_df2 = get_rule(\n",
    "        rxn_dict, molsig_r1, molsig_r2, novel_decomposed_r1, novel_decomposed_r2)\n",
    "\n",
    "    X = rule_comb\n",
    "\n",
    "    ymean, ystd = loaded_model.predict(X, return_std=True)\n",
    "\n",
    "    result = {}\n",
    "    # result['dG0'] = ymean[0] + get_ddG0(rxn_dict, pH, I)\n",
    "    # result['standard deviation'] = ystd[0]\n",
    "\n",
    "    # result_df = pd.DataFrame([result])\n",
    "    # result_df.style.hide_index()\n",
    "    # return result_df\n",
    "    return ymean[0] + get_ddG0(rxn_dict, pH, I, novel_mets), ystd[0], rule_df1, rule_df2\n",
    "    # return ymean[0],ystd[0]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_novel_molecule(add_info):\n",
    "    result = {}\n",
    "    for cid, InChI in add_info.items():\n",
    "        c = Compound.from_inchi('Test', cid, InChI)\n",
    "        result[cid] = c\n",
    "    return result\n",
    "\n",
    "\n",
    "def parse_novel_smiles(result):\n",
    "    novel_smiles = {}\n",
    "    for cid, c in result.items():\n",
    "        smiles = c.smiles_pH7\n",
    "        novel_smiles[cid] = smiles\n",
    "    return novel_smiles\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "db_smiles = load_smiles()\n",
    "molsig_r1 = load_molsig_rad1()\n",
    "molsig_r2 = load_molsig_rad2()\n",
    "\n",
    "loaded_model = load_model()\n",
    "ccache = load_compound_cache()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimating dG for reaction with novel metabolite"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "rxn_str = 'C01745 + C00004 <=> N00001 + C00003 + C00001'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'C01745 + C00004 <=> N00001 + C00003 + C00001'"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rxn_str"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "add_info = {\"N00001\":\"InChI=1S/C14H12O/c15-14-8-4-7-13(11-14)10-9-12-5-2-1-3-6-12/h1-11,15H/b10-9+\"}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'InChI=1S/C14H12O/c15-14-8-4-7-13(11-14)10-9-12-5-2-1-3-6-12/h1-11,15H/b10-9+'"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "add_info['N00001']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "pH = 7 # any number between 0-14 \n",
    "I = 0.1    #min_value=0.0, max_value=0.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'N00001': 'Oc1cccc(/C=C/c2ccccc2)c1'}\n"
     ]
    }
   ],
   "source": [
    "try:\n",
    "    novel_mets = parse_novel_molecule(add_info)\n",
    "    novel_smiles = parse_novel_smiles(novel_mets)\n",
    "    novel_decomposed_r1 = decompse_novel_mets_rad1(novel_smiles)\n",
    "    novel_decomposed_r2 = decompse_novel_mets_rad2(novel_smiles)\n",
    "\n",
    "except Exception as e:\n",
    "    novel_mets = None\n",
    "    novel_smiles = None\n",
    "    novel_decomposed_r1 = None\n",
    "    novel_decomposed_r2 = None\n",
    "\n",
    "print(novel_smiles)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "rxn_dict = parse_formula(rxn_str)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABLAAAADICAIAAAB3fY8nAAB2lElEQVR4nO3dd1xTZxcH8F/ClOEWRQXrqopb0Kq490DrwC1q3ZO6tS7c0rqwbn2tolYtamtx1y0qLlzgFhVEFJQ9ZeS8fzwxIi4gNwTI+X76xzUkz32SQu49zzhHRkRgjDHGGGOMMaZ75NruAGOMMcYYY4wx7eCAkDHGGGOMMcZ0FAeEjDHGGGOMMaajOCBkjDHGGGOMMR3FASFjjDHGGGOM6SgOCBljjDHGGGNMR3FAyBhjjDHGGGM6igNCxhhjjDHGGNNRHBAyxhhjjDHGmI7igJAxxhhjjDHGdBQHhIwxxhhjjDGmozggZIwxxhhjjDEdxQEhY4wxxhhjjOkoDggZY4wxxhhjTEdxQMgYY4wxxhhjOooDQsYYY4wxxhjTURwQMsYYY4wxxpiO4oCQMcYYY4wxxnQUB4SMMcYYY4wxpqM4IGSMMcYYY4wxHcUBIWOMMcYYY4zpKA4IGWOMMcYYY0xHcUDIGGOMMcYYYzqKA0LGGGOMMcYY01EcEDLGGGOMMcaYjuKAkDHGGGOMMcZ0FAeEjDHGGGOMMaajOCBkjDHGGGOMMR3FASFjjDHGGGOM6SgOCBljjDHGGGNMR3FAyBhjjDHGGGM6igNCxhhjjDHGGNNRHBAyxhhjjDHGmI7igJAxxhhjjDHGdBQHhIwxxhhjjDGmozggZIwxxhhjjDEdxQEhY4wxxhhjjOkoDggZY4wxxhhjTEdxQMgYY4wxxhhjOooDQsYYY4wxxhjTURwQMsYYY4wxxpiO4oCQMcYYY4wxxnQUB4SMMcYYY4wxpqM4IGSMMcYYY4wxHcUBIWOMMcYYY4zpKA4IGWOMMcYYY0xHcUDIGGOMMcYYYzqKA0LGGGOMMcYY01H62u4AYyyHCQ3FiRMwNkbHjjA21nZvGGOM6bZTpxAYiAYNULmytrvCWN7EM4SMsTRev0bfvrC0BIBu3ZCaqu0OpRERgdGjMWwYBg3C06cAQKTtPjHGGNOkiRPh64vKlTFvHk6e1HZvPvbnnxgwAIMGYfNmgC9JLBfjGULGWBp//olx49CiBQDcuIFLl9C4sbb79N6sWejXD/b2ePECI0bg558RGIjevWFuru2eMcYY04CkJNy+jRUrAKBiRQwbhlattN2n9x48wMGD2L0bMhkGD0bp0ggKQrVqaNBA2z1jLNM4IGSMpREejmLFlMcWFnj7Vqu9+di9e7C3BwArKyQloVUr6Olpu0+MMcY0JjoahQopjwsXRnS0VnvzMW9vtG8PmQwAfvwRDx9i/Hgtd4mxrOIlo4yxNGrUwJUryuNr11CjhlZ78zG5HAqF8jg1FXL++mKMsTytaFGEhiIlBQCuX4eNjbY7lIaJCeLjlcdxcTAx0WpvGFOLjHjFM2NMJTUVo0ZBLkdMDOrUwU8/Yf58WFtj4kRt9wxYuRIpKRg3DseP4/BhbNqk7Q4xxhjTsP/+w7p1KF0ar19jzRpcvow1a3DwIPLl03LHwsPh6Ijt22FggP794e6OkiW13CXGsooDQqZzvLxw+TKmTAERfv0VZmYYNAhmZggIwMWL6NtX2/3LCZKSIJPBwABnzqBFC5ib49EjlCihzS49eQKFAteu4fJlVKqEESNgYKDN/jDGmESGD8eCBSheHNu2oXJlBAaiZ08AWLIEv/yi7c7lEImJMDCATIZ69eDjg0WLMGOGlrt07BjKlcPWrVAoMGAAqlbVcn8YUwOvuWI6JyAAW7bg2jUQ4eJFXLyId+8AICICt25puW85xYkTqFYNZ8+ieXM4OCAmBvPmZbGpxETcuSPBXsQRI1C9OvLnx+rVGDuWo0HGWJ5x44Yyurl9Gy9f4vZt5eNnzmixUzlJaChGjUL37pDLsWwZACxZglevstja8+d48EDdHNoHD6J9ezg7Y8kS/PorR4Mst+OAkOmiSZMwZ86Hy8GePdi2DZ6eWu1TjnLrFh49woQJUCiwbBkMDLB5M3x9M93O/fvo3h3nzmHKFGzYkPX+eHjg9GmYm6Nhw6w3whhjOVKpUiheHKdOKf955w62bcO2bXj9WqvdylH++Qf//osTJ9CsGTp3RmwsXFyy0s7o0di4EQcP4scfEReXxc68e4fJkwGgY8cstsBYDsMBIdNFBQqgd2+sX6/8Z7lyqFgRZcpotU85yqRJKFMGt25h+3ZUqoSRI5GampX8ab/+Cjc3jBuHrVvx559ZHJGNj8fUqQDg6ooiRbLSAmOM5WwzZ+LXX5GUBACFC6NiRVSsqP1dcjmFhQWmTQOAKVOQmooVK2BoiD/+wJ07mWvn1i3IZFiyBFOmoGdP7NiRxf4sXYpHj1C1KkaNymILjOUwHBAyHTVgAE6fVh7Xqwd7e9SsCQBhYXjzRov9yhmMjbFoEQDMnIm4OMyZg3JlUTEFT45lrp3AQFSooDwuVSqLn+yiRQgIgK0tBg/OyssZYyzHMzXFqFHYswcASpeGvT3s7VGgAFJT8eQJFzwHJkzAd9/h9m1s24by5TFyJJrb4Zpr5hp5/hwVKyqPK1fGs2dZ6UlQEFxdAWDNGuhz8TaWR3BAqNuuX8fIkRg2DMePA1DGAAC8vHD+vBb7pTm+vihRAsWLQybDwoWoVQu1a8PQEAAKFECFCvD0xNy5iI3Vdke1rm9fNGyI4GAsXYqiRbFjLCzP47/JUKRk6OX//IOtW2FlhadPlY8EB8PAAJs2ZW6eMOwJjnpCLsfatVxngjGWx7x9i5cvlcORXbuiTx9YWKBSJeVPf/gBR4/Czw+rV2uxjzmDsTEWLwaAmTMRHY3Fc9DqEYJ24/HhDL08KAgTJqBkSTx+rHzk4UOULYvduxEQkLme7FqBd+/QuzeaNcvcCxnLwTjLqA4LDUW/fvjnHxgZoU8fzJmDyZPx338AsHOnMmtW3nL7NurXR9u22Lv3a0lJZs3CnDnKKFGnXbwIl0loGo3xJ2BSDOuqIfwxOq6D3VcXydy4gYkTce4czM3h6Ynly9G1q7J+1OPH+P131KqF1avRqFGG+rC7Ex4eRelpGLro209mjLHcgwidOuHqVfz999e+Ee/exfXrGDgwG3uWMxGhiwO+i8aPzdFiPi6vxPGJKFoZo+5A/uUrenw8Vq/GokWIicGSJfD3h7U1ihTB4cNYsgT16kEux9SpmD4dxsbf7sOz09jeErIaGHIEpUpJ+OYY0y4ecddhFy6gc2eYmcHAAIMH49gxxMdjyxZs2YJz57TdOenFxaFPHyQmwtLyi9GgQoENG9CrF0eDAAB7e4ywhuI+Ts2AniFauQLAGRe8i/r884ODMWIE6tXDuXMoXBgzZ6JhQ+zahSpVMGMGxo5FixYoWxa3bqFJEwwY8O0ccf7H8egQjE3Re6zEb40xxrRt2TIcPoyUFFhZffE5N2/i4sW8NzybJTIZ1s5G4YvwXoaoQNQbiyLf4+0D+HyhJi0R9u5FlSqYPh0xMXBwQM+e2LwZ3bujXj0cOAALCzg6IjER8+ahatVvZ5ZTpODYeABo0YejQZbHcECow9JODstkyrpzVauialWULq29bmnK6NG4fx/VqmHFii8+JywMcjm8vREeno09y8la/wp9Y9zegeBrqNINZZog/g28Fqd/WnIcDqxAxYrYtAn6+pgyBU+fYto0GBrC3BwNGih/o378EXfvYu5cGBtjxw5Urgz3tVAkf/7UqUk46gwAzRfAzFKTbzL3uX379qZNmwoWLNi6deu5c+cePHgwIiJC251ijGXCtWuYNQsyGbZu/VpKsxcvAOTVPRyZV7o+qvZESgJOz4Lc4P0w5RwkfHLNfn4GHdqgZ08EBqJePVy4gIMHUa4cANjYwM4OBgYoUQI7d+LcOdSsiadP8eOP6NABoU++eParqxHqi8IV0GCCxt5hrpSamrp48eL+/ftXrVp1wIABmzZtunv3Lq9AzF14yagOCwnBgAH4918YGMDJCdOmYcqUvLpk9I8/MGQITE1x7RqqVNF2b3KXk9Nw8TdYNcRPF/D6Bna0RfP5qDta+VMi3N+H/6Yg6gX+qQSr8nBzQ/ny32gzKAgzZmDnTqyoCf14tHNDhfbpn+O1GKdnolhVjLz5teVAusff379Ro0YhISFpv7319PSqV6/esGHD+vUbNGzY/5v/BxhjWhQZiTp18OwZJk7E8uXa7k3uEv0CayojOQFDL6NUPWxrhnyF0GEtzEsqnxD+GKdm4t5ehNvj7wDMno2hQ7+xBV2hwM6dmDQJ7cuh4k3UHYXmC2CU/6PnxIViTSUkRqLvYVTsoKl3lwsR0ciRIzdt2mRqahqXppJHiRIl6tev37Bhw0aNHGvVKsspc3M4Dgh12/nz2L4dADp1wo8/wsVFWX/83DkoFGjeXLu9k8qjR7CzQ0wM3N3zUpCbXZJisPp7xL5Gz/2o0g0pCdB//73+/CyOT8TrmwBQqh5auqFsg0y0fOUUro1B2EMAqNINbVegQJpx8h1t8PQEBp7Bd80keid5QWhoaOPGjR89etSyZcv//e9/vr6+Fy9evHDhwvXr19+9ewfgu+++f/78YYECqFsX9vZo1AgNG8LEBABWrECbNqhWDRcuwMwMtWpp960wpru6d8fff8PODhcv8g6FzDv1Cy64onQDDL6I1MQPl6SEcJybj2vroEiGoTnsp6PuhEzU7nj7Ft5LcNMNpIB5KbRZiqq9IZMpf/rgH+zvi3Kt0YdrFn9k5syZixcvzpcv37///psvXz4fH5+LFy+eO3cuNDRUPKFu3cM3b3b4/ns0agR7e9jaompVALh9G2fOKAtarV6NceO09hYYOCDUabduoWNHjByJ2bO13RUNio+P79Qp9PTp7wYNwtat2u5NLnV9PZ4cR5ul8N2NUD/I9VCxA+7/gwf/AEABa7RcjGp9P1w4M06RjGvrcGYO3kVDPx/sp8K0KKKCYFoMdYbh1Q2OBtOKiopq1qzZrVu36tWrd+rUKTMzM9WP4uPjfXx8Ll26FBRksnfvuJCQD6/S10etWmjWDHfvwsgI+/djyxYULYquXbXwFhhjW7feGzzYpkAB3LihXMPIMuddNLY1xQ8/o2JH7GiNUvWQ+g4V2uHIWCSEQ6aH2j+h+QKYlchK469u4Og4vLgEAGWaoP54BF2BXB+VOsGkKOQGKGAt7bvJ1dauXTt27FgDA4MDBw506PDRxOnjx48vX77s7e3t57fg0qUiafOLlyqFBg3QsCFWr8aOHbC3R8eOOJyxfLFMQzgg1GHTp+PXXzF2bN5OaD106NB9+/a3aLF7x452pqba7k1ud2oGbLqjRG1sa4r8VnjoiXpj0WQWDM2+/dqviH2FE9NwZyf0jVHSFn088fQUXl5F698k6ndekJCQ0LZtWy8vr6pVq547d65IkSIA3N3dL1682KBBgwYNGlSqVEn2PiYPDsbFi7hwAT4+uHYNSUlo0ADFiqFzZ6SkQC7ngJAx7bh582aDBg1sbbtOnLize3c9bXcnl4t/g2MT0G0nvBbD0BzHx+O7ZmizHCVqqdUsEe7swImpiAuBjSNK/YAa/eDhiAGnoJ+BTKQ6488//xwwYAAR/fHHH4MGDQKQkJDQq1evH374oWHDhnXr1lWNWsbF4eZN+Pjg4kWcOYO3bwHg119hYIBjx3DoELp04YBQyzgg1FVEKFcOz5/DyyujBQCymUKBjRtx4wasrTFxIhISsGcPypSBlRWsrFCkSEba+Ouvv3r37m1sbOzt7V2LV8ip79QMGBeAIgXRQWg6B5BlcQj2swIvIPwx/Pag/3E8OIAXF9F6qWSN53LJydSt24+HDh0sU6bMhQsXSr9P+9S1a9cDBw6I48KFC9evX79ZM/c6dYrWqwdzc+VrY2Nx7RpSU7F6Nf75B506oXFjVKrEASFjmXTnDjZvBhH690f9+ti3D4mJyqtSqVJfq2X0XmxsrJ2d3cOHD8eMGbNmzZps6HIeF/8Gf3VDnWG4tQ1dtyP+DUrUlqzxd1G4tAwWVREfjup9sbsTnE5wQKhy+PDxrl07JScnr1q1ytnZWTx4/vz5pk2bimOxs71z54kVKzo1aPAhtwARHj7E5cswN0d4OMzN8fIlTp/mgFDLOCDUVZcvK3M/BgTk0HrfixahcGGMHInz57FhA4YPR4sWH35qYoIyZWBtHV69+hpz8zJlyli9Z/y+lNCTJ09sbW2jo6M3btw4fPhw7byLPOOhJ267o8j3KGAN60YoVi0rC0Qzwr0FilREwbKoOxqRz1CoHAzNv/2qPI0IQ4bg8WP/J0+anTlzonLlyqof+fj4nD9//tKlS97e3i9fvjQ3zx8fH56aqqenh2rV0LAhOnZEx47KJ//4I/79FzdvomNHrF3LASFjmRESgn79sHcv9PXRqxdWr4ajI27dUv5ULoelpQgOt9vYRBYsqLoqFStWTNWGk5PTzp07q1evfuXKlXycZEMdMcE4PQsV2+P6BjSbi+I1YFRAIyd6cADey2Fpi6o9UdIWIXdQsq5GTpSrXL6MPn1SDA379O5tM0/kngAAhIeHHz9+XKwUvXXrVnJycuPGG728hgOwsECDBmjQACNGoGBBAPD0REgIhg2DoyMiI3HypJbeDAPAAaHumjABbm6YNAnLlmm7K1/QogVOnVJGHS1aYOVKbNyIFy8QEIAXLxAZKZ4V1KCBlbd32tcVL17cysqqU6dOBw4cuHnzZs+ePf/6669s730e8voW/puEZ6cB4PtOaOYCS1vpz5IUq1x3urMt+h8HgJQErLVByju0ckUNJ03Fn7nB+PFYtQr58+Ps2bjatb+47jkwMNDHx//cuebe3rh5E8nJADBiBDZsUD5hzhzMnw8AixejcWM0bpwNfWcsr9izB9HREGOL+/YhJATh4bh/H4GBCAjAq1d4v0eqnqXltTRFVvPly2dtbW1lZdW+fftJkyaZmZldu3Yt7bAOy5zkeFxahou/ITkOxWvCohq67ZT+LERIiYeBKR4cQEywMrH2hSU4NRM1+qPNUpgWl/6kucTt22jWDJGR+PnnJDe3L+ZEEjvbb96sePp0icuXIXa26+khPBz58wOAtzciItChAx4+xIYNWLkyu94A+xx9bXeAaYNCgb17AaBXL2135atUMYBcjmrVsG7dhx/FxIjLcGJ4+MwWLQIDAwMDA1+8eBEUFBQSEhISEpKYmOjn51e+fPnNmzdrpe95QewrnJ2LG1tAqchXGPZTUdIW5qUAgBTw24OqPaQpCEGp2NYUhcqh43oUrqh8MO4NzCwR5I0DA3Hjf+iwGsVrSnCu3GbBAqxaBUND7N2Lr0SDAKytra2trcW8X0ICfHzg7Q07uw9PENHgzp3Yvh1FinBAyFhmpKRA//0tk74+UlI+yseWkoLgYAQGIjCwV1BQvcBA1VUpPDz84cOHjx49unTpEoD169dzNJhFosrRiamIfA4A3zug5WLlMYC3D5AQDquG0pzrznacnoWO62BcEIqUDx3QM8SdHXh0EM3mou4YyHXuLvrpU7Rvj8hIdO2K5cu/liHXxMSkcePGjRtDrCf194e3N54/V0aDABo0AIDQUEyapKy3ybSIZwh10vnzaNoU5crhyZOcO/EyZQpatED79vD3x6RJaN4cvr6wsvqwjdDaGsbpV/MT0evXr0+fPu3k5GRoaHjp0qU6depopfu5W3w8VrtBbw3iXkHPUJk5xrjQhyfs64W7Hmjnhh9+luB019fj8GgUsMaYe1CkQs9QuU+DCHe24+R0xL6GTA92I9F0PkwLS3DGXGLDBowaBT097NkDR0dp2ty6FYMHo317HDkiTYOM6YTnzzF+PPbtg54e+veHkxO2bUPp0rC2Vl6VSpeGhcWnr4uLiwsICBgwYICPj8/QoUN5jDKLzp/HvY0I2QUAJe3QZjnKNPnw00AvuLdEwe8w2g96atfxSIrB6kqIfYVuf6J6X8SFwvT9/9kIfxwbj0eHAMCiGtqtQdmm6p4u9wgNRePGePQILVrgyBEYGUnQZmoqLC3x5g3u3eMy0drEAaEu2jl7du3z58u1aJHPxUXbffmyuDjMm4eICBgYYNYsDBmCY8fSP8fCYnPz5kfevRNzI2LDRpkyZfz9/fv37//dd9+dO3dOG13PzYiwbx+mTsXz55jRFNXM0XYlCldI/7RHh7C7E4wLwfkx8mUowc8XJYRjTSXEv0Wvv1G5Kw4Ox7NT6Lrjw0BvchwuLsWFJUgqgS3JmDoNY8dCL+9n5/vnH/ToAYUCmzZh6FDJmg0LQ4kSkMkQGqrcyMEYy5BDh7B/P2QytG2LQoXQtm36Jxgbhzdr1jM5WVyJrK2tS5cuLS5MQ4YMOXDgwJMnT6ysrLTR9dwsMBCzZmHnTlgWxC+WaPAz6gyB7ONLAKVifU28uYu2K1B/grpnPD4Rl1fCyh4/eeHNXWyui7qj0Wb5hyc8PYmjznh7H5fqw7Ao1qxBmTJfbi6PiIpC8+a4eRP16uHUKZipl1w8rcGDsXUrFi/GL79I1ibLNGI6JiUlpXjx4gBu3bql7b5khpcXbdxIs2bRgAHUrBmVL09GRgQ4fW5ASV9fXy6X58uX7+3bt9rud65y/jzZ2RFAANna0tmzX3vyjjY0F3RsvLonPTSS5oK2tyIiCvaheXJaYEhvH6R/Wuhd+nmQsm92duTtre55c7aYGCpWjABydZW+8aZNCaDdu6VvmWmUQqEIDg7Wdi8YERG9ekV//kmurjR2LHXqRDVrUuHCBDxr+Pkli6VKlQKwaNEibfc7VwkPp4kTydCQADI1pXnzKDb2i09+dJjmglwLUdwbtU4aepfmG9A8PXp9i4hoW3OaCzr6c/qnpSTS8ZVkaqrs2+LFlJio1nlzvJkzCSAbG5L8xurffwmgevUkbpZlCs8Q6pz//vuvbdu2lStXvn//vrb7oh4ivH59KyjIPzDwxYsXAQEBL9579epVs2bNzp49u2TJkunTp2u7oznYnj04exbFi2P8eCxahOXLAaBUKSxejP79v5F+9s09bKgJyDDaD0W+z2IHXtzCVjvI9TDyNopUwh8NEXQZjaaj5ZLPP//gQTg74/lzyGTo3x9Ll6J43tnZn5KCJk1w4AAsLDB9Orp2xYkTmDVL+hO5uWHCBPTujd27pW+cqe/du3cvX74MDg5+9erV0/eCg4MDAgLi4uIWLlzYo0eP77/P6h8d05z4+LeBgT4BAUFBQS9evHj+/Lk4CAwMbNeu3YEDB0qWLPn8+XODDBSo0FEBAVizBvHx6N4dRYqgZUuEhUEux8CBWLgQJUt+4+U728H/OH5wRrtVWe+Dewc8P4q6Y9BhDfz2YH8fmFpg7EMYF/zMk4ODMX06du4EEcqXh5sbHByyfuqc57ffYGqKMWPg6wtvbzx6hPHj8b7skWQSElCsGOLjERgofeMso7QdkbLsNnjwYABz587Vdkc0KDEx8ejRowBKlSqVlJSk7e7kVO7uNGsWJSfTjRvUrh0dOEAmJjRtGkVHZ7QFz6E0F7R+TBY7oFBQw4bUpj4dmkNEdPMPmgtaXorexXztVXFxNGsWGRsTQH360LFj1LMnDR1Ko0dTQkIWe5IzJCVRrVo0eDARUY8eXxsKV9Pz5ySTkbl5nh/RzulevXp15cqVffv2rVix4ueff+7SpYutrW3aKgWfMjc3B1CrVq2EXP7brmsSEhKqV68OYDdPzX9JfDw1a0aBgRQfTz170qVLZGNDzZvTjRsZbSHkDs3To3nW9ORhFvuwdy+VtKSZrSk+jJLiaGUZmgu6seUbrzp9mqpVI4AMDMjPjwYPpmHDqEcPungxi93IMcaOpWbN6PlzunSJ5s/X4Im6dSOA1q7V4CnY13FAqFuSkpIKFy4M4O7du9rui2YpFIoqVaoA8PDw0HZfcqouXSgiQnns6Ehv3lBISOZaiHlFzeoTQGfOZKUD7u4EUPHiFBlJURG01Jrmgnx3Zei1T55Q9+7k50ctWlByMhHR9u20cmVWupFjJCVR9+40dSqdOqXZgJCIatUigI4e1eAp2Fc8ePCgefPmX4r6jIyMKlas2KJFi0GDBrm4uGzZsuXEiRMPHjxISEiIiYkROSpHjRql7TfBMmfDhg0AGjRooO2O5FTnz9Ps2crjixdpxoxMX5KIaP1sypePOnXKSgfi4+m77wigTZuIiI7Op7mgzfVIkfrt1yYnk5sbzZtHv/xCx48TEcXFUaNGpFBkpSc5xtixdOYM9eih8YBQ3A60bq3BU7Cv07mEuTru2LFj4eHhtWrVsrGx0XZfNEsmk40dO3bMmDG///57jx49tN2dHCklBaqVSwYGSElBiRKZa8GsBFo54OxlTJgAH59vLDFNJyZGuX982TIUKIAJE3CAMH4AqvbO0MvLl8e+fbhzB1WrKnPBt2qFSZMy1/8cafZsdO8OExPNnqV79xR9/YhLl261a9das2din4iPj2/YsGFycnLBggVLlixZsmTJcuXKlStXztLSUhx/99138k/+mpKTk8VSQw8Pjx9++GH9+vWNGjXq27evNt4BywonJ6eZM2d6e3tfvXq1Xr162u5OzpOc/NElKTn5s1lbv6HbWEz/HQcP4sQJtM7kl9vixXj+HHXqYPBgPH4Mx1/RqRmG/gpZBi5t+vr4+WcAaNMGolC7iQlKl8br17C0zPS7yEnKlkX16jh06NO07lJycECFCkmpqT6RkVUKcrozbcjMDRzL/Xbv3g2gVw4vPyiRQYMGFS5c+MKFC9euXdN2X3Kk5s2V5SjfvsWrV1ncjDdpEsqUwa1bcHfP3AvnzkVwMBo2RL9+uHcPa9fiRTCaTcxcHZQCBRAVpTyOiMgbSTPNzDB8OI4f1+xZOnW6e/26xaZNTgqFQrNnYp/YsWNHeHh4jRo1IiIi7t69e+LEiY0bN06bNm3AgAGtWrUqVKjQzZs39+7du2rVqunTp/fs2dPOzq5kyZJ9+vQRL69evfry5csBjBgx4uHDh1p9KywTTExMxJaN1atXa7svOVKdOvDyQmIiAHh4oGXLrDRiYYFp0wBgyhSkpmbihU+fYtkyyGRYuxZ6epgwAXFxMKuAUpkM3c3NP1yVIiNRoEDmXp4jTZ2KQ4c0e4rChVGmTIfTpxse4YJI2qLtKUqNCw0NPXbs2OLFix0dHa9evRqShRUIuVlUVJSXl5ebm5uTk5ONjY2pqalcLr9+/bq2+5VNJk2aBGDAgAHa7kiO9O4dzZ5NQ4fSkCH04JOsnhm3cycBVLJkJtY4+vuTgQHp6dHNm0REbdoQQOPGZeXsP/5Ix47R8+fUs2duTz2akkI7diiPp0/X+Aa/cuXKAbh06ZJmT8M+IfaS7dmzR/zz6NGjTk5OTZo0KVu27FfSjTRq1ChtI2JusEaNGryZMBcJCAjQ19c3MDAICgrSdl9ypHPnaOBAGjyYVq/OeiMJCcqVn//7XyZe1aMHATRoEBGRpycBlD8/vXqV6bN7etLw4RQQQPv2Ue6/9zhwQLmz5PRp0vT+m99//x1Ajx49NHsa9gV5MMuoGHP1ee/+/fuq91i4cOGKFSuePn3aRNPrsbTnzZs3N9J4+vRp2p/q6+unpKTY29ufPHnSWKPT/zlDQEBA+fLl9fT0AgICSmR2PaQuuHQJbm5o2BDjx2e9ESLUr4+rV7FtGwYOzOir9u7F/fuYMwceHujVC0WK4OFDFMl8ScOEBOzYgdBQdOiAOnUy/fIcqV07HD8OHx/lG0pNxbt30i8inThx4sqVK6dNm+bq6ipx0+zLTp482bp167TZJn///fefxWIzAEChQoXSLh9VHZctW1aWZv48NjbWzs7u4cOH48aNEzdSeZKo6h6YRkBAQNeuXe3t7evWravt3mWFo6Pj/v3758yZM08sLGTpjByJ8HBs2IDChbPeyO7d6NsXNWrg9u2MviQoCLNm4ddfUbAgatTAo0f4/XeMG5eVs1+8iFOnYG2N/v2VOxpyuX/+Qf/+6NYNO3YoH4mJgbm5xGd5+fKllZWVqanpmzdvdOEGNcfRdkQqgYCAgH/++Wf27NkdO3a0/GSttrm5eePGjcePH7927dqyZcsCcHBwSElJ0XavJRMeHp52DlD28Yo7AwMDGxsbJycnNzc3Ly+vZ8+elSlTBoCjo2Nqagb2Sed+Xbp0ATBv3jxtdyRH2rVLmatTTVev0pEjdPYsrVlDGZlxevCATp2isDBKTaWKFQmgzZuzfvY6dQggH5+st5DDjBlDAM2ZQ0S0ZQsVLUpLlkh/lrNnzwKoWLGi9E2zL+vcuTOAhQsXqh65f//+tm3bTp069fjx48TMzAvfvn1b3Dbt379fAz3NVuHh4devX/f09BSrZ3v06GFra2tpaSn73BryokWLmpubP3yY1UySWnXu3DkAxYoV46ndzytRgoCsTM2lpVCQmxs9fEibN9POnRTz1czVRBQVRWfO0J07RESrVxNA1asr05VlgVg1069fFl+e8zx5QgAVLEjv3lFYGNnZUYkSpIlbSFtbWwCHDh2Svmn2LbkyIHz9+rWHh8f06dPbtGlTtGjRdJeKggULtmjRYvLkybt27Xrw4EHasOfx48cio/fIkSO12H9JrFu37ksBcJMmTcaPH799+3Y/P79PQ9+7d+8WKlQIwOTJk7XS82x2+vRpABYWFpm609IVO3YQQP37S9DUb7+RiwvduEETJtDGjbRvH23c+PkEawsX0pQptGsXOTjQzZv05AlNmqTWtaVmTQLo1q2st5DD/Pef8oaEiA4cIIB++EH6s6SkpIjvw3v37knfOvucZ8+e6enpGRkZvX79WpIGxdxgwYIFnz59KkmD2UahUPj4+DRr1qx8+fKGhoZfGrPOly9fpUqVWrduPXjw4Llz527btu306dPdu3cHUL169fj4eG2/j6yws7MDsG3bNm13JEcqVowACg1Vt51Xr6hlS7pwgQ4epFat6O5d2riRbt/+zDP9/aldO9q5k+bMoYkTKTmZVq6k8+ezfmqRMTP3rxdNS5TVOHGCiJSjuBcuSH+WBQsWABg2bJj0TbNvyX0B4f37938RyQnTRID29vbOzs7u7u5+fn6Kryb5PX/+vBhSXZmbM9TPmTNHlaMsf/78ad/+1+f9RFG+s2fPGhkZAfj999+zq8vaVLNmTQA7VNuzmMq2bQTQwIESNNW4sfIgJYWaNv3i0xITqUUL5fGTJ+TkJMGpxZXK11eCpnKGpCQqVIgAevSI4uPJ1JRkMnrxQvoTiRQXixcvlr5p9jkTJ04E8NNPP0nYZrdu3QDUrVv33bt3Ejarae7u7mIvpVCoUCEbG5tWrVoNHz7c1dXVw8PDy8vL39//s1c0Ve2NXDq2u3XrVgC1a9fWdkdypMKFCaCwMHXbWbmS/v5beTx3Lp08+cVnTphAV64oj7t0yUqti3T++IMAkvTPXOtmzSKAxowhIpo0iQCaMkX6s/j6+ooR/Ly0ji+3yH0B4dChQwFUqVJlzpw5Bw4cCAwMzOAL//e///37779E9Ndff8nlcrlcvm/fPk32VFOio6Pz5csnl8s3btzo7++f8Rc+fvy4UqVK4kPYvXu3TCaTy+X//POPpjqaY/zvf/8DYGtrq+2O5DwSXrfSBoFfCQiDgqhvX+VxSgq1bCnBqatUIYDy1jRXv34E0LJlRERduxJA69ZJ0OzLly8bNGhw+/0w+V9//QWgSJEikyZN2r9/f3BwsATnYF8QFxcnysBeu3ZNwmYjIiLEbogJEyZI2KymiVmyKVOmPHz48LOLJ5OTkwMDA728vHbu3Ll48eJRo0Z17Njx5cuX4qd37tzJly8fgJ07d2ZvxyWQmJhYvHhxAF5eXtruS85ToAABFBmpbjtz5tDZs8rjNWto794vPrNXrw9BoLOzBItNNm8mgIYOVbednOTaNQKoVClSKMjLiwCqUEGalqdMmTJt2jTVP0uWLCmTydq3b7969err168nZ3nhLsukXBYQJiUlFSlSBJmvq37mzBkApqamV69eJaIlS5YAMDY2vnjxomZ6qkE7duwA0PQr99xfsGjRIgBmZmY+Pj5ENH/+fLEgxzuX52b8psTERAsLC3BCxU9JeN1q357eviUiCgwkR8cvPi05mZo0UR7fvEmSFNf+/nsC1EqUmvN4eBCgnHYV64/atlW3zbCwsGrVqgFwcHAgosTExBYtWujp6aVdcFGmTJk+ffr8/vvv165dS07WiW3G2Wbt2rUAGqvm0qVz9epVQ0NDmUz2t2pKJGfz8vJKu48uJSXl0KFD69ev/+WXX/r169eoUSNra2v9z2XjSBtBrV+/XlzRHuTCv/3Zs2eDEyp+lpkZARQdrW47//zzYe91v3706NEXnzlrFp06pTxu316CWHTDBgJoxAh128lJFAqytiaArl6l1FTlTk8/P3WbXbp0KQBDQ0M/Pz8i2rdvn1wuL5wmn5CpqWnTpk1/+eUXT0/P0NC3ErwT9gW5LCA8ePAggJo1a2bhtaNGjQJQtGjRx48fE9GYMWPE0Hiu25ju4OAAYF3m5wsUCsXAgQMBWFpaPn/+nD75TPIwscy4d+/e2u5IDiPhdev2berWjYYNo+7d6eu/Tlu3Ur9+5OpKHTrQs2cSnLpCBQK+cdLcJiaGjI1JLqdXrygsjFq0cKlcuVakGncqcXFxjRo1AlCtWrWwsLDU1NSePXuKb4OtW7e6uLi0adMmf/78qsuwoaFp/vxka0vOzuThQWm3vJ05o/ywX7wgX1/y9FQ+fvOmRta15hkiGvfQTO72FStWqLOZMDEx8fHjx2JPQTYQv3uzZ88W/0xNTf10G6FcLi9ZsmTDhg179eo1ZcqU1atXe3p6hoeHp22nX79+AGrUqJHrNhMGBwcbGhrq6ek9k+Q7MC8xMSGA4uLUbUehoHHjaPBg6tuXvr5FKCSEOnakxYtp8GByc1P3vES0di0B0gx35iTjxhFAM2cSEc2Y4fXDDz+6ui5Xp8Ht27eLpWq7d+8motOnT4stXePHj9+wYcPAgQO///77tN8JTZtes7SkHj3IzY28vEi1Rv7t2w+rg48cIR8fEisJFAo6eFCdDuqWXBYQ9u/fX7Xj5a+//rqVmZn9pKSktm3bAqhcuXJYWFhKSsqPP/4IoHz58rmoOGFERISRkZGent7r168jIiJat279v8xU2klKSmrdujUAGxub8PDwdJ+J5rqtdS9fvjQwMNDX13/Bd6xprVsn8XUrg7lhoqPp8WOSapNA2bIEUG5LqvFNQ4f6NGkyxt19LxG1aNECwJ9//pm1ppKSktq3bw/AyspKLLMfN24cgAIFCtwUpSCJiCg1NdXX13fjxo2DBg1q23aaTEbAh//Kl6f+/WntWnJyovbtKTWVjh0jNzdq3Vr58t9+o//+U/NN51nHjx8HUKpUKQ0FXQqFQmRUrlev3lc2E4aHh/v5+Z04cUIk83RycmrVqlW5cuXERPHUqVNPnz6tie6lJb6NDQwM0n4bDxw4cNiwYQsWLHB3dz979uzTp08zsiUyJiamUqVKAMaOHavJLmuEqCSZdrEcIyIyMiJAsjKsGbwkpaaSv78Ec4OCyFOaC38nv+706bgmTVY4OAyi99MzdnZ2WW7N09NTrAJwc3MjomvXrpmbm3/6t/zmzZuDBw/OmDGjWbNmdnZRaS9J5ubUogXNmkUeHlS4MInVAx060KJFdOYMEVFKCrVrl/X3q2tyU0CYkJBQoEABAI8fP05ISBCD2ZkaDY2Ojq5Vq5ZYtJOYmBgXF/fDDz+IK2ic+sNR2WLLli0A2rRpozpurbody5ioqKgaNWqIRaeJiYnpPhPN9DpH6NWrF4DOnTt7eXk9f/4828bCc7LdGzc2tbL6VRN7w7PRsx9+8LeySn7+XNsdkdjmzZsBdOrUiYhWrVqV5QVmqqUBxYoVE4vrZsyYIZaLn/9qJr2oKDpxglxcyMFBmeRGZB4fMICWL6cNG5QBYZMm9PgxPX5M06ZxQPhFYmWHRvP3hIeHf/fdd2Izob+//5kzZ9zd3efPnz906NC2bdtWqVJFbLr7LAMDA5Fy1tLSUqoMqF8ifv36qF/thoiIbt++Ld7Xrl27JGkw21y/fh2Aubn5vn37fH19o6KitN2jHKFF2bJNraySc/MFOnzdOn8rq7fTp2u7IxJLTk4Wif0fPnyYmJhobm4uk8kynsgjrTNnzojJQFES7NGjR2JXbb9+/b6SGTElhfz8yN2dhg8nGxtSDVnOnUsTJlC7dpSUpAwId+6kx4/p4UMOCDMhNwWE+/fvF8Gb6rhu3bqZbeTly5fW1tYAevbsmZqaGhwcLEldvsjIyHQ5kTS0gqVNmzYAtmzZojrO1AyhEBQUVLp0aXFJVigU6T4TDfQ6R1i+fLkYUFApVKiQra2tg4NDRvLa5Ulubm4Afv75Z213RC2lSpUCEBQUpO2OSCwkJESUKIiOjn7x4oVMJjMzM8tC7bLx48cDyJ8/v9g8vGbNGhEAZKrWU0oK3bpFa9fSb7/RwIEUFkYdOtCOHeTmRjVqkKurcgkwB4Sf9eTJE7lcbmRkpOnVKN7e3vr6+mZmZl8K/IyNjcuVKyeSebq4uGzcuPHEiRP+/v4pKSmpqamtWrUC0Lx5c82l+NPEju7Vq1cDKFiwYKayrGldYmJikSJFrKysPv2/4+TkNG3aNPF/x8/PLzY2VtudzT6i8mSuvgovX74cwMSJE7XdEekNGDAAwG+//UZEjo6OAFavXp3ZRm7fvl2wYEEAo0ePJqKgoCBxH+7g4JCpFDKvX9OBAzR1Kh06RLNn019/0a+/KgPC4cPJ1ZWWLOGAMBNyU0Aodh0sX75cdbxM5ODLJD8/P/G7+Msvv1CaunyTJk3KeCORkZFpy8HL5fJZs2aJBJ5E5Ovr+91330mexTQ0NFRfX9/AwCAsLOzNmzeq4yw0defOHREdzZkzhz75TCR36tSpdu3atWvX7oImKtdkwIsXL8TIVp06dezt7UuXLp0ul0ZaovJVq1atfvrpJxcXly1btuTG5EMZkTeuW6IaZ57MkCl2/YlvkqxV7HVxcRE3mmfPniWinTt3yuVymUy2devWLPdq4EAKD6ebN8nGhpeMZoiIyYcMGaKJxmNjY1U761JSUiwsLIyNjUuWLNmoUaP+/fvPmDFj/fr1hw8f9vPzi/5Woo7Xr1+XKFECwMKFCzXRVXq/sKVOnTrSNtunTx8AdnZ2Ga+90bt373bt2mXqui8t1R7+Nm3a2NjYmJqafumSJJPJLC0t69ev7+joOGHCBDc3t/3790dERGir55qjUCjE+9V2R9Ty22+/AZiSy5fefJaYjLG3t6f3CQ5btWqVqRaePHkivmT69OmTmpr69u1bGxsbAPXr18/ywMetWyT2I3fvTvXrq7Vk9N9//xU3q5nakpY35JqAMDY21tTUVCaTPX/+PC4uTnWctdZOnz4ttrCvWbOGiM6dO/fNunyvXr06fPjwggULunbtKgYz0hIvNzExEVlMly1bJuIKadNarlu3Du+XkIljkS0wa44cOSIWcK9fv54++UyktX37dvFBaaXUR0pKitiC1a5du7TjjuHh4devX/fw8HBzc5s2bVqPHj1sbW0tLS3FCGVajRo12rhxY/b3XNPEdWvy5Mna7ohaxIRDLtoJnHEiA5uTkxMRrV27dvTo0X6ZSesmviX09PTE393BgwfFn7wYVsuyUaNI3ItOmkQbNlCXLsrHV62iLVto1aoPuf0YEcXExIjRNw3dYbi5uZmamoqLl7hdq1ixYhYmWN68eUNEp0+f1tPTk8vlp1R5FyUlxjW2b98ubbMRERHlypUDMH78+Ay+RAwRNmzYUNqeZNDevXvFnYOYtxfi4+P9/f3FDk8XF5fhw4eLHZ4GBgafRolz587NewlpUlJSxFeWtjuiFldXV+TR3aFxcXEmJiZyufzVq1fh4eHt27fP1CK1ly9fiho5rVq1Evu2GjZsCKB69erp8kVliq8vLVpERPTwIdWpQ8uXK/cTpqZS9+60di1NmkQZXNkqdmcA0NAXYE6WawLCXbt2ifvydMdZ9scff4jvnXR1+VRpu1++fOnp6eni4tKjRw8xgJGWmZmZra2tk5OTm5ubl5dXYmJiNmQxbdq0Kd6XXRLHahZbF5uUDAwMjh8/Tp98JhLSbkA4c+ZMAMWLF3/16pV45Pbt2z4+PqGhoZ99fmxs7N27d48dO7Z58+Y5c+b07t1bX19fX18/ywMQOVbeuG6JUjRv3+bBhNRPnjwBULBgwSzUHN+9e7eYDBQXbG9vbzEFoUrtqDkiDR0TxILGLBQKygiFQiES8R04cICImjVr9vWRzS/ZtGmTqanpyZMn6f0evxIlSqi+MKVy7tw5ABYWFllY+fxN165dy1TtDS0GhAEBASKxvmq5XWRk5Llz5/z9/T/7l56SkhIUFHTx4sU9e/YsXbrU2dlZbBbN7Ys7PvXu3TsAhoaG2u6IWkSJrxkzZmi7IxrRuXNnAJs2bcrsCyMjI0XGih9++CEmJiYpKaldu3YAypUrp+kFPh4edONGhp7JAWEuIDKCim/PtMfqEIWATExMrly5QkQLFiwQI3aNGjUScw5pFShQoFmzZhMnTty5c+e9e/c+HX9NTk4WGTsrVKgQGhoqeRbT4OBgPT09Y2PjqKiotMdqNjt9+nQA+fPnF9km5syZI/osbc4VLQaEqgFvca8jdOrUSfTHyMioXLly9vb2PXr0EHs2PD09r1+/nu6DFcnNp06dms2d17S8cd0SS77VGV/MyUStghMnTmTqVaIwHYClS5cS0Z07d8SnNHLkSM10U0mhoPXrJShOlWcoFIrKlSun++rLQnj/JSLX33fffZeSkuLr6yuTyczNzbNwXZg7d64I1YKDg7+0pEJ9YtORi4uLhG2mtXLlSmS49oa2AsKkpKQGDRoAcHBwUCgU4sFTp06pbjYysrn9xo0b4sKdx1LRxMfHAzA2NtZ2R9Qiijxnw9CbVohV3x06dMjsC1u2bAmgatWqou6RyPNnYWGh6dpvd+7Qhg0ZfTIHhDldVFSUsbGxXC5/+fJl2mM1m1UoFAMGDDAyMtq/f794pHXr1iVLllRFgPb29s7Ozu7u7n5+fhm5Lmo0i6m41HXv3j3dsZoUCoWTk9OAAQPEPYpCobC0tKxSpcqjr1RxzTxtBYQhISFig5nIZKUyfvz4mjVrpi1++qmiRYuqBh1EOrhChQpl3+b+iAg6doyuXqX3dwya0Lhx4zxw3RIJh7t06ZLraopmxJAhQ4yNjfPly9eqVSsXF5cTJ05k5MskJSVl2LBh06dPJyJ/f3/xV9ClSxfNJQsRvLxo4UJSb+FCnqJQKGbPnp0vXz6xMlChULi6utasWVOqrGOijJBYAzx06FBkNUdUamqqaKpZs2ZiPkrES66urpL0k4iCgoJEtQnN5X9SKBRdu3bFt2pvCNoKCCdMmADAysoq7f5/Ly+vRo0aWVtbi0Xdn2VkZJS2CLPYYLx27drs67qPDx09Sppci3H58mUApqammjtFNhDDK9WqVVu9erVCk1dwrXj9+rVcLjczM7OxsRk+fLi4Q87IC8+dO1etWrWAgAAicnZ2FrfZNzI4c5dVKSn0yy+0cSM9eZKh52clIAwIUFY3iY+ntKsqUlPpwQN68ECjd3ESyh0B4datWwG0bNky3bH63r17l3YFvxgWHTt2bJYXB2ooiykR1a9fH++LGqc9Vl9SUpLqO+vJkycymSx//vzSLunRSkCousVp2rTpl+6DExISxJ4Nd3d3V1fX4cOHOzg42Nraino4G9IMK4kx3Q0ZH2hSR2AgtWlDHh7k6kpDh2ruPCITg0gslHupSlqbmZlJvjdJu5KTk6tWrSqXy9PeFxoYGNStW9fZ2XnXrl1f+aZSKBQKhSIkJEQsKWzRokXeriuTY/30008AbGxsYmNjVXXzRowYoX7LDx8+lMlkJiYmYWFh4eHhJiYmMplMrPXIgtevX4uBg/nz5xPRkSNH5HK5vr6+l9iOozaxGqVfv36StPYl4eHh4pr7zZxzWgkIDx06JJPJ9PX1v5KoTGxu9/T0FLUie/ToYW9vX65cOblcXrlyZdXTPDw8kNX9olkxfjwtWEB791K7dnT/voZO8vbtW/ErraH2s8fkyZNVX9ddunTJY0We//e//wFId1UqVaqUo6Pj8uXLL1269JULjfhdnTVrFjJQ90grshIQOjoqK1heu0aqvIwpKeToSIsW0YIF5OgoWdVlTcodAaFYZyyWLKc9lpbI4WloaKjmX68qY6cYnldlMVUndUdAQIBMJjM1NY2NjU17rE4/P0usmx04cKA6jTx+/Hj3x0aPHi3+xiZMmJDuR5pL5TRv3jwAFhYWWZtMfv36ddrVOLt37xZ3ddkx4DdrljJJFhH170+S5lJPTk5WfV8rFAp7e/tt27ZJ2H42CwkJUQWEwsCBA2NiYqQ8R0wMjR5NQ4dS79507ZqULX/L77//DqB8+fJPnz718PBwdna2t7dP935LlCjh4OAgJg/TzTtFRUXVrl0bQN26dSX+TPI0aYfDEhISatasCeCnn34iIl9fX1E3T+wGV4fYuC5Sty9ZsgRAx44d1Wkw3QL7adOmAShdurTIN6OO+Ph4sdf38uXLajb1TStWrMifP//mzZvTPrh37950lx4x6vf999+ne1yqkdZPvXjxQnwIYiF3ZsXHx7948UL1z5SUFLGT8OjRo9L18QvevKEff1Qe37ghecn1tDczHh4etWrVytVlJ0RJMBVra2upRlU+2L+f+valQYPIxSU7J6DCw8NFwdLt27eLTPs9evQQYysq+vr6tra2Ynndp4u3165dC0BPTy+De3016t69e+n+/EXNXgCzZs1K96P7XxoH+WxA+O+/pFpesXAhZTI9uFbkgoDwzZs3Yp3J27dv0x5LfiKx9f9H1beeGs6cOZPZLKZfJ5J/9O3bN92x5MSGpSNHjqjTiMhtmEEa2ph37tw5cWcj8uWoLzk5WRRvzOxurqwYMuTD+oYZM+jCBZLubn7KlClppwQvXbqUbj1t7jJixIhPf6kqVap08+ZNyc7h4kKenkRE0dHUuHG2XX3fvn0rbh8PHjyY9vHY2FjVlVhcm9NeiVXLeO7evdukSRMAFStWzMkpWJOTk2NzkkePHhUtWnTixIkShtB+fn4mJiYA3N3diWjDhg0AzMzMsjybR0QRERGi3qCvr68qNjh27JiaXRVb60UKruTkZLEusUOHDmoOhIkEZnZ2dmp277MUCkWzZs1mzZoVGxtLRGIzf7qgS1TBzggDAwNNdDI5Odne3h5A+/btpRpVFDcD7du3l6S1r/H1pTFjlMdv31K3bhJekgIDA62trVV/bgqFYujQoervCdKWu3fvflrUSl9ff+HChZJFua9eUceOymWKc+dSNkZW48aNA9C4ceN0v8P+/v7u7u7Ozs62trbpJg8tLS0dHBxcXV29vLy2bdsmUp398ccf2dbnrxBpFDLoi+vnHR1p4EAaOpS6dPkQEK5aRe83o5GHB2kge7/kckFAKK6dYgOrOFZzEPRLxJVv165dkrT2559/ymQyPT09kfxNlcX0n3/+yUJrYphfJP9Meyyt+/fvAyhUqJCaOQ82b95s8jERDwMwMjJK96NZs2ZJ1X+VsLAwsXB3pqTpDhctWlSgQJkxY05L2Obn/fYbqX5PunShCxcoXz5yciI1bh9VxCi1WMefB4itCJ8yNjZWP++UUps2pPqL6N2bsqvg4ZAhQwC0adPm60978ODBtm3bRowYUb169XQ3IjKZzNraOjCD+ba1RBTpyYHSpqFSn8jhbGpqeu/ePXqfp6pGjRpZ3kwoPjfx6yFqGHz//ffq33GmpqaK3A8tWrRISUlRVXDNWtVfFTFH+ueff6rZvc86fvw4gFKlSiUlJT169Egul4tltGmfU6RIkXSXHvF/WTw5rYIFC2qik2LFrCTTrSphYWEmJiZ16kx++FDKJHCfER1Nqi+iU6do5kxq147q1CF3d0nWwjk5OWniTkArnj179qWvlBYtWkiTTvPIkQ+zT1euUHYVPPTz8xMZ12/fvv2Vp0VGRh4/fnzu3Lnt2rUTFXdURKbrlStXZk+Hv2np0qXp/vxVC3CMjY3T/WjFihWfb+WzM4Senh/+Hy1YQNkwja+2XBAQNm/eHO8HVsUeP3EsrRcvXogLg4SjwiJjZ7ospqampnfv3s1UO6rU84mJiREREdbW1gUKFNDEdiAxNjx8+HDJW87OPYQKhUKkRa5fv760uVJDQ6OLFFHI5fT4sYStfiIlhSIjqVMnmjOHfvqJfvuNVq4kmYwA0tOjnj2zsGoxODg4bZagtWvXanond7YJCgoaPnz4l4b/1dq/ERJC8+fTypXUrduHPAodO0o4NP4V169fl8vlhoaGmZpEioqK+u+//+bNm9e+fXsDAwN9fX1NjBxJa9WqVSY5iarmm+Qr8QYMGACgWrVqcXFxqs2EY1QTL5mRmpoqyu4dOnSI3meHkiq/iKo8/eLFi4no8OHD39z29hXh4eEiLaGlpaWE6VXTcnBwUPVWFHzKSDbd7NxDeOrUKblcrqenJ3nqwokTYwDJl3B+TIR8a9dS//40dy45ONDDh1SsGAEEUMWKtGEDZXJcQ6FQpN0/FhgYmKsXqqTzxx9/iL/uT1lYWGT9iyU1lQ4epO7d6cSJD7V9Tp6k+fOl6vnXiao2GS/1SUSpqal+fn6bN2/+6aefSpUqZWRkVL9+fc31UH3q7iGcOpVWr6alSyk5mbp3p99+o19/pR49KDcsgc7pAeGrV6/09PSMjIwiIyOJ6MWLF8uXLxfH0hKjrT179pSwTYVCIZYjW1paitwPo0aNGjJkSMajlKioqHPnzoltk6p9fQqF4kkG8yVlkkg+oYlku9kZEP76668AChUqpImygYMHE0BZSuOXMevWUdOmyi+Xly8pOlr5+KNH5OxMxsbKa7C9PXl6Znzt4h9//NG4cWPN9DhHCAkJcXFx+WzaWCsrq0zvXL95k4YPp3z5CKBixWjPHho7lt6+pWPHqFcvzbyDj6SmporEUb+ohhszT8zz/PfffxJ2TBeEhYVVrFjRzc1N8oyssbGxVapUATBs2DAiun37tthMmIVlKQcOHABQoUKF1NTUmzdvQuoKBKdOndLT09PX1xd/O5MmTRJ/Sl/ZrCFSoXh4eLi6ujo7O/fo0cPW1lZspwdQtGjRLl26iGfGx8ePGTNGzJSq78mTJ3K53MjIKCQkJDo6WuQcvnPnzjdfmG0BoSrAXrhwoeSN+/mRTEZmZqSB2yIiInr5kqpWpcOHiYhiY0m1iTExkdzdqVIl5SWpWDFyccl4AtLY2NjSpUuLsfI8KTU11dPTU+SiS0cmkzk7O2ducCQ6mjZupCpVlJ/2hg3UvDndvEkvX5KDgyRLh75JFAC3sLCIiIjIWguHDh1C9qxwVoO6WUZfvqSEBFq0iKKjKTWV7twhX99cEQ1Szg8I7969a2dnV65cOU1n8qhbty4AVf0JqSQlJYlElzY2NuHh4d+8w4iKihK7g5ycnGxsbFRLsS0sLFq2bKnRD8HHxwdAiRIlNJGYPtsCwitXrojaxGKlruR8fUkmI3NzzVx9d+8muZxkMvrS72FQEE2ZQvnzE5Csr9++eXN3d/eMXFdSU1NdXV01NDyfc0RFRS1durRUqVLpLsAZ3b+RmkqentSqlfKiK5dTq1bKwPu//2jaNFqzhhISKDKSNDzttmUL2doG29o6qrNgQSxKzGNpV7OH5hJa+Pr6mpiYANixYwcRrVmzBoC5uXlmK6aIxTJubm70Poup5DXKf/nlF7xf4picnNywYUMA7dq1u3PnzsGDB9euXTt9+vS+ffva29tbWVl9umlKxdzcXOxvzJcvn6+vLxHNnDkTQPXq1SWpvTF+/HgAQ4YMIaIVK1YAaNWqVUZemD0BYWpqaqtWrQA0b95cQ0VfWrYkgL60nE0t4eFUrRoB1KTJ58cfU1Lor7/Izk58Z65p1WrcuHEZKQJJRIcOHbqWvTm6tOLs2bPt27eXyWTp/i7q1auXoQ/qyROaNo0KFVJelcqWJVdXCguj0FBasoRmziSxevPvvzW6dCUujurUeVe//rotW7K+90+U76pdu7aEHZOcBHUInz0jNzdJO5VNcnpAGBsbK7ImiIydGuLv7y/q+UpVGyqtqKioGjVqAGjatOmn6zxDQkKOHj26aNGi7t27ly1bNt1XhrGxcd26dR0dHcXCa43Wi5s6dSqAcePGaaLx7AkIIyIixJ3HpEmTNHeWZs0IoFWrpG73xAkyMiKAvpmALiqKfvttd8eO4iMtXbr00qVLP50cuHnzpvppDHOjd+/eubu729jYpPtrat68+ZcSFYSFhf36669h4sYKoAIFaOLEzyd3jY+n6tVJX58+TvQioagosrQkgHbvVqsdMavz66+/StQvJo1NmzYBMDMzEznrROkXW1vbjO8CePDggbhgRUZGvnnzRhTmlXzZSHJycuPGjY2NjQ8fPkxEz58/L1SokJh/+yxRTr1Hjx7Ozs6inPr169dVf3FiQ+z3338fHR2dbqZUHTExMWKT0q1btxQKhVjkksFl0tkTELq4uACwsLCQZvPY5/z7LwH03XdSZ7aPj6dGjQigqlXpmwvvT59O7tixtKWlGIDr1avX9evX0z0lJSVlzpw50qbwzS3u3Lnj5OSkWo4u5M+f/0urAxQKxX///bd95EiSy5VXpebN6e+/P///eO1aksmofXtKTtZQ/6dPJ4Ds7NSa6woKCgJgaWkpXb+kp25AGBxMXbrQxo2U1XlULcrpASF9krFTE0SiIScnJw21HxQUJBJU9u7dOygo6MSJE6o5wHTjRoaGhjY2Nk5OTm5ubl5eXqpbhKNHj4p6tevXr9dEDxUKhQilLly4oIn2syEgVJUkrlu3rkanwv7+mwCqUEHSVQCXL5OZGQGU4YGPpKQkd3d3kRUWgLm5ubOzc9r0IU+ePJE2K0buolAoPD09xbSGSrFixcTdrcqjR4+cnZ3FgMvS+vWpQgVydf3GV7mLCwGULx9duqSJnjs7E0CNGqmbzXTp0qXQwMQRU1///v3xfoosOjpahDHOzs4Zb+H8+fMiTd/ChQsBdO7cWRP9fPHiRdqyQBcuXFi1apWNjU27du2GDx++cOHC7du3nzt37tmzZ9/cB6GqvdGnTx/6ZKY0y0Ru8KZNmxKRp6cngO+++y6DE3HZEBCePXtWJLvW6Mrt1FQqX54ASRcuJCVRhw4EkJUVZTgr1Z07d9Lu6La3t/f09Ey7smnHjh1xcXHS9TKXef78uepyo+Lk5BSbpupGYmKiu7t79erVAejL5UkVKpCTE329Ope/PxUvTgANGaKJbj9+TEZGJJeTmit8k5OTxU5aDU2VS0KCGcJcKxcEhPRJxk7JiRm8Q5qsE3Ljxg1zc3MR1KUbImratOmECRN27Nhx9+7dr/ydiGKgBgYGUtVRSOvSpUsArKysNLQqNRsCwvj4+BYtWpiYmGRwvUqWpaRQuXIESDdF5OdHhQsTQAMGZDYIUCgUhw8fFlu9xYDCJc1EKbmXl5eXg4ODauRF7N9ISEg4ePBg69atxeMymaxt27bHjx7N6Oc/ahQBVLQopUnVIwk/PzIwID29b9wAZMSOHTugsfo0TB0xMTGVK1fG+9wn165dMzQ0LFy4cGarg2RrLRy1PXz4UJT+E4VP082UZoFCoRAfo7isiN0Zy5cvz+DLsyEg/PPPP42NjdUpQZxBK1YQQC1bStScQkEDByp3BmZ+f9qLFy8mT56smkzmr6B0IiMj3dzcLC0tVfeBlStXvn37dkBAwNSpU1U74UuVKrVo0aLwDOakvXaNTE0JIA3k5hEjA2pP5xMRiRV/r1+/lqAtzeCAMBdIl7FTQlLVWvimo0ePrlq1Kn/+/Pb29qJkp5+fX6Y2q4hNHebm5lLWWCMiop9//hmA5q5bT58+3bFjx44dOzSXAf/Zs2cymSxfvnzh4eEaOoXKsmUkl0uT2evZs2cXO3cmgLp2VWfFz40bN/r27WtoaNitWzcJupXn+Pr6pl2xo8qFbWZmNnr06EzfkqakkPi/Vr48SVrir3lzAigzc0VfdOLECQDNmzeXoC0mtTt37qQtT+/h4ZGFJFh79uwBYGNjo+k99lIReSlUqbbTzpRmobWjR48CsLa2Tk5OvnfvnkwmM/mk2sRX7N27d8eOHZoYYFXp3r07gEWLFol/hoSEaGjcOTKSzM2pXDlpNrdvmDEj1cKC8ucnNZJRR0dHu7m5FShQoFq1ahr9kHOpuLi433//XazMEmP9qo249evX3717d6ZzpB86RPr6JJORpCX+DhwggAoVotBQCVoTM5+S38FK6MGDB+JmNSdHrRqSawLCTzN2SkUs8R86dKiEbX6JmrkKFAqFk5MTgJIlS0oYWaWmpoo8HLl9h3ebNm3wST1iTYiKolmzKCiIiOjMGXr2LIvthISEVKxYEcDJoUNJ7Z0Vjx8/BlC2bFk128nDAgICxo8fb2RkVLRo0SJFiri4uGS9KEV0NNWuTQA1bJgk0d7j3bsJIAsLaXYf+Pr6FixYqEULBwnaYhqwbt06ZLI8vUKhCA4O9vb29vDwWLZsmai2qqF9BBoiruNVq1YVtTfEFN+oUaOy0FSHDh3wvlr0qFGjkNUaHppz7NgxvC+Q+Pr1a2NjY1NTUw0NWe7erdx1nJxMIpPUuXNZaWfevHkA2lWsmJq1139MVATJWvllXZCcnLxz584qVaoULFhQT0+vR48eam3bWbeOADIwSJJoiXJiIlWsSIBkZdU7dHAoWdLqxInL0jTHJJVrAkL6JGOnVM2K3e25YskNEb17905kk69atWqWk/+KdlTHZ86cAZANqVw17fDhwwDKlCmTrLGt1SqNG5PYc7p4Mf36K/n6ZrqFqKgoW1tbADVr1lTnf6XKqVOn8H47DfuKkSNHZmpp2RcFB6eULz+mSRMHBwdJNkV06kQA/e9/6rdERPTmjXJkl+VYffv2xefK0ycmJvr7+3t5eYkqDsOHD2/VqpWNjU263UcmJiYjRoxIuwEp54uNjRUJnwYPHkxpZkozW7P+8ePHcrk8X758b9++jYiIMDMzk8lkUpWykIpCoRA3GB4eHqThIcstW6hcOXr4kBITqUMHGjuWBgygTNY8Vg5S6OnpSbW5o1atWgA+TTDD0nr69CmA0qVLS9DWL79crlq1XJkykpQaPn6cDAyoZk3J8hU5ORFAW7dK0xqTVm4KCOlbGTu/7unTp/v27ZsxY0bv3r1VD4oKTsWKFcuGEEIqkZGRYtq9efPmGV/m+vLlS09PTxcXlx49etjY2PTo0UP1I3F/PFNV5zTXUigUohqs5OVDPtW5M82bR4cO0eLF1KgR/f238vH9+zNUN/7du3didKNChQpSrUzYtm0bNJkbKc8QSaQkSVz89N49seUjC+l5w8Ppf/+je/do+XJlcrgDB2jHDsmSFSkUZGhIgPpzz0xTVOXpmzZt+vPPP3ft2tXW1lZss/mS4sWL29nZdevWbfz48StXrsxFVy4VPz8/kVHG3d2dsjRTSkTOzs54n6dU5E9q27atpnqshrVr1wKwt7en93XYNDRkuWULrVlDXbpQYiI1bkyNG5OI6UJDM1qR4u+//9bT05PJZP+TalCKSHw9hkqy3DDvSkhIAGBoaCjBoLxCMWTQILGOLLOL6ZKT6coV2ryZTp1S3sa8ekXbt6uzcDi9KVMIoCVLJGuQSSiXBYT0ccbOr//xqEIgBweHdJfYV69eiedMnz49B64z+aZnz56JQrdikPWz/P399+7d+8svv7Rt2/bTO4xatWqJpyUnJ1tYWCBjlXxzvt9//z17Zsk6d6aEBGrdmlxc6PvvSTUA2q0bqYZWd+/e7eXl9elrU1JSHB0dxVKiZ1leb/qJ+fPnA5gxY4ZUDeZVW7ZsATBw4EBJWjt//rzIqvfNgf/ERPL2plWrqF8/Kl1amUvc2JgKFCAxW9m3r8SpqsVZJF1izyR28+ZNY2PjQoUKpf1+NjQ0LFeunL29fY8ePaZNm7Zx40ZPT8/r169HR0dru7/SEDnSTE1NxZyeaqY04zUJDh482LRp09u3b6emppYrVw4aTguXZXFxcSIounr1qkaHLLdsob//pmXLaNs2srOjfv2Uj3t5kSpvzps3b1xcXD778pMnTxoZGQH47bffpOpSXFwcAGNj49y++CgbFCxYEEDWtzCkkanFdM+e0V9/0aRJVL++cgARoGrVqH59Skiga9dozhz1e/TB8uUE0PjxUrbJpJL7AkIiunPnjsgJMefjX9X4+PgdO3aMHz++SZMmn5ZLKlGiRIcOHWbNmvX333+r1tiUL18ewPnz57XxPtRy/fp1sXxowYIFaR8/f/58ixYtxJdLWkWKFGnduvW0adP++uuvx48fq76gxSaHypUra+NNSC86Olr8bmh6jYpI837kCJUsSQUKfNhvbWf3ITVzjx49/vrrL3Hs7u6u2lgvZqiKFClyN7MLer5q2LBhyG0birRCLC2WcD7Bw8NDLpfLZLJ0VeBTU1Pv3bu3bdu2/v1dypZ9Y2CgvNym+69aNXJ0pMBA6QPCunUJoMu8XyMH27dvH4BixYotW7Zs7969ly9f1lzBuhxF7IevVq1aXFxcVFSUuBbPmjUrs+38888/ACpUqKDmFn3NERVBBwwYQJocshQBYVIStWlDNjb0yy/Kx3ftItWiqEuXLv3www/i+O3bt6pKJ6GhoSIB7NSpUyXs0r179wB8//33EraZV4nNtH5+fpK0FhUVJaq8NGnSJN0gS3h4+LFjx+bPn1+z5qmiRVM/e0kyM6O1a2n2bOkDwj//JIB69ZKyTSaVXBkQEtGRI0c+rcuXmJiYtu6npaWlg4ODi4uLh4dHuj+ziIiI06dPizUnJUuWzLEXkq87dOiQvr6+TCYTWbyFc+fOibdfqFChtOlM0w7RRUdHnz9/3s3NbcCAARYWFjKZbPTo0dp4Bxoxfvx4AD/99JOG2n/3jkJCaMIE5T9HjqRVqz5UKzh9mmJilMdXrlxRLQf19fX1f1/ofP369YaGhpLX1Wzbti2AdHX22Kd8fHwA1KxZU8I2ly9fLiZ2/vrrrwMHDsyYMaNly5aqXKbA7s9ed8V/FSrQw4fUt6/0AaHYlKiZYj1MGk2bNoUmq+zmWKry9MOHDyei69evOzo6ZmGGpEWLFgDc3Nw00EdpPH/+XE9Pz9DQ8NWrV6ohSzUTLT55QosWfdjZFRRE//5Lp08TEXl50ZAhH6pFvHhBqlO9efNGVZcoNjZWlVj/2bNn1apVs7a2lnYqTww3t2rVSsI28ypRO0rC0sEvX74USaccHR0vXry4atWqfv36ff/99+8rMNX7yiVJT48ePaI+fWjXLokDwlOnCCBOdJAz5daAkL5Ql2/ixImLFi06evRouoJOkZGRXl5eqnLwcrlc3Ka1bdt2a27e37px40bxIaiS4sTGxh46dCjdGHNUVNRn375qeVLp0qXzzLD0s2fP9PT0jIyMNJQ1eNIksrTMYgI3YfHixSIgkbY8q7i78s1Cfhsd8/LlSwDFixeXtllRu0WVOvxjW75y9bWyIiKaNYvKlpU4IBw2jADasEHKNpmEfH19ZTKZubl5VFSUtvuiBeqUp09JSXnx4sWOHTvEBxgpSb0FjenatSuAuXPn0vshy69s9/imd+/Izo4AZemjhw/J3JzGjMl65o/Q0FARpkpbH2Lz5s1qvlPd0bt3b7wvQiOVO3fumJmZiZXAn6j2lUsSQNeuUWAgVa8ucUB47x4BVKmSlG0yqeTigJDe7wD8bF2+4ODgQ4cOzZ8/v0uXLmKYJC0TE5P69euPHj06Z+46yBSxHCV//vy3b99WPfilAFgwMDCwsbFxcnJyc3M7depUgwYNANSpUydGNbeVy3Xq1AmfLKaVxOHDJJORvj6pkxo6ISFBVB+ScO8+EYk1P7p5Z5kpycnJenp6crlc2tQOx48fByD2E36igUwW96VLr9jhEx9PP/xAz56Rp6dkXZo9mwCaO1eyBpm0hgwZAmCCar2B7hFjml8pT5+QkODv73/ixAl3d3dVwtVy5cqpVgPVqFHDWZLCnZokUnlbWFgkJiY+ffpUzSHLn38mgKytKSyMEhOV5W/UXIbn6uoKwMbGRsJvRVE++ku7FllaYphg2bJlEraZkpJSoUIFsZjuE8bA4S9dkiwsSOTrdXOjlSvJw4Pi4qTpUng4AVSggDStMWnl7oAwbV2+a9euqVLIiC3maZmZmdnb2w8fPtzNzc3Ly0vTNeizk0KhEDvyixcvPnny5K8HwP/73/9u3LiRruDp27dvxU739u3b58acdZ86efKk+K2Q9n/0ixdUtCgBpP6ue1GguXjx4lLFb2FhYWJcQJLW8rzixYsDkHBWPDk5uVq1anifkiEhISE8PDw8PDwkJMTf39/f3//Bg4jr10n8d+oUnThBJ07QgQN07NiHRmJiqHJlMjCgo0el6dWePdS0KW3ZIk1rTFphYWEmJiYymezhw4fa7os2ifL0NjY2Z86c2bNnz2+//TZu3LhOnTrVrFlTpGP5LJlMZmlpKUp454qcJWJPl5gLFUOWCxcuzEI7hw6RTEYGBiTWfo4apVx2ruaV5N27d2Ib56ZNm9RqKI1BgwZJPu6ZV4mAfMqUKRK2uWbNGgDW1tZxcXFEFP7e8+fP/f39Hz16duuW8pJ08aLyknT0KHl40Pusi0Tv08B07ixN5QmFguztqVMnyhN3mnlN7g4IiSgxMVHswUgnf/78aXfQ5dJdghmUkJBgb2+fNg4UAbCzs/PGjRszEgD7+/uLW+ShQ4dmT581TZQn2bVrl1QNJidTo0YEUPv2JEVqaIW9vT2A2bNnS9E7unXrFoDq1atL0lqeJ27OJKnUJIjE9xUrVsxsOZx0Zs4kgMzNyccn0699+pTS5u1/+ZJatSJR3y7Hz6DoIrF0vFOnTtruiJZFR0dXrFgxzYbbjxgZGaVLuHrixAk/P7/cVX2R3ic3rlOnDhGdOHEia0OWAQEB9etfBUjMJO3dSwAZGUlTGGDPnj0ALCwspBqmFNs7/5OoSHreJnnVqLCwsCJFigA4oN4m8idPyMKCABo2LNOvffeO3m9TVZo8WZmG/eRJ+vdfdfrFpJfrA0IiCg8P37ZtW5UqVdq0aTN9+nQPDw9V9g7d8fbt27t3706ePHnXrl0PHjzIQgDs7e0tCgT/6uqqiR5mM7ESyc7OTqoGZ8wggEqVojdvpGnQ29tbJpPly5cvs8WCPsvT0xNAhw4d1G9KF4gEPEeOHJGktdevX4vbWfUbVCho0CACyNKSMl6RxMuLevQgfX3q0uXDg8+eka0tifKi/HuR4yQnb+vcuby1tbS7tnKpe/fu7d27t169eo6OjhMmTHBzc/vnn3+uX7+el+rXJSYmigpPFy9epCwNWSYnJ9vb2xsaGo4Z869CQf7+VKAAAbRunWSdbNSoEaSrXVSxYkUAmSovqbMkT8AzfPhwAK1bt1a/qStXyNSUAFq8OKMvCQkhV1cqXZrkckq7AKJbN3JwoMhI2ruXl67kOHkhIGRS8fT0LFqkyJnatenj7Pm5UXx8vBgeu5wm6f6TJ0/8/PyyMPx5+vQruZz09NTKJfMpsY+8f//+6jclFoeMHDlS/aZ0wcCBAwFskeiKJBaud+3aVZLWROJ4gKpUoa/nXIyPj9+yJblGDeXGDyMj+umnD3Xtnz0jZ2dycqJ79zggzHn27iUgtWrVXLHckUlixowZAHr16kXvhywz9Y0tkiaULl36zZs3iYmJzZuHAeToKGUPr1y5IpPJjI2N1R+mVCgUYog5183laoVY41OtWjVJWvPx8RGJbaWKxj09SU+PZDJyd//GM2/ejBk4kIyMlFelGjXI2/vDT7t1Iy8vGjuWA8KciANC9pGoNWuUt5bShj7aMG3aNAB9+/ZVPTJixAixDMnY2LhcuXKtWrVycnL65jKkV69eFS9evFGjgYsXZ7RocgYFBgaKTURXVLULs0q82UWLFknSsTxPwo/r4sWLYqb36dOn6rcmREdTrVoEUOPG9NlK3a9evXJxcSlatGijRo8BKl6cpk2joCAiotRUOnyY2rShgwfJ2ZmCgqhLFw4Icx6xAF3CyR2W4718+dLAwEBfXz8wMDAuLu7q1asZf+2xY8fkcrm+vr6XlxcR/fzzz6VKWbdu7S15dlWRkiDtdTNrQkJCABQtWlSSXuV5r1+/lurjUigUIlPgtGnT1G9NZf16AsjAgD67BDg1NdXT07NVq1YVKnSRyUgup1atyNNTub8mMJCmT6dx46hbN0pOpuHDaeZMDghzHA4I2SdECrPChekLad9yi6CgIAMDAwMDgxcvXohHZs2aZWNjY2pq+pVEBSVLlpwvknkTEVFqamrr1q0BNG/eXNoqEYIY9G3QoIGaEwXiKp6F7O26aeXKlQDUT06Ymppat25dvE8oL6GgILKyotq1Tw4YMCTtCvBLly716tVLlWLxxx9/2bmTxEak6GhavZq+/145NNuzp3Lr4IoVVLKktL1j6rlxgwAqWJDySmJnlkFiVcgvqrLxGZOcnCwSUy9evJiIDh48KJPJDAwMVEUFJfTixQsxTHlBnVTaRNeuXQNQu3ZtqTqWt6Wmpoqy0upnwtu6dSuAEiVKSJ5yfPJkMjWNadKk1507d1QPhoWFubq6qnJYFCxYcO7cN6rR0QsXlNsZADI0pI4dKTmZ3r4lKysOCHMcDgjZJ1JTqVs3Aui770gzpfyyTY8ePQDMmjUr3ePx8fEilfnGjRtdXFzSpTJPmyZ73rx5Yp+9huo0RkdHlyhRAsDevXvVaUfs/Th79qxUHcvbdu/eDaBnz55qtrNu3TqkSeMmLT+/2GLFLABMnDhRPPLu3Tvx22JgYNCrVy/V7aC/P02bRoUKKUPBkiXJxYUePaI9e4iIkpNJoj1BTCJin+ikSdruB8tuly5dAlC4cOHMfmPcvXt37NixqampgYGBYjfE8uXLNdRJsbS1fv366gxT7t+/H8CPP/4oXb/yuJIlSwJQjV9nTVRUlKWlJYA///xTqo6ppKbS6NGzxLplVT/nz58vQsFKlSqtXbtWVC979448PKhBA+UlycCAevSgEydo/XrlpoYDB9Sq3cU0gQNC9jnx8co/ZVtbys0bALy8vACYm5t36NBhxIgRCxcu3L59+9mzZ/39/T87DieKHYeEhIh/njt3ThSs02ietPXr1wMoW7asOgkqxfichKsW8zZRFqxJkybqNBIWFla0aFEAf//9t1QdS+fs2bOirPDKlSvFIytXrpwxY0ZQUBARKRSKkydPjh69Sy5XXnebNqX9+6XJD840JTSUjI1JLifdS37GiKhevXoymczGxqZv377Tpk1bs2aNp6fnrVu3wr6+Y5iI3ueVEfnDNLf7NCYmRgQVf/31V5YbcXNzAzBu3DgJO5a31alTB8C1a9fUaUTUM7S3t9fQr8e7d+9atmwJoGrVqhEREUQUEhLSqVOnY8eOiTOK7Qy1ayeKS1KxYjRzJr18qYm+MIlxQMi+IDSUypfPA7tcdu/eLe7aP1WoUCFbW1sHB4fhw4e7urp6eHh4eXn5+/uLFXqhoaGlSpWCdJUhviQlJaV69eoAli5dmuUWDAwM5HJ5XiqwqVH3798HULFiRXUaEVtSJcwL91l79uyRy+VyuTztHHJiYqK7u7v4tdHTMyxfPqlHD0qTPonlYAsWEEA8c6KrvL29xV37p760uT3m/dLiqVOn4n1eGY12ctOmTWKYMuGzm5gzYNKkSXhfl5VlRIcOHQAcPHgwyy3cvXvXwMBAT0/v1q1bEnYsncjISHHpadasWdpRbB8fn+HDhxsbGwNo2vTfWrVo40bJKtqzbMABIfuyhw9JhCj+/rR9O50+TUQUFka3byufcOaMtrqWKXfu3Dl48ODatWunT5/er1+/Ro0alSlTRrUL61NGRkYVK1YUQ6RNmjTRxNbBdE6dOiVmMl9nZo1ucHDwwYMH586d27p1a5lMZmpqqprbZF8XEREhPvAst3Djxg2Rxu2+5rfaiprFxsbGXl5eAQEB06ZNUxXsLlmy5IIFC0JCvj23wHKKhQspf346eVLb/WBaExoaeu7cue3bty9atGjkyJEdO3asVq3alyoxCiVKlLCxsZHJZKq8MhqVkpIiCmO4ZqYMVWJi4tWrV9etWzdkyJBChQoBmMTrojNs8ODBADZt2pTlFtq0aQNg7NixEvbqs168eFG6dGkAvXv3fvfu3Z49e0QaGwB6enrdunU7d+6ipvvAJCcjoq98BzGGixfh5oaxY+Hjg+fPMWAA/v4bixcDQJs2+O8/bfcv6yIiIp4+fRocHPzq1aunT5+qjp8/f65QKEqUKNGkSZMVK1aIeUJN69Chw9GjR8eMGSMKSHxWZGSkn5+fz3v37t1T/cjQ0DApKalEiRJbt25t165dNnQ4tzMxMUlISBg5cmSTJk0aNmxYpkyZjL+WiJo1a3b+/PkpU6b89ttvmuukypgxY9atW5cvX77k5OSUlBQAP/zww88//+zo6PiVoQ2mTRcuICUFzZoBwMKFmDED69bh9m2YmWHMGJQvD5lMyz1kOUxiYmJwcPCnV6WAgIC4uDgAQ4YMqVKliph807TTp0+3bNnS3Nz80aNHYuvyp1JSUh4+fOiTRmJiouqnMpmMiEaMGLFixQoTE5Ns6HOuNnPmzMWLF1etWnXw4MENGjSoU6eO2C+QQR4eHr169SpSpMjDhw/FLlONunnzZtOmTWNiYgoUKBAVFQWgUKFCQ4cOHTNmTKYupiwH0XJAynK+vn0pMFB53KYNnT9PEyZQWBiFhVHz5lrtmabEx8ffv39f/VIQmXLv3j19fX09PT0/Pz/Vg5GRkV5eXm5ubk5OTjY2Nun+eM3Nze3t7Z2dnd3d3Tdv3qyvrw9AJpMNHz5cEzlO8pibN2+KMllCiRIlHBwcXFxcTpw4ER8f//XXuru7i5dESp70/QtSU1O7d+/+448/GhgY9OjRQ80EgCw7/Pknbd2qPG7dmlavJjc3IqJHj6hNG+11i+U+CoXi5cuX3t7eb9++zc7zduzYEcCoUaNUjyQnJ/v5+bm7uzs7O9vb26f9ChXKlSvn5OTk5uZ26tSpH374QTxYqVKl69evZ2fPc6NXr14NHTpU9Unq6+vb2tqK6/s3swPEx8eLMGzz5s3Z01siOnXqVMOGDQsVKlSxYkU3NzcuOJnb8Qwh+5a2bXHwIAwNAcDJCb17Y+FCNGkCAPv2wd9fu73LS0aPHr1+/fr69et36tRJjLYGBASkfUL+/Plr165tZ2dna2trZ2dXoUIFWZpJhi1btqguJzY2Nn/++WetWrWys/+5S0pKipeX16VLly5fvuzt7R0WFqb6kaGhYZ06derXr9+gQYOGDRuKtTEqMTExlStXDg4O3rlzZ79+/bKtw0SUkJAQGRkpktGxnG7XLpw8iYYNAWDrVhQsiH37IG6gu3fH5s14v+6XsZzpwYMHNWrUUCgUc+bMef36tY+Pz+3bt9+9e6d6glwur1Spkq2trbgk1apVy8zMTPXTt2/f1q9f39/fH4C+vv7MmTNnz56tp6enhXeSSzx79uzs2bPe3t7e3t737t1TKBSqH1lbWzds2FBclWrXrp1uYcisWbMWLVpka2t79epVuVyebR0moqdPn5YrV07G6x1yPw4I2bdMnoxu3dCwIRQKNGuGZctw4EDeWDKa07x588bKysrIyCg6Olo8YmZmVrNmTdv3qlSp8vXv+vHjx69atUoctzcwmOnsbL90Ka9MUyGiwMDAzy5oCQ4Ovnjx4oULF3x8fK5evZqcnKz6kaWlpa2tbaNGjezt7e3s7GbOnLlixYqGDRteuHCBr4Lsi3btwoMH6NgRAKZNg6kpB4Qs1+nateuJEyfEglVBfB8K9vb2hb/6a+zr69uoUSNxRbMGRtjYTD9yRM5LCtN4/vx5mTJlPr2UxMbG3rp1S1yV0g1ZGhgY1KhRw97e3tbWtmnTpikpKdWqVXv37t2lS5fq16+fvd1neQcHhOxbQkIwejSsrPD8OX76CaVK5Zk9hDnN/v37HR0dZTLZkCFDmjRpYmdnV6lSpUyN9qWmpnbu3PnIkSPTgSXioTZtsG0bLC010+VcxsfHp0+fPiIV21eeFh0dfeXKFW9v78uXL1++fFlkoBGMjIzEFr7r16/zBCz7ml27kJSEQYMAoE0bdOwIQ0OMGoWAAAwfjuPHtdw9xr7l3bt31apVe/LkSeXKlYcOHWpra1unTp38+fNnqpEjR4507tzZOjX1ImAJKAoUkG/ahJ49NdTnXKdly5ajRo1ydHT8ynMUCsX9+/e9vb3FkpYHDx6kvXU3MzOLjY396aef/vjjD833l+VZHBCyjImLg6mp8phIOemkUCAbFyfkbTExMTY2NkFBQfXq1bty5Yo67TS2t7/s62useqhoUfzvf/jxRym6metFRESI9HcZ9/TpUzFzePHixZs3bw4YMKBq1aqTJ0/WUA9ZHnH+PFJS0KIFAMydi1mzsHIlHj2CsTGmT0e2pKpiTB3z5s2bO3cugAsXLoj6h1mzatWqpPHjp6R9aNAg/P47zM3V7GEeEBUVlT9//kwtNomJibl9+7aYPLx48WJKSsqIESMmT55cvHhxzfWT5XkcEDKWI0ycOHHlypUAlixZMn36dHWaCnj2zLJyZcOkpI8eHToUq1ZBB1O9vX6N0FDUqCFJY5GRkSkpKV+qbMkYY3nDkydPqlevnpiYWLx48ZcvX6q59+9+69ZVTp786CFra/z1F3RzieN//6FNG0laSk1NDQwMLFu2rCStMV3G0zuMaZ+vr6+q2kSXLl3UbK1M2bKG584p8wCp/O9/sLLCwYNqNp77+PigQwc8fy5JYwULFuRokDGW5/3888+ihkSXLl3UzwRT5ehRZQkWlcBANGyIPn2Qmqpm47lMYiKmT8eSJZI0pqenx9EgkwQHhIxpmUKhGDFihMhiYmNjU7lyZQkarV8f7u7p08mEh8PVVYLGc5eOHTF9OrZs0XY/GGMsd9izZ8+RI0fEcdeuXSVoUV8f//yDKlU+epAIe/YgTboUnWBsjEOH4OmJqChtd4WxDzggZEzLtq9fX87b+zsAUkwPftC7N1xcAHwUFurIHoPDh3Hnzod/jh2LBQskbD4xMXHevHmvXr2SsE3GGMsJYsLC/jdu3I+AHChQoEDz5s2labdgQRw8iHQrLAwMEBkpTfs5WWQk1q//8M+SJeHtjQIFJDzDmTNntm3bJmGDTNdwQMiYVoWGOv38807gCTBMqrFYlTlz0KcPVPuEa9RApUoYNgwDB8LHBy9f4tdflT/65RekSSyei6WkgAjx8ejYEYGBGjqJoaFh2roUjDGWZ5gMGnTy7dsDgA/QuUMHw3S7D9RRvjz+/htGRsp/6ulh3jwsWoShQ7F4MRQKTJ4Msfvd1xebN0t2Xu1KToZcjo0bsXy55k4SGRn59RIgjH2dvrY7wJhumzxZLzUVgB4wrXDhcra2UjYuk+GPP/D8Oby9Ub48nJxgbIwlSxATAwcHbNz4YWfd48e5byPHs2dwdcXGjQDQsyfWr8eMGSBCbCxq1cKECWjfHhcuIJM5RTNCLpcvXLhQ8mYZY0zLfHz03i8WrQVMbdxY4vYbN8amTRg4EDIZ3Nzwzz84eBAmJnB1xa5dePQIohp7TAxevpT41Nmga1fs3AlTU1y+jCNH0KgR1q+HpSVCQrB1K7p2RfHi6N9fM2eWdDSZ6R4OCBnTnvPnsXOn6l/l+/eXvoi8sTGOHsWTJ7C2xvjxWLYMAMzNUbUqAgJw7RpmzwaAe/ckPm/2+/dfVKmC8eMBwMkJs2cjLg6hoZoICBljLA9SKDBmjDIkA2BsXM3JSfqzDBiAOnUgkyEhAQEBytzXjo7KFSvz50NPD0FByAP16xcswJkz0NeHjw82b4anJ44e1XafGPs8DggZ05KkJIwcibR1X7p108iJChSAmHg0NkZCgvLBxEQYG6NaNYwZAwC3bmnk1Jp2/jyGDQOA69dhZfXhA6xXD76+ylhX08LCUKRIdpyIMcY0av16pK2C27YtzMw0cqJq1QDAz+/DJSkhAcbGADBiBIyMcP06rl/XyKk1bexY6OsjNBTff4/ixaGvDwB16mDKFNSoIVUBpK+JiYGRUfo044x9C+8hZExLli3D/fsf/lm0KNSo/Jsh3btj2TIkJeHuXbx8iVKlkC8fSpRAiRIfNnXkLk2aYPNmbN4MOzsULYrXr5WPv34NCwuNn/3UKTRrhvbtNX4ixhjTtNevMXPmR49oehVilSp48ACPHiExEStWwNERAIoXR4kSyL3b4daswebN+OUXGBkhPFz54Nu3KFhQ46eOiMD48fjuO1y8qPFzsTyHZwgZ04aAACxe/NEjP/6oHErUnA4dQIQJE1C0KHbtAgBV+rgOHXJrTKji6IjRo1GkCCIicPs25s/X+BkTEzF6tKbmdRljLDtNnPhRIQS5HB07avaMenr480+sXo2ICPTrh6ZN8eiR8jpoaYkfftDs2TVNLkedOli1CvXrY/Vq5WIcjcqXDwUK4NYtWFlp/Fwsz5FR2hVrjLHsMWoUNmz46JFDhzR+9c1j4uJw9y7q1QMALy/Ur4+ICBw7BhMTODgoVx8xxhj7pkePUKnSR4/Uq/fR8lGWEefPw94eenp4+xbBwahRAydOIDAQ9vaQpMIwYxrDASFj2rB5M4YP//BPfX3ExHAMwxhjTAuSk2Fh8VFJwAULMGuW1vrDGMtevIeQMc1bvx4REQAQHIxt2xAZCX9/tGuHihWVT6hVi6PBXCwqCm5ueJ+rnTHGcrTAwA8Jrl1dQYSDB9GhA2rXhrk5AOjro08fLXaQqevsWUybpu1OsNyEA0LGNO/KFWUuteho3LiBYcPQsyeOHsWsWRg2DIcO4do1bXeRZUlUFMaNQ7lyuHwZxYsrH1QlbScCL8FgjOU0kZEfMkufO4czZ3D0KHbswIkTqFkTGzfizRuUL6/NHrIs27cPtWph9GiUKqW8GKW9EqkuT4x9jANCxrLFnj344w/8/TeIEBGBOnUAYMAAPH7MWwdzMRMTlCiBn36CoSE2bsTw4VAo0K6d8qd//okdO7TaP8YY+xw/P/zxB/74Ay9f4uhRDB8OuRxFiqB5c1Spkh0pMZmGyGQYNQoVKuDJE3TujBs3sG0b9uxR/lR1eWLsY5xllLFsUbEiihaFiQlevfrocckr0bPsZGCAzp2xejW2bwcAV1f8+6+2+8QYY99SpAiqVAEAU1MgzZVIJuN1Dblb9+7o1Al//IFixRAejn790LOntvvEcgEOCBnLFra2KFkShQrh0iXkz49792Bjg337NF57kGmanx/s7JTH9erhwgVERWHYMAB48gQ//aTFrjHG2OdZWqJBAwAoWBCtW8PdHXZ2iI6GlxemTNF255h6YmNRrBgAFC6s3K6yeTNOnwaAmBhtdozlYBwQMqZ5tWsrc8aYm6NaNcyejUWLEBcHKytO45brFSuGR4+Ux8HBKF4cBQpg82YA2LmTN2wwxnKc/PlRtaryuH59tGqF168xZAj09LBiBczMtNo5pjaZDCkp0NeHQqGc7x02TJklqE0b7XaN5VhcdoIxxtSQkoJOnTB8OExN8dtv2LcPPXviv/+A9wHhgAHa7iJjjDGdsXUrbt1Cr17Yvx+VKsHAAMbGHwJCcXli7GMcEDLGmHoSEnD4MJKT0bYtChfGuXNo2hQAAgJAhO++03L3GGOM6ZQbN3D7NqpXh50dnj2Dnh6srQF8uDwx9jEOCBljjDHGGGNMR3HZCcYYY4wxxhjTURwQMsYYY4wxxpiO4oCQMcYYY4wxxnQUB4SMMcYYY4wxpqM4IGSMMcYYY4wxHcUBIWOMMcYYY4zpKA4IGWOMMcYYY0xHcUDIGGOMMcYYYzqKA0LGGGOMMcYY01EcEDLGGGOMMcaYjuKAkDHGGGOMMcZ0FAeEjDHGGGOMMaajOCBkjDHGGGOMMR3FASFjjDHGGGOM6SgOCBljjDHGGGNMR3FAyBhjjDHGGGM6igNCxhhjjDHGGNNRHBAyxhhjjDHGmI7igJAxxhhjjDHGdBQHhIwxxhhjjDGmozggZIwxxhhjjDEdxQEhY4wxxhhjjOkoDggZY4wxxhhjTEdxQMgYY4wxxhhjOooDQsYYY4wxxhjTURwQMsYYY4wxxpiO4oCQMcYYY4wxxnQUB4SMMcYYY4wxpqM4IGSMMcYYY4wxHcUBIWOMMcYYY4zpKA4IGWOMMcYYY0xHcUDIGGOMMcYYYzqKA0LGGGOMMcYY01EcEDLGGGOMMcaYjuKAkDHGGGOMMcZ0FAeEjDHGGGOMMaajOCBkjDHGGGOMMR3FASFjjDHGGGOM6SgOCBljjDHGGGNMR3FAyBhjjDHGGGM6igNCxhhjjDHGGNNRHBAyxhhjjDHGmI7igJAxxhhjjDHGdBQHhIwxxhhjjDGmozggZIwxxhhjjDEdxQEhY4wxxhhjjOkoDggZY4wxxhhjTEdxQMgYY4wxxhhjOur/BxdCVdX3XLcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<PIL.PngImagePlugin.PngImageFile image mode=RGB size=1200x200 at 0x1711E6902E0>"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "draw_rxn_figure(rxn_dict, db_smiles,novel_smiles)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "dG = -121.79 ± 100.57 kJ/mol\n"
     ]
    }
   ],
   "source": [
    "mu, std, rule_df1, rule_df2 = get_dG0(rxn_dict, 'R00801', pH, I, loaded_model, molsig_r1, molsig_r2, novel_decomposed_r1, novel_decomposed_r2, novel_mets)\n",
    "\n",
    "print(\"dG = %.2f ± %.2f kJ/mol\" % (mu, std))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bulk estimation of dG for a list of KEGG reactions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "KEGG_rxn_list = {\"R00010\" : \"C01083 + C00001 <=> 2 C00031\",\n",
    "                 \"R00303\" : \"C00092 + C00001 <=> C00031 + C00009\",\n",
    "                 \"R00304\" : \"C00103 + C00001 <=> C00031 + C00009\",\n",
    "                 \"R07294\" : \"C15524 + C00001 <=> C02137 + C00010\",\n",
    "                 \"R01252\" : \"C00148 + C00026 + C00007 <=> C01157 + C00042 + C00011\",\n",
    "                 \"R00406\" : \"C00091 + C00149 <=> C00042 + C04348\"\n",
    "                }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "R00010\n",
      "C01083 + C00001 <=> 2 C00031\n",
      "dG = -12.45 ± 3.49 kJ/mol\n",
      "R00303\n",
      "C00092 + C00001 <=> C00031 + C00009\n",
      "dG = -12.40 ± 3.30 kJ/mol\n",
      "R00304\n",
      "C00103 + C00001 <=> C00031 + C00009\n",
      "dG = -18.78 ± 3.37 kJ/mol\n",
      "R07294\n",
      "C15524 + C00001 <=> C02137 + C00010\n",
      "dG = -14.46 ± 31.43 kJ/mol\n",
      "R01252\n",
      "C00148 + C00026 + C00007 <=> C01157 + C00042 + C00011\n",
      "dG = -427.04 ± 41.12 kJ/mol\n",
      "R00406\n",
      "C00091 + C00149 <=> C00042 + C04348\n",
      "dG = -3.27 ± 4.37 kJ/mol\n"
     ]
    }
   ],
   "source": [
    "pH = 7 # any number between 0-14 \n",
    "I = 0.1    #min_value=0.0, max_value=0.5)\n",
    "\n",
    "for keys in KEGG_rxn_list:\n",
    "    kegg_rxn_string = KEGG_rxn_list[keys]\n",
    "    kegg_rxn_dict = parse_formula(kegg_rxn_string)\n",
    "    mu, std, rule_df1, rule_df2 = get_dG0(kegg_rxn_dict, keys, pH, I, loaded_model, molsig_r1, molsig_r2, [], [], [])\n",
    "    print(keys)\n",
    "    print(kegg_rxn_string)\n",
    "    print(\"dG = %.2f ± %.2f kJ/mol\" % (mu, std))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}