response.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "f70ab501-b03b-4ca6-a043-4e231d8008de",
   "metadata": {},
   "source": [
    "# Remove response"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "560555ea-4484-40b1-b890-0ad3ea30b167",
   "metadata": {},
   "source": [
    "## import function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "33efc07b-2cfa-48ab-868d-e97605845620",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as  mpl\n",
    "import csv\n",
    "import obspy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "id": "b9ebd26c-c99a-4367-a116-b457ade43660",
   "metadata": {},
   "outputs": [],
   "source": [
    "def find_response_file(path, idb = None, ifs = None, iFilter = None, iLC = None, mode = 'COUNTS',inverted=True, plot='True') :\n",
    "    '''\n",
    "    function to find metadata in the logfile\n",
    "    \n",
    "    path : path of the logfile - is a string\n",
    "    idb : force the gain - From 0 to 36 B - is a integer, float or string\n",
    "    ifs : force the sampling frequency - 50, 100, 125, 250, 500, 1000, 2000 or 4000 Hz - is a integer, float or string\n",
    "    iFilter : force the Final Filter Phase - could be LP or MP - is a string\n",
    "    iLC : force the DC Filter - could be 'OFF', '1' or 'DC' - is a string\n",
    "    \n",
    "    mode = 'COUNTS' or 'mV'\n",
    "    interted = True or False - if True, the polarization has been already inverted\n",
    "\n",
    "    '''\n",
    "    \n",
    "    \n",
    "    with open(path) as file :\n",
    "        readcsv = csv.reader(file, delimiter=' ')\n",
    "        line = 0\n",
    "        db = np.zeros(4)\n",
    "        for row in readcsv :\n",
    "            #print(row)\n",
    "            if row != [] :\n",
    "                if row[0] == 'Project' :\n",
    "                    pj = row[-1]\n",
    "                if row[0] == \"Serial\" :\n",
    "                    sr = row[-1]\n",
    "                if row[0] == \"Channel\" and row[2] == \"Gain\":\n",
    "                    c = int(row[1]) - 1\n",
    "                    db[c] = round(np.exp(0.1155*int(row[-1])))\n",
    "                if row[0] == \"Anti-alias\" :\n",
    "                    if row[-1] == 'Linear' :\n",
    "                        Filter = 'LP'\n",
    "                    else :\n",
    "                        Filter = 'MP'\n",
    "                if row[0] == \"Sample\" and row[1] == \"Rate\" :\n",
    "                    fs = int(row[-1])\n",
    "                if row[0] == \"Low\" and row[1] == \"Cutter\" :\n",
    "                    if int(row[-1]) == 0 :\n",
    "                        LC = 'Off'\n",
    "                    elif int(row[-1]) == 1 :\n",
    "                        LC = '1'\n",
    "                    else :\n",
    "                        LC = 'DC'\n",
    "                        \n",
    "    if idb != None :\n",
    "        db[0] = round(np.exp(0.1155*int(idb)))\n",
    "    elif idb == None:\n",
    "        idb = int(np.log(db[0])/0.1155)\n",
    "    if ifs != None :\n",
    "        fs = int(ifs)\n",
    "    if iFilter != None :\n",
    "        Filter = str(iFilter)  \n",
    "    if iLC != None :\n",
    "        LC = str(iLC)\n",
    "\n",
    "    from obspy.clients.nrl import NRL\n",
    "\n",
    "    nrl = NRL()\n",
    "\n",
    "    response_node = nrl.get_response(\n",
    "    datalogger_keys=['DTCC (manufacturers of SmartSolo','SmartSolo IGU-16HR3C', '%s dB (%i)'%(idb, db[0]), '%s'%fs, 'Linear Phase', '%s' %LC],\n",
    "        sensor_keys=['DTCC (manuafacturers of SmartSolo)', 'DT-SOLO' , '5 Hz', 'Rc=1850, Rs=430000'])\n",
    "    print(response_node)\n",
    "    response_node.response_stages[0].stage_gain = 76.7\n",
    "    response_node.recalculate_overall_sensitivity()\n",
    "\n",
    "    if plot == True :\n",
    "        response_node.plot(min_freq=0.001, output='VEL')\n",
    "\n",
    "    if inverted == True:\n",
    "        inv = 1\n",
    "    else :\n",
    "        inv = -1\n",
    "\n",
    "    if mode == 'COUNTS' :\n",
    "        paz_node = {'poles': response_node.get_paz().poles,\n",
    "        'zeros': response_node.get_paz().zeros,\n",
    "        'gain': 3355.4428*inv, # need to invert the trace for some reason\n",
    "        'sensitivity': response_node.instrument_sensitivity.value}\n",
    "\n",
    "    elif mode == 'mV' :\n",
    "        paz_node = {'poles': response_node.get_paz().poles,\n",
    "        'zeros': response_node.get_paz().zeros,\n",
    "        'gain': -1.,\n",
    "        'sensitivity': response_node.instrument_sensitivity.value}\n",
    "    else :\n",
    "        raise ValueError('mode should be either COUNTS or mV')\n",
    "    \n",
    "    \n",
    "    return paz_node "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "id": "534c76d0",
   "metadata": {},
   "outputs": [],
   "source": [
    "from obspy import read_inventory\n",
    "\n",
    "def simulate_mseed(stream, logfile_path, pre_filt = None, plot=True) :\n",
    "    \n",
    "    '''\n",
    "    Function to remove the response of minions based on the obspy remove_response function\n",
    "    \n",
    "    stream --> obspy stream object\n",
    "    \n",
    "    logfile_path : path of the logfile - is a string\n",
    "    pre_filt : prefilter of obspy remove_response function\n",
    "    output : output parameter of obspy remove_response function - could be 'VEL', 'DISP', 'ACC' or 'DEF'\n",
    "    zero_mean (bool) : If True, the mean of the waveform data is subtracted in time domain prior to deconvolution.\n",
    "    water_level (float) : Water level for deconvolution.\n",
    "    taper (bool) : If True, a cosine taper is applied to the waveform data in time domain prior to deconvolution.\n",
    "    taper_fraction (float) : Taper fraction of cosine taper to use.\n",
    "    \n",
    "    Parameters for the find_response_file function :\n",
    "    \n",
    "    idb : force the gain - From 0 to 36 B - is a integer, float or string\n",
    "    ifs : force the sampling frequency - 50, 100, 125, 250, 500, 1000, 2000 or 4000 Hz - is a integer, float or string\n",
    "    iFilter : force the Final Filter Phase - could be LP or MP - is a string\n",
    "    iLC : force the DC Filter - could be 'OFF', '1' or 'DC' - is a string\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    trace = stream[0]\n",
    "\n",
    "    #net = trace.stats.network\n",
    "    #station =  trace.stats.station\n",
    "    #ch =  trace.stats.channel\n",
    "    #loc = trace.stats.location\n",
    "\n",
    "    ifs = trace.stats.sampling_rate\n",
    "    \n",
    "    if trace.stats.mseed.encoding in ['STEIM2','STEIM1','INT32','INT16'] :\n",
    "        mode = 'COUNTS'\n",
    "    elif trace.stats.mseed.encoding == 'FLOAT32' :\n",
    "        mode = 'mV'\n",
    "    else :\n",
    "        raise ValueError('Encoding not recognized')\n",
    "    \n",
    "    if plot not in [True, False] :\n",
    "        raise ValueError('plot should be either True or False')\n",
    "\n",
    "    paz_node = find_response_file(logfile_path, ifs = ifs, plot=plot)\n",
    "\n",
    "    \n",
    "    if pre_filt == None :\n",
    "        pre_filt = (0.05, 0.06, ifs//2 - 5, ifs//2)\n",
    "\n",
    "    stream.simulate(\n",
    "        paz_remove=paz_node,\n",
    "        remove_sensitivity=True,\n",
    "        pre_filt=pre_filt)\n",
    "\n",
    "\n",
    "    return stream"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "id": "754fd72b-bc11-4844-9729-9a0792cb95e6",
   "metadata": {},
   "outputs": [],
   "source": [
    "def remove_resp_mseed(trace, logfile_path, pre_filt = None, output = \"VEL\", zero_mean=False, water_level=60,plot=False,\n",
    "                idb = None, ifs = None, iFilter = None, iLC = None) :    \n",
    "    '''\n",
    "    Function to remove the response of minions based on the obspy remove_response function\n",
    "    \n",
    "    --> Needs improvements, currently not working\n",
    "\n",
    "    trace --> obspy trace\n",
    "    \n",
    "    logfile_path : path of the logfile - is a string\n",
    "    pre_filt : prefilter of obspy remove_response function\n",
    "    output : output parameter of obspy remove_response function - could be 'VEL', 'DISP', 'ACC' or 'DEF'\n",
    "    zero_mean (bool) : If True, the mean of the waveform data is subtracted in time domain prior to deconvolution.\n",
    "    water_level (float) : Water level for deconvolution.\n",
    "    taper (bool) : If True, a cosine taper is applied to the waveform data in time domain prior to deconvolution.\n",
    "    taper_fraction (float) : Taper fraction of cosine taper to use.\n",
    "    \n",
    "    Parameters for the find_response_file function :\n",
    "    \n",
    "    idb : force the gain - From 0 to 36 B - is a integer, float or string\n",
    "    ifs : force the sampling frequency - 50, 100, 125, 250, 500, 1000, 2000 or 4000 Hz - is a integer, float or string\n",
    "    iFilter : force the Final Filter Phase - could be LP or MP - is a string\n",
    "    iLC : force the DC Filter - could be 'OFF', '1' or 'DC' - is a string\n",
    "    \n",
    "    '''\n",
    "    from obspy.core.inventory.response import Response\n",
    "\n",
    "\n",
    "    #net = trace.stats.network\n",
    "    #station =  trace.stats.station\n",
    "    #ch =  trace.stats.channel\n",
    "    #loc = trace.stats.location\n",
    "\n",
    "    ifs = trace.stats.sampling_rate\n",
    "    \n",
    "    if trace.stats.mseed.encoding in ['STEIM2','STEIM1','INT32','INT16'] :\n",
    "        mode = 'COUNTS'\n",
    "    elif trace.stats.mseed.encoding == 'FLOAT32' :\n",
    "        mode = 'mV'\n",
    "    else :\n",
    "        raise ValueError('Encoding not recognized')\n",
    "\n",
    "    paz_node = find_response_file(logfile_path, ifs = ifs)\n",
    "    resp = Response.from_paz(paz_node['poles'], paz_node['zeros'], paz_node['gain'], paz_node['sensitivity'], input_units=mode, output_units='M/S')\n",
    "    print(resp)\n",
    "    pre_filt = (0.05, 0.06, ifs//2 - 5, ifs//2)\n",
    "\n",
    "    trace.remove_response(resp, output=output, pre_filt=pre_filt, zero_mean=zero_mean, water_level=water_level, plot=plot)\n",
    "\n",
    "    return trace"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a8101d92-d5f0-4773-bbbb-4cba476e9c06",
   "metadata": {},
   "source": [
    "## Get metadata information from log file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "id": "905c28f7-24b5-4a0f-abc1-4f08d00fc280",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Channel Response\n",
      "\tFrom M/S (Velocity in Meters per Second) to COUNTS (Digital Counts)\n",
      "\tOverall Sensitivity: 4.9001e+06 defined at 25.000 Hz\n",
      "\t7 stages:\n",
      "\t\tStage 1: PolesZerosResponseStage from M/S to V, gain: 76.6\n",
      "\t\tStage 2: ResponseStage from V to V, gain: 64\n",
      "\t\tStage 3: CoefficientsTypeResponseStage from V to COUNTS, gain: 1000\n",
      "\t\tStage 4: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 5: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 6: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 7: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "{'poles': [(-22.2111-22.2178j), (-22.2111+22.2178j)], 'zeros': [0j, 0j], 'gain': 3355.4428, 'sensitivity': 4906497.352699268}\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\julie\\AppData\\Local\\Temp\\ipykernel_13596\\751132783.py:59: ObsPyDeprecationWarning: \n",
      "        DEPRECATED\n",
      "\n",
      "        Direct access to online NRL is deprecated as it will stop working when\n",
      "        the original NRLv1 gets taken offline (Spring 2023), please consider\n",
      "        working locally with a downloaded full copy of the old NRLv1 or new\n",
      "        NRLv2 following instructions on the\n",
      "        `NRL landing page <https://ds.iris.edu/ds/nrl/>`_.\n",
      "        \n",
      "  nrl = NRL()\n"
     ]
    }
   ],
   "source": [
    "path = 'D:/GitHub/ULB-SeismoLab/Minions-Toolbox/ExampleData/Log/453001250/20231110162230/DigiSolo.LOG'\n",
    "paz_node = find_response_file(path)\n",
    "print(paz_node)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "797b7ce4-f75a-4e83-aed0-0dec67e51f84",
   "metadata": {},
   "source": [
    "## Remove response"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "id": "edf86ed1",
   "metadata": {},
   "outputs": [],
   "source": [
    "file_path = 'D:/GitHub/ULB-SeismoLab/Minions-Toolbox/ExampleData/Miniseed/'\n",
    "file = '453004444..0.0001.2023.11.10.11.09.45.000.N.miniseed'\n",
    "st = obspy.read(file_path+file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "id": "d7378d30-40a0-4ced-88c6-fdb26f04ec89",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "         network: BE\n",
       "         station: 04444\n",
       "        location: BE\n",
       "         channel: GPN\n",
       "       starttime: 2023-11-10T11:09:45.000000Z\n",
       "         endtime: 2023-11-10T12:56:33.000000Z\n",
       "   sampling_rate: 1000.0\n",
       "           delta: 0.001\n",
       "            npts: 6408001\n",
       "           calib: 1.0\n",
       "         _format: MSEED\n",
       "           mseed: AttribDict({'dataquality': 'D', 'number_of_records': 3749, 'encoding': 'STEIM2', 'byteorder': '>', 'record_length': 4096, 'filesize': 15355904})"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "file_path = 'D:/GitHub/ULB-SeismoLab/Minions-Toolbox/ExampleData/Miniseed/'\n",
    "file = '453004444..0.0001.2023.11.10.11.09.45.000.N.miniseed'\n",
    "st = obspy.read(file_path+file)\n",
    "trace = st[0]\n",
    "trace.stats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "id": "a682df5f-58f9-4fd1-b061-94f30165dd44",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\julie\\AppData\\Local\\Temp\\ipykernel_13596\\751132783.py:59: ObsPyDeprecationWarning: \n",
      "        DEPRECATED\n",
      "\n",
      "        Direct access to online NRL is deprecated as it will stop working when\n",
      "        the original NRLv1 gets taken offline (Spring 2023), please consider\n",
      "        working locally with a downloaded full copy of the old NRLv1 or new\n",
      "        NRLv2 following instructions on the\n",
      "        `NRL landing page <https://ds.iris.edu/ds/nrl/>`_.\n",
      "        \n",
      "  nrl = NRL()\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Channel Response\n",
      "\tFrom M/S (Velocity in Meters per Second) to COUNTS (Digital Counts)\n",
      "\tOverall Sensitivity: 4.90185e+06 defined at 50.000 Hz\n",
      "\t7 stages:\n",
      "\t\tStage 1: PolesZerosResponseStage from M/S to V, gain: 76.6\n",
      "\t\tStage 2: ResponseStage from V to V, gain: 64\n",
      "\t\tStage 3: CoefficientsTypeResponseStage from V to COUNTS, gain: 1000\n",
      "\t\tStage 4: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 5: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 6: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n",
      "\t\tStage 7: CoefficientsTypeResponseStage from COUNTS to COUNTS, gain: 1\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "No artists with labels found to put in legend.  Note that artists whose label start with an underscore are ignored when legend() is called with no argument.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAHQCAYAAACSr3/KAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAAB5l0lEQVR4nO3dd3gUVdsG8Ht7eiCVdEJLoRN6B+kdRVCUqr4oKiAqir5YEPFVFLGBgBT1syAqiIhApPceIIQaEtIJSSA9W+f7I2QlhIWU3Z3s5v5d18LuzOzMs082u0/OnDlHIgiCACIiIiIyK6nYARARERHZIxZZRERERBbAIouIiIjIAuRiB0BERLVTYmIifv/9d1y9ehVarVbscKrN1dUVPXr0wMCBA+Hg4CB2OFSHSNjxnYiI7vbuu+/inXfegUqlQlhYGFQqldghVYsgCMjKykJiYiL8/Pywc+dOhIeHix0W1REssoiIqJwffvgBTz75JN566y288sorcHV1FTukGjt37hwee+wx5OTkICEhAUqlUuyQqA5gkUVEROX0798fer0eO3fuFDsUszpz5gxat26Nv/76C0OGDBE7HKoD2PGdiIiMBEHAvn37MHToULFDMbuWLVsiMDAQ+/btEzsUqiNYZBERkZFer4darYanp6fYoZidRCKBp6cn8vPzxQ6F6ggWWUREVIFEIhE7BIuw19dFtROLLCIiqpR33nkHEomk3K1BgwYWPeZvv/2GyMhIqFQqREZGYsOGDeXW7927F8OHD4e/vz8kEgk2btxo0XiIqoJFFhERVVrz5s2Rnp5uvJ09e7ba+1q7di169+5tcv2hQ4cwbtw4TJgwAadPn8aECRMwduxYHDlyxLhNYWEhWrdujS+//LLacRBZCossIiKqNLlcjgYNGhhv3t7exnUajQZz5sxBQEAAnJ2d0alTJ+zevbvax1qyZAn69++PuXPnIjw8HHPnzsVDDz2EJUuWGLcZPHgwFixYgIcffrgGr4rIMlhkERFRpV2+fBn+/v4IDQ3FY489hqtXrxrXTZkyBQcOHMDPP/+MM2fO4NFHH8WgQYNw+fLlah3r0KFDGDBgQLllAwcOxMGDB2v0GoishUUWERFVSqdOnfDdd99h27ZtWLlyJTIyMtC1a1dkZ2cjPj4eP/30E9avX48ePXqgcePGeOWVV9C9e3esWbOmWsfLyMiAr69vuWW+vr7IyMgwx8shsjjOXUhERJUyePBg4/2WLVuiS5cuaNy4Mb799lsEBQVBEAQ0a9as3HPuHA4iKSkJkZGRxnU6nQ5arRYuLi7GZU8++SS+/vpr4+O7rwYUBIFXCJLNYJFFRETV4uzsjJYtW+Ly5csICAiATCbDiRMnIJPJym1XVkT5+/sjJibGuPz333/Hb7/9hh9++MG4zM3NzXi/QYMGFVqtMjMzK7RuEdVWLLKIiKha1Go1zp8/jx49eqBt27bQ6/XIzMxEjx497rm9XC5HkyZNjI99fHzg6OhYbtmdunTpgujoaLz00kvGZdu3b0fXrl3N+0KILIRFFhERVcorr7yC4cOHIzg4GJmZmViwYAHy8vIwadIkhISE4IknnsDEiRPxySefoG3btsjKysLOnTvRsmXLas0VOHPmTPTs2RMffvghRo4ciT/++AP//PMP9u/fb9ymoKAAV65cMT5OSEhATEwMPDw8EBwcbJbXTVRd7PhORESVkpKSgscffxxhYWF4+OGHoVQqcfjwYYSEhAAA1qxZg4kTJ+Lll19GWFgYRowYgSNHjiAoKKhax+vatSt+/vlnrFmzBq1atcLatWuxbt06dOrUybjN8ePH0bZtW7Rt2xYAMHv2bLRt2xZvvfVWzV8wUQ1JBEEQxA6CiIhqB51OB4VCgbVr12LSpElih2N2bdu2Rbdu3Th4KVkFW7KIiIiILIBFFhERGUmlpV8LWq1W5EgsQ6PRQC5nd2SyDhZZRERkJJVKERgYiNjYWLFDMbvCwkIkJCRUu48YUVWxyCIionJGjRqFn3/+GWlpaWKHYlbLly9HcXExRo8eLXYoVEew4zsREZWTmJiIHj16QCqV4qmnnkJUVBRUKpXYYVWLIAjIysrC5s2b8eOPP+KFF17AF198IXZYVEewyCIiogri4+Px9ttvY9OmTcjPzxc7nBpr3rw5pkyZgtmzZ3NaHrIaFllERGSSRqNBVlYWNBqN2KFUm5ubGzw8PMQOg+ogFllEREREFsCO70REZHY5OTmo7+EJqVRqttudcxgS2QIOFkJERGaXlpaGWzdz4NpuGBTeDWu8v/zjG3H6zJmaB0ZkRSyyiIjI7MLDwxHSMBRZhTlw7f9sjfalvZmOm9HLMGL4y2aKjsg6eLqQiIjMTi6X4+235qHw4kHkn9yM6nb/1Rfews0tn8LT0xPTpk0zc5RElsWWLCIisoiJEyfi9OnT+Oyzz6C5Hg+P/tMhkSsq/Xx1xhXc/GMhXBXApj82wtHR0YLREpkfry4kIiKL+vbbb/HMf6ZB5uIJ587j4Ny8NyRSmcntdXmZyDu8HoVnotGmTRv8sXEDAgMDrRgxkXmwyCIiIouLjY3Fm//9Lzb98QccPANKi63IXuWKLV1eFnIP/4KiM9vh5uaG11+bg5kzZ8LBwUHEyImqj0UWERFZzalTp/DW229j859/QuXuBYfIh6D0bYyiuF0ounLUWFy98MILcHFxETtcohphkUVERFYXExOD5cuX47v/+z8UFRSgectWeG7afzBx4kS4urqKHR6RWbDIIiIi0RQWFiI1NRVNmzblnIJkd+p8kWUwGJCWlgZXV1f+ghMRkd0SBAH5+fnw9/eHVMoRnKyhzg/hkJaWhqCgILHDICIisork5GRerWkldb7IKjv3n5ycDDc3txrvT6vVYvv27RgwYAAUisqPB1MXMDemMTemMTemMTem2XpuijQ6dHx/BwDg6JsPwUlZ86/rvLw8BAUFsc+bFdX5IqvsFKGbm5vZiiwnJye4ubnZ5C+2JTE3pjE3pjE3pjE3ptl6buQaHaQqJwCl30/mKLLKsGuM9fCkLBEREZEF1PmWLCIiotpGKZPiq/HtjPfJNrHIIiIiqmXkMimGtvITOwyqIZbHRERERBbAliwiIqJaRqc3YNu56wCAgc19IecpQ5vEIouIiKiW0egNeP7HkwCAuPkDWWTZKP7UiIiIiCyALVlERERUZTdv3sQPP/yAP/74AykpKdDr9WKHZFXOzs7o1KkTnnjiCfTo0eOe27DIIiIiontKTU3F5s2bUVRUhJYtW6Jfv34AgBs3bqBv3764cOEC+vXrh8GDB0MurzslhSAIyM3NxdatW7F8+XJ8+eWXeP755ytsV3cyQkRERA9kMBiwfPlyrFn7LY4dPQKJTAaJRAZXFyfcunkTAPDKK6/g+vXrOHPmDCIiIkSOWDyCIOCll17CCy+8gIEDB6JJkybl1rNPFhEREQEA8vPzMfrhhzH9+ecRd0sCz2EvI/CFH+DeYwK0Wh0AQKPRYOPGjZg+fXqdLrCA0imKFi5cCCcnJ6xfv77CerZkEREREa5cuYJhw0cgPjEJ3g/Pg1OTjvfcLjk5GXl5eSb7IdU1Tk5O6NChA86ePVthHYssIiKiWkYhk2LRmFbG+5Z27do1dOjYCSUyJ/g8+TEUnkEmty0pKQFQWlxQKUdHR2Ne7sQii4iIqJZRyKR4tL3pQsec1Go1Rj/8CIqhgPf4jyBzdLPKce2JRCK553L2ySIiIqrDZs6cidNnzqD+iLlmLbA++OADSCQSzJo1677bffXVV4iIiICjoyPCwsLw3XffmS0GU86ePYtevXrB0dERAQEBmD9/PgRBKLeNWq3Gm2++iZCQEKhUKjRu3BirV6+u0nHYkkVERFTL6PQG7L18AwDQs6m3xUZ8//7777F8+XJ4DHwBqgZNHvyESjp27BhWrFiBVq1a3Xe7ZcuWYe7cuVi5ciU6dOiAo0eP4plnnkH9+vUxfPjwah07MTERoaGhFYqmMnl5eejfvz/69OmDY8eO4dKlS5g8eTKcnZ3x8ssvG7cbO3Ysrl+/jlWrVqFJkybIzMyETqerUixsySIiIqplNHoDpq49jqlrj0OjN1jkGMnJyZj27HNwafEQXFoPNNt+CwoK8MQTT2DlypWoX7/+fbf9/vvvMW3aNIwbNw6NGjXCY489hqeeegoffvhhue3WrFmDiIgIODg4IDw8HEuXLq12fD/88ANKSkqwdu1atGjRAg8//DDeeOMNLF682FiYbd26FXv27MGWLVvQr18/NGzYEB07dkTXrl2rdCwWWUREVI4gCDAYBOgNAnR6A7R6AzS6mt20+urddNW86Q0CDAKgv/06KnszVPMmCNW7iWnmrFnQyx1Qv980k32KquP555/H0KFDjQOX3o9arYaDg0O5ZY6Ojjh69Ci0Wi0AYOXKlXjzzTfx/vvv4/z581i4cCHmzZuHb7/9tlrxHTp0CL169YJKpTIuGzhwINLS0pCYmAgA2LRpE9q3b4+PPvoIAQEBaNasGV555RUUFxdX6Vg8XUhEVIvoDQKu55Ugq0CN7EINsgs0yClU42aRFkVqHYo0+ts3HQrUOmTckGFZwiFjkaA1GKDTC9DdLpDK7gsQIAiAAAACyj0WBOH2/+K+dvOT46XD0WIHcV+B9R0xq18zjIkKtOpx9+3bhw2//w6v4a9CqjLfVYI///wzTp48iWPHjlVq+4EDB+Kbb77BqFGj0K5dO5w4cQKrV6+GVqtFVlYW/Pz88N577+GTTz7Bww8/DAAIDQ1FXFwcli9fjkmTJlU5xoyMDDRs2LDcMl9fX+O60NBQXL16Ffv374eDgwM2bNiArKwsTJ8+HTk5OVXql8Uii4hIBMUaPc5n5OFcWh4uZuThWnYRknOKkHqrGFp9VaodCVCYb7E4ybJSbhbjlfWn4eYgx4DmDax2XIPh9ilIM7ZgJScnY+bMmdi+fXuF1ilT5s2bh4yMDHTu3BmCIMDX1xeTJ0/GRx99BJlMhhs3biA5ORlPPfUUnnnmGePzdDod3N3djY+bN2+Oa9euAYCxhdDFxcW4PiQkBOfOnTM+vrvlruw5ZcsNBgMkEgl++OEH43EWL16MMWPG4KuvvoKjo2OlXh+LLCIiK8gp1ODw1Wzj7UpmAQwmaim5VAIvFxU8nJXwdFHC01mJek5KuKjkcFLJ4KSQwUklh4MMiD19Cp07dYCDUgGFTAqZVAKFVAq5TAK5VAK5TAq5tPSLQyIp/RKRlN2H5Pb/AO56fPd2kFT/+7jaLWTVfJ4AAVqtDtHR0ejfvz8Uisp91VU3zuq+PJ3BgCX/XMaPR5Lw4dYL6BfhC6nUfEXP/fTq1QsjRo7E1l2r4diovVlas06cOIHMzExERUUZl+n1euzduxdffvkl1Go1ZDJZuec4Ojpi9erVWL58Oa5fvw4/Pz+sWLECrq6u8PLywo0bpZ3/V65ciU6dOpV77p372rJli/H0YmpqKnr37o2YmBjjeoVCYbzfoEEDZGRklNtXZmYmgH9btPz8/BAQEFCukIuIiIAgCEhJSUHTpk0rlRMWWUREFhJ/owBbYzOw/VwGTqfkVljv5aJCc383RPi5oZGXM4I8nBDs6YQGbg6QVeLLVqvVQkgS0KOJV7kvEQK0WgmcFUA9J0Wtzs3cweH4MyYN8TcKcTLpJto39LDasT//7DOEhUfg5s6V8Bw8s8b7e+ihhyqMej5lyhSEh4fjtddeq1Bg3UmhUCAwsPSU6c8//4xhw4ZBKpXC19cXAQEBuHr1Kp544gmTzw8JCTHeL5uo+u55BMt06dIFb7zxBjQaDZRKJQBg+/bt8Pf3N55G7NatG9avX4+CggJji9ilS5cglUqNcVYGiywiIjPKzCvBbydT8fvJFFzOLCi3LszXFV0ae6JzIw+0C64PH7fKnVIh++XqoEDfCB/8EZOG6LjrVi2yQkJC8PWypZgyZQpUARFwaTWgRvtzdXVFixYtyi1zdnaGp6encfncuXORmppqHAvr0qVLOHr0KDp16oSbN29i8eLFiI2NLdep/Z133sGMGTPg5uaGwYMHQ61W4/jx47h58yZmz55d5TjHjx+Pd999F5MnT8Ybb7yBy5cvY+HChXjrrbeMpwvHjx+P9957D1OmTMG7776LrKwsvPrqq5g6dWqlTxUCLLKIiGpMbxCw4/x1/HI8Gbsu3oD+9nlAhUyCLo29MKh5A/SL9IGPK4sqqqhveGmRdTghx7hMIZNi/sjmxvuWMnnyZOzfvx9rv1sOpW9jKH0bW+xYAJCeno6kpCTjY71ej08++QQXL16EQqFAnz59cPDgwXId059++mk4OTlh0aJFmDNnDpydndGyZcsHDnJqiru7O6Kjo/H888+jffv2qF+/PmbPnl2uYHNxcUF0dDRefPFFtG/fHp6enhg7diwWLFhQpWOxyCIiqqYCtQ7rjydj9YEEJOf8e2l3VEh9PBoViMEt/eDuWHtPVVHt0DaodCyp82l5UOv0UMllUMikmNiloVWO/+WXX+L4iZO4sOl/8J7wKWQOLg9+UiXt3r273OO1a9eWexwREYFTp049cD/jx4/H+PHjK3XMhg0bPnB4jJYtW2Lv3r333SY8PBzR0TW7OpVFFhFRFaXnFmPNgUT8dDQJ+SWlI0DXc1JgbPsgjG0fiCY+riJHSLYkyMMRHs5K5BRqcCE9H62D6ln1+A4ODtjw+29o2y4KWT/OgceoN6HwCLBqDLbOVFHHIouIqJJSbxVj6a4r+OV4snGYhUbezpjaLRSPtAuEo9J0x14iUyQSCSL8XHHgSjYuZxagdVA96A0Cjt4+fdgx1KNSF0LURGhoKI4cPoRhw0fg2v+9jPrDXoVjo6h7bls2iGdVB+a0Z8XFxfcc3d4uRnxPSEhAnz59EBkZiZYtW6KwsFDskIjIjqTcLMIbG86i96Jd+OFIErR6AR1DPbB6cnv881IvPNk5hAUW1Ugjr9JTdFdvlF4sodbp8fjKw3h85WGodXqrxBAWFobjx46if99euPHbu8j64wMUxu2BQV1UbrugoCC4urpi//79VomrtisuLsbx48fRvHnzCuvsoiVr8uTJWLBgAXr06IGcnJxyQ+UTEVVX6q1ifLnzCn498W/LVZdGnpjZryk6N/IUOTqyJ428nQGUDvshJnd3d/y5aRO+/PJLrPn2O8T8uQhSuQISiRTOTqUXbqhUKowcORJfffUVxo0bh7CwMFFjFpMgCJg3bx4KCwvx6KOPVlhv80XWuXPnoFAo0KNHDwCAh4f1Ln8lIvt0s1CDr3ZdwXeHrhkn5+3a2BMzH2qKTiyuyAIaeZe1ZIl/JkYqlWLGjBmYMWMGrl27hs2bN6OoqAgtW7Y0bvPxxx+jT58+aNWqFQYOHIjw8PBaPR6ZuQmCgNzcXGzbtg3x8fFYsmQJmjVrVmE70YusvXv3YtGiRThx4gTS09OxYcMGjBo1qtw2S5cuxaJFi5Ceno7mzZtjyZIlxqLq8uXLcHFxwYgRI5CSkoIxY8bgjTfeEOGVEJGtK9bosfpAAr7eHY98dWmH9k6hHnh5QBg6hvIPOLKcYI/SEddTbxWLPnH0nUJCQvD8889XWO7r64t9+/bh+++/xx9//IENGzZAp9OJEKF4nJ2d0bt3b6xcuRJ9+vS55zaiF1mFhYVo3bo1pkyZgkceeaTC+nXr1mHWrFlYunQpunXrhuXLl2Pw4MGIi4tDcHAwtFot9u3bh5iYGPj4+GDQoEHo0KED+vfvL8KrISJbpNMbsP5ECpb8cwnX89QAgPAGrnh9cDh6NfOuMM8Zkbn5uZeeiivS6JFXrINCXvvfc56enpg1a1a1x6uqC0QvsgYPHozBgwebXL948WI89dRTePrppwEAS5YswbZt27Bs2TJ88MEHCAwMRIcOHRAUFAQAGDJkCGJiYkwWWWq1Gmq12vg4Ly8PQOn0FGXzHtVE2T7MsS97w9yYxtyYZsncCIKA6POZ+CT6Mq5mlXbuDajngJceaoLhrfwglUpq9V/nfN+YZmu5kQGo76TAzSItkrPzEeTx76jiWq0WWknNW7dsJRf2RPQi6340Gg1OnDiB119/vdzyAQMG4ODBgwCADh064Pr167h58ybc3d2xd+9eTJs2zeQ+P/jgA7z77rsVlm/fvh1OTjWfILNMTQcws2fMjWnMjWnmzk18HrDpmgyJBaUtBs5yAQMCDejuWwB5Wgy2psWY9XiWxPeNabaUGyfIcBMSbNqxH03cBJR9RW/bth0qM1y8WlRU9OCNyKxqdZGVlZUFvV5vnBW7jK+vr3EGbblcjoULF6Jnz54QBAEDBgzAsGHDTO5z7ty55YbOz8vLQ1BQEAYMGAA3N7cax6zVau+Y+b3udAKsDObGNObGNHPnJvlmERZtu4y/z10HADgqpJjStSGe7h4CVwfbyj3fN6bZYm7+yDmF1Is3ENisBYa2DUBW/WsAgKGdQ6CU13zEpbIzN2Q9tbrIKnN3fwhBEMote9ApxzupVKp7DvGgUJh3pnZz78+eMDemMTem1TQ3+SVaLN0dj1X7E6DRGSCVAOM6BOGlfs1sfqJmvm9Ms6Xc+NcvPUWYma+Fs6MK0/tUvFqtJmwlD/akVhdZXl5ekMlkxlarMpmZmRVat4iI7kVvELD+eDI+3n4JWQWl/TG7NfHEf4dGIsKv5q3XRObi7VJa7GcXqh+wJdmKWl1kKZVKREVFITo6GqNHjzYuj46OxsiRI0WMjIhswcH4LLy3+TzOp5eeJgn1csYbQyLQL8KHVwxSrePhogQAZBdooDcIiE3NBQC0CHC3+LQ6ZBmiF1kFBQW4cuWK8XFCQgJiYmLg4eGB4OBgzJ49GxMmTED79u3RpUsXrFixAklJSXj22WdFjJqIarPErEIs3HIe2+NK+125Ocgxs18zTDBT3xYiS/Byvl1kFWqg1ukx8qsDAIC4+QPhpBT965qqQfSf2vHjx8sN4lXWKX3SpElYu3Ytxo0bh+zsbMyfPx/p6elo0aIFtmzZgpCQELFCJqJaKrdYiy93Xsbag4nQ6gXIpBI80SkYs/o1g8ftLzCi2qrsPZpTqBE5EjIX0Yus3r17P3B02+nTp2P69OlWioiIbI1Ob8BPx5LxafQl4xdUr2be+O/QCDT1dRU5OqLK8XQpvSirrO8g2T7RiywioprYc+kG3v8rDpeul06s28THBW8OjUCfMB+RIyOqGs/bLVn5JTpodAaRoyFzYJFFRDbpSmYB3v8rDrsu3gAA1HNSYHb/Zni8YzAUMva7Itvj7qiATCqB3iDgZiFHZ7cHLLKIyKbcLNTgsx2X8f3ha9AbBMilEkzq2hAz+jaFuxPHASLbJZVKUN9JiawCNXKKeMrQHrDIIiKboNUb8P2RBHy24zJyi0v/yu8X4Ys3hoSjkbeLyNERmUd9JwWyCtTG9zjZNhZZRFSrCYKA2JsSfPblQeMkzuENXDFvWCS6NfESOToi83JzLG2NLVLrMfOhpgAAuZSnv20ViywiqrUuZuRj/p/ncCBeBqAIns5KvDwgDOM6BHFwRrJLbg6lX8tFWj1e6m/eaXXI+lhkEVGtk12gxuLoS/jpaBIMAiCTCJjaLRQv9msGNxubxJmoKspasvJ4utAusMgiolpDrdPj24OJ+GLHFeSrdQCAgZE+6KhKw8SBzTjBLdk919stWbnFWly6ng8AaOLtAilbbm0SiywiEp0gCNh27jo++Ps8rmWX9rtq7u+Gt4ZFol2QG7ZsSRM5QiLrKGupzSnUYMCnewFwWh1bxp8aEYkqNjUX722Ow5GEHACAt6sKcwaG4ZF2gZBKJdBqedqE6o6y04X5JXzf2wMWWUQkisz8Eny87SLWn0iBIAAquRT/6dkIz/ZqDGcVP5qobiprycov0YscCZkDP8mIyKpKtHqs2p+ApbuuoFBT+kUyorU/XhscjoB6jiJHRyQuN8fSr+U8tmTZBRZZRGQVgiBg85l0/O/vC0i9VQwAaBNUD/OGRSIqpL7I0RHVDsaWLF5daBdYZBGRxZ1OvoX3Nsfh+LWbAAA/dwe8NigcI1r786opojsYh3BgS5ZdqFaRdevWLfz666+Ij4/Hq6++Cg8PD5w8eRK+vr4ICAgwd4xEZKPSc4uxaOtF/H4qFQDgqJDh2V6N8Z+ejeColIkcHVHt43K7P2KRhn2y7EGVi6wzZ86gX79+cHd3R2JiIp555hl4eHhgw4YNuHbtGr777jtLxElENqRYo8fyvfH4ek88SrQGAMDD7QIwZ2A4Grg7iBwdUe3lrCr946NIo8N/ejYCwGl1bFmVi6zZs2dj8uTJ+Oijj+Dq6mpcPnjwYIwfP96swRGRbTEYBPxxOhUf/n0RGXklAID2IfUxb1gkWgfVEzc4IhtQNh6WzgC8MiAMSjkLLFtW5SLr2LFjWL58eYXlAQEByMjIMEtQRGR7TlzLwfw/43A6JRcAEFDPEXOHhGNoSz9IJOx3RVQZTnecRi/S6KCUK0WMhmqqykWWg4MD8vLyKiy/ePEivL29zRIUEdmOlJtF+HDrRfx5unRUdmelDM/3bYKp3ULhoGC/K6KqUMikUMqk0OgNuJJZAF83BwTUc+QFIjaqykXWyJEjMX/+fPzyyy8AAIlEgqSkJLz++ut45JFHzB4gEdVOBWodlu2+gm/2JUCtM0AiAcZGBeHlgc3g48p+V0TV5aSSQVNkwJivDwHgtDq2rMo/tY8//hhDhgyBj48PiouL0atXL2RkZKBLly54//33LREjEdUiBoOAX0+mYNG2i7iRrwYAdG7kgXnDItHc313k6Ihsn7NSjltFHMLBHlS5yHJzc8P+/fuxc+dOnDx5EgaDAe3atUO/fv0sER8R1SKHr2bjvc1xOJdW2mUgxNMJbwyJwIBIX/a7IjITJw5vYjeq3f7Yt29f9O3b15yxEFEtdS27EB9suYCt50ovbnF1kGNG36aY2DUEKjm/EIjMyYlzd9qNSv0kP//880rvcMaMGdUOhohql7wSLb7aeQVrDiRCozdAKgHGdwrGS/2awdNFJXZ4RHbJmS1ZdqNSRdann35a7vGNGzdQVFSEevXqASgdAd7JyQk+Pj4ssojsgE5vwLrjyVi8/RKyCzUAgB5NvfDfoZEIa+D6gGcTUU2wk7v9qNRPMiEhwXj/xx9/xNKlS7Fq1SqEhYUBKB2+4ZlnnsG0adMsEyURWc3+y1lY8FccLmTkAwAaeTvjv0Mj0CfMh/2uiKygbNR3sn1VLpfnzZuHX3/91VhgAUBYWBg+/fRTjBkzBk888YRZAyQi64i/UYCFf53HjguZAAB3RwVe6tcUT3QOgULGUaeJrKWsJat1oDtaBdaDjGNk2awqF1np6enQaiteWqrX63H9+nWzBFUdRUVFiIiIwKOPPoqPP/5YtDiIbE1ukRaf7biM7w4lQmcQIJdK8GTnEMzq1xT1nDjaNJG1lfXJ6tzYE3MHR4gcDdVElYushx56CM888wxWrVqFqKgoSCQSHD9+HNOmTRN1GIf3338fnTp1Eu34RLZGqzfgxyNJ+PSfS8YxeR4K98EbQyPQ2NtF5OiI6q6yIRyK1HqRI6GaqnKRtXr1akyaNAkdO3aEQqEAAOh0OgwcOBDffPON2QOsjMuXL+PChQsYPnw4YmNjRYmByJbsupiJBZvjEH+jEADQzNcF/x0aiZ7NODUWkdhUt6ejyi3WILtADQ9nJftD2qgqd7Tw9vbGli1bcOHCBaxfvx6//PILzp8/jy1btsDHx6fKAezduxfDhw+Hv78/JBIJNm7cWGGbpUuXIjQ0FA4ODoiKisK+ffvKrX/llVfwwQcfVPnYRHXNpev5mLj6KKasOYb4G4XwcFZiwagW2DKjBwssolpCJS/9at50Oh1RC/5BsZYtWraq2teJNmvWDM2aNatxAIWFhWjdujWmTJlyz7kP161bh1mzZmHp0qXo1q0bli9fjsGDByMuLg7BwcH4448/jLEcPHjwgcdTq9VQq9XGx2WTXWu12nv2Nauqsn2YY1/2hrkxzdK5ySnU4POd8fj5eAr0BgEKmQQTOwdjeq9GcHNUQDDooTXUzg9yvm9MY25Ms+XcKO5q/tBqtdBKhBrv1xZzYeskgiBU6Sc3derU+65fvXp19YORSLBhwwaMGjXKuKxTp05o164dli1bZlwWERGBUaNG4YMPPsDcuXPxf//3f5DJZCgoKIBWq8XLL7+Mt956657HeOedd/Duu+9WWP7jjz/Cycmp2rET1UY6A7AvQ4JtKVIU60tPN7TyMGBEsAHejiIHR0T3dCRTgh/j/x3G4aOOOphjVIeioiKMHz8eubm5cHNzq/kO6YGqXGSNHj263GOtVovY2FjcunULffv2xe+//179YO4qsjQaDZycnLB+/fpyx505cyZiYmKwZ8+ecs9fu3YtYmNj73t14b1asoKCgpCVlWWWN51Wq0V0dDT69+9v7LNGpZgb08ydG0EQsOPCDfxv6yVcyykCAEQ0cMUbg8PQuZFHjfdvTXzfmMbcmGbLufnrbAZm/XLG+Pj0vL5mGaA0Ly8PXl5eLLKsqMo/tQ0bNlRYZjAYMH36dDRq1MgsQZXJysqCXq+Hr69vueW+vr7IyMio1j5VKhVUqorTgSgUCrP+Ipp7f/aEuTHNHLmJS8vDgr/icDA+GwDg5aLCqwObYUxUkE2Pt8P3jWnMjWm2mBtnh/JDp5S+hpoXWbaWB3tglrH7pVIpXnrpJfTu3Rtz5swxxy7LufuqCkEQ7nmlxeTJk81+bCJbcSNfjcXRF/HzsWQIAqCUS/F091BM79MELpxwlshmlHV8J9tntk/e+Ph46HQ6c+0OAODl5QWZTFah1SozM7NC6xZRXVWi1WPNgUR8tesKCtSlv4NDW/nh9UHhCPJgP0MiW+Og4LQ69qLKRdbs2bPLPRYEAenp6fjrr78wadIkswUGAEqlElFRUYiOji7XJys6OhojR44067GIbI0gCPg7NgMf/H0eyTnFAIBWge6YNywSHRraVr8rIvpXWUuWo0KKIS39bfo0f11X5SLr1KlT5R5LpVJ4e3vjk08+eeCVh/dSUFCAK1euGB8nJCQgJiYGHh4eCA4OxuzZszFhwgS0b98eXbp0wYoVK5CUlIRnn322yscishdnU3Lx3uY4HE3MAQD4uqkwZ2A4RrcNgJQfyEQ2rawly0kpxydjW4scDdVElYusXbt2mTWA48ePo0+fPsbHZS1lkyZNwtq1azFu3DhkZ2dj/vz5SE9PR4sWLbBlyxaEhISYNQ4iW3A9rwSLtl3EbydTIAiAg0KKaT0bY1qvRma5+oiIxOdwe6CsEg5CavOq/KlcNkxDvXr1yi3Py8vDqFGjsHPnzirtr3fv3njQKBLTp0/H9OnTqxoqkd0o1uixct9VfL0nHkWa0g/e0W0D8OrAMPjX44BXRPZEJS9tySrR6lGk0cFRIeO0OjaqykXW7t27odFoKiwvKSmpMN0NEdWMIAjYdDoNH/59AWm5JQCAdsH18Nbw5mgTVE/c4IjIIspasvQCEPnWNsTNH8iWahtV6Z/amTP/DowWFxdX7oo/vV6PrVu3IiAgwLzREdVhJ5Nu4r3NcTiVdAsAEFDPEa8NDsfwVn78q5bIjvHqQvtR6SKrTZs2kEgkkEgk6Nu3b4X1jo6O+OKLL8waHFFdlHarGB9uvYA/YtIAAE5KGab3boynezTihy9RHaCUcZwse1HpIishIQGCIKBRo0Y4evQovL29jeuUSiV8fHwgk/ELgKi61HpgyY4rWHUgESVaAyQSYEy7QLw6MAw+bg5ih0dEViKVSqCQSaDV13xSaBJXpYussqv5DAaDxYIhqosMBgG/n0rFwlMy5GqvAgA6hnrgrWGRaBHgLnJ0RCQGlVwKrZ5XF9q6ShVZmzZtwuDBg6FQKLBp06b7bjtixAizBEZUFxxLzMH8P+NwNjUXgASB9R3x5pAIDGrRgP2uiOowB4UMBWoWWbauUkXWqFGjkJGRAR8fH4waNcrkdhKJBHpW3kQPlJxThP/9fQF/nU0HADirZHjIV4OFU7rBxbHiBOZEVLco2C/LLlSqyLrzFCFPFxJVX36JFkt3x2PV/gRodAZIJcC4DsGY0ScUR/fu4MSwRASgdIJ3AOjcyANStmrbLA68QWQFeoOA9ceT8fH2i8gqKB1nrlsTT8wbFonwBm7QarUiR0hEtUnZH1wzHmrKq4ptWKWKrM8//7zSO5wxY0a1gyGyRwevZOG9v87jfHoeACDUyxlvDonAQxE+7HdFRPdUdrqQVxjatkoVWZ9++mmldiaRSFhkEd2WkFWIhVvOIzruOgDAzUGOmf2aYULnEOOpACKiezEWWTp20bFllSqyEhISLB0Hkd3ILdbiix2X8e2hRGj1AmRSCZ7sFIxZ/ZqhvrNS7PCIyAbIbrdyP/3dcU6rY8Nq9FMrm9iZpzyIAJ3egJ+OJuHTfy4jp7C031XvMG+8OSQCTX1dRY6OiGyJQs7vVXtQrXMWq1atQosWLeDg4AAHBwe0aNEC33zzjbljI7IZey7dwODP9mHeH+eQU6hBEx8XrJ3SAWundGSBRURVxiEc7EOVW7LmzZuHTz/9FC+++CK6dOkCADh06BBeeuklJCYmYsGCBWYPkqi2upKZj/f/Oo9dF28AAOo7KfBS/2YY3zEYcn5IElE1sciyD1UuspYtW4aVK1fi8ccfNy4bMWIEWrVqhRdffJFFFtUJNws1+GzHZXx/+Br0BgFyqQSTujbEjL5N4e6kEDs8IrJxChlPF9qDKhdZer0e7du3r7A8KioKOp3OLEER1VZavQHfH7qGz3ZcRm5x6dhW/SJ88caQcDTydhE5OiKyFwopiyx7UOUi68knn8SyZcuwePHicstXrFiBJ554wmyBEdUmgiBg54VMvP/XeVzNKgQAhDdwxbxhkejWxEvk6IjI3vB0oX2o1tWFq1atwvbt29G5c2cAwOHDh5GcnIyJEydi9uzZxu3uLsSIbNGFjDws2Hwe+69kAQC8XJR4eUAYxrYPgox/bRKRBShvj/Ie6uXMaXVsWJWLrNjYWLRr1w4AEB8fDwDw9vaGt7c3YmNjjdtxWAeyddkFaiyOvoSfjibBIABKmRRTu4fi+T6N4erAfldEZDmOt4usYa38OK2ODatykbVr1y5LxEFUa6h1eqw9kIgvd15Bvrq0n+GQlg3w+qAIBHs6iRwdEdUFZR3fNXqO+G7LOIQs0W2CIGDbuetYuOU8knKKAAAtAtwwb2gkOjXyFDk6IqpL/p1Wh3MX2rIqF1klJSX44osvsGvXLmRmZsJgKF9lnzx50mzBEVlLbGou3tschyMJOQAAH1cVXh0YhkfaBULKfldEJJJvDyXilYHNOK2OjaryT23q1KmIjo7GmDFj0LFjR/a9IpuWmVeCj7dfxPoTKRAEQCWX4j89G+HZXo3hrOKHGhGJo6wlS29gS5Ytq/K3yF9//YUtW7agW7duloiHyCpKtHqs2p+ApbuuoFCjBwCMaO2P1waHI6Ceo8jREVFdx8FI7UOVi6yAgAC4unIuNrJNgiBg85l0/O/vC0i9VQwAaBNUD/OGRSIqpL7I0RERleI4Wfahyj/FTz75BK+99hquXbtmiXiqLDk5Gb1790ZkZCRatWqF9evXix0S1VIxybcw5utDePGnU0i9VQw/dwd89lgb/P5cVxZYRFSrsMiyD1VuyWrfvj1KSkrQqFEjODk5QaEoP15QTk6O2YKrDLlcjiVLlqBNmzbIzMxEu3btMGTIEDg7O1s1Dqq90nOLsWjrRfx+KhVA6fgzz/VujGd6NIKjkuPPEFHto2SRZReqXGQ9/vjjSE1NxcKFC+Hr6yt6x3c/Pz/4+fkBAHx8fODh4YGcnBwWWYQijQ4r9l7F13viUaItvQr2kXaBeHVgGBq4O4gcHRGRaVLWWHahykXWwYMHcejQIbRu3dosAezduxeLFi3CiRMnkJ6ejg0bNmDUqFHltlm6dCkWLVqE9PR0NG/eHEuWLEGPHj0q7Ov48eMwGAwICgoyS2xkmwwGAX+cTsWHf19ERl4JAKBDw/qYNywSrQLriRscEVEllJ0udHOQc1odG1blWjk8PBzFxcVmC6CwsBCtW7fGl19+ec/169atw6xZs/Dmm2/i1KlT6NGjBwYPHoykpKRy22VnZ2PixIlYsWKF2WIj23PiWg5GLz2Al9adRkZeCQLrO+Kr8e3wy7QuLLCIyGaUdWWI9HfjtDo2rMotWf/73//w8ssv4/3330fLli0r9Mlyc3Or0v4GDx6MwYMHm1y/ePFiPPXUU3j66acBAEuWLMG2bduwbNkyfPDBBwAAtVqN0aNHY+7cuejatet9j6dWq6FWq42P8/LyAABarRZarbZKsd9L2T7MsS97Y8ncpN4qxqJtl/FXbAYAwFkpw3O9GmFyl2CoFDLodDqzH9Oc+L4xjbkxjbkxzdZzI9we6FunN5jtNdhqLmyZRBCEKo10Jr19ovjuvliCIEAikUCv11c/GImk3OlCjUYDJycnrF+/HqNHjzZuN3PmTMTExGDPnj0QBAHjx49HWFgY3nnnnQce45133sG7775bYfmPP/4IJyfOS2drSvTAP6lS7EqTQCdIIIGAzj4ChgQZ4KYUOzoiouo5myPBNxdlaOgi4KWW1f9evVNRURHGjx+P3NzcKjeIUPWYdYLoU6dO1SiYu2VlZUGv18PX17fccl9fX2RklLZYHDhwAOvWrUOrVq2wceNGAMD333+Pli1b3nOfc+fOxezZs42P8/LyEBQUhAEDBpjlTafVahEdHY3+/ftXaOWr68yZG71BwO+n0vDpP5dxo0ADAOgcWh9zB4ch0s/2Pjz4vjGNuTGNuTHN1nMjPZcBXDyD5CIpevd7yCzT6pSduSHrqfJPrVevXuUe5+bm4ocffsA333yD06dPY9asWeaKzchUqxkAdO/evcL8ifejUqmgUqkqLFcoFGb9RTT3/uxJTXNz+Go23tsch3NppR8YIZ5OeHNIBPpHin+1a03xfWMac2Mac2OareZGdTtmvUG4/RpqXmTZYh5sXbV/ajt37sTq1avx+++/IyQkBI888ghWrVplztjg5eUFmUxmbLUqk5mZWaF1i+zftexCfLDlAraeK30/uDrIMaNvU0zsGgKVnB1Dich+8IpC+1ClIislJQVr167F6tWrUVhYiLFjx0Kr1eK3335DZGSk2YNTKpWIiopCdHR0uT5Z0dHRGDlypNmPR7VTXokWX+28gjUHEqHRGyCVAOM7BeOlfs3g6VKxVZKIyNbJpCyy7EGli6whQ4Zg//79GDZsGL744gsMGjQIMpkMX3/9dY0CKCgowJUrV4yPExISEBMTAw8PDwQHB2P27NmYMGEC2rdvjy5dumDFihVISkrCs88+W6PjUu2n0xuw7ngyFm+/hOzC0n5XPZp64b9DIxHWgPNnEpH9knOCaLtQ6SJr+/btmDFjBp577jk0bdrUbAEcP34cffr0MT4u65Q+adIkrF27FuPGjUN2djbmz5+P9PR0tGjRAlu2bEFISIjZYqDaZ9/lG1iw+TwuXs8HADTydsa8oZHoHeZt8/2uiIgehKcL7UOli6x9+/Zh9erVaN++PcLDwzFhwgSMGzeuxgH07t0bDxpFYvr06Zg+fXqNj0W1X/yNAiz86zx2XMgEALg7KvBSv6Z4onMIJ0wlojpDztOFdqHS31pdunTBypUrkZ6ejmnTpuHnn39GQEAADAYDoqOjkZ+fb8k4yc7lFmkx/884DPx0L3ZcyIRcKsGUbg2x59XemNwtlAUWEdUpZacL5VIJW7VsWJW/uZycnDB16lTs378fZ8+excsvv4z//e9/8PHxwYgRIywRI9kxrd6AtQcS0OvjXVh9IAE6g4CHwn2w7aWeeHt4c9Rz4oiiRFT3OKtKTzTVc1JyWh0bVqPmgbCwMHz00UdISUnBTz/9ZK6YqI7YdTETg5bsxTt/xuFWkRZhvq74/qmOWDW5Axp7u4gdHhGRaGS3W6/0VRgHkmqfmo9uBkAmk2HUqFHG6XCI7ufy9QL8b/tl7L10AwDg4azE7P7N8FiHIMh5WpCIyDiEg95QpZnvqJYxS5FFVBk5hRqsvyrFoSOHSkcxlkkwtVsonu/bBG4OHImYiKiMTl9aXOWX6FCs0cNRyVOGtohFFlmcRmfAd4cS8dmOy8gvkQIQMLC5L+YOjkBDL2exwyMiqnWkt1uyBAAC2Jplq1hkkcUIgoDouOtYuOU8ErOLAAABTgI+fKwDujfjtEhERKbIeEWhXWCRRRYRl5aH9zbH4dDVbACAl4sKs/s1gWPGaXQK9RA5OiKi2o3T6tgHFllkVjfy1VgcfRE/H0uGIABKuRTP9AjFc72bQCUVsGXLabFDJCKq9e4cjNTAzu82i0UWmUWJVo81BxLx1a4rKFDrAADDWvnhtUHhCPJwAgBotVoxQyQishmSOy601rHIslkssqhGBEHA37EZWLjlPFJuFgMAWge6Y96wSLRvyNOCRETVcWdLFodxsF0ssqjazqbk4r3NcTiamAMAaODmgDmDwjCqTYDxyhgiIqo6ufTfpizWWLaLRRZV2fW8Eny09SJ+O5kCAHBQSDGtZ2NM69UITkq+pYiIasrV4d/PUiUHabZZ/EakSivW6LFy31Us2x2PYq0eADC6bQDmDAqDn7ujyNEREdmPO68u1HFqHZvFIoseSBAEbDqdhg//voC03BIAQLvgenhreHO0CaonbnBERHZIcsc4WTxbaLtYZNF9nUy6ifc2x+FU0i0AQEA9R7w+OBzDWvmV+xAgIiLzKdbo73mfbAuLLLqntFvF+HDrBfwRkwYAcFLK8HyfJniqeygcFJxDi4jIku6cSofjZNkuFllUTqFah+V74rF871WodQZIJMCjUYF4ZUAYfNwcxA6PiKjOYY1lu1hkEYDSv5R+P5WKj7ZeQGa+GgDQMdQDbw2LRIsAd5GjIyKquwwCqyxbxSKLcDQhB+9tjsPZ1FwAQLCHE94YEo6BzRuw3xURkchYYtkuFll1WHJOET74+zy2nM0AALio5HixbxNM7tYQKjn7XRER1Qbsk2W7WGTVQfklWny1Kx6r9ydAozdAKgEe6xiM2f2bwctFJXZ4RER0B54ttF0ssuoQvUHAL8eT8cn2i8gq0AAAujfxwn+HRSC8gZvI0RERURkJJJDg9qlCCassW8Uiq444eCUL8zfH4UJGPgCgkZcz3hwagb7hPux3RURUyzgqZXBRyZGv1kHJ7hs2i0WWnUvIKsTCLecRHXcdAODmIMfMfs0woXMIlHLOh0VEVFuV/f3LqwttF4ssO5VbrMUXOy7j20OJ0OoFyKQSTOgcgpkPNUV9Z6XY4RER0QNIb89fyBrLdrHIsjM6vQE/HU3C4uhLuFmkBQD0DvPGf4dGoImPq8jRERFRZZRo9cgv0d2+rxM5GqouuzhftHnzZoSFhaFp06b45ptvxA5HNHsu3cDgz/Zh3h/ncLNIi6Y+Llg7pQPWTunIAouIyIYYBAH620M36A0iB0PVZvMtWTqdDrNnz8auXbvg5uaGdu3a4eGHH4aHh4fYoVnNlcx8LPjrPHZfvAEAqO+kwOz+zfB4x2DIZXZRRxMR1VkChyO1WTZfZB09ehTNmzdHQEAAAGDIkCHYtm0bHn/8cZEjs7ybhRp8tuMyvj98DXqDAIVMgkldGuLFvk3h7qQQOzwiIjIDA1uybJbozRx79+7F8OHD4e/vD4lEgo0bN1bYZunSpQgNDYWDgwOioqKwb98+47q0tDRjgQUAgYGBSE1NtUbootHoDFi1PwG9P96NtQcToTcI6B/pi+0v9cJ/h0WywCIisiO8utB2id6SVVhYiNatW2PKlCl45JFHKqxft24dZs2ahaVLl6Jbt25Yvnw5Bg8ejLi4OAQHB0O4x5vvfuM+qdVqqNVq4+O8vDwAgFarhVarrfHrKduHOfZ1N0EQsOtSFv7390UkZBcBAMJ9XfDGkDB0aeRpseOaiyVzY+uYG9OYG9OYG9NsPTfaOzq763Q6s34/kfVIhHtVKSKRSCTYsGEDRo0aZVzWqVMntGvXDsuWLTMui4iIwKhRo/DBBx/g4MGDWLRoETZs2AAAmDlzJjp16oTx48ff8xjvvPMO3n333QrLf/zxRzg5OZn3BZlRWiGw4ZoUl3JLGx9dFAKGBhnQ2UeAlGOJEhHZFbUemHO0tB1kRnMdGpthUo6ioiKMHz8eubm5cHPjLB/WUKuLLI1GAycnJ6xfvx6jR482bjdz5kzExMRgz5490Ol0iIiIwO7du40d3w8fPgxPT897HuNeLVlBQUHIysoyy5tOq9UiOjoa/fv3h0JR89N22QVqLNkZj1+Op8AgAAqZBFO6huDZno3g6iB6Q2SVmDs39oS5MY25MY25Mc3Wc1Ok0aHtgp0wCMD/TW2PTqE1v5grLy8PXl5eLLKsqFZ/S2dlZUGv18PX17fccl9fX2RkZAAA5HI5PvnkE/Tp0wcGgwFz5swxWWABgEqlgkpVcRJkhUJh1l/Emu5PrdNj7YFEfLnzCvLVpc3GQ1o2wOuDIhDsWXtb3CrD3Lm2J8yNacyNacyNabaaG3eFAv71HJFysxjODkqzvAZbzIOtq9VFVpm7+1gJglBu2YgRIzBixAhrh2URgiBg27kMLNxyAUk5pf2uWgS44a1hzdHRDH/JEBGRbZDe/p4z1JrzTVRVtbrI8vLygkwmM7ZalcnMzKzQumUPYlNz8d7mOBxJyAEA+LiqMGdQOB5uG2CcXoGIiOqGso/9WtSrh6qoVhdZSqUSUVFRiI6OLtcnKzo6GiNHjhQxMvPKzCvBx9svYv2JFAgCoJJLMa1nI0zr1RjOqlr9IyIiIgso0eqRkVtivE+2SfRv8IKCAly5csX4OCEhATExMfDw8EBwcDBmz56NCRMmoH379ujSpQtWrFiBpKQkPPvssyJGbR4lWj1W7U/AV7uuoEhT+ks0so0/5gwKR0A9R5GjIyIisRgEASW60lFI9TxfaLNEL7KOHz+OPn36GB/Pnj0bADBp0iSsXbsW48aNQ3Z2NubPn4/09HS0aNECW7ZsQUhIiFgh15ggCNh8Jh3/+/sCUm8VAwDaBNXDvGGRiAqpL3J0RERUm7DGsl2iF1m9e/d+4Pnm6dOnY/r06VaKyLJikm/hvc1xOHHtJgDA390Brw0Ox/BW/ux3RXZNEATodDro9RVPfchkMsjl8vsOJExUV7FPlu0SvciqK9Jzi/HR1ovYcKp0yh9HhQzP9W6MZ3o0gqNSJnJ0RJal0WiQnp6OoqIik9s4OTnBz88PSqXSipER1X4ssWwXiywLK9LosHzPVSzfG48Sben59UfaBeLVgWFo4O4gcnRElmcwGJCQkACZTAZ/f38olcpyLVaCIECj0eDGjRtISEhA06ZNIZWKPq0qEVGNsciyEINBwO8nU/DR1ovIyCu9QqRDw/qYNywSrQLriRsckRVpNBoYDAYEBQWZnLrK0dERCoUC165dg0ajgYMD/wAhKsOzhbaLRZYFJOQDj644gjOppZNPB9Z3xBtDIjC4RQP2OaE660GtU2y9IipPImGBZetYZJlJVoEaey/dwJYzafjnghxAHlxUckzv0xhTu4XCQcF+V0REVDlOSjnCfF1xISMfKgX/ALFVLLKqSW8QEJN8C3suZmL3pRs4m5pr/ItDAgFj2wfilYER8HatOE8iERER2T8WWVWQVaDGnos3sPvSDey7fAO3irTl1kf6uaFnU0+4517G0yObczJOIiKqtrLuJTxlaLtYZN1HaWvVTey+eAO7L5a2Vt3JzUGOHk290SvMG72becPHzQFarRZbtlwWKWKi2stgMNRoPVFdUqLVIzmndMgTjY6/G7aKRdZdbuSrsefSDey+mIl9l7OQW1y+taq5vxt6h3mjd5gP2gbVg1zGc+VE96NUKiGVSpGWlgZvb+/7DuEglUo5ThYRSqfVKVDrAACNvJ1Fjoaqi0XWbZ//cxmHU4sQe/uKwDLujgr0aOqF3mE+6NnMCz6uvLScqCqkUilCQ0ORnp6OtLQ0k9s5OTkhODiYVxkS3YVjKtouFlm3rdh3FVJV6Rg+LQPcb7dWeaN1IFuriGpKqVQiODiY0+oQUZ3CIuu2wS0aYECbUPRs5s0rAoksQCKRQKFQ8IIQIqozWGTdtujR1nBzcxM7DCIiIrITPA9GREREZAF1viVLuD0ASV5e3gO2rBytVouioiLk5eXxtMhdmBvTmBvTmBvTmBvTbD03RRodDOrSIRzy8vKgU9b867rse07gwFtWIxHqeLZTUlIQFBQkdhhERERWkZycjMDAQLHDqBPqfJFlMBiQlpYGV1dX45VNHTp0wLFjx8ptd/cyU4/z8vIQFBSE5ORks/bxuldM5niOqW2qspy5Mb3cHnNzv/Xmys2d95kb5uZ+y2p7bmryWXO/ddXJzdGjR5Gfnw9/f38OlWIldf50oVQqrVDRy2SyCr+Udy970GM3Nzez/mLfKyZzPMfUNlVZztyYXm6PubnfenPl5l7bMzfMzb2W1fbc1OSz5n7rqpMbd3d3uLu7Vzp2qrk6X2TdqyVr6tSpFfpo3b3M1OOyZebq42XqeOZ6jqltqrKcuTG93B5zc7/15srNnfeZG+bmfstqe25q8llzv3XVyU1ubi5bsqyszp8uZJ8sIiKqS9gny3rqfEuWq6srAJjtvL1Wq8X27dsxYMAAm7yixZKYG9OYG9OYG9OYG9NsPTdFGh06vr8DAHD0zYfgZKarC4OCgozfe2R5db7IKjtFaK7z9lqtFk5OTnBzc7PJX2xLYm5MY25MY25MY25Ms/XcyDU641Rvbm5uZimyynD6Kuup80UWERFRbaOUSfHV+HbG+2SbWGQRERHVMnKZFENb+YkdBtUQy2MiIiIiC7DZlqzw8HBcvHjxnus+++wzzJgxo0r7i029Bedcg/HxnZdc3nkBZvnlFZdqtTpczQNOXLsJmVx+j23vvb9y61HuwQOO/e/297pOVCIBJJCg7BS85PY/ZcskKD0//+/90g3KrTO1rzv3U2GfxqMZl+n1OtwoBq7lFEEpV5Tfz+3jlO3nzudBAkglEsgkEkilEsikZfcBmaT0MfsYEJE90ekN2HbuOgBgYHNfyHnK0CbZbJG1YcMGREZGYseOHYiMjIRWq0WjRo3w008/Yfjw4Safp1aroVarjY/LxhB5bMURYyfDmpPjs3NVGy257pBjQcx+s+9VIsG/RZgEt/8vLcCkEgnk0orryv1/u2CTy6RQyCRQlPtfCuUd9xUyCZRyabnHd2+jkkvhqJTBUSGDg0Ja+r9SBieFDA4KGRwV0nIfmlqtttz/9C/mxjTmxjRbz02RRofnfzwJADg9r69ZOr7bai5smc0WWRkZGZDL5ejWrRtUKhViYmKg0+nQo0cPqFQqk8/74IMP8O6771ZYXk8pQKYUYKpBxFQ7ieQeG0jutf7u591jxQOPUY3jCLdvEO5qCSv7XzC93Z3blFt2j33d87l3Pb7fNhViK7fdg1upBAHQCQJgsJ1h32QSAUopoJQCChmglMqwJHYnHGUCHOWAgwxwlAEOcqH0fxluLy9d73z7Jq8jf+BGR0eLHUKtxdyYZqu5UeuBsq/obdu2QyWr+T6LiopqvhOqEpstss6ePYtmzZoZC6qYmBh4e3vD19f3vs+bO3cuZs+ebXxcNm7IrjkPmW0Ih+joaPTv398mLxu2pJrkRhAE6A0C9AJgMAjQC8Jd/wN6gwDD7e1K//93W72h7Pl3Pg/lHusNAnR6ATqDAK3eAK3eAI3+3/taXdn9O9ff3kZnMK7T6A3Q6Awo0epRpNGjRKtHsfb2Y63eWNjqBQmK9UCxHkC5PzCrdurT1UEOT2clPJyV8HBSlP5/++bprEQDdxUauDnA180BKhusyPg7ZRpzY5qt56ZIo8OcozsBAAMHDjDbOFlkXTZbZJ05cwYtW7Y0Po6JiUGrVq0e+DyVSnXPli6FQmHWX0Rz78+e1OXcCEJpEVaiMaBIq0OxRo9irR4FxRrs3n8Ika3boVgrIK9Ei/wSnfH/fOP/pffzSnS4VaSBQYBxeWL2g/9KLS26HODn7nD7f0eEeDqhoaczQjyd4OpQe38udfl98yDMjWm2mhuF8O8fW6WvoeZf17aYB1tn00XWqFGjjI9jYmLQrl078QIiqgSJRAKVXAaVXAZ3/PuBp9VqkVFPwKDmvpX+IDQYBOQWa5FdqEF2gRo5hRpkF2qQc/uWXajBjfwSXM9TI+1WMdQ6Q+m2hRqcS7v3X7Sezso7ii5nNPN1QVgDV4R4OkMm5cUFRERVYZNFlsFgwLlz5/DWW28Zl129ehWjR48WMSoi65JKJajvrER9ZyWa+Ljcd1tBKC3I0nNLkJFbgrTcYmTkliD1ZjGu5RThWnYhsgo0xiLsZNKtcs9XyaVo6uuCMF83hDdwRXN/N7QKqgcXlU1+hBARWYVNfkLGx8ejqKio3OnB1q1b45133kHbtm3Rs2dPEaMjqn0kEgnqOSlRz0mJCL979z3ML9HiWnZR6S2nEFdvFOLS9Xxcup6PEq0Bsal5iE3Nu2OfQFMfF7QJqofWQfXQJqgewhu4scWLiOg2myyymjZtWm6sKQD4888/RYqGyD64OijQIsAdLQLcyy3XGwQk5RThYkYeLmTk40J6Ps6m5iL1VjEuXS/ApesF+OV4CgDAzUGOTo080aWRJ7o28UQzH1dIWXQRVZlCJsWiMa2M98k22WSRRUTWI5NKEOrljFAvZwxq8e80H5n5JTidnIuY5Js4nZyL08m3kFeiQ3TcdUTHlQ6i6OGsRPcmXngowge9w3zg7siOt0SVoZBJ8Wj7ILHDoBpikUVE1eLj6oD+kQ7oH1k6bIpOb8C5tDwcupqNg/HZOJaQg5xCDTadTsOm02mQSyXoGOqBhyJ8MbhFA/jXcxT5FRARWRaLLCIyC7lMita3+2c926sxNDoDYpJvYeeFTOw4fx2XMwtwML60AHtvcxw6hnpgVJsADGnZAPWclGKHT1Sr6PQG7L18AwDQs6k3p9WxUSyyiMgilHIpOoZ6oGOoB14fHI7ErEL8c/46tsddx9GEHOPt7U2x6BPmg8c7BaNXU2/24SICoNEbMHXtcQBA3PyBLLJsFIssIrKKhl7OeLpHIzzdoxHSbhXjz9Np2BiThvPpedgeV1p8BXk4YnzHEIxtHwhPF9PTYxER2QIWWURkdf71HDGtV2NM69UYl67n4+ejyfj1RDKSc4rx4dYL+DT6Eh6JCsC0no0R4M5TiURkm1hkEZGomvm64q3hkXh1YBj+PJ2G/ztyDWdScvHT0WSsO5aMwc0bIJJnSojIBrHIIqJawVEpw9gOQRjbIQjHEnOwdNcV7Lp4A3/FZuAvyHHOcBqvDopAqJez2KESEVUK/z4kolqnQ0MPrJnSEVtm9MDQlg0ggYAtsdfRf/Ee/HfjWWTml4gdIhHRA7HIIqJaK9LfDUvGtsKrrfTo3cwLOoOA/zuchL4f78Hq/QnQ6Q1ih0hEZBJPFxJRrRfgDKx8tB1OJOfhgy3ncTolF/M3x2H9iRQsGNUCUSH1xQ6RyKwUMinmj2xuvE+2iT85IrIZnRt5YsP0blg4uiXcHRU4n56HR5YdxDubzqFEqxc7PCKzUcikmNilISZ2acgiy4bxJ0dENkUqlWB8p2DsfLkXxkQFAgDWHkzE0M/34UzKLXGDIyK6A4ssIrJJni4qfPxoa6yZ0gE+rirE3yjE6KUH8dWuKzAYBLHDI6oRvUHAofhsHIrPhp7vZ5vFIouIbFqfMB9sm9UTQ1v5QW8QsGjbRTzz3XHkFmnFDo2o2tQ6PR5feRiPrzwMtY6nwm0Viywisnn1nZX4anw7fPRIK6jkUuy4kIlhX+5DbGqu2KERUR3GIouI7MbYDkH47bmuCPJwRHJOMR79+hD+ibsudlhEVEexyCIiu9IiwB2bX+iBns28UazV4z/fH8e3BxPFDouI6iAWWURkd9ydFFg1qT0e7xgEgwC8vekcFm45D0FgB2Iish4WWURklxQyKRaObok5g8IAACv2XsW8P2J55SERWQ2LLCKyWxKJBNN7N8FHj7SCRAL83+EkzP39LC+JJyKr4LQ6RGT3xnYIgkIuwcu/nMa648nQGQQsGtMKUqlE7NCI7kkulWLu4HDjfbJNLLKIqE4Y3TYQSpkMM34+hd9OpsDdUYF5wyIgkbDQotpHKZdiWq/GYodBNcTymIjqjKGt/LBoTCsAwOoDCfhq1xWRIyIie8Yii4jqlIfbBeKtYZEAgI+3X8JPR5NEjoioIr1BwOnkWzidfIt9CG0YiywiqnOmdg/Fi32bAADmbYzFwfgskSMiKk+t02PkVwcw8qsDnFbHhrHIIqI6aXb/ZhjR2h86g4DpP5xEYlah2CERkZ1hkUVEdZJEIsFHY1qhdaA7bhVp8dS3x5BXwkmlich8WGQRUZ3loJBh5cT2aODmgPgbhZj721mOCk9EZsMii4jqNB83Byx7sh0UMgn+OpuO/zt8TeyQiMhOsMgiojqvbXB9vD44AgDw3ubziE3NFTkiIrIHLLKIiABM7dYQ/SN9odEb8PyPJ1Go1okdEhHZOBZZREQo7Qj/8ZjWCKjniGvZRfjf3xfEDonqMLlUipkPNcXMh5pyWh0bxp8cEdFt7k4KfHR7RPjvD1/Dvss3RI6I6iqlXIqX+jfDS/2bQSnnV7Wt4k+OiOgO3Zp4YWKXEADAnF/PcFgHIqo2FllERHd5fXA4QjydkJ5bgoV/nRc7HKqDDAYBl67n49L1fBg4rY7NYpFFRHQXJ6UcHz/aGgDw87FkHE/METkiqmtKdHoM+HQvBny6FyWcVsdmscgiIrqHDg09MK59EADgvxtjodUbRI6IiGyN3NIH8PDwqNL2EokEJ0+eREhIiIUiIiKqnNcHh2N7XAYuZORjzYEE/KdnY7FDIiIbYvEi69atW1iyZAnc3d0fuK0gCJg+fTr0ejaNEpH46jsrMXdIBOb8egZL/rmMYa384V/PUeywiMhGWLzIAoDHHnsMPj4+ldr2xRdftHA0RESVN6ZdINYfT8axxJtYtO0iPh3XRuyQiMhGWLxPlsFgqHSBBQD5+flo1KiRBSMiIqo8qVSCt4Y1BwBsOJWKMym3xA2IiGwGO74TET1Ay0B3PNw2AACwYPN5CAIvqSeiB7P46cJNmzZVetsRI0ZYMBIioup7ZWAY/jqbjqOJOdh27joGtWggdkhkx+RSKf7Ts5HxPtkmixdZo0aNKvdYIpGU+ytQIpEY77PDOxHVVv71HPGfno3wxc4r+N/f5/FQhA8UMn75kWUo5VK8MSRC7DCohqzSJ6vstn37drRp0wZ///03bt26hdzcXGzZsgXt2rXD1q1bLR0KEVGNTOvVGF4uSiRmF+G3Eylih0NEtZxV/wybNWsWPvvsMwwcOBBubm5wdXXFwIEDsXjxYsyYMcOaoRARVZmLSo5ne5WOlfXFzitQcyRushCDQUByThGSc4o4rY4Ns2qRFR8ff8/xstzd3ZGYmGjNUIiIquXJziHwdVMh9VYx1h1LFjscslMlOj16fLQLPT7axWl1bJhVi6wOHTpg1qxZSE9PNy7LyMjAyy+/jI4dO1ozFCKianFQyPBCnyYAgC93XkGJll+ARHRvVi2yVq9ejczMTISEhKBJkyZo0qQJgoODkZ6ejlWrVlkzFCKiahvbIQgB9RyRma/G/x2+JnY4RFRLWWXE9zJNmjTBmTNnEB0djQsXLkAQBERGRqJfv37lrjIkIqrNVHIZZjzUBK/9dhZf77mKJzuHwEEhEzssIqplrFpkAaVDNgwYMAADBgyw9qGJiMzm4XaB+HzHFaTeKsb6EymY0JmT2hNReVYvsgoLC7Fnzx4kJSVBo9GUW8crDInIVihkpYNFvr3pHJbvicfjHYIg57hZRHQHqxZZp06dwpAhQ1BUVITCwkJ4eHggKysLTk5O8PHxYZFFRDZlbPsgfL7jMlJuFuPPM2kY3TZQ7JCIqBax6p9dL730EoYPH46cnBw4Ojri8OHDuHbtGqKiovDxxx9bMxQiohpzVMowtXsoAGDZ7niOZ0RmI5NKMKFzCCZ0DoFMyj7LtsqqRVZMTAxefvllyGQyyGQyqNVqBAUF4aOPPsIbb7xhzVCIiMziyc4hcFHJcel6AXZcyBQ7HLITKrkM741qgfdGtYBKzosqbJVViyyFQmG8itDX1xdJSUkASgcjLbtPRGRL3B0VePJ2p/dlu6+IHA0R1SZWLbLatm2L48ePAwD69OmDt956Cz/88ANmzZqFli1bWjMUIiKzmdq9IZQyKU4m3cLJpJtih0N2QBAEZBeokV2ghiDwNLStsmqRtXDhQvj5+QEA3nvvPXh6euK5555DZmYmVqxYYc1QiIjMxsfVASPa+AMAVu1PEDkasgfFWj2iFvyDqAX/oJizCtgsq11dKAgCvL290bx5cwCAt7c3tmzZYq3DExFZ1NRuofj1RAq2xmYg9VYxAuo5ih0SEYnMai1ZgiCgadOmSElJsdYhiYisJtLfDV0be0JvEPDtwUSxwyGiWsBqRZZUKkXTpk2RnZ1trUMSEVnVU7eHc/jpaBIK1TqRoyEisVm1T9ZHH32EV199FbGxsdY8LBGRVfQJ80EjL2fkl+iw/niy2OEQkcisWmQ9+eSTOHr0KFq3bg1HR0d4eHiUuxER2TKpVIIp3RoCANYcTISeg5MS1WlWnVZnyZIl1jwcEZHVPRIViI+3X8K17CLsOH8dA5o3EDskIhKJVYusSZMmWfNwRERW56SU4/GOwfh6TzzWHEhkkUXVIpNK8Ei7QON9sk0WP12Yl5dXpe3z8/MtFAkRkXVM7FI639yhq9k4n161z0AioHRanU/GtsYnY1tzWh0bZvEiq379+sjMrPx8XgEBAbh69aoFIyIisiz/eo4YdLsFa+2BRHGDISLRWPx0oSAI+Oabb+Di4lKp7bVarYUjIiKyvCndGuKvs+nYGJOK1waHw8NZKXZIZEMEQTCO9O6okBnn/SXbYvEiKzg4GCtXrqz09g0aNIBCobBgRERElhcVUh8tA9xxNjUXPx1NwvN9mogdEtmQYq0ekW9tAwDEzR8IJ6VVu1CTmVj8p5aYmGjR/Z8/fx7dunVD/fr14erqiitXrqB9+/bYvXu3RY9LRHQ/EknpcA6zfzmN7w9dw396NoJCZtVRc4hIZDb/Gx8REYEOHTpg48aNiImJQYMGDbB161axwyIiwtBWfvByUSEjrwR/x2aIHQ4RWZnNF1kAcOnSJYSFhSE/Px8ymQwODg5ih0REBJVchic7BwMA1hxIEDkaIrI2my+y8vLy4OjoCKVSidjYWERGRoodEhGR0ROdQqCUSXEq6RZikm+JHQ4RWZHN96Q7d+6csbA6f/48mjVrdt/t1Wo11Gq18XHZOF5ardYsVzaW7YNXSVbE3JjG3Jhm67mp5yDF0Ja+2BCTjlX74rH40VZm27et58aSbD03Wq3ujvtaaCU1n6LJVnNhy2y+yIqLi0Pz5s0BAM7Ozvj777/x2muvmZwL8YMPPsC7775bYfn27dvh5ORktriio6PNti97w9yYxtyYZsu5aawHADn+OpuODooUuJt5NAdbzo2l2Wpu1LffMwCwbdt2qMwwHmlRUVHNd0JVIhEEwaozmO7btw/Lly9HfHw8fv31VwQEBOD7779HaGgounfvbvHj36slKygoCFlZWXBzc6vx/rVaLaKjo9G/f38ORXEX5sY05sY0e8nNYyuP4kTSLUzv1Qgv9TPPcA72khtLsPXcqLV6vPJbLADg40daQKWoeZWVl5cHLy8v5ObmmuX7jh7Mqi1Zv/32GyZMmIAnnngCp06dMhY7+fn5WLhwIbZs2WLxGFQqFVQqVYXlCoXCrL+I5t6fPWFuTGNuTLP13Ezt3ggnfjyJdcdTMKNfMziY4UuzjK3nxpJsNTcKhQJfT2hv9n2SdVm14/uCBQvw9ddfY+XKleV+2F27dsXJkyertU+JRHLfGxFRbTCwuS/83R2QXajBptNpYodDRFZg1SLr4sWL6NmzZ4Xlbm5uuHXrVrX2mZSUhF69eiEiIgItW7bEL7/8AkEQjDciotpALpNiQpeGAIA1BxL5+URUB1i1yPLz88OVK1cqLN+/fz8aNWpUrX3K5XIsWbIEcXFx+Oeff/DSSy+hsLCwpqESEZnd4x2D4KCQ4nx6Ho4k5IgdDtViRRodGr7+Fxq+/heKNLoHP4FqJasWWdOmTcPMmTNx5MgRSCQSpKWl4YcffsArr7yC6dOnV2uffn5+aNOmDQDAx8cHHh4eyMnhhxcR1T71nJQY3TYQAAcnJaoLrNrxfc6cOcjNzUWfPn1QUlKCnj17QqVS4ZVXXsELL7xQ4/0fP34cBoMBQUFBZoiWiMj8pnRriJ+OJiE67jqSc4oQ5GG+oWOIqHax+ojv77//PrKysnD06FEcPnwYN27cwHvvvVfj/WZnZ2PixIlYsWKFGaIkIrKMZr6u6N7ECwYB+O5QotjhEJEFiTKtjpOTE9q3b4/w8HD8888/OH/+fI32p1arMXr0aMydOxddu3Y1U5RERJYxpVtDAMDPx5JRqGZ/GyJ7ZdUia+zYsfjyyy8BAMXFxejQoQPGjh2LVq1a4bfffqvWPgVBwOTJk9G3b19MmDDBnOESEVlEnzAfNPR0Qn6JDr+fShU7HCKyEKsWWXv37kWPHj0AABs2bIDBYMCtW7fw+eefY8GCBdXa54EDB7Bu3Tps3LgRbdq0QZs2bXD27Flzhk1EZFZSqQSTujYEAKw9kACDgcM5ENkjq3Z8z83NNc4puHXrVjzyyCNwcnLC0KFD8eqrr1Zrn927d4fBYDBnmEREFjcmKhCfbL+E+BuF2HclC72aeYsdEtUiUokEfcK8jffJNlm1JSsoKAiHDh1CYWEhtm7digEDBgAAbt68CQcHB2uGQkQkKlcHBR5tz+Ec6N4cFDKsmdIRa6Z0NOsUTGRdVi2yZs2ahSeeeAKBgYHw9/dH7969AZSeRmzZsqU1QyEiEt3krg0hkQC7L95A/I0CscMhIjOzapE1ffp0HD58GKtXr8b+/fshlZYevlGjRtXuk0VEZKtCPJ3xULgPAODbg4niBkNEZmf1IRyioqIwevRouLi4GJcNHToU3bp1s3YoRESim9ItFADw64kU5BZrRY6GaosijQ4R87YiYt5WTqtjw6za8R0AUlJSsGnTJiQlJUGj0ZRbt3jxYmuHQ0Qkqq6NPRHm64qL1/Ox/ngynu5RvXlcyf4Ua/Vih0A1ZNUia8eOHRgxYgRCQ0Nx8eJFtGjRAomJpbPRt2vXzpqhEBHVChKJBJO7NcTc389i7cFETOkWCpmUV5MR2QOrni6cO3cuXn75ZcTGxsLBwQG//fYbkpOT0atXLzz66KPWDIWIqNYY1SYA9ZwUSLlZjOi462KHQ0RmYtUi6/z585g0aRIAQC6Xo7i4GC4uLpg/fz4+/PBDa4ZCRFRrOCpleKJTMADgy12XIQgcnJTIHli1yHJ2doZarQYA+Pv7Iz4+3rguKyvLmqEQEdUqU7uFwlEhQ2xqHnZeyBQ7HCIyA6sWWZ07d8aBAwcAlF5R+PLLL+P999/H1KlT0blzZ2uGQkRUq3i6qDCxSwgA4LMdbM0isgdW7fi+ePFiFBSUDrj3zjvvoKCgAOvWrUOTJk3w6aefWjMUIqJa55mejfDdoWs4k5KL3RdvoM/tMbSo7pFKJOgU6mG8T7bJqkVWo0b/Xprs5OSEpUuXWvPwRES1mpeLCk92DsbKfQlYsuMyeod5Q8Iv2DrJQSHDumldxA6Dasjqg5ECgEajQUpKCpKSksrdiIjquv/0bAwHhRSnk29h96UbYodDRDVg1SLr0qVL6NGjBxwdHRESEoLQ0FCEhoaiYcOGCA0NtWYoRES1krerCk92Ku2btWjrRRgM7JtFZKuserpwypQpkMvl2Lx5M/z8/NgMTkR0D9P7NMG6Y8mIS8/DH6dTMbptoNghkZUVaXTo/uEuAMD+1/rASWn1CVrIDKz6U4uJicGJEycQHh5uzcMSEdkUD2clnu3dGIu2XcTH2y5hcAs/OChkYodFVpZTqHnwRlSrWfV0YWRkJMfDIiKqhKndQtHAzQGpt4rxf4eviR0OWdmdp4l/O5EiYiRUExYvsvLy8oy3Dz/8EHPmzMHu3buRnZ1dbl1eXp6lQyEishmOShle6t8UAPDFzivILdKKHBFZ0/o7Cquujb1EjIRqwuKnC+vVq1eu75UgCHjooYfKbSMIAiQSCfR6zjhORFTmkXaB+GZfAi5nFuDTfy7hnRHNxQ6JrCAjtwSfbLtofOxXz0HEaKgmLF5k7dq1y9KHICKyS3KZFG8Pb44nVx3Bd4cSMbZ9ECL93cQOiyxIEAT8d+NZFGjY6GAPLF5k9erVy9KHICKyW92bemFIywbYcjYDb2+KxS/TuvDKbDu2+Uw6/jmfCbkU0BnEjoZqyiod34uKivD8888jICAAPj4+GD9+PDvAExFV0n+HRsJRIcOxxJvYcCpV7HDIQjLzS/D2pnMAgGk9G6NVoDtaBbpzWh0bZpUi6+2338batWsxdOhQPPbYY4iOjsZzzz1njUMTEdk8/3qOePGhJgCAhVvO41YRL+23N4Ig4NX1Z5BTqEF4A1fM6NcUm17ojk0vdOfwHTbMKuNk/f7771i1ahUee+wxAMCTTz6Jbt26Qa/XQybjm4eI6EGe7t4Iv59MxZXMAsz/Mw6Lx7UROyQyo28PJmLPpRtQyaX4/PG2UMn53WgPrNKSlZycjB49ehgfd+zYEXK5HGlpadY4PBGRzVPKpfhoTCtIJcDvp1Kx4/x1sUMiM7l0PR8L/74AAJg7OBzNfF1FjojMxSpFll6vh1KpLLdMLpdDp9NZ4/BERHahXXB9PNW9dJ7XNzacRV4xx86ydYVqHV748SQ0OgN6NfPGpK4NAQDFGj26/W8nuv1vJ4p5paHNssrpQkEQMHnyZKhUKuOykpISPPvss3B2djYu+/33360RDhGRzXp5QBh2nM/E1axCzP/rAvo6iR0RVZcgCHjttzO4dL0APq4qLHq0lfHKUQECUm8VG++TbbJKkTVp0qQKy5588klrHJqIyK44KGRY9GgrPPr1IfxxOh0uTSQYInZQVC1rDiRi85l0yKUSfPVEO/i4ctBRe2OVImvNmjXWOAwRUZ0QFeKBWf2aYXH0Jay/KsWkrEI086sndlhUBUcTcrBwy3kAwBtDItChoYfIEZElWHWCaCIiMo/n+zRBp9D60BgkmLnuDEq07LdjKxKzCjHt++PQGQQMa+WHKd0aih0SWQiLLCIiGySTSvDJmJZwlgs4n5GPeRtjIQjsu1Pb3SzUYMraY7hZpEWrQHd8NKYVR/C3YyyyiIhslK+bAyY2NUAqAdafSMGaA4lih0T3odbpMe3/TiAhqxAB9RzxzaT2cFJapdcOiYRFFhGRDQuvJ+D1QWEAgPe3nMf+y5yyrDbS6Q2Y9XMMjibkwFUlx+rJHe7b0V0CCZr6uKCpjwskYEuXrWKRRURk4yZ3CcYj7QKhNwiY/sMJXMzIFzskuoPBIGDOr2fwd2wGlDIplj0ZhbAG9x9w1FEpQ/TsXoie3QuOSo7+bqtYZBER2TiJRIL3R7dA+5D6yCvRYeLqI0i5WSR2WITSsbDm/RGL30+lQiaV4MvxbdG9qZfYYZGVsMgiIrIDDgoZvpnUHs18XXA9T42Jq44iu0Atdlh1msFQWmD9cCQJEgmweGxrDGjeQOywyIpYZBER2Yl6Tkp8O7Uj/N0dcDWrEBNWHcXNQo3YYdVJOr0Br6w/jf87XFpgffhwK4xsE1Dp5xdr9Oi/eA/6L97DaXVsGIssIiI74ufuiO+e6gQvFyXi0vPw+MrDbNGyMrVOjxd+PGU8RbhkXBuM7RBUpX0IEHA5swCXMws4rY4NY5FFRGRnmvi44Of/dIa3qwoXMvLx2IrDyMwvETusOiGnUIMnvzmCredKO7l//WRUlVqwyL6wyCIiskNNfFyx7j+d0cDNAZczC/Dw0oO4ksmrDi3pSmYBRi89gGOJN+GqkmPNlA7oH+krdlgkIhZZRER2qpG3C36Z1gUNPZ2QcrMYDy89iMNXs8UOyy7tuXQDDy89gGvZRQjycMTv07uiWxNeRVjXscgiIrJjwZ5O+H16N0TdHt5hwqoj+OloEqfgMRO9QcAn2y9i8pqjyCvRoV1wPWyY3g1Nfe8/DhbVDSyyiIjsnIezEj883QlDW/lBqxcw9/ezeGX9GV61VkOZ+SV48psj+GLnFQgC8GTnYPz4TGd4uajEDo1qCU6aRERUBzgoZPjisbZo4e+ORdsu4LeTKTiXlosvHm/LVpcqEgQBf55Jx1t/xOJWkRZOShk+eLilWTu4SyBBQD1H432yTSyyiIjqCKlUgud6N0aboHp48aeTuJCRj6Ff7MerA8IwtXsoZFJ+mT9IdoEa8/6IxZazGQCASD83fP54WzTxcTHrcRyVMhx4va9Z90nWx9OFRER1TJfGntgyowd6h3lDozPg/S3n8diKQ7z68D70BgE/HLmGhxbvwZazGZBLJZj5UFP88UI3sxdYZD/YkkVEVAf5uDlgzeQOWHcsGe9tjsOxxJsYtGQfnuoeihcfagoXFb8eypy4dhNvb4pFbGoeACC8gSs+frQ1WgS4ixwZ1Xb8LSIiqqMkEgke6xiMbk288O6f5/DP+Uws33sVG06l4uUBzfBIu0DIZXX3hMfFjHx8sv0itsddBwC4Osjxcv9meLJziMXzUqLVY+zyQwCAX6Z1gYNCZtHjkWWwyCIiquOCPJzwzaQO2HnhOt79Mw7Xsovw2m9n8fWeq5jVrymGtfKvU/21Ll/Pxxc7r+DPM2kQBEAqAcZEBWLOoHCrXTloEAScSck13ifbxCKLiIgAAH3DfdG1sRe+P3QNy/bEIyGrEDN/jsFnOy5jardQPNIuEI5K+2xREQQBB+OzsXLfVey+eMO4fGhLP7zUvyma+PAKTKo6FllERGTkoJDhmZ6N8HinYKw9kIAVe6/i6o1C/HdjLD7efhHjOwZjXIcghHg6ix2qWdws1GBjTCrWHUvGhYzSjv8SCTAg0hcv9m3KfldUIyyyiIioAheVHC/0bYrJ3ULxy7FkrDmYgOScYizdHY+lu+PRoWF9jIkKxOCWfnBzUIgdbpWUaPXYe+kGNsak4p+4TGj0BgCAo0KGse0DMaVbKBp62UcRSeJikUVERCa5qOSY2j0Uk7o2RHTcdfx0NAn7Lt/AscSbOJZ4E29uiEXnRp7oH+mLfpG+xgE0a5usAjX2Xb6BbbHXsefSDRRr/x3tvrm/G8a2D8LINv6o56QUMUqyNyyyiIjogWRSCQa1aIBBLRogI7cEG06l4veTKbicWYD9V7Kw/0oW3t50Do28nNGpkQc6N/JEx1APNHBzgERi3U7zgiAgI68Ep7MlOPHXBRxJuImL18uPARZQzxGDWjTAI+0CEenvZtX4qO5gkUVERFXSwN0Bz/VujOd6N8bVGwX45/x1RMddx/FrN3E1qxBXswrx09FkAICnsxKR/m6I9HNDWANXBHs4IdjDCd6uqhoXX4IgIKdQg8TsIlzLLsTVG4U4l5aLs6l5yCpQA5ABSDJuH97AFf0jfTGweQM093ezevFXVR7ObFWzdSyyiIio2hp5u+A/3i74T8/GyC3S4lhiDo4kZOPw1RzEpechu1CDfZezsO9yVrnnOSik8Hd3RH1nJeo7KeHhrICbgwIKuRQKmRRKmQQSiQQanQFqnQEanQHFWh2yCzTILtQgu0CNG/lqFJqY5FoqAXwdBPRpGYTuTX3QKdQDnjY0cbOTUo6T8/qLHQbVEIssIiIyC3cnBfrd7psFlHYwv3Q9H+fS8hCXlofLmflIzilGem4xSrQGXM0qBLIKa3RMiQTwd3dEiKcTQjydEOHnhhYB7mji6Yhd/2zDkCGRUChsq2M+2Q8WWUREZBEOChlaBdZDq8B65ZZrdAak5xYj7VYJbhVpcLNIi5tFGuQVa6HVC9DqDdDqDdAbBKgUUqjkMqjkUjgoZPBwVsLLRQlPFxU8nZXwr+d4z9HQtVqtlV4lkWkssoiIyKqUcilCPJ3tZqwtSyjR6jFp9VEAwLdTO3JaHRvFIouIiKiWMQgCjiTkGO+Tbaq7M38SERERWRCLLCIiIiILYJFFREREZAEssoiIiIgsgEUWERERkQXYfJG1dOlShIaGwsHBAVFRUdi3b5/YIREREdWYo0IGRw7dYNNsegiHdevWYdasWVi6dCm6deuG5cuXY/DgwYiLi0NwcLDY4REREVWLk1KO8+8NEjsMqiGbbslavHgxnnrqKTz99NOIiIjAkiVLEBQUhGXLlokdGhEREdVxNtuSpdFocOLECbz++uvllg8YMAAHDx40+Ty1Wg21Wm18nJeXB6B0CgZzTMNQtg9O6VARc2Mac2Mac2Mac2Mac1MRc2F9NltkZWVlQa/Xw9fXt9xyX19fZGRkmHzeBx98gHfffbfC8u3bt8PJycls8UVHR5ttX/aGuTGNuTGNuTGNuTHNVnOjNQCrL5aebJoaZoDCDOedioqKar4TqhKbLbLKSCSSco8FQaiw7E5z587F7NmzjY/z8vIQFBSEAQMGwM3NrcbxaLVaREdHo3///pz5/S7MjWnMjWnMjWnMjWm2npsijQ6vHNkJAOg/oB+clDX/ui47c0PWY7NFlpeXF2QyWYVWq8zMzAqtW3dSqVRQqVQVlisUCrP+Ipp7f/aEuTGNuTGNuTGNuTHNVnOjEP5tLCh9DTX/urbFPNg6m+34rlQqERUVVaEpODo6Gl27dhUpKiIiIqJSNtuSBQCzZ8/GhAkT0L59e3Tp0gUrVqxAUlISnn32WbFDIyIiojrOpouscePGITs7G/Pnz0d6ejpatGiBLVu2ICQkROzQiIiIqI6z6SILAKZPn47p06eLHQYRERFROTZfZNWUIAgAzHfVhVarRVFREfLy8tjJ8C7MjWnMjWnMjWnMjWm2npsijQ4GdemQC3l5edCZ8erCsu89sjyJUMeznZKSgqCgILHDICIisork5GQEBgaKHUadUOeLLIPBgLS0NLi6uhrH1+rQoQOOHTtWbru7l5l6XDbuVnJyslnG3TJ1PHM9x9Q2VVnO3Jhebo+5ud96c+XmzvvMDXNzv2W1PTc1+ay537rq5Obo0aPIz8+Hv78/pFKbHVzAptT504VSqbRCRS+TySr8Ut697EGP3dzczPqLfa+YzPEcU9tUZTlzY3q5PebmfuvNlZt7bc/cMDf3Wlbbc1OTz5r7ratObtzd3eHu7l7p2KnmWMrew/PPP//AZQ96bI2YzPEcU9tUZTlzY3q5PebmfuvNlRtL56W6x2Buqv+cupqbmnzW3G9dbcwNVVTnTxeaW15eHtzd3ZGbm2vWv57sAXNjGnNjGnNjGnNjGnNDtQFbssxMpVLh7bffvufUPXUdc2Mac2Mac2Mac2Mac0O1AVuyiIiIiCyALVlEREREFsAii4iIiMgCWGQRERERWQCLLCIiIiILYJFFREREZAEsskSSn5+PDh06oE2bNmjZsiVWrlwpdki1RnJyMnr37o3IyEi0atUK69evFzukWmX06NGoX78+xowZI3Yootu8eTPCwsLQtGlTfPPNN2KHU6vwfXJv/Hwha+IQDiLR6/VQq9VwcnJCUVERWrRogWPHjsHT01Ps0ESXnp6O69evo02bNsjMzES7du1w8eJFODs7ix1arbBr1y4UFBTg22+/xa+//ip2OKLR6XSIjIzErl274Obmhnbt2uHIkSPw8PAQO7Rage+Te+PnC1kTW7JEIpPJ4OTkBAAoKSmBXq8H691Sfn5+aNOmDQDAx8cHHh4eyMnJETeoWqRPnz5wdXUVOwzRHT16FM2bN0dAQABcXV0xZMgQbNu2Teywag2+T+6Nny9kTSyyTNi7dy+GDx8Of39/SCQSbNy4scI2S5cuRWhoKBwcHBAVFYV9+/ZV6Ri3bt1C69atERgYiDlz5sDLy8tM0VuWNXJT5vjx4zAYDAgKCqph1NZhzdzYuprmKi0tDQEBAcbHgYGBSE1NtUboFsf3kWnmzI2tfb6Q7WGRZUJhYSFat26NL7/88p7r161bh1mzZuHNN9/EqVOn0KNHDwwePBhJSUnGbaKiotCiRYsKt7S0NABAvXr1cPr0aSQkJODHH3/E9evXrfLaasoauQGA7OxsTJw4EStWrLD4azIXa+XGHtQ0V/dq+ZVIJBaN2VrM8T6yV+bKjS1+vpANEuiBAAgbNmwot6xjx47Cs88+W25ZeHi48Prrr1frGM8++6zwyy+/VDdE0VgqNyUlJUKPHj2E7777zhxhisKS75tdu3YJjzzySE1DrDWqk6sDBw4Io0aNMq6bMWOG8MMPP1g8VmuryfvI3t4nd6tubuzh84VsA1uyqkGj0eDEiRMYMGBAueUDBgzAwYMHK7WP69evIy8vD0DpbPF79+5FWFiY2WO1NnPkRhAETJ48GX379sWECRMsEaYozJGbuqIyuerYsSNiY2ORmpqK/Px8bNmyBQMHDhQjXKvi+8i0yuTGXj9fqHaSix2ALcrKyoJer4evr2+55b6+vsjIyKjUPlJSUvDUU09BEAQIgoAXXngBrVq1skS4VmWO3Bw4cADr1q1Dq1atjP0tvv/+e7Rs2dLc4VqVOXIDAAMHDsTJkydRWFiIwMBAbNiwAR06dDB3uKKqTK7kcjk++eQT9OnTBwaDAXPmzKkTV+dW9n1UF94nd6tMbuz184VqJxZZNXB3/w9BECrdJyQqKgoxMTEWiKp2qEluunfvDoPBYImwaoWa5AZAnbqC7kG5GjFiBEaMGGHtsGqFB+WmLr1P7na/3Nj75wvVLjxdWA1eXl6QyWQVWh8yMzMr/AVV1zA3pjE3lcdcmcbcmMbcUG3DIqsalEoloqKiEB0dXW55dHQ0unbtKlJUtQNzYxpzU3nMlWnMjWnMDdU2PF1oQkFBAa5cuWJ8nJCQgJiYGHh4eCA4OBizZ8/GhAkT0L59e3Tp0gUrVqxAUlISnn32WRGjtg7mxjTmpvKYK9OYG9OYG7IpYl3WWNvt2rVLAFDhNmnSJOM2X331lRASEiIolUqhXbt2wp49e8QL2IqYG9OYm8pjrkxjbkxjbsiWcO5CIiIiIgtgnywiIiIiC2CRRURERGQBLLKIiIiILIBFFhEREZEFsMgiIiIisgAWWUREREQWwCKLiIiIyAJYZBERERFZAIssIiIiIgtgkUVEdmvy5MmQSCSQSCTYuHGjWfe9e/du475HjRpl1n0TkX1gkUVkQ+4sGu683TlhLpU3aNAgpKenY/DgwcZlpoquyZMnV7pg6tq1K9LT0zF27FgzRUpE9kYudgBEVDWDBg3CmjVryi3z9vausJ1Go4FSqbRWWLWWSqVCgwYNzL5fpVKJBg0awNHREWq12uz7JyLbx5YsIhtTVjTceZPJZOjduzdeeOEFzJ49G15eXujfvz8AIC4uDkOGDIGLiwt8fX0xYcIEZGVlGfdXWFiIiRMnwsXFBX5+fvjkk0/Qu3dvzJo1y7jNvVp+6tWrh7Vr1xofp6amYty4cahfvz48PT0xcuRIJCYmGteXtRJ9/PHH8PPzg6enJ55//nlotVrjNmq1GnPmzEFQUBBUKhWaNm2KVatWQRAENGnSBB9//HG5GGJjYyGVShEfH1/zxN4lMTHxnq2GvXv3NvuxiMg+scgisiPffvst5HI5Dhw4gOXLlyM9PR29evVCmzZtcPz4cWzduhXXr18vd4rr1Vdfxa5du7BhwwZs374du3fvxokTJ6p03KKiIvTp0wcuLi7Yu3cv9u/fDxcXFwwaNAgajca43a5duxAfH49du3bh22+/xdq1a8sVahMnTsTPP/+Mzz//HOfPn8fXX38NFxcXSCQSTJ06tUIL3urVq9GjRw80bty4egm7j6CgIKSnpxtvp06dgqenJ3r27Gn2YxGRnRKIyGZMmjRJkMlkgrOzs/E2ZswYQRAEoVevXkKbNm3KbT9v3jxhwIAB5ZYlJycLAISLFy8K+fn5glKpFH7++Wfj+uzsbMHR0VGYOXOmcRkAYcOGDeX24+7uLqxZs0YQBEFYtWqVEBYWJhgMBuN6tVotODo6Ctu2bTPGHhISIuh0OuM2jz76qDBu3DhBEATh4sWLAgAhOjr6nq89LS1NkMlkwpEjRwRBEASNRiN4e3sLa9euvW++Ro4cWWE5AMHBwaFcHp2dnQW5XH7P7YuLi4VOnToJw4YNE/R6faWOQUTEPllENqZPnz5YtmyZ8bGzs7Pxfvv27ctte+LECezatQsuLi4V9hMfH4/i4mJoNBp06dLFuNzDwwNhYWFViunEiRO4cuUKXF1dyy0vKSkpdyqvefPmkMlkxsd+fn44e/YsACAmJgYymQy9evW65zH8/PwwdOhQrF69Gh07dsTmzZtRUlKCRx99tEqxlvn000/Rr1+/cstee+016PX6Cts+9dRTyM/PR3R0NKRSngAgosphkUVkY5ydndGkSROT6+5kMBgwfPhwfPjhhxW29fPzw+XLlyt1TIlEAkEQyi27sy+VwWBAVFQUfvjhhwrPvbNTvkKhqLBfg8EAAHB0dHxgHE8//TQmTJiATz/9FGvWrMG4cePg5ORUqddwtwYNGlTIo6urK27dulVu2YIFC7B161YcPXq0QhFJRHQ/LLKI7Fi7du3w22+/oWHDhpDLK/66N2nSBAqFAocPH0ZwcDAA4ObNm7h06VK5FiVvb2+kp6cbH1++fBlFRUXljrNu3Tr4+PjAzc2tWrG2bNkSBoMBe/bsqdDCVGbIkCFwdnbGsmXL8Pfff2Pv3r3VOlZl/fbbb5g/fz7+/vtvi/T7IiL7xnZvIjv2/PPPIycnB48//jiOHj2Kq1evYvv27Zg6dSr0ej1cXFzw1FNP4dVXX8WOHTsQGxuLyZMnVzgl1rdvX3z55Zc4efIkjh8/jmeffbZcq9QTTzwBLy8vjBw5Evv27UNCQgL27NmDmTNnIiUlpVKxNmzYEJMmTcLUqVOxceNGJCQkYPfu3fjll1+M28hkMkyePBlz585FkyZNyp3mNLfY2FhMnDgRr732Gpo3b46MjAxkZGQgJyfHYsckIvvCIovIjvn7++PAgQPQ6/UYOHAgWrRogZkzZ8Ld3d1YSC1atAg9e/bEiBEj0K9fP3Tv3h1RUVHl9vPJJ58gKCgIPXv2xPjx4/HKK6+UO03n5OSEvXv3Ijg4GA8//DAiIiIwdepUFBcXV6lla9myZRgzZgymT5+O8PBwPPPMMygsLCy3zVNPPQWNRoOpU6fWIDMPdvz4cRQVFWHBggXw8/Mz3h5++GGLHpeI7IdEuLujBRHVeb1790abNm2wZMkSsUOp4MCBA+jduzdSUlLg6+t7320nT56MW7dumX1KHWsfg4hsE1uyiMgmqNVqXLlyBfPmzcPYsWMfWGCV2bx5M1xcXLB582azxrNv3z64uLjcs7M/ERHAju9EZCN++uknPPXUU2jTpg2+//77Sj3no48+wn//+18ApVdTmlP79u0RExMDAPccIoOIiKcLiYiIiCyApwuJiIiILIBFFhEREZEFsMgiIiIisgAWWUREREQWwCKLiIiIyAJYZBERERFZAIssIiIiIgtgkUVERERkAf8PdGlwZvhP9t8AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "path = 'D:/GitHub/ULB-SeismoLab/Minions-Toolbox/ExampleData/Log/453004444/20231110141508/DigiSolo.LOG'\n",
    "rst = simulate_mseed(st, path)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b7a84b89-e3a6-47db-9bb1-06708d9d3cae",
   "metadata": {},
   "source": [
    "## Save data with the removed response"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "57d546c0-3cb3-4cdb-9649-ec389faf419b",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\julie\\anaconda3\\envs\\bas2\\lib\\site-packages\\obspy\\io\\mseed\\core.py:770: UserWarning: The encoding specified in trace.stats.mseed.encoding does not match the dtype of the data.\n",
      "A suitable encoding will be chosen.\n",
      "  warnings.warn(msg, UserWarning)\n"
     ]
    }
   ],
   "source": [
    "output_st = obspy.Stream(rtrace)\n",
    "output_st.write('./ResponseRemoved/CR_'+file, format=\"MSEED\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "07f1fbc2-b077-4a0e-94c7-0d04dfd060f1",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "bas2",
   "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}