ToOSimulations_with_redback.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "ca02b38c",
   "metadata": {},
   "source": [
    "## Notebook to show how to use redback to generate toO style observations for any model. \n",
    "\n",
    "You will need to install Redback. Instructions available at https://redback.readthedocs.io/en/latest/. I suggest installing from source via GitHub."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d27dcd12",
   "metadata": {},
   "source": [
    "## Incase you do it on google colab "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01e209f3",
   "metadata": {},
   "outputs": [],
   "source": [
    "!git clone https://github.com/nikhil-sarin/redback"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97f1c4eb",
   "metadata": {},
   "outputs": [],
   "source": [
    "cd redback/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8d388819",
   "metadata": {},
   "outputs": [],
   "source": [
    "pip install -r requirements.txt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eb48bea8",
   "metadata": {},
   "outputs": [],
   "source": [
    "!pip install lalsuite nestle"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "64bde9af",
   "metadata": {},
   "source": [
    "### You will also need to download latex for google colab"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "057c026d",
   "metadata": {},
   "outputs": [],
   "source": [
    "!sudo apt-get install -y texlive texlive-latex-extra texlive-fonts-recommended dvipng cm-super"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "ba8c5fc6",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:837: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:899: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:962: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1029: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1097: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1164: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1232: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1300: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1369: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1441: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1515: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n",
      "c:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:1588: SyntaxWarning: invalid escape sequence '\\l'\n",
      "  \"\"\"\n"
     ]
    },
    {
     "ename": "ModuleNotFoundError",
     "evalue": "No module named 'sncosmo'",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mModuleNotFoundError\u001b[0m                       Traceback (most recent call last)",
      "Cell \u001b[1;32mIn[1], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\n\u001b[0;32m      2\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mpandas\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mpd\u001b[39;00m\n\u001b[0;32m      3\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01msimulate_transients\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m SimulateOpticalTransient\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\__init__.py:1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m analysis, constants, get_data, redback_errors, priors, result, sampler, transient, \\\n\u001b[0;32m      2\u001b[0m     transient_models, utils, photosphere, sed, interaction_processes, constraints, plotting, model_library, simulate_transients\n\u001b[0;32m      3\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mtransient\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m afterglow, kilonova, prompt, supernova, tde\n\u001b[0;32m      4\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01msampler\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m fit_model\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\analysis.py:6\u001b[0m\n\u001b[0;32m      2\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mnumpy\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mnp\u001b[39;00m\n\u001b[0;32m      4\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mpandas\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mpd\u001b[39;00m\n\u001b[1;32m----> 6\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mmodel_library\u001b[39;00m\n\u001b[0;32m      7\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mutils\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m logger, find_nearest\n\u001b[0;32m      8\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mresult\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m RedbackResult\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\model_library.py:1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mtransient_models\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m afterglow_models, \\\n\u001b[0;32m      2\u001b[0m     extinction_models, kilonova_models, fireball_models, \\\n\u001b[0;32m      3\u001b[0m     gaussianprocess_models, magnetar_models, magnetar_driven_ejecta_models, phase_models, phenomenological_models, \\\n\u001b[0;32m      4\u001b[0m     prompt_models, shock_powered_models, supernova_models, tde_models, integrated_flux_afterglow_models, combined_models\n\u001b[0;32m      5\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mutils\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m get_functions_dict\n\u001b[0;32m      7\u001b[0m modules \u001b[38;5;241m=\u001b[39m [afterglow_models, extinction_models, fireball_models,\n\u001b[0;32m      8\u001b[0m            gaussianprocess_models,  integrated_flux_afterglow_models, kilonova_models,\n\u001b[0;32m      9\u001b[0m            magnetar_models, magnetar_driven_ejecta_models,\n\u001b[0;32m     10\u001b[0m            phase_models, phenomenological_models, prompt_models, shock_powered_models, supernova_models, tde_models, combined_models]\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\transient_models\\__init__.py:1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mtransient_models\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m afterglow_models, \\\n\u001b[0;32m      2\u001b[0m     extinction_models, kilonova_models, fireball_models, \\\n\u001b[0;32m      3\u001b[0m     gaussianprocess_models, magnetar_models, magnetar_driven_ejecta_models, phase_models, prompt_models, \\\n\u001b[0;32m      4\u001b[0m     shock_powered_models, supernova_models, tde_models, integrated_flux_afterglow_models, phenomenological_models, combined_models\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\transient_models\\afterglow_models.py:5\u001b[0m\n\u001b[0;32m      3\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mutils\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m logger, citation_wrapper, calc_ABmag_from_flux_density, lambda_to_nu\n\u001b[0;32m      4\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mconstants\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m day_to_s, speed_of_light, solar_mass, proton_mass, electron_mass, sigma_T\n\u001b[1;32m----> 5\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01msed\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m get_correct_output_format_from_spectra\n\u001b[0;32m      6\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mastropy\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01munits\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01muu\u001b[39;00m\n\u001b[0;32m      7\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mnumpy\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mnp\u001b[39;00m\n",
      "File \u001b[1;32mc:\\Python312\\Lib\\site-packages\\redback\\sed.py:4\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mtyping\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m Union\n\u001b[0;32m      3\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mnumpy\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mnp\u001b[39;00m\n\u001b[1;32m----> 4\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01msncosmo\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m TimeSeriesSource\n\u001b[0;32m      6\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mconstants\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;241m*\u001b[39m\n\u001b[0;32m      7\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mredback\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mutils\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m nu_to_lambda, bandpass_magnitude_to_flux\n",
      "\u001b[1;31mModuleNotFoundError\u001b[0m: No module named 'sncosmo'"
     ]
    }
   ],
   "source": [
    "import redback\n",
    "import pandas as pd\n",
    "from redback.simulate_transients import SimulateOpticalTransient\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f675ed98",
   "metadata": {},
   "source": [
    "We first design a strategy. This takes the form of a dataframe which specifies sky pointings, time, limiting mags, and the bands. Here now we have updated the pointings file manually for the Golden Stratergy for the GW case"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "6bcc6760",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "         expMJD  _ra  _dec filter  fiveSigmaDepth\n",
      "0  59582.180147  1.0   1.5  lsstg            25.0\n",
      "0  59582.492505  1.0   1.5  lsstu            25.0\n",
      "0  59583.492485  1.0   1.5  lsstz            23.0\n",
      "1  59583.658149  1.0   1.5  lsstg            25.0\n",
      "0  59583.669905  1.0   1.5  lssti            23.0\n",
      "1  59583.850379  1.0   1.5  lsstu            25.0\n",
      "1  59584.937064  1.0   1.5  lsstz            23.0\n",
      "2  59585.042648  1.0   1.5  lsstg            25.0\n",
      "2  59585.607560  1.0   1.5  lsstu            25.0\n",
      "0  59585.872426  1.0   1.5  lsstr            24.5\n",
      "3  59586.028325  1.0   1.5  lsstu            25.0\n",
      "2  59586.245133  1.0   1.5  lsstz            23.0\n",
      "1  59586.562660  1.0   1.5  lssti            23.0\n",
      "3  59586.893282  1.0   1.5  lsstg            25.0\n",
      "4  59587.355231  1.0   1.5  lsstu            25.0\n",
      "3  59587.823244  1.0   1.5  lsstz            23.0\n",
      "4  59587.895677  1.0   1.5  lsstz            23.0\n",
      "5  59588.444323  1.0   1.5  lsstz            23.0\n",
      "4  59588.597719  1.0   1.5  lsstg            25.0\n",
      "8  59589.617223  1.0   1.5  lsstz            23.0\n",
      "5  59589.624317  1.0   1.5  lsstu            25.0\n",
      "6  59589.715564  1.0   1.5  lsstz            23.0\n",
      "2  59589.832316  1.0   1.5  lssti            23.0\n",
      "7  59590.351273  1.0   1.5  lsstz            23.0\n",
      "7  59590.486153  1.0   1.5  lsstu            25.0\n",
      "5  59590.605219  1.0   1.5  lsstg            25.0\n",
      "9  59590.647158  1.0   1.5  lsstu            25.0\n",
      "9  59590.805927  1.0   1.5  lsstz            23.0\n",
      "6  59591.183543  1.0   1.5  lsstu            25.0\n",
      "3  59591.351637  1.0   1.5  lssti            23.0\n",
      "8  59591.410339  1.0   1.5  lsstu            25.0\n",
      "1  59591.703325  1.0   1.5  lsstr            24.5\n",
      "6  59591.870588  1.0   1.5  lsstg            25.0\n",
      "7  59592.569476  1.0   1.5  lsstg            25.0\n",
      "4  59593.718029  1.0   1.5  lssti            23.0\n",
      "8  59593.932699  1.0   1.5  lsstg            25.0\n",
      "9  59594.934338  1.0   1.5  lsstg            25.0\n",
      "5  59596.395199  1.0   1.5  lssti            23.0\n",
      "2  59596.906750  1.0   1.5  lsstr            24.5\n",
      "6  59598.434857  1.0   1.5  lssti            23.0\n",
      "7  59600.912609  1.0   1.5  lssti            23.0\n",
      "3  59601.071891  1.0   1.5  lsstr            24.5\n",
      "8  59603.302651  1.0   1.5  lssti            23.0\n",
      "9  59605.926214  1.0   1.5  lssti            23.0\n",
      "4  59606.209735  1.0   1.5  lsstr            24.5\n",
      "5  59611.315049  1.0   1.5  lsstr            24.5\n",
      "6  59616.778593  1.0   1.5  lsstr            24.5\n",
      "7  59622.079047  1.0   1.5  lsstr            24.5\n",
      "8  59627.068726  1.0   1.5  lsstr            24.5\n",
      "9  59631.340752  1.0   1.5  lsstr            24.5\n"
     ]
    }
   ],
   "source": [
    "pointings = pd.read_csv(\"pointings_data.csv\")\n",
    "print(pointings)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0da052b",
   "metadata": {},
   "source": [
    "We now define the function to convert the distance for Megaparsec to redshift"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a876273c",
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.cosmology import Planck18 as cosmo\n",
    "import astropy.units as u\n",
    "from astropy.cosmology import z_at_value\n",
    "\n",
    "# Given distance in Mpc\n",
    "distance_mpc = 350 * u.Mpc\n",
    "\n",
    "# Calculate the redshift corresponding to the distance of 350 Mpc\n",
    "redshift = z_at_value(cosmo.luminosity_distance, distance_mpc)\n",
    "redshift"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e200dfbe",
   "metadata": {},
   "source": [
    "We now specify a redback model (or a user implemented model) and the parameters we want to simulate given the above cadences. Here we use the Metzger 2017 Kilonova model at a luminosity distance of 350 MPc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "9bee8aac",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "16:28 redback WARNING : [Errno 2] No such file or directory: '/Users/nikhil/Documents/postdoc/redback/redback/priors/one_component_kilonova.prior'\n",
      "16:28 redback WARNING : Returning empty PriorDict.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'mej': 0.05, 't0_mjd_transient': 59582.0, 'redshift': 0.01, 't0': 59582.0, 'temperature_floor': 3000, 'kappa': 1, 'vej': 0.2, 'ra': 1.0, 'dec': 1.5}\n"
     ]
    }
   ],
   "source": [
    "model_kwargs = {'output_format': 'magnitude'}\n",
    "# Any redback model can be referred to as a string.\n",
    "# If the user has their own model, they can pass a function here instead.\n",
    "# There are over a 100 models implemented in redback, lots of models for kilonovae, GRB afterglows,\n",
    "# supernovae, TDEs and other things\n",
    "model = 'metzger_kilonova_model'\n",
    "# Load the default prior for this model in redback and sample from it to get 1 set of parameters.\n",
    "# We can sample from the default prior for this model for a random kilonova.\n",
    "parameters = redback.priors.get_priors(model=model).sample()\n",
    "\n",
    "# We fix a few parameters here to create a nice looking kilonova.\n",
    "# You can change any of the parameters here or add additional keyword arguments\n",
    "# to change some physical assumptions. Please refer to the documentation for this and units etc\n",
    "parameters['mej'] = 0.05\n",
    "parameters['t0_mjd_transient'] = 0\n",
    "parameters['redshift'] = 0.074816848\n",
    "parameters['t0'] = parameters['t0_mjd_transient']\n",
    "parameters['temperature_floor'] = 3000\n",
    "parameters['kappa'] = 1\n",
    "parameters['vej'] = 0.2\n",
    "parameters['ra'] = 1.0\n",
    "parameters['dec'] = 1.5\n",
    "print(parameters)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e2aa0d32",
   "metadata": {},
   "source": [
    "Now we use the Metzger 2017 Kilonova model at a luminosity distance of 700 MPc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "43b21755",
   "metadata": {},
   "outputs": [],
   "source": [
    "model_kwargs = {'output_format': 'magnitude'}\n",
    "# Any redback model can be referred to as a string.\n",
    "# If the user has their own model, they can pass a function here instead.\n",
    "# There are over a 100 models implemented in redback, lots of models for kilonovae, GRB afterglows,\n",
    "# supernovae, TDEs and other things\n",
    "model = 'metzger_kilonova_model'\n",
    "# Load the default prior for this model in redback and sample from it to get 1 set of parameters.\n",
    "# We can sample from the default prior for this model for a random kilonova.\n",
    "parameters = redback.priors.get_priors(model=model).sample()\n",
    "\n",
    "# We fix a few parameters here to create a nice looking kilonova.\n",
    "# You can change any of the parameters here or add additional keyword arguments\n",
    "# to change some physical assumptions. Please refer to the documentation for this and units etc\n",
    "parameters['mej'] = 0.05\n",
    "parameters['t0_mjd_transient'] = 0\n",
    "parameters['redshift'] = 0.14309633\n",
    "parameters['t0'] = parameters['t0_mjd_transient']\n",
    "parameters['temperature_floor'] = 3000\n",
    "parameters['kappa'] = 1\n",
    "parameters['vej'] = 0.2\n",
    "parameters['ra'] = 1.0\n",
    "parameters['dec'] = 1.5\n",
    "print(parameters)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ca22933",
   "metadata": {},
   "source": [
    "We now simulate a kilonova with the above parameters and the strategy designed above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "d955414c",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "16:28 redback INFO    : Using the supplied as the pointing database.\n"
     ]
    }
   ],
   "source": [
    "kn_sim = SimulateOpticalTransient.simulate_transient(model='metzger_kilonova_model',\n",
    "                                       parameters=parameters, pointings_database=pointings,\n",
    "                                       survey=None, model_kwargs=model_kwargs,\n",
    "                                        end_transient_time=20., snr_threshold=5.0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "1b3eae74",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "            time  magnitude  e_magnitude   band system  flux_density(mjy)   \n",
      "0   59582.180147  20.409959     0.003165  lsstg     AB           0.024909  \\\n",
      "1   59582.492505  20.239631     0.002705  lsstu     AB           0.029142   \n",
      "2   59583.492485  18.429289     0.003219  lsstz     AB           0.154529   \n",
      "3   59583.658149  19.876845     0.001942  lsstg     AB           0.040607   \n",
      "4   59583.669905  18.601859     0.003752  lssti     AB           0.132599   \n",
      "5   59583.850379  21.623112     0.009777  lsstu     AB           0.008064   \n",
      "6   59584.937064  18.750489     0.004338  lsstz     AB           0.114662   \n",
      "7   59585.042648  21.581913     0.009428  lsstg     AB           0.008362   \n",
      "8   59585.607560  24.393896     0.127071  lsstu     AB           0.000620   \n",
      "9   59585.872426  20.468957     0.005278  lsstr     AB           0.023675   \n",
      "10  59586.028325  24.567203     0.150044  lsstu     AB           0.000525   \n",
      "11  59586.245133  19.330835     0.007331  lsstz     AB           0.067853   \n",
      "12  59586.562660  19.930000     0.012938  lssti     AB           0.038450   \n",
      "13  59586.893282  22.524063     0.021727  lsstg     AB           0.003629   \n",
      "14  59587.355231  25.340487     0.277944  lsstu     AB           0.000284   \n",
      "15  59587.823244  20.150582     0.015586  lsstz     AB           0.031918   \n",
      "16  59587.895677  20.167839     0.016147  lsstz     AB           0.030808   \n",
      "17  59588.444323  20.469292     0.020944  lsstz     AB           0.023751   \n",
      "18  59588.597719  23.356844     0.049400  lsstg     AB           0.001596   \n",
      "19  59589.617223  20.990594     0.033668  lsstz     AB           0.014775   \n",
      "20  59589.624317  25.760812     0.758183  lsstu     AB           0.000104   \n",
      "21  59589.715564  20.994217     0.034814  lsstz     AB           0.014289   \n",
      "22  59589.832316  21.473338     0.055395  lssti     AB           0.008980   \n",
      "23  59590.351273  21.219385     0.042220  lsstz     AB           0.011782   \n",
      "24  59590.486153  25.789167     0.983188  lsstu     AB           0.000080   \n",
      "25  59590.605219  24.252614     0.099297  lsstg     AB           0.000794   \n",
      "26  59590.647158  26.373819     1.024206  lsstu     AB           0.000077   \n",
      "27  59590.805927  21.376891     0.047366  lsstz     AB           0.010502   \n",
      "28  59591.183543  25.785831     1.156168  lsstu     AB           0.000068   \n",
      "29  59591.351637  21.997492     0.081512  lssti     AB           0.006103   \n",
      "30  59591.410339  25.458105     1.209595  lsstu     AB           0.000065   \n",
      "31  59591.703325  22.916330     0.047650  lsstr     AB           0.002622   \n",
      "32  59591.870588  24.558999     0.128641  lsstg     AB           0.000613   \n",
      "33  59592.569476  24.273943     0.142803  lsstg     AB           0.000552   \n",
      "34  59593.718029  22.348474     0.113220  lssti     AB           0.004394   \n",
      "35  59593.932699  25.304951     0.165987  lsstg     AB           0.000475   \n",
      "36  59594.934338  24.498299     0.179152  lsstg     AB           0.000440   \n",
      "37  59596.395199  22.347813     0.133996  lssti     AB           0.003712   \n",
      "38  59596.906750  23.260227     0.074631  lsstr     AB           0.001674   \n",
      "39  59598.434857  22.577039     0.140353  lssti     AB           0.003544   \n",
      "40  59600.912609  22.561513     0.141376  lssti     AB           0.003519   \n",
      "41  59601.071891  23.590757     0.077399  lsstr     AB           0.001614   \n",
      "\n",
      "    flux_density_error  flux(erg/cm2/s)    flux_error  time (days)  detected   \n",
      "0             0.000073     4.163814e-14  1.214800e-16     0.180147       1.0  \\\n",
      "1             0.000073     3.065054e-14  7.644000e-17     0.492505       1.0   \n",
      "2             0.000458     6.139773e-14  1.823467e-16     1.492485       1.0   \n",
      "3             0.000073     6.803572e-14  1.214800e-16     1.658149       1.0   \n",
      "4             0.000458     8.224144e-14  2.863284e-16     1.669905       1.0   \n",
      "5             0.000073     8.571283e-15  7.644000e-17     1.850379       1.0   \n",
      "6             0.000458     4.567434e-14  1.823467e-16     2.937064       1.0   \n",
      "7             0.000073     1.414848e-14  1.214800e-16     3.042648       1.0   \n",
      "8             0.000073     6.679318e-16  7.644000e-17     3.607560       1.0   \n",
      "9             0.000115     2.228902e-14  1.088188e-16     3.872426       1.0   \n",
      "10            0.000073     5.693895e-16  7.644000e-17     4.028325       1.0   \n",
      "11            0.000458     2.676294e-14  1.823467e-16     4.245133       1.0   \n",
      "12            0.000458     2.420106e-14  2.863284e-16     4.562660       1.0   \n",
      "13            0.000073     5.940861e-15  1.214800e-16     4.893282       1.0   \n",
      "14            0.000073     2.793160e-16  7.644000e-17     5.355231       0.0   \n",
      "15            0.000458     1.257868e-14  1.823467e-16     5.823244       1.0   \n",
      "16            0.000458     1.238033e-14  1.823467e-16     5.895677       1.0   \n",
      "17            0.000458     9.378885e-15  1.823467e-16     6.444323       1.0   \n",
      "18            0.000073     2.758908e-15  1.214800e-16     6.597719       1.0   \n",
      "19            0.000458     5.802702e-15  1.823467e-16     7.617223       1.0   \n",
      "20            0.000073     1.896557e-16  7.644000e-17     7.624317       0.0   \n",
      "21            0.000458     5.783369e-15  1.823467e-16     7.715564       1.0   \n",
      "22            0.000458     5.841163e-15  2.863284e-16     7.832316       1.0   \n",
      "23            0.000458     4.700175e-15  1.823467e-16     8.351273       1.0   \n",
      "24            0.000073     1.847668e-16  7.644000e-17     8.486153       0.0   \n",
      "25            0.000073     1.209009e-15  1.214800e-16     8.605219       1.0   \n",
      "26            0.000073     1.078359e-16  7.644000e-17     8.647158       0.0   \n",
      "27            0.000458     4.065477e-15  1.823467e-16     8.805927       1.0   \n",
      "28            0.000073     1.853353e-16  7.644000e-17     9.183543       0.0   \n",
      "29            0.000458     3.604439e-15  2.863284e-16     9.351637       1.0   \n",
      "30            0.000073     2.506390e-16  7.644000e-17     9.410339       0.0   \n",
      "31            0.000115     2.339600e-15  1.088188e-16     9.703325       1.0   \n",
      "32            0.000073     9.117490e-16  1.214800e-16     9.870588       1.0   \n",
      "33            0.000073     1.185489e-15  1.214800e-16    10.569476       1.0   \n",
      "34            0.000458     2.608824e-15  2.863284e-16    11.718029       1.0   \n",
      "35            0.000073     4.586639e-16  1.214800e-16    11.932699       0.0   \n",
      "36            0.000073     9.641736e-16  1.214800e-16    12.934338       1.0   \n",
      "37            0.000458     2.610414e-15  2.863284e-16    14.395199       1.0   \n",
      "38            0.000115     1.704445e-15  1.088188e-16    14.906750       1.0   \n",
      "39            0.000458     2.113581e-15  2.863284e-16    16.434857       1.0   \n",
      "40            0.000458     2.144023e-15  2.863284e-16    18.912609       1.0   \n",
      "41            0.000115     1.257103e-15  1.088188e-16    19.071891       1.0   \n",
      "\n",
      "    limiting_magnitude  \n",
      "0                 25.0  \n",
      "1                 25.0  \n",
      "2                 23.0  \n",
      "3                 25.0  \n",
      "4                 23.0  \n",
      "5                 25.0  \n",
      "6                 23.0  \n",
      "7                 25.0  \n",
      "8                 25.0  \n",
      "9                 24.5  \n",
      "10                25.0  \n",
      "11                23.0  \n",
      "12                23.0  \n",
      "13                25.0  \n",
      "14                25.0  \n",
      "15                23.0  \n",
      "16                23.0  \n",
      "17                23.0  \n",
      "18                25.0  \n",
      "19                23.0  \n",
      "20                25.0  \n",
      "21                23.0  \n",
      "22                23.0  \n",
      "23                23.0  \n",
      "24                25.0  \n",
      "25                25.0  \n",
      "26                25.0  \n",
      "27                23.0  \n",
      "28                25.0  \n",
      "29                23.0  \n",
      "30                25.0  \n",
      "31                24.5  \n",
      "32                25.0  \n",
      "33                25.0  \n",
      "34                23.0  \n",
      "35                25.0  \n",
      "36                25.0  \n",
      "37                23.0  \n",
      "38                24.5  \n",
      "39                23.0  \n",
      "40                23.0  \n",
      "41                24.5  \n"
     ]
    }
   ],
   "source": [
    "# We can print the observations that were simulated to see what the data looks like. \n",
    "# This will include extra stuff like non-detections etc\n",
    "print(kn_sim.observations)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1de0cce4",
   "metadata": {},
   "source": [
    "One can also save the simulated model to a dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "7c535c0e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# We can also save the observations to a file using the save_transient method.\n",
    "# This will save the observations to a csv file in a 'simulated' directory alongside the csv file\n",
    "# specifying the injection parameters.\n",
    "kn_sim.save_transient(name='my_kilonova')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f26e070d",
   "metadata": {},
   "source": [
    "Now we plot the detection on a plot using the redback format"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "90532c0b",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "16:28 redback WARNING : [Errno 2] No such file or directory: 'kilonova//my_kilonova_metadata.csv'\n",
      "16:28 redback WARNING : Setting metadata to None. This is not an error, but a warning that no metadata could be found online.\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(24.0, 18.0)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnYAAAHdCAYAAACQZzRJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA8MklEQVR4nO3db2wj54Hn+Z/a7cj2tMVqCovkMumBu+gAxuEOE5Pii/jVnLuYHDA5TGKT0ty8OHgBk4Rx96onUbVezfrNyBQyebfIkHJwvgNuDi1W7AwuWZzNau/smwwQibRfLG59QLPagHdwCTaiSuoe24od8V70kBZFUiLF/6XvBzDcKj6sekgWW79+/s7V6/W6AAAAMPMuTboCAAAAGA6CHQAAQEAQ7AAAAAKCYAcAABAQBDsAAICAINgBAAAEBMEOAAAgIAh2AAAAAUGwAwAACIjLk67AWXzfVzqd1srKipLJZNdyGxsb2t3dbT4nFospk8mMq5oAAAATN7XBLpVKKRwOS5Icx9HKykrXstlsVrZtyzTN5rFCoaBsNqt8Pj/yugIAAEyDqQ12xWJRkuR5ngqFQtdyruvKMIyWUCdJmUxGsVhMvu/LMIxRVhUAAGAqzPwYu0qlIs/zOj5mmmbXxwAAAIJm5oOdaZpyHKdjq16lUlE0Gp1ArQAAAMZv5oNdMpmUaZrKZrNKJBLyfV++7yuVSjW7cwEAAC6CqR1j149yuaxUKiXXdXX16lVFo1HduXNnaGPrDg8PdXh4KEk6OjpSrVbT4uKi5ubmhnJ+AACAbur1uu7fv6+vfvWrunTp9Da5QAQ7wzCUSqVkGIYcx1GlUlE6nR5ai936+rpeffXVoZwLAADgPD766CN97WtfO7VMIIJdIpFQNptVsViU53nKZrNyHEeRSETlcnnglru1tTXdvHlTkrS/v68/+qM/0kcffaSFhYUh1B4AAKC7g4MDXbt2TU8++eSZZWc+2GWzWaVSqebixaZpqlQqNdexs2174LXs5ufnNT8/33JsYWGBYAcAAMamlyFgMx/sCoWC6vV62/FMJiPf93X79u0J1AoAAGD8Zn5W7GlO24IMAAAgaGY+2EWjUbmu2/Ex13VP3YoMAAAgSKY+2Pm+L0mq1WodHy8Wi8pms207TFQqFZVKJa2uro66igAAAFNhasfY2bYtz/NUqVSaP5dKJYXD4ZbJEI3JErlcruX5kUiEBYoBAMCFMlfvNPMAXR0cHCgUCml/f59ZsQAAYOT6yR5T3xULAACA3hDsAAAAAoJgBwAAEBBTO3kCAM7r6OhIn3/+ecfFywFg3Obm5nT58mVdujT69jSCHYDAODw81KeffqpLly6N7S9RADjL73//ex0eHuro6EiPPfZY2zalw0SwAxAIh4eH+t3vfqdQKDTpqgBAV/fv35ekkYU7/jkLIBA+/fRTPfnkk5OuBgCc6sknn9Snn346svMT7ADMvKOjI7pdAcyMS5cu6ejoaDTnHslZAWCMPv/8c12+zMgSALPh8uXL+vzzz0dyboIdgJlXr9dpsQMwMy5dujSyWfv8ExcAhuhnd9/Sx599rCcefULfffp7k64OoL/75Yf658PP9Qfzl/UXzz016epgxAh2ADBEf3/3Le1+uqvFxxYJdpgK/+c/fqj/cnCof7UwT7C7AOi7AAAACAiCHQAAQEAQ7AAAAAKCYAcAkO/7SiQSikQiisVik64OwD15TgQ7AIAMw1CpVFI0Gh35tSqVysivgdnHPXk+BDsAQFM8Hh/5NfL5/MivgeDgnuwPy53MqA9+/oY+++SBHn38ip75zkuTrg4A9MzzvElXAWgRpHuSFjsAwNg4jjPpKgAtgnZP0mI3o2ilAzBrPM9TOp3W0tLSpKsCSArmPUmwAwCcyvM8OY4j0zRVq9Xk+75M09T29rZyuVxPZQqFgkqlksLhsHZ2dpRKpZrnz+VyMk2z5ZqFQkG+72t3d1eRSESZTEaFQkGSVCwWVSqVxvcGYOqM6570fV83btyQ7/syDEPlcrlZZmNjQ/l8XrVaTcViUZZljf196IRgBwDoyvd92batYrHYcrxQKDTHJfVSJpPJKJPJKJvNyvO8trLHpVIpxeNxra6uNs9z9epV7e3tDfOlYUaN855shLlUKtU2Dm91dVXJZFKRSGSYL29gjLEDAHTluq7C4XDb8Uwm0zzeS5l+ruc4TjPUNc7j+35zLFQmk+nrnAiWcd+TUveZuSdbmqcBwQ4A0JVpmtra2uq4zlej66qXMr0qlUodf1lGo1Hdvn27r3MhmMZ9T84agh0AoKtoNCrLshSLxRSJRJTNZpstZ40xRb2U6VUkElGtVms73hgfBYz7npw1BDsAwKkakxUsy5LrukqlUopEIi2tIb2U6YVlWfJ9v2U8U+PnlZWVob0mzLZx3pOzhmAHAOiqEbAsy1I+n1e1WtXe3p6i0ahu3LjRc5nTHF9HzDTNlgHtjeUo8vn8WLaWwvQb9z05awh2AICuXNeV67otxwzDULFYVDgclud5PZU57mRX68nWuUQioWKxqEqlokqlos3NTSZMoGnc9+RpfN/v/wWMGMEuKH75I+nf/5uH/weAIeq2NMnxFrReykgPx9Cd9stwZ2dHtVpNhmEomUwqmUzKMIy+64xgG+c9KT0MhZ3KnAyP04BgBwBocfIXWKfZhY0yjQkNvZSRHnaNeZ7XfKxSqbT8ol1aWlIul1OlUpHv+1PZIoLxm+Q9KUnLy8stZRrnarT0TdNesyxQHBTP3Zx0DQDMMN/3lU6n5bqufN9XKpVSNptVOBzWnTt3muPdGr/IfN9vtoj0UqYhGo2qWCzKtm3FYjFJrevSGYahbDbbfOz48UboY6zdxTAt96T08P4rlUpKp9OKx+PNVuTGeFDbtlUqlU5deHtc5ur1en3SlZglBwcHCoVC2t/f18LCwqSrA0DS4eGhJGl+fn7CNZH+9f/9P2n3010tPrao//W//98nXZ2ZY9t2cwux4xrjpmzb1ubmppLJ5IRqOHv+h7/5B/2Xg0P9q4V5/V9/+SeTrg7U/99Z/WQPumIBAFPBcRx5ntdxokRjtuzm5qbW19cnUDtgNhDsAABToVarnbndk2EY59oSCrgoGGMHAEP0Z09/Tx9/9rGeePSJSVdl5mQyGSUSCTmO07Gr1fM85XI55XK5CdRudv2P33xK/3z4uf5gnl/5FwFj7PrEGDtg+kzTGDsMznEcbW9va3FxsTlIvTEbMZPJsPwJZt4ox9gR7PpEsAOmD8EOwCwZZbCjXRYD+eDnb+izTx7o0cev6JnvvDTp6gAAcKExeQIAACAgaLHDQGilAwBgetBiBwAAEBAEOwAAgIAg2AEAAAQEwQ4AACAgmDwBAMP0yx9JhwfS/IL03M1J1wZgWaoLhmAHAMP0jz+S7v+T9OQfEuwwFT74xRv6pPYbPR7+MsHuAqArFgAAICCmvsXO932l02mtrKx03BS6YWNjQ9VqVYZhyPO8M8sDAAAEzdQGu1QqpXA4LOnhhtArKytdyyYSCZmmqXw+3zwWiUQkiXAHAAAujKkNdsViUZLkeZ4KhULXcoVCQa7ram9vr+V4LpdTOp0m2AFAD3zfVyqVkud5MgxD5XJ50lXCBcc9eT4zP8Yul8spGo3KMIyW49FoVL7vy3GcyVQMAGaIYRgqlUqKRqMjv1alUjn18Ugkomw2O/J6YLpN0z05S2Y+2Hme1+yyPc40TUlSqVQad5UAYGbF4/GRX+P4sJluTv5jHRfXtNyTs2Jqu2KHZWdnZ9JVAAAc43neqY9Xq9Ux1QR46Kx7cpbMfLAzTVO1Wq3tuO/7Lf8fxOHhoQ4PDyVJBwcHA58PAC4qhsdg2gTtnpz5rthsNtuxb9x13aFdY319XaFQSKFQSNeuXRvaeQHgIvE8T+l0etLVAJqCeE/OfIvd6uqqbt++rUKhoEwm0/KYYRhDGaextrammzcfriB/cHBAuANwoXieJ8dxmj0kvu/LNE1tb28rl8v1VKZQKKhUKikcDmtnZ0epVKp5/lwu1xwXnUgkmAWJM43znpybm5Npmkomk82l1KSH4/I8z2tblWPSZj7YSVK5XJZt29rY2GgeW11dVSqV0tLS0sDnn5+f1/z8/MDnAYBZ4/u+bNtuLkHVUCgUmuOSeimTyWSUyWSUzWbleV5b2YZSqaRsNsv4aHQ17nvSNM22cZ+VSkXZbLbrcyYpEMFOUjOhNzTG1h1P1wCA/riu23HlgUwm02xR66VMPyKRCMEOXY3znvR9v+PSO6lUSplMZirXyp35MXbdNMbYneyeBQD0zjRNbW1tdRzL3Oi66qUMMCzjvCdrtVrbOnq2batWq03tEikzH+wcx9HVq1fbZr/m83nlcjnWQgKAAUSjUVmWpVgs1lw4uDGL0LKsnssAwzLOe9I0zZbyrutqY2NDd+7cGeIrGq6pD3aNwNZpSZPGcdM0WwKcbdsyDEOrq6tjqCEABFuxWFSpVJJlWXJdV6lUSpFIpKU1pJcywLBM6p5MpVJaXV0dy24Y5zW1wc62baVSqWaTaePnk33dmUxGlmUpm802/4tEIlM5oBEAZk1joLllWcrn86pWq9rb21M0GtWNGzd6LnOaoK0jhtGa1D2ZSqUUDodbxvR7njd1ixtP7eSJk5MhhlUWANA713XbuqMMw1CxWFQkEpHneT2VaSwdIbX3wEzbL0ZMt0nck47jyHGctokXlUpl6lrvprbFDgAwHbr1gBz/hdZLGenhjNdh7AiEi22c96Tv+0qn08rlcm3P3d7e7jj7dpIIdgCAFid/yXWaXdgo02j16KWM9LBrzPO85mPdWjwIfzhukvdkKpWSaZpt4/Z931ehUJi6SZpT2xULABifRquE67ryfb85pjkcDuvOnTvNsUSNLivf95stIr2UaYhGoyoWi7JtW7FYTNIXy1KdrEMikejYSoKLYRruScdx5LquMpmMCoVC8xzValVbW1tjeR/6NVev1+uTrsQsOTg4UCgU0v7+vhYWFiZdHQCSDg8PJWk6doj5m69J9/9JevIPpb/8z5OuDaCfvfIn+qT2Gz0e/rK+++N/mHR1oP7/zuone9AVCwAAEBAEOwAAgIBgjB0ADNM3b0qHB9I8QzUwHZ7505f02ScP9OjjVyZdFYwBwQ4Ahum5m5OuAdDime+8NOkqYIzoigUAAAgIgh0AAEBAEOwAAAACgjF2mFkf/PyN5oBgxpAAAECLHQAAQGDQYoeZRSsdAACtaLEDAAAICFrsAGCI7ucLqj94oLkrV/RkNjPp6gC4YAh2ADBEDwqbOvr1r3XpK18h2AEYO7piAQAAAoJgBwAAEBAEOwAAgIAg2AEA5Pu+EomEIpGIYrHYpKsD4JwIdgAAGYahUqmkaDQ68mtVKpWRXwO4qAh2AICmeDw+8mvk8/mRXwO4qAh2AICx8jxv0lUAAotgBwAYG8dxJl0FINAIdgCAsfA8T+l0etLVAAKNnScAAKfyPE+O48g0TdVqNfm+L9M0tb29rVwu11OZQqGgUqmkcDisnZ0dpVKp5vlzuZxM05TneUqlUvI8T8vLy83nSVKpVFKpVJrI6wdmCcEOANCV7/uybVvFYrHleKFQaI6V66VMJpNRJpNRNpuV53ltZSXJNE2Vy2UlEgnVajUVCgWtrq7KdV3Zti3P82Sa5oheKRAMdMUCALpyXVfhcLjteCaTaR7vpUw/otGoKpWKLMuSJFmWpb29PUId0AOCHQCgK9M0tbW11XHtuUZ3ai9l+lWr1VrW1DMM41znAS4agh0AoKtoNCrLshSLxRSJRJTNZpszWxstar2U6Retc8D5EOwAAKcqFosqlUqyLEuu6yqVSikSibS00PVSph/n6cIFQLADAJyiMfnBsizl83lVq1Xt7e0pGo3qxo0bPZc5DWvbAcNDsAMAdOW6rlzXbTlmGIaKxaLC4bA8z+upzHG1Wq3lZ3aiAIaHYAcAOFWnpUkktUxu6KWMJEUiEfm+P7S6AWhFsAMAtDgZvDrNeG2UaUxy6KWM9LC71vO85mOVSqUt/Pm+39aqB6A3LFAMAJDv+0qn03JdV77vK5VKKZvNKhwO686dO/I8T57nNQOX7/vNVrpeyjREo1EVi0XZtq1YLCbp4Xp30sMuWdu2W+pgmmZzdwsAZ5ur1+v1SVdilhwcHCgUCml/f18LCwuTrg4ASYeHh5Kk+fn5CddE+v9icR39+te69JWv6L8qb0+6OgCmUL9/Z/WTPeiKBQAACAiCHQAAQEAwxg4AhuhKJq36gweau3Jl0lUBcAER7ABgiJ7MZiZdBQAXGF2xAAAAAUGwAwAACAiCHQAAQEAQ7AAAAAKCyRMB8bO7b+njzz7WE48+oe8+/b1JVwcAAEwALXYAZt7c3JyOjo4mXQ0A6MnR0ZHm5uZGcu6pbrHzfV/r6+vyfb+5/+Da2pqSyWTH8hsbG9rd3dXi4qKq1aoSiUTXskFDKx0ussuXLze36AGAaff555+PbAvEqQ12vu/Ltm3lcjkZhiFJqlQqisViSiaTbRtLZ7NZRSKRls2iE4mEarVac4NpXDC//JF0eCDNL0jP3Zx0bTBCly5dosUOwMw4OjrSpUuj6TSdq9fr9ZGceUC2bWttba0Z6ho2NjZk27ZKpZIsy5L0ReA7+VK6HR9EPxvxYsL+5mvS/X+SnvxD6S//86RrgxE7PDzU7373Oz355JOTrgoAdHX//n196Utf6qvFrp/sMbUtdo7jyHEcVavVluONMFcsFpt/zufzikajbedoHHMc58J0yeKYb978osUOgdf4S3J/f1+XLl3S5cuXR/YvYgDox9HRkT7//HMdHR3pscceG1k3rDTFwc40Te3s7LQdb7Tg1Wq15jHXdTsGu0b5UqlEsLuI6H69cObn5zU/P9/8S3RKOyQAXDCPPPKI5ufnx/KPzakNdqVSqePxSqUiSYrH481jnuc1W+9OCofDHQNiPw4PD5sDsw8ODgY6F4DRu3Tpkr70pS9NuhoAMHYz10+Rz+dlGEZzQoTv+6eWNwzjzDJnWV9fVygUUigU0rVr1wY6FwAAwKjMVLBzXVeu62pzc7NtUsUora2taX9/X/v7+/roo4/Gdl0AAIB+TG1XbCepVEr5fL5lvNxZAW/Q1jrpi3E7AAAA02xmWuxSqZTW1tb6XpOuVquNtXUPAABgUmYi2Nm2rXg8rtXV1Y6PG4bRMkv2ON/3tbS0NMrqAQAATIWpD3aFQkGLi4ttoa5QKDT/vLy8LM/zup4jkUiMrH4AAADTYqqDneu68n2/Y0vd8bFzqVRKlUqlbTyd67qS1HUpFAAAgCCZ2i3FPM9TIpGQZVnNMXKN4OZ5nrLZbMskilQqJdM02/aKPVluUGwpBgAAxqmf7DG1wS4SiZzavVoul9t2m9jY2NDu7q4WFxdVrVaVSCSGvuMEwQ4AAIxTIILdtCLY4Swf/PwNffbJAz36+BU9852XJl0dAMCM6yd7zNQ6dsAs+OAXb+iT2m90+fErBDwAwFhN9eQJYBY986cv6b9J/s+am5P+o/Nv9cEv3ph0lQAAF8TQWuw+/PBDua6rarWq9fV1SdL+/r7K5bKef/75YV0GmHqN1rlHj7XYAQAwDkMZY3fr1i1tbGzIMAzNzc1pd3e3+dh7772nO3fu6Pvf//6gl5kKjLEDAADj1E/2GLgrdnNzU57nqVqtqlartSwcLEnPPvus0um0Xn/99UEvBQAAgFMM3BVbrVa1tbXV/Hlubq6tTCgU0tWrVwe91Ez5u19+qH8+/Fx/MH9Zf/HcU5OuDgAAuAAGbrGLRCItP3fr2b13796glwIAAMApBm6x69RC10m1Wh30UjOFVjoAADBuA7fY7e3t6c0332z+3CnovfLKK4rFYoNeCgAAAKcYyqzYpaUlzc3NaWVlRdVqValUSr7va3t7W4VCQZZl6fbt28Oo78QxKxYAAIzTRLYU29jY0K1btx6edG5O9XpdhmEol8spnU4P4xJTgWAHAADGaaJ7xd67d0+e58k0TV2/fn2Yp54KBDsAADBOE90r9vr164EMdAAAANNubHvFvvLKK+O6FAAAwIXUc4vd+++/f+6L1Go1bW1t6cc//vG5zwEAAIDT9Rzsnn/+ee3v73dcgLixxMlpjwEAAGC0eu6KDYfD2tnZ0d7eno6Ojpr/7e3tKZ1ONx87+d/bb7+tdDqtvb29Ub6OifjZ3bf0d//p/9DP7r416aoAAAD03mKXTCb17LPPth0vFova2NjoOkvDsizF43FtbW3p5ZdfPn9NAQAAcKqeg91rr73W8Xi9Xj9z6m0oFOq6h+ws++7T35t0FQAAAJrGtlcsY+2As33w8zf02ScP9OjjV/TMd16adHUAADNm4GB39+7dnspVq9VBLwUE3ge/eEOf1H6jx8NfJtgBAPo2cLBLJBJ65ZVXTl3KZG1tTfF4fNBLAYH3zJ++1GyxAwCgXwMHuxs3buidd97R4uJic6KEYRjyfV/ValVbW1taXl7WCy+8MIz6AoFGKx0AYBBD2yvWdV3dunVLlUqlecw0TeVyOb344ovDuMRUYK9YAAAwTv1kj6EFu+Pu3bsX2P1iCXYAAGCc+skeI9krNqihDgAAYJqNJNh1srKyMq5LAQAAXEgDT5549913zyzjeZ5c1x30UgAAADjFwMEumUzK9/2uj8/NzQVy1wkAAIBpM3BXbDgcVrVa1dHRUdt/e3t7un37tvL5vI6OjoZRXwAAAHQxcLDLZrNdJ0uEQiElk0ktLy/r9ddfH/RSAAAAOMXAwe4HP/jBmWVCoRDdsQAAACM2tlmx+/v747oUAADAhTSWYHdwcKBf/epX47gUAADAhTXwrNh4PH7q477vy/M8lUqlQS8FAACAUwwc7Pb29hSNRrsGPMMwtLy8rFAoNOilAAAAcIqBg51hGNra2hpGXQAAADCAgcfY3blzZxj1AAAAwIAGDnZ0sQIAAEyHgbtiP/zwQz311FNtx3/605/K8zzNzc0pGo3q+eefH/RSAAAAOMVcfcCVg1dWVnT79u1Ty7z33nsql8t6+eWXB7nUVDg4OFAoFNL+/r4WFhYmXR0AABBw/WSPgbtie8mFzz777KCXAQAAwBkGDnZzc3Nnljk4OFC5XB70UgAAADhFX2Ps3nvvPeXz+ebYOUna2dnRt7/97Y7la7Vac4HiYrE4eG0BAADQ1bnH2DmOo0wmo8XFRd24caNruUgkIsuyAtMdyxg7AAAwTv1kj4EmT3iep1u3bl2oBYoJdgAAYJzGNnnCNE0lEolBTgEAAIAhGXjyRDqd7qnc+++/3/e5fd+XbdvKZrNKJBKKxWJyHOfM56RSqTPLAQAABM3ACxT3yrZtvf322z2Xb4S6XC4nwzAkSZVKRbFYTMlksm0yRiqVUjgclvRw/N/KysrQ6g4AADAL+gp2r7zyiq5evaq//uu/bh6Lx+NnPq8xM7Yf6+vrLaFOkqLRqHK5nGzbluu6siyr+Vgj6Hmep0Kh0Ne1AAAAgqCvYPfOO+/o0qVLLcFub29PlmUpEol0fV69Xlcul+urYo7jyHEcVavVluONMFcsFluCHQAAwEXXV7A7GbKkhxMo/vZv//bM5+7s7PRzKZmm2fE5jRa8Wq3W1/kAtPvg52/os08eaO/e/6Or1/9rPfr4FT3znZcmXS0AwDkNPMbunXfe6anc5uZmX+ctlUodj1cqFUm9dQEDON0Hv3hDn9R+I81d0j+V/70eD3+ZYAcAM2xskydCodBQzpPP52UYhjKZzFDO14vDw0MdHh5KeriWDAAAwDQaeLmTXg1jlqrrunJdV5ubmy2TKkZtfX1doVBIoVBI165dG9t1AQAA+jG0Frv333+/67g33/ebXaiDSKVSyufzSiaTA5+rH2tra7p586akhy12hDsExTN/+pI+++SB/t9/97/ps48fTLo6AIABDRzs7t27p1gsJt/3Ty03Nzc30HVSqZTW1tbG2gXbMD8/r/n5+bFfFxi1xni6Rx+/os8+eaBHH78y2QoBAAYycLCzbVubm5uyLOvUcXTf+ta3BrpGPB7X6urquc8BoDsmTABAMAw8xi4ej+vFF188c3LEefeULRQKWlxcbAt1LEIMAADQauBgd/Xq1Z7K/eAHP+j73K7ryvf9ji11Z3X9AgAAXDQDd8XW63UdHBxoYWHh1HJvvvmmXnjhhZ7P63mestmsLMuSbduSvghzjcc6aZRhAWMAAHDRzNXr9fqgJ3n99de1tLSkb3zjG13LrKys6Pbt2z2fMxKJnLq/bLlcVjQabf5s27Y8z1OlUpHneTIMQ5ZlKRwOK5/P93zdsxwcHCgUCml/f//MMAsAADCofrLHwMHu29/+tqSHW4b5vi/TNNvWmPN9X57n6fe///0gl5oKBDsAADBO/WSPgbtit7e3ZVmW0um0FhcXO5b57W9/q9dff33QSwEAAOAUAwc70zS1tbV1Zrl79+4NeikAAACcYuBg5zhOT+U2NzcHvRSAKfTBz99oLm7MengAMFkDL3eSyWR66mYNhUJ68803tbS0pG9/+9t6//33B700gCnwwS/e0H90/q0++MUbk64KAFx4A7fYpVIpvfzyy2eWe++99/Tyyy9rc3NToVBItm2rWCwyAQEAAGBIBg52DW+++aZu374t0zT153/+5/rjP/7jlsdv3bqlXC6nF198UdLD5UwKhYK+//3vD6sKAAAAF9rAXbHpdFrhcFjJZFLValXFYlHRaFQ/+clPWsrt7Oy0bCt2/fp1DWEJPQAAAPyLgYPd66+/rlwup6OjI+3s7Oju3bva3d3V22+/rYODg2Y53/f11FNPtTz35Hp3AAAAOL+Bu2Lv3r2r1157reWYYRjK5XJyXVcvvPCC9vf3Oz53bm5u0MsDmLBn/vSl5qxYAMBkDRzsui1KPDc3d2ZXa2NfVwCziyVOAGB6DNwVaxiGfvjDH7Yc+/DDD5XNZptj6jzP6xjyfvvb3w56eQAAAPyLoUye+NWvfqVHHnlEX//617W4uKhIJKLr16/LdV1tbm4qlUopl8u1rHf35ptvKh6PD3p5AAAA/Iu5+pCmprquK9d1FYlEZFmWrl+/rv39fe3s7GhpaUmhUEi3bt3SvXv3FA6H5Xme3n777WFceqz62YgXQCt2qQCA/vWTPYYW7Dq5c+eOPM9TPB7XN77xDUnST3/6U9VqNaXT6VFddqQIdsD5/eyVP9Entd/o8fCX9d0f/8OkqwMAM6Gf7DFwV+xpbty4oXQ6rXq93uyGffHFF2c21AEAAEyzkQa7BtM0VSwWx3EpAACAC2sowe7dd99VPB7X4uJi23+PPPKIwuFwy64TAAAAGL6B17F77733lEwmlclklMlkVC6XFYvFFA6HVavVVC6XlUgkmnvEAgAAYDQGDnaFQkH37t1TKBSSJC0tLenq1avN7cPS6bTu3bund999V88///yglwMwar/8kXR4IM0vSM/dHOqp2aUCAEZr4K7YaDTaDHXSwwWLXddtKXP9+nV5njfopQCMwz/+SPoPrz78/5A9852X9N+m/heWOgGAERk42J3c7/X69esql8uDnhYAAAB9GjjY7e3tSXq4jdi7774rSQqFQvrJT37SUq5UKg16KQAAAJxi4DF2mUxGt27dkuM42tvb0+7urm7duiXTNJXP53Xjxg25rqulpaVh1BcAAABdDNxiFwqF9Nprr6lYLGpnZ0fSw3F2Ozs7+v3vf69cLqdwOKwf//jHA1cWAAAA3Q3cYtfw7LPPtvxsmiZj7QAAAMZoaMEOQEB88+YXy50AAGZKz12xa2trA13ohz/84UDPBzAmz92U/rt/M/Q17AAAo9dzi12lUhnoQtvb2wM9H9PvZ3ff0seffawnHn1C3336e5OuDgAAF07Pwa5UKunrX/+6DMPo+yK+77NA8QXw93ff0u6nu1p8bJFgBwDABPQc7EKhkK5fv36uixwdHZ3reZgtf/b095otdgAAYPx6DnYffvihtra2JEmWZfUd8paXl/urGQbyd7/8UP98+Ln+YP6y/uK5p8ZyTVrpAACYrL5a7NLptCTpvffe0+uvv65wOCzLsrSwcPbsuXg8fv5aAgAA4Exz9Xq9PsgJ7ty5I8/zFIlE9Pzzzw+rXlPr4OBAoVBI+/v7PQVaAACAQfSTPQYOdg37+/tyXVd7e3taWlrSN77xjWGcduoQ7AAAwDhNJNgdd+/ePbmuq7m5OVmWpaeeemrYl5gYgh0AABiniQe741ZWVuT7vrLZrF544YVRXmosCHYAAGCc+skePe880Y93331Xr7zyih555BEVi0Xt7u6O4jIAAAA4Zmh7xX744YfK5/MqFAra29uTaZp67bXXlEwmz73+HQAAAHo3ULA7ODjQ1taW8vm8KpWKQqGQlpeXlc1m9eyzzw6rjgAAAOjBuYLdm2++qdu3b8txHNXrdVmWpa2tLb344ovDrh8AAAB61PMYu/fff1+vvPKKFhcXlUwmVa1W9dprr2lvb0/vvPPOmaFubW1t4MoCAACgu55nxV66dElXr15VJpNRJpPpa9zc/v6+TNMMxCQKZsUCAIBx6id79NUVm06nFQ6H9dOf/rR5rJEL5+bmOj7nt7/9rQqFgvb39/u5FAAAAPrUc4vdt771Lb3zzjvnvtDTTz+tu3fvnvv504IWOwAAME4jWccukUgMVKlsNjvQ8wEAAHC6ke88ETS02AEAgHGa+M4TAAAAGD+CHQAAQEAMbUuxUfB9X+vr6/J9X57nqVaraW1tTclkcqCyAAAAQTS1wc73fdm2rVwuJ8MwJEmVSkWxWEzJZFLFYvFcZQHMmF/+SDo8kOYXpOduTro2ADDVprYrdn19vSWoSVI0GlUul5PjOHJd91xlAcyYf/yR9B9effh/AMCppjbYOY6jWCzWdtyyLElqaYXrpywAAEBQTW2wM01TtVqt7XijVe74Y/2UBQAACKqpHWNXKpU6Hq9UKpKkeDx+rrIAAABBNbXBrpt8Pi/DMJTJZIZa9jSHh4c6PDyU9HCRQAAAgGk0tV2xnbiuK9d1tbm52TJRYtCyZ1lfX1coFFIoFNK1a9cGOhcAAMCozFSwS6VSyufzPa1N10/Zs6ytrWl/f1/7+/v66KOPBj4fAADAKMxMV2wqldLa2lpP3ar9lO3F/Py85ufnh3IuAACAUZmJFjvbthWPx7W6ujrUsgAAAEEy9cGuUChocXGxLagVCoWBygIAAATNVHfFuq4r3/c7tr75vn/usgBmyDdvfrGlGADgVHP1er0+6Up04nmeEomELMtqzmptBDTP85TNZpsTI/opO6iDgwOFQiHt7+9rYYFfNAAAYLT6yR5TG+wikYg8z+v6eLlcVjQa7bvsoAh2AABgnPrJHlPbFVutVkdSFgAAIKimfvIEAAAAekOwAwAACAiCHQAAQEAQ7AAAAAJiaidPAAAugF/+6It1Cp+7OenaADOPYAcAmJx//JF0/5+kJ/+QYAcMAcEOACaJFisAQ0SwA4BJosUKwBAxeQIAACAgCHYAAAABQbADAAAICIIdAABAQBDsAAAAAoJgBwAAEBAEOwAAgIAg2AEAAAQEwQ4AACAg2HkCAACcje3vZgLBDgCA87pIYYft72YCwQ4AgPMi7GDKEOwAYJK+efOLFh8AGBDBDgAm6aK38hBsgaEi2AEAJueiB1tgyAh2uLB+dvctffzZx3ri0Sf03ae/N+nqAAAwMIIdLqy/v/uWdj/d1eJjiwQ7AEAgEOxwYf3Z099rttgBABAEBDtcWLTSAQCChi3FAAAAAoJgBwAAEBAEOwAAgIAg2AEAAAQEwQ4AACAgCHYAAAABQbADAAAICIIdAABAQLBAMQAAONs3b0qHB9L8wqRrglMQ7AAAOK+LFHaeuznpGqAHBDsAAM6LsIMpwxg7AACAgCDYAQAABATBDgAAICAIdgAAAAFBsAMAAAgIgh0AAEBAEOwAAAACgmAHAAAQECxQDGCm3c8XVH/wQHNXrujJbGbS1QGAiSLYAZhpDwqbOvr1r3XpK18h2AG48KY62Pm+r/X1dfm+L8/zVKvVtLa2pmQy2bXs4uKidnd3ValUlEgktLq6OoGaAwAAjN/UBjvf92XbtnK5nAzDkCRVKhXFYjElk0kVi8W2svl8vuXY9evXtb293VIWAAAgqKZ28sT6+npLqJOkaDSqXC4nx3Hkum7z+NbWlgqFQssxwzBkWZYcx5HneeOsOgAAwERMbbBzHEexWKztuGVZktTSCre0tNQSABvC4bCkh613AAAAQTe1XbGmaWpnZ6fteCPA1Wq15rFoNKq9vb22sq7ryjRNRaPRkdUTwGRdyaSbs2IB4KKb2mBXKpU6Hq9UKpKkeDx+6vMdx1GtVtOdO3cGrsvh4aEODw8lSQcHBwOfD8DwMBMWAL4wtcGum3w+L8MwlMm0/2XueZ4cx9H29rZ839e9e/c6dtH2a319Xa+++urA5wEAABiluXq9Xp90JXrluq4SiYSKxWLHJU8afN9XoVBQPp9XsVgcuCv2ZIvdtWvXtL+/r4WFhYHOCwAAcJaDgwOFQqGessdMtdilUinl8/lTQ530cBze6uqqqtWqYrGYqtWqTNM893Xn5+c1Pz9/7ucDAACMw9TOij0plUppbW2tYxfsac+RJNu2R1UtAACAqTETwc62bcXj8a67SMRiMSUSibbjjeVOGhMuAAAAgmzqg12hUNDi4mJbqCsUCs0/VyqVjosQN5ZEYbkTAABwEUx1sHNdV77vd2ypO77osGVZKpfLbWUaS6asrKyMrI4AAADTYmpnxXqep0QiIcuymkuWNMKc53nKZrPNSRSe58m2bW1ubjbLep6nSCSiTCbTsofsoPqZmQIAADCofrLH1Aa7SCRy6h6v5XK5pYvV933Zti3DMOT7flv4GxaCHQAAGKdABLtpRbADAADj1E/2mOoxdgAAAOgdwQ4AACAgZmrnCQDAxXY/X1D9wQPNXbmiJ7O9L1gPXBQEOwDAzHhQ2NTRr3+tS1/5CsEO6ICuWAAAgIAg2AEAAAQEwQ4AACAgCHYAAAABQbADAAAICIIdAABAQLDcCTClfnb3LX382cd64tEn9N2nvzfp6gAAZgDBDphSf3/3Le1+uqvFxxYJdgCAnhDsAAAz40om3dx5AkA7gh0AYGaw2wRwOiZPAAAABATBDgAAICDoigWm1J89/b3mrFgAAHpBsAOmFDNhAQD9oisWAAAgIAh2AAAAAUGwAwAACAiCHQAAQEAQ7AAAAAKCWbEAAKBv9/OF5vZu7AgyPQh2AACgbw8Kmzr69a916StfIdhNEbpiAQAAAoJgBwAAEBAEOwAAgIAg2AEAAAQEwQ4AACAgCHYAAAABQbADAAAICIIdAABAQLBAMQAA6NuVTLq58wSmB8EOAAD0jd0mphNdsQAAAAFBsAMAAAgIgh0AAEBAEOwAAAACgmAHAAAQEAQ7AACAgCDYAQAABATBDgAAICAIdgAAAAFBsAMAAAgIgh0AAEBAEOwAAAACYqqDne/7sm1b2WxWiURCsVhMjuP0/PxYLDbC2gEAAEyXy5OuQDeNUJfL5WQYhiSpUqkoFospmUyqWCye+vxsNqtKpTKGmgIAAEyHqW2xW19fbwl1khSNRpXL5eQ4jlzX7frcSqWinZ2dMdQSAABgekxtsHMcp2NXqmVZknRqi93t27e1srIysroBAABMo6kNdqZpqlartR1vtOB1ekySNjY2tLa2NsqqAQAATKWpHWNXKpU6Hm+Mm4vH4x0fM02zpfsWAADgopjaYNdNPp+XYRjKZDJtj92+fVu5XG7o1zw8PNTh4aEkaX9/X5J0cHAw9OsAAACc1Mgc9Xr9zLIzFexc15XruioWi22tcqPsgl1fX9err77acuzatWsjuRYAAEAn9+/fVygUOrXMTAW7VCqlfD6vZDLZcnzUXbBra2u6efOmJOno6Ei1Wk2Li4uam5sbyfUwXgcHB7p27Zo++ugjLSwsTLo6GCI+2+Disw0uPtt29Xpd9+/f11e/+tUzy85MsEulUlpbWxtrF2zD/Py85ufnmz8zhi+YFhYW+EskoPhsg4vPNrj4bFud1VLXMBPBzrZtxeNxra6utj3mOI4qlYqy2WzL8cY6do3jJ9fEAwAACJqpD3aFQkGLi4ttoa5QKCiTySiZTLZ1zUpf7DyRz+fHVVXMqPn5ef3VX/1VS6ssgoHPNrj4bIOLz3Ywc/VeplhMiOu6qlQqHVvqNjY2Oh5vyGazKhQKPc0gAQAACIKpbbHzPE/ZbFaWZcm2bUkP9489/thpui1gDAAAEFRT22IXiUTkeV7Xx8vlsqLRaNvxQqGgUqkkx3EkPdxfdmlpiS5ZAAAQeFMb7AAAANCfqd0rFgAAAP0h2AEAAAQEwQ4AACAgCHYAAAABQbADAAAICIIdAABAQBDsAAAAAoJgBwAAEBAEOwAAgIAg2AEAAATE5UlXAJiERCKhaDSqlZUVRaNReZ6nfD4v3/fZV3iG+L6vdDqtlZUVJZPJruU2Nja0u7urxcVFVatVJRKJU8tj8nr5bPkezx7f97W+vi7f9+V5nmq1mtbW1rp+xnx3+0eww4VUq9W0sbGhjY2N5jHLslQqlSZYK/QqlUopHA5LkhzH0crKStey2WxWkUhEuVyueSyRSKhWqymTyYy8ruhPP58t3+PZ4vu+bNtWLpeTYRiSpEqlolgspmQyqWKx2FKe7+75zNXr9fqkKwGMWyqVUjwe1/b2tkzTVCKRkGVZk64W+uR5niKRiIrFYsd/xTd+aZz8a67bcUyPsz5bie/xrLFtW2tra81Q17CxsSHbtlUqlZqfH9/d86PFDhdSOBzW6urqpKuBEcvn84pGo23HG8ccx6FbZ4bxPZ4tjuPIcRxVq9WW440wVywWm3/mu3t+TJ4AEFiu68o0zY6PGYZBlx0wRqZpqlartR1vtOAdf4zv7vkR7HChVSoVFQoFVSqVSVcFI+B5XnO81knhcFg7OztjrhFGge/xbCiVStrb22s73vjc4vF48xjf3fMj2OFCqtVqsm27OQi3VqspFovJ87xJVw1D4vv+qY8bhnFmGUw3vsfBkM/nZRhGc0IE393BMMYOF1IikWiZVWVZllZWVpRIJNrGfwCYTnyPZ5/runJdV8VisW1SBc6HFjtcSJ2myluWJc/z5DjOBGqEYTvrlwT/4p99fI9nXyqVUj6fb5kIwXd3MAQ74F80BuoyKPdiqNVqtBAEEN/j2ZFKpbS2ttb3mnR8d09HsMOFk0qlFIvFuj7eadYWZpNhGF0/T9/3tbS0NOYaYVj4Hs8227YVj8e7LlfDd/f8CHa4cCqVSsfZVo2/RI7PzMJsW15ePnUgfSKRGGNtMEx8j2dXoVDQ4uJiW6grFArNP/PdPT+CHS6cZDLZsZumMSaHrWqCI5VKqVKptI3JcV1XktilYIbxPZ5NruvK9/2OLXXHv6d8d8+PYIcLZ21tTdlstuVYpVLR+vo6M7NmTOMv/W5dNpZlKZlMan19veV4Lpfjs55yZ322fI9nj+d5ymazqlarsm1btm0rm80qm80qkUi0LEjMd/f82CsWF1JjM2rDMJrN/Wtrax23sMH0sW1bnuepUqnI8zwZhiHLshQOh5XP59vKb2xsaHd3V4uLi6pWq0okEmxHNKX6+Wz5Hs+WSCRyavdquVxu++z47vaPYAcAABAQdMUCAAAEBMEOAAAgIAh2AAAAAUGwAwAACAiCHQAAQEAQ7AAAAAKCYAcAABAQBDsAAICAINgBAAAEBMEOADBSjuMoEokoEokokUhMujp9y2azzfpvbGxMujrAqS5PugIAgODLZrNaXV2ddDXOpbFHreM42t7ennBtgNMR7IAZEIlEVKvVFA6HZZqmDMOQJLmuK0kyTVOmaUp6uDH6zs6OJCmXyymTyUykzmjneZ5SqZQ8z9PS0pJKpdKkqzR1bNuW4ziq1Wq6c+dO26bwAE5HsANmgOd5yufzbSHNtm1tbGyoWCw2g530MNylUilVq9WW8pFIRJZlNVsgZtksvhbTNFUul2eyO7JSqYwlZOVyOa2srCgWi438WkAQMcYOmHK+7yuZTJ7a8hYOh1t+NgxD+Xxevu+3lW209gXBrL6WWWyFGmeAPv6PFAD9ocUOmHK1Wk0rKyt9P+94l23DyRa8WRak1zILPM8b27VmNbAD04AWO2DK+b5/7haMxcXFIdcGF5HjOJOuAoAeEeyAKWcYxrmD3Sx2+WG6eJ6ndDo96WoA6BFdscCUG2S8kWVZzT8nEgl5nifDMFQulyV1nqVZKBQkPezqrFQqyuVyikajqlQqzdm25XJZkUik6/IVGxsbze60arWqSCTS8+xcz/PkOI5M01StVmu2WG5vbyuXy53rtfi+33x+t/ezsZTF4uKidnd3u9Z5kNd2ku/72traav5cLpdl23bbRJgbN27I9/2W19uoSz6fV61WU7FYbH7ex9+L5eVl5XK55ntRKpWas3HPeq8LhYJKpZLC4bB2dnaUSqWa1z7tveyH67oqlUqKRCLNY8vLy13L+76vQqHQ/AzK5bKy2WzzHzGNWbWN+2Nzc1PJZLLlPbNtW4ZhKJPJKJfL9XTPATOjDmBmra6u1iXV9/b2eiqfyWTq0Wi07XgymaxbllXP5/Mt5yoWi3XDMOrlcrleKpVanmOaZj2fz7edKxqN1ovFYtv5V1dXz6zf3t5ePZlMth3P5/Ntx8/7WjpZXV1tO//e3l5bnQd5bSevd9r7d/K9blyn0+utVqt1SR2fY1lWPZlM1nO5XL1er9dLpVJdUr1arfb9XluW1fPrO6lYLDbrcNzq6mo9k8m0HNvb26tnMpm6pHq5XO74nJPlDcNoef17e3t1SR2vWa8/fF+q1WqzbK/vQ7FY7PuzBsaNrljgAjneKnKcaZra2dnR0tJSy8B1y7Lk+77y+XxL61/jsWKx2HKssSr/8RYS6WHrzsbGxpkD8F3XbZvhK0mZTKbt+HlfS6VSabvmxsaGNjc3W443WnEadR70tXXSqaUvl8u1tIw1xOPxjuc4rdWs0dLa+Owsy9Le3p5M0+zrvR4F13VVKBTaZtsahqFsNtvxOZVKpeUzaZTPZDKybbvt2O3btzueJ5FINN+3Sb8PwLAR7ABocXFRvu+3jclrBKNOIcowDNVqtZZjtm13nMHbmKHbWFC5G9M0tbW11Ra+JHUMO52c9Vo61TmZTLbNxKzVai1lB31tJ3Wb+dkIoMPauqpWq7W8F43rDuO9HoRt223/WGjoFlYb99zJEB2JRNqOZbNZVSqVjkH+eDif9PsADBtj7ABIOrvl5yyNX6zVarXrLMqzliiJRqOyLEuxWEymacqyLCUSCSWTya4hoJN+xn4db9E6rtG6JQ3ntfWjMb5rWOfqZFjv9XlVKpW+l/ExTbP5mUgPPxff91vGHTZEo1FFo1Hl8/mWVsGTn/ek3wdg2Ah2ACQNvnZYI/ykUqmOvxBPdmF2UywW5bpu8/+FQkGmaapYLPY8y7fX19Ko81nLwgzrtfXKMIyhrRt3WnfiMN7r82i8tvPcc77va319Xb7vKxaLNUPZ8UkoDdlsVtlsthnsui0dNKn3ARgFumIBDMXxvWrPq/ELv7FVWLVa1d7enqLRqG7cuDGMarZo1Hl3d7encoO8tn4MsnZhrwZ9rwdZ2+6876fnebp+/boikUhziz3TNLuG18YYxsaM4EKh0BbCx33PAaNGsAMwFI2xZqd1IZ71i9x13baxaoZhqFgsKhwOj2T3A9M0Tz1vI2QN+tr64Xle18kSw7puv+/1ybGJg34Wpmn23X2dSqVkmmbbpJOTdTv+ujKZzKnboU3ingNGiWAHYGiOr5d20vF18E5zcqZtw6i6xHK5nBzH6RiQPM9r/tIfxms7rlsgcxxHhmG0rRFoGEbH5/Q7aeO4Xt/rSCQy9NbKXC7XsftUUtf3stt4yHK53FK/42GsMYliY2Oja5f5uO85YJQIdsAMa/wyO9li0ctzjtvd3T31F3evv9QzmYwsy2pbrsL3fbmu29Ng9E4zFBvXP9k9OYzXkkwmlclkOs6AdBynGQaG8dpOXvdkd2Zj/NidO3fayi8vLzcnCxwv3/jsO7UsHX+8k17fa8uyWq5dqVQGDj3JZFLLy8sty5Q0rt8IWifrbllWW5D1PE+JRKL53EqloqWlpebjjUkUt2/f7tq93c89B0y7uXq9Xp90JQD0rlKpNAeP7+zsNLsKo9GowuFwx24n3/eVTqfluq5835dlWc2dA44fTyaTWllZUTKZ1MbGhkqlklzXbZ6/WCzK8zzZtt3ynGw22xJsCoVCc3eKxmr+vezO0Fj9vxFSGr/Yfd9vtmAN87WcvHZjB4TT6nze13byHI3dFU7uPJHL5bpOKnBdV/l8XvF4vFkmk8lobm5OhmE01xbs9BmZptmyi0Iv73Wn9ycWizWv26vG2nPdzru9vd2288TVq1dlGIaWl5db7ulsNqtardYMc42ZrI3POJVKtdWt0dLaqc79vA+NurIbBaYZwQ4AMFKnBbtx2NjYGMq1CXaYBXTFAgAABATBDgAQGI1txxpO7jQBBB3BDgAQGPl8Xuvr682fPc9jAgQuFHaeAAAERmPSSKFQ6DoRBAgyJk8AAEbKcRyl02mFw2GZpqlSqTTpKvUlm83Kdd3mDGgmT2CaEewAAAACgjF2AAAAAUGwAwAACAiCHQAAQEAQ7AAAAAKCYAcAABAQBDsAAICAINgBAAAEBMEOAAAgIP5/ZNRn2uxBPbgAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "kn_object = redback.transient.Kilonova.from_simulated_optical_data(name='my_kilonova', data_mode='magnitude')\n",
    "\n",
    "# Make a dictionary for colors on the plot\n",
    "band_colors = {'lsstg':'#4daf4a', 'lsstu':'#377eb8', 'lsstr':'#e41a1c', \n",
    "               'lsstz':'#a65628', 'lssti':'#ff7f00'}\n",
    "ax = kn_object.plot_data(show=False, band_colors=band_colors)\n",
    "ax.set_ylim(24, 18)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "821caf8f",
   "metadata": {},
   "source": [
    "The above plot only shows the detection and not the input lightcurve or non-detections. Let's add those in. As the axes is returned with can use the attributes stored in the kn_sim object directly in the original redback format."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "77448ec0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.1, 20.0)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAokAAAHdCAYAAABxIONdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAACiS0lEQVR4nOz9aXRUaZ7neX7vYvtupl0IkIkdfGPx8IjKpSoCMrunsqdyAffqdbImA5iY7jeTleW4nzNv8pyZ8sA7K2fO6e6pgohez0x1u0NEVk2vmchjcqmMCHdAvrCDZBIIENpMtslk673z4kpCAgFazGRm0v8Tx8JM167ZfSTk0k/P8n8U0zRNhBBCCCGEWECtdQOEEEIIIUT9kZAohBBCCCGeIyFRCCGEEEI8R0KiEEIIIYR4joREIYQQQgjxHAmJQgghhBDiORIShRBCCCHEcyQkCiGEEEKI50hIFEIIIYQQz9Fr3YBXSSQSnDx5kvfee4/jx4+/8LyPP/6YycnJ+dccOnSIU6dOrVczhRBCCCE2lLoNiSdOnCAcDgNw8eJF3nvvvReee/r0ac6cOUM0Gp0/dv78eU6fPs25c+eq3lYhhBBCiI2mbkPihQsXAIjFYpw/f/6F5/X29hIMBhcFRIBTp05x6NAhEokEwWCwmk0VQgghhNhwGn5OYl9fH7FYbMnnotHoC58TQgghhBAv1vAhMRqNcvHixSV7G/v6+jh48GANWiWEEEII0dgaPiQeP36caDTK6dOnOXbsGIlEgkQiwYkTJ+aHrIUQQgghxMrU7ZzElbh69SonTpygt7eXUCjEwYMH+eyzzyo2FzGfz5PP5wEwDIN4PE4kEkFRlIq8vxBCCCHEi5imSTqdpqOjA1Vdv/69DRESg8EgJ06cIBgMcvHiRfr6+jh58mTFehI/+ugj/viP/7gi7yWEEEIIsRrDw8Ns2bJl3a63IULisWPHOH36NBcuXCAWi3H69GkuXrxIT08PV69eXXOP4ocffsgf/uEfApBMJtm6dSvDw8P4/f4KtF4IIYQQ4sVSqRRdXV34fL51vW7Dh8TTp09z4sSJ+ULb0WiUS5cuzddJPHPmzJprJTocDhwOx6Jjfr9fQqIQQggh1s16T3Nr+JB4/vx5TNN87vipU6dIJBJ88sknNWiVEEIIIURja/jVzS/zsm38hBBCCCHEizV8SDx48CC9vb1LPtfb2/vS7fyEEEIIIcTS6j4kJhIJAOLx+JLPX7hwgdOnTz+3s0pfXx+XLl3i/fffr3YThRBCCCE2nLqdk3jmzBlisRh9fX3zH1+6dIlwOLxoIcrcQpWzZ88uen1PT48U0xZCCCGEWCXFXGrVh3ihVCpFIBAgmUzK6mYhhBBCVF2tskfdDzcLIYQQQoj1JyFRCCGEEEI8R0KiEEIIIYR4Tt0uXBFCiNUyDINSqbRkoX0hhFhviqKg6zqq2lh9cxIShRAbRj6fJ5fLoapqQ/5AFkJsTOVymXw+j2EYOJ3O57b6rVcSEoUQG0I+n6dQKBAIBGrdFCGEeKF0Og3QEEFR/swWQmwIuVwOn89X62YIIcRL+Xw+crlcrZuxLNKTuAGViwUK00kKmRSF6RSlfBaWmJql2R3Y3D7sbi82tw+by4siw3OiARmGIUPLQoiGoapqQ/zckpC4QRRnpsmMDpNPT1EuWH+hqDY7dk8Aj68DRVEWnW+aJqX8DIVMkuzE4/njNrcPd7gVV7gVzWZf189BiNUqlUrouvw4E0I0Bl3XKZVK2O31/XtWfqo2uMJ0ivTIfXLJCVSbHXeoBZsngN3rR7c7l/UehlGmNDNNMZsml5oi+WiA5MN+HP4w7kg7zmAEVdWq/JkIsXqmadb9X+RCCDFHVdWGqL4gIbFB5VNTpJ/cJ5+eQnO4CG7bgzvcuqrhYlXVsHv82D1+PM2dGKUiM1NjZCefMDV4A0XT8TZ34m3tQtVtVfhshNg4/lX/n5EtZnHb3Pz2jt+pdXOE4F/+fIjpfAmPQ+ff+872WjdHNBAJiQ3GMMokhm4zMzWGzeUl1L0fV6j5ueHktVB1G57mTjzNnZRyWaYnHpMZGyYz/ghf2za8LVtk7qIQL/Cv+/+MydwkEWdEQqKoC//9L4YYT+Vp9jskJIoVkZDYQEr5GeID1yjlZwh178Mdbq36NXWnm8CWHXhbukiNDJF6NMD02EN8Hd24I20VDadCCCGEqB8SEhtELjlJfPAmqm6jefchbG7vul5fszsIbduNr7WL1ONBEvdvMz02TGj7vnVvixBCCCGqT8YM65xpmqRHhpjs/wa7x0/LnvUPiAvpTjfh6H6a9xwCYOz2FdJP7jfEBFwhhBBCLJ/0JNYx0zSZGrzJzNQYvvbt+Nq3Pz+8a5pQmIZcAvJJMA1QNFB1UGfvNQe4wqBXbqm93eOnec9hUo8HST2KkUtMEurei+5wVewaQoj1k0gkOHHiBLFYjGAwyNWrV2vdJLHJyfdk7UlIrGOpRzFmpsYIR/fjCrU8fSKfgdSwFQxzCTBK1nGbxwqGRhnMMlz/FIrToLtg72+DIwDuJvA0W6FxjWVtFFUlsKUHZzDC1OAtxm5eJtC1A09Tx5reVwix/oLBIJcuXZr/pVxNfX19HDx4sKrXEI1PvidrT4ab61R28gmZ0QcEtux4GhAL0zDyJQz9JSSHrV7C8A7Y8g7s+Lcg+l3Y/uvWfc8xaN4LkV3Wre1NcPgg/Qge/hL6/zd43Ae55Jrb6vAGadl3BFeohcT9O8QHb2IaxprfVwix/o4cOVL1a5w7d67q1xAbh3xP1o70JNahQibJ1P3buCNteFu7oDgDk/cg+QB0B7QcgMBWeFUZmu/84eKPA13WfT4D06OQGIL7fw3uZgj3WD2Mq6RqOqHte3D4Q0wN3aJcyBPu2Y9WwSFuIcTGUO1eISFWSr4nlyYhsc6UCjkmB65jd/sJdu2C8VswFbN6DZv3QnD7moeJcXitWygK6RGI91u9iw4/hHeCf/XDxe5wK5rdQbz/OuO3+2ja8Tq607229gohNoyLFy/WuglCLCLfky8mw811xDDKxPuvoSgK4eg+lNGvID5gBbfu71m9fZXcHk9RrEC4/ddgy7etBS4jV2H4l1DIrvptHd4gzXsPoaAwdvsq+fRU5doshGhYsViMkydP1roZQsyT78mXk57EOpIYuk0pn6V55xtoo1/BzCR0HAZfW/Uv7mmybpkxGP3GmvfYtAdC3VaYXCHd4aJ5z0HisetM3Pua0LY9uCPr8HkIISouFotx8eJFotEo8XicRCJBNBrl8uXLnD17dlnnnD9/nkuXLhEOh7ly5QonTpyYf/+zZ88SjUYXXfP8+fMkEgkmJyfp6enh1KlTnD9/HoALFy5w6dKl9fsCiLqzXt+TiUSC733veyQSiedWWH/88cecO3eOeDzOhQsXOHr06Lp/HapNQmKdmB5/ZK1k3rYT2+Q1KGSsBSnuyPo2xNsCrr8LE7dh/Ia10KX1DXD6V/xWqm4jsuMNEg/uMjV0C9M08TS1V77NQoiqSSQSnDlzhgsXLiw6fv78+fl5XMs559SpU5w6dYrTp08Ti8WeO3ehEydOcOTIEd5///359wmFQkxNyaiEWN/vyblguNQK6/fff5/jx4/T09NTyU+vrshwcx0wyiVSj4dwBYK40vegmIWu76x/QJyj6dB6ALb+ilVe5/5fQ3x1k3oVVSW0fQ/upg4S92+TjY9WuLFCiGrq7e0lHA4/d/zUqVPzx5dzzkqud/HixfmAOPc+iURifu7YqVOnVvSeYmNZ7+9JePEK62d7wDcaCYl1IDM6jFHIECiPWDUOt/4KOAO1bha4QrDt160FLuM3rPI7qyxtE9y6C3ekjanBm2TjYxVuqBCiWqLRKJ9++il9fX3PPTc3PLecc5br0qVLS/7iPXjwIJ988smK3ktsTOv9PbmZSUissXIhT+bJfQJqyioXs/VXwO6pdbOeUlVo2QftByH9GIb/Foq5Fb+NoigEt+3BFWphavAGM4nxKjRWCFFpBw8e5OjRoxw6dIienh5Onz4936M3NwdrOecsV09PD/F4/Lnjc/PJhFjv78nNTEJijaUeD6LPjON26VbBa5uz1k1amr8Tuv4OlHLw4G9gZuVzgxRFIdS9D1eohXjsBrnkZBUaKoSotLmFIkePHqW3t5cTJ07Q09OzqJdmOecsx9GjR0kkEovmf819/N5771XscxKNbT2/JzczCYk1VJyZZmZ0EJ8tixrZYa0urmeuIGz7NWubv+GfQ/Lhit9CURRC2/fi9EeYHLhGZuwRpmlWvq1CiIqYC2tHjx7l3LlzDAwMMDU1xcGDB/ne97637HNeZmGdumg0umgxwVyJknPnzsm2aQJY/+/JzUxCYg0lh+/hmBnGGW6Hpr21bs7y6A5rUY2vE558CVNDK34LRVUJR/fjjrSTHL7L+K3L5DOJijdVCLF2vb299Pb2LjoWDAa5cOEC4XCYWCy2rHMWenY4+dlew2PHjnHhwgX6+vro6+vjRz/6kSxWEfPW+3vyZRKJxMo/gQYiIbFG8qkpSqO38HpdKO0HrRXFK1DT3jdVhfY3rQUtY9dWtfJZUVVC23bTvOcQKCoTd74kHrtBuZCvfHuFEGvyonI1C3v2lnMOWHMOX/aL9cqVK8TjcYLBIMePH+f48eMEg8EVt1lsbOv5PQlWwFzqnGeD6EYjIbEGTNMkOXQdV3EMx5Y3XjjMbGSzFO/codD3JfkvviD/tz8n91d/Ra73M/Kf/Yz8F19QvHuX8tgYZr4G4aplP4R3WCuf4wOregu7x0/znkMEt+0mn5pi9MbnpEeGMMqlCjdWCLFcz/4yXGqV6Nw5c4tJlnMOWMN/sVhs/rm+vr5Fv7QPHz7M2bNn6evrI5FIbPieGrE8tfyeBHj33XcXnTP3XnM9kBt172cppl0DM5NPUEa/xtWxDZr3Pfd8OT5F+f4Q5fEJFIcdNRBAcXtQbDqKzQ42HQwDI5Wi/GSU0tB9AFS3G7WtDX1rF4rdvj6fTPNeQIHxm2AaENm54rdQFAVPUweuYDPpkSFSI0OkR4fxtXbhae5E1W2Vb7cQYpFEIsHJkyfp7e0lkUhw4sQJTp8+TTgc5rPPPpufHzj3SzGRSMz31CznnDkHDx7kwoULnDlzhkOHDgGL6x4Gg0FOnz49/9zC43MBUuYmbg718j0J1vffpUuXOHnyJEeOHJnv3Z6bP3vmzBkuXbr00iLxjUgxZdXAiqRSKQKBAMlkEr9/5buQAIx9/mc4Zx7iP/IPwdMMWL2LxsgTSg/uY6TSqF4v2rZtaO1tKOrLO3zNmRmMZBIjPkV5ZAQAbUsn+rZtKM51Wi09cRcm70BkNzTtWtNblQo5Mk8eMD3xGEXV8LZswduyRcKieKH8bE+6w+GocUvgH/1v/xGTuUkizgj/9b/139W6OQ3nzJkz89vwLTQ3z+zMmTP86Ec/4vjx4zVqYeP5d/7ZXzKeytPsd/A//uO/W+vmCFb+M6sS2WM1pCdxnRXSU6jxezh3fXtRQCxeu0b5yShaJIL90EG0yPJ3W1FcLjSXC62tDX1HD6UHw5SHH1AefojW3o7WvR3V7a7Wp2Rp2mXt8TxxGxQVIjtW/Va63Ulw6y58bdvIjD6YvQ3j7+jG07IFZRV7SQsh6t/FixeJxWLze+8uNLfqORwO89FHH0lIFGIdSEhcZzMPr6OqJrYtbwJgGoYVEMfGsb/xOlpr65reX7Hbse3oQd++jfLDh5Tu36c8MoLeE0Xbvr26ASuy09oxZuIW2Nzg71jT22l2B4GunXjbtpEeGSL5sJ/s5BOC2/dgd/sq1GghRL2Ix+Ov3DItGAyuals1IcTKSUhcR6ZpUhy5jjPSjeLwWQHxm28wJiatgNjSUrFrKbqOvn07WlcXpViM4r1+yhMT2A4cQHW5Knad5zTvsfaefvIl6E5wr/2HuWazW9v6hVtJPLjD+K0reFu78HV0o6paBRotROX8gx2/Q7aYxW2rcu/9BnTq1CmOHTvGxYsXl+wpnOtlXKqnUbzYv/vt7UznS3gc8itfrIzMSVyhtcwLyI0PMf3F/wff4RPYmnfMB0TbG6+jNTdXqcUWY2qKwvUbUCyi796F3tlZxYsZ8PCXUEhXfJtB0zDIjD4gNTKEZnMQ3LYbp196FTa7epqTKNbu4sWLXL58mUgkMr9AYG5V6alTp6Qkjmh4jTInUULiCq3lHyrZ968x4vcJ/r3/M8VvrmFMxrG9+QZa0/rstGKWSpTu3KX06BFacxO2/furtwq6XIQH/wZM0wqKemWvU8plSTy4Qz6dwBVuJdC1w9r7WmxKEhKFEI2kUUKi1ElcJ0apQGnsDvaO/RSv37AC4ltvrltABGsI2rZ/H/Y338BIJil8cRkjm63OxTQbdH4LjBI8vmzNVawg3emmaddbBLftJpecZOz6F2TjYxW9hhBCCLGZSUhcJ/lH16BcRPdtpzw6hm3/vhWtYK4kraUF+7e+BUDhiy8wksnqXMjuhs63IZeEJ19V5RKepg5a972N3RdkavAGk/3fyK4tQgghRAVISFwnhUfXUINb4PEEasCP1t5e0/aoLhf2t4+guN0UrlylPDFRnQu5gtD+FqQfr2qf5+XQ7A4iPQcI9xygkE0zeuNzpiceV+VaQgghxGYhIXEdlDNjlKYeoTk7MZIp9J0r35WkGhS7HfuhQ6jhEIUvv6L0qErBytcOwW5r+75cqjrXAFzBZlr3vY0r1Ezi/h2rV7FYqNr1hBBCiI1MQuI6yD/4GkO1o6d1tEgErY5qfCmahu3NN9E7OijeuEEpNlidCzXvA7sXHl+BKu7LrOo2Qtv3Wr2KmRRjN75gZkrmKgohhBArJSGx2spFik9uoejNKIUi+s7V70RSLYqiYNu/D72nh2J/P6WhocpfRFWh4zCU8zD6TeXf/xmuYDMt+9/G7g0Qj90gPngTo1Ss+nWFEEKIjUJCYpWVJgYoTqfR8h60tlbUdVy6vlK2nih6tJvi3XvVGXq2e6D1dUg/guRw5d//GZrNTmTHawS37SGXmGDs5mXymUTVryuEEEJsBFJ+vcoKD7+inFdw+fzoO+qvF/FZth07oFCkePMmik2v6C4wAPg7ITsBo9fAGQRH9bfX8zS14/CHmIrdZOLOl/g7uvG2bZM9oEV1/PxPIZ8Chx++84e1bo0Q3P6f/huKMxlsLi97fuv3a90c0UAkJFbTzBT58ftQ8GHfuhXV3RjbdOl792CWihS/uQYH36r8HMqWAzAzBY+vWoW2tep/G+p2J0273yL9eJDU40Hy6SlC2/eh2aX4sqiwX/yp1Vvu65SQKOrC7f/5v2EmPoor3CohUayIDDdXUXHyPoWJSZy+LvTu7lo3Z9kURbH2eA6FKH71deXrKKoadByy9ni+/1cwPV7Z938BRVHwd0aZGrzFvb/4H7j63/zfmUlUqfSPEEII0eAkJFZRYXQAI1PCtXsfitNZ6+asiKKq2N58A8XrofDllxjT05W9gMMH238dbG5rn+eRr6C0PuVqdKcbd3MHmt1JfOAayYf9mIaxLtcWQgghGoUMN1dR8dFdFHsIWwP1Ii6kaBr2t96i8MVlil99jf1bb6PoFfyWsXug69uQeADjN2F61BqK9ndW7hpLWDjckhkdJvmwn0ImRTi6X4afhRBCiFnSk1gt+Qyl+Dha01YUm63WrVk1xWbD9tabmPk8xevXMU2z8hcJboXuvwfuJhjpg+FfwvT6DAN7W7to3n2QciHH2M3L5JKT63JdIepNIpHg2LFj9PT0cOjQoVo3Rwj5nqwDEhKrpJx+QimTxtFRH7urrIXqdmM7sJ/y2DjlwaHqXER3WPMUO45YtRQf/gKG/hpSj6DKQ8F2b4DmfYexeXxM9n9D6lGsOmFYiDoWDAa5dOkSBw8erPq1+vr6Xvp8T08Pp0+frno7RH2rp+/JzUpCYpUUhu9iKHYcWxtzqPlZWksLejRKsb+f8mQVe9t8bdZcxS3vgGa3ehYHfwbxWFXDoqbbiex4HX9nlPST+0ze+0q29BOb0pEjR6p+jXPnzr3ynGAwWPV2iMZQL9+Tm5HMSayS0sg9TLsfW1NrrZtSMXpPFDOVovjNNyjvvIPqclXvYp5m65ZLwdSANWcxeR9a3wB3dbY1VBQFX9s27B4/8dhNxm5dJhI9gN0bqMr1hNisYrHYS58fGBhYp5YIYXnV9+RmJT2J1VAuUpp4hN68DVV9+ZfYMA3ypRypfJKp3BRlo7xOjVw5RVGwvXYAxWan+NVXmOV1aKvTD+1vWb2Lqg2G/9YqxF3F/Z8dvhAtew+j252M3/2S6fFHVbuWEJvNxYsXa90EIRaR78kXk57EKjDTYxQzSWy7n5+PaJomj6cfk8onyZfzFI3iovlvqqIScASJuCIEHAE0RVvPpr+SYrNhe+N1a8XzrVvYDxxYnws7fLD170BiCMZvQWbU2uLPW+EdYWZpdgdNu94i+bCfxIO7FKZTBLbuQlXr699DiEYSi8U4efIkhw8frnVThADke/JVJCRWQfFRP2XDxNu1OCSWzTKxxACJfIKQM4TP7sOhObBrDhyaVXolnosTz8Xpn7qHqqiEnCFa3W147d5afCpLUn0+bPv2Urh2nXIkgtbevj4XVhQIdYOnFUa/gUefg38LNO8H3V75y6kqwa27sHv8JB7coZjNEO45gO6o4jC7EHUoFotx8eJFotEo8XicRCJBNBrl8uXLnD17dlnnnD9/nkuXLhEOh7ly5QonTpyYf/+zZ88SjUYBOHbsGLFYjGAwyNWrV2vy+Yr6t57fk4qiEI1GOX78OD09PfPnnDt3jlgsxtTU1Pp+8utIQmIVFB/dxXQEsUeanx4zitybuku2mGVnaBdBR3DJ13Z4O+jwdpAr5YjnJpmYmeBW/CbtnnY6vJ2oSn3MENDa29EmJijevo0aDKJUc37is+xu6HoHksMwdgOmx6paX9EdacPm8jI5cI3xW1cJR/fj8Ieqci0h6k0ikeDMmTNcuHBh0fHz58/Pz+NazjmnTp3i1KlTnD59mlgs9ty5cy5dusTp06e5cuVKFT4bsRGs9/dkNBp9bp5sX18fp0+ffuFrNor6SBwbiWlSGhtCCbTO9zjlSjluTt4kX86zJ7z3hQFxIafupMPbyWtNr9Pp7WRkeoQbE9eZLlZ455M1sO3ZA5pG4cbN2pSMCXQtrq/48HMozlTlUja3l+a9h7C5PUzc+4rM2MOqXEeIetPb20t4if3bT506NX98OeesxMLeGiGetZ7fk4lEYslyTCdOnODUqVMcP358BS1vPBISK8xIjlGaTqG37wAgU0hzc/IGCrAvvG/Fw8aKotDh7WR/5ACKonJz8gYP0w8xzNpvI6fYbNj278eIxykPD9emEQvrK+ZTMPSXMDUIVQitmm4nsvNNvK1dJIfvMTV0S7bzExteNBrl008/XbKO3Nzw3HLOEaJS1vN7Mh6PP1en8cyZM8Tj8U1RNkdCYoUZIzHKpQL2LbvIFrPcjt/GqbvYG9mHQ1/9/s1um5t9kX10eDsZmX7MzckbFMq1r+OnRSLoW7so3b1X+f2dV8LXBtv/Lvg6Yew6DP8CCtmKX0ZRFAJbdhDavpdsfJTxO32UC/mKX0eIenHw4EGOHj3KoUOH5otcz60GPXr06LLPEaJS1vN7MhqNLjq/t7eXjz/+mM8++6yCn1H9kpBYYaWRAcouP/ZQEyPTI9g0G3vCe7Cpa9+aT1VUOr2d7Ivsp2iUuBW/Sb6Uq0Cr10bfuRPF5aJ47Vpte9Y0G7S9Dlu+DcUs3P8ra1/oKnBH2mjecwijWGDs1hUKmWRVriNEPbhw4QKXLl3i6NGj9Pb2cuLECXp6ehb10iznHCEqpVbfkydOnOD9999fl11g6oGExAoyTZPSaAzT24zqdBPPTdLsaqn4YhOPzcO+yD5A4Wb8Ftli5XvMVkLRNGwH9mOkM5TqoSCpp8mqq+hth9Gv4eEXUKp8b5/d7aN572F0h8uqpzgxUvFrCFFrc5P8jx49yrlz5xgYGGBqaoqDBw/yve99b9nnvIzUqRMrUavvyRMnThAOh+dXT89dZyMX4paQWEFmfJxSNoEa2cJUOQVAs6v5Fa9aHYfmYF94HzZV51b8JplCpirXWS41EMDW00NpcAgjkahpWwCrV7H9TWuuYi5hzVVMVz7EaTY7TbvexB1uI3H/NsmH/bLvs9hQent76e3tXXQsGAxy4cIFwuEwsVhsWecsFI/HF328kX/JisqrxffkxYsXuXjx4nOrmTd6T7mExAoyngxRNkronTsYzY4RdoaxaWsfZn4Rayh7Ly7dxZ2p26TytR3y1Lq3o/q8FG/dqp+gNDdX0RWGx1fgyTdQ4V1tFFUltH0Pga6dZEaHmez/BqNUrOg1hKilF5X5WDjktpxzwFq5nKiHPyRFQ1vP78lEIsHJkyc5e/bsc6+9fPnyqlbwNwoJiRVkjA5RsjsoOF3kSzma3dXZDWQhXdXZHd6Dx+blztQdEvlE1a/5IoqiYNu7FyOdofygOnMBV0W3Q+cRa4eW1DDc/xvIpyt+GW/LFiI736AwnWL8dh+lXG2nAQixWs/+wlxqlejcOXMFh5dzDljDf7FYbP65vr6+Jed3SZAUC9Xye/LEiRNEo1Hef//9597r/PnzBIPBVX5W9U+KaVeIWSpRmniA4YmQVHM4dSd+u39drq0pGrtCu+hP3GMg0c/e8D7cNve6XPtZaiCA3rWF0kAMrbUVxbn6Fd0VF9w226N41QqKLfutYxXk9Idp2XOIyf5rjN22Cm87/Rv3r0yxccz1lvT29pJIJDhx4gSnT58mHA7z2Wefzc+9mhuWSyQS8z01yzlnzsGDB7lw4QJnzpzh0KFDgFW7bqk2HDt2bMneG7E51MP35MWLF+nt7eXUqVOcP39+/j0GBgb49NNP1+XrUEuKWTfjgo0hlUoRCARIJpP4/U9DYHlslOz/+v8k2d7D43372Brqps2zTtvVzbXBKHMrfouiUWR/ZD92rfJb1S2HWSyS//nPUYNB7G+8UZM2vJRRtnZqSd4HX4fVw1jhaQFGqUh88Cb5VJxA1y68LdXZDUZY8nlrYZLD4ahxS4B/tgXSj6xyTP9Yiq6L2vtXP/i7zMRHcYVb+e1//pe1bo5g5T+zXpQ9qk2GmyvEGL1P2SyT9nhQdRtNVVqw8jKaavUoAtyduku5wnPvlkux2bDt2kV5dIzyxERN2vBSqmaVyuk4DNPjcP+vIVfZ+ZyqbiOy43W8LVtIDt8l8eCuFN4WQgjRUOp6uDmRSPDRRx+RSCTmu4s//PDDF26D8/HHHzM5OUkkEmFgYIBjx46t25Y5xsRDyg6dhF2n1RlGV2vzpbVrdnaHdnMrfpOBZD87g7tQFGXd26G1t6M9HqF46zbqd76Nomnr3oZX8rWDI2AtaHnwb6z9nys4/KwoCoGuneguD4n7dyjlpglHD6Dq1VvMJIQQQlRK3YbEuc25z549Oz8ptK+vj0OHDnH8+PHn5hScPn2anp6eRfWLjh07Rjwen59bUE1mZpyMAYbTSYu7terXexm3zU1PcAf3pu4ynH7AVn9l590tl753D4Wf/4LS4CC2HTtq0oZXsrth66/A+A0Y/Qayk7PDz5X7T8PT1IHucBEfuMHY7atEdryGzemp2PuLOvPtP7S2iHSs35CQEC+z5+//PsWZDDbXyraFFaJu5ySeOXOGDz/88LlVQx9//DFnzpyZr6IOT8Pjs5/Ki46vxVLzAkzTJPs//JBBM0XpO/973tz+7Ypdby1Gp0e5nxpiq38bbZ62mrSh2D9AeWgI+7ffQfXUeTBKPbKCou60hqIdvoq+fSk/w2T/NcrFvCxoqbC6mpMohBCvIHMS1+jixYvzq4wWmguGC3sSz507t+Tqt7ljVa/mXyhQyCWYxqAtUpteu6W0elpp87QxnH5AplD5ki/LoUe7UZxOirdu1+T6K+LvhK2/CooK/9sfWref/2nF3l53uGjecxC7x8/kva/JjMmiBiGEEPWrbkNiNBp9rgI6MN+zuPC53t7eRTWPnj3/0qVLVWnjHDOXYzo3ieoL0eyp7VDzs7b4uvDYPPQn+ikZpXW/vqKq6Lt3YcTj9bmI5VkOrzX87PBZO7SkR6CCC05UTV+woOWeLGgRQghRt+o2JF66dImpqannjs8Vxjxy5Mj8sVgs9sKK5+FwmCtXrlSnkbPMXIZcPos31Imm1tcCDVVR6Qn0UDbLDKUGa9IGrbkZNRSidO9e/ezE8jKqBv/Wn8LRf2rtAf3wF1DMVezt5xa0BLftZnr8EZP9X8sOLUIIIepO3YbEFzl37hzBYHBR8dWXCQaDa67cn8/nSaVS87dnmdk0pXIRe6BpTdepFofuZLu/m/hMnPHsWE3aYNu5AyOdwRh5UpPrr0qoG7q+A8WsVSYn+3zP9lp4mjpo2vUmhWyGsdtXKeamK/r+QgghxFo0VEic27D7Rz/60bpug/PRRx8RCAQIBAJ0dXU993w5m8RQTBzu+l2IEHFFaHY3cz91n5nSzLpfXw0G0VqaKfb3N9bwqjsM234N7F4Y/jlMVbY31uEL0bLnEAoK47f7yKUqG0SFEEKI1WqokHjixAnOnTu3qPbhq8JiJfb//PDDD0kmkySTSYaHh597vpCOY6oKdtfL21JrW/3bsGt2BhL9GOb6BzV9507MfJ7yEl/DuqY7YMs7Vs/i2HV48nVF5ynqTre1oMXtm13Q8qhi7y2EEEKsVsOExBMnTvDhhx+uuOZhPB5fc6+jw+HA7/fP356VT09gaDZcrsCarlNtmqKxI7iTmdIMD9PrH9RUjwe9o4PS4CBmscHm4Kmqtddz25uQemj1KlZwnqKq24jsfEN2aBFCCFE3GiIknjlzhiNHjvD+++8v+XwwGFxyJTRYPYmHDx+uZvMoTicwbXbsdndVr1MJbpubLt9Wnkw/IZFPrPv19Z4olA1K9++v+7UrItAFXX8HSjPWPMWZ5xdXrdb8gpatu2RBSwNLnztP6p/9Kelz52vdFCGEWJO6D4nnz58nEok8FxDPn3/6A/jdd98lFou98D2OHTtWtfYBFLNTKDYXuq0xCvm2edrwOwIMJQfXfX9nxelE27aV8v0HmLnK9cStK1dwdp6ix+pRTDyo6Nt7mjtlQUsDy5z/Eek//X+QOf+jWjdFCCHWpK5DYm9vL4lEYskexIVzDU+cOEFfX99z8w97e3uBpwW4q6WcS6PY3A21J2+3fzslo8SjzPoXdNa3bQNVoRSrTUmeitAdsOXb4O+C0a9h9DpUsLyPwxeiZe9hFEVl/NZVcsnJir23EEIIsRx1u3dzLBbj9OnTHD16lDNnzgBPg+Hcc3OOHj3K8ePH+eijjxbt3Xz27FkuXLhQ1ZXQZrGIUZhG8zWjqHWduRdx6E46fVsYTj0g7Izgta/fnp6KzYYejVK8ew+tvQ01FFq3a1eUqkLb69YevWPXoZCBjkOgVeaPBWu/5+ukRga5duE/p2nnG3iaO9nzW79fkfcXQgghXqZuQ+KxY8eIxWKLhpUXWhgGwdqmb25f50gkwsDAAKdPn160EroazHyecjmHVueLVpbS5m5jcmaSodQg+yL7UZX1C7laVxfl0VHyV65i29GDtn07iqKs2/UrKrTdKpHz+Arc/xvofNvauaUCFFXFFWrhwc//F+ID13D4w+z+3/1HDfUHiRBCiMZUtyFxYGBgxa950cKWairnZjBLOWyexusNUxSF7kA3Nydv8GT6CR3ejvW7tqpiP3yYUixG8V4/5ck49gP7UZzOdWtDRXmarHmKj76AB38D7QfBu/YtGud6De1uH9nJJxRnppm4+yXhntfQbPY1v78QcxKJBCdOnCAWixEMBrl69WqtmySEqDHpjlijwvQUpmlg9z1fSDuZLTAwmiaTq98Vqh6bh1Z3G48yD8mV1nchiaKq2HbswH74EGZ2mvwvf0l5rDY7wlSE3W3t++xussJifOV/6LzInt/6fQ7+Hz7grf/gjygVcozdukxh+vndf4RYrbl97g8ePFj1a81tryqEqG8SEteomJ7CxMTuDs4fMwyTwfEM3wwnGE3l+PL+FHdGUuSK67uSeLk6fZ3YNXvt9nYOh3G88w5qMEjhq68pXL+OscT2hw1B06HjMIR3wvhNGPmqooW37d4ALXsOo9udjN/pIzvZQNscioZw5MiRql/j3LlzVb+GEGLt6na4uVEUslOLdluZzpW48yTFTKHM9iYPnSE3o6kcDyammUjn6Qi66Iq40bX6yeeaorHd382d+G3Gs+M0u5vXvQ2K3Y79zTcpPXxIKRYj/3gE1e9D27IFra0NRW+gb1VFgeY94PDBk6+gOG0FR70yJZI0u4OmXW+ReHCHqaFbFGcy+Dt7GndOp9h0XlayTAhRPxroN299Kk0nUFQNze5leHKaB5NZXHaNN7eG8DitL2970EWL38mjeJaHU1meJGfY0eqj2V8/8+8CjgARV4Th9ANCzhC6WptvDX3LFrTOToyJCcoPH1K8dZvSnbtobW3o0W4Ul6sm7VoVfyfY3PDo8tMFLc7nd+xZDUVVCW3fi83lJfmwn2I2Qzi6v6HKMInN6eLFi7VughBimeqnO6tBlaYTqLrO/ZTJ0MQ0HUHXooA4R1MVtjZ5ONwdIeSxc2ckxdR0oUatXtpW3zYM0+BxprZ7ByuKgtbcjP2tt3D+6q+gbd9OeWKc/C9+SenR45q2bcVcIdj2q1ZZnAf/BtKVHR72tnYR2fkGhWyasVtXKGYzFX1/ISopFotx8uTJWjdDCLFM0pO4RqVsAsXmJFkw6epws7355aVP7LrK7nY/ZSPJrcdJXu8K4nXWR++PTbPR4e3gUeYRLe5WnHrtezoVpxNbTxR9axfFO3cp3riBMTaKbe/exlkJbXNZC1pGvoTHl6FpL0R2VOztnf4wLXsPEx+4xvidqwS37cEdXvvKaiHmxGIxLl68SDQaJR6Pk0gkiEajXL58eb4c2avOOX/+PJcuXSIcDnPlyhVOnDgx//5nz54lGo0Si8XmV1i/++67868DuHTpEpcuXarJ5y/EZiUhcY3K2SSm7qaEjm+ZYU9RFPZ0BLg2nODGwyRvbA3htGtVbunytHraGMuOMZx+wM7Qrlo3Z55is2E/sJ9yawvFmzfJ/+IX2HbvQetor3XTlkfVrELbE3dg4hYU0tD6hlWQuwJ0h4umPYdIDN1mavAmxWxa5imKikgkEpw5c4YLFy4sOn7+/Pn5uYXLOefUqVOcOnWK06dPE4vFnjsXIBqNcvXqVY4dO0Y8Huf8+fO8//779Pb2cubMGWKxGNFotEqfqRDiWTLcvAZGsUi5ME1Jc6FoOh7H8jO3pirs6wygqQrXHyYoliq3AnYtNEWjy9fFVG6KVKH+Vhhrzc04vv1t1KYmCtevU/jqq8bZA3puQUv7QUg/hoe/gFK+Ym+vqhrh6H4CW3aQGR1m8t5XlEv1NaVBNJ7e3l7C4edLfJ06dWr++HLOWYmDBw/S19c3v6Xq0aNHmZqakoAoxDqTkLgGxdw0lAuUVBeazY7DtrIvp11X2b8lQMkwufEoSdmo3N6/axF2RvDYPDxIPcCs4H7ElaLY7dhfew37G69jJJPkf/4LSsPDddnWJfk7oes7UJi2FrTkKhvGva1dNO16k2J2mrGblylkkhV9f7G5RKNRPv300yVrG84NGS/nnJWKx+OLajZWc3tVIcTSJCSuQS6bRjGKlHQfHqd9VUN7LrvO/s4A2XyJOyOpugg6iqKw1b+NbHGaiZmJWjfnhbTWVhzf+Q5aayvFW7cpXLmCkWmQhRsLF7QM/23FF7Q4fCFa9h2x6ine/ZLMWG0XI4nGdfDgQY4ePcqhQ4fo6enh9OnT8yuU53r6lnPOSkmvoRC1JyFxDYrZNBgFynbfioaan+Vz2djd4Wcyk2ckMVPBFq6ez+4j7ArzMPOQslmfRcDBmqto278P++HDkC9Q+OXnFPsHMCtYwLpq5ha0uJusBS2T/RV9+7l6imM3L/P1f/+n9P23P8Qw6vffUtSvCxcucOnSJY4ePUpvby8nTpygp6dnUc/hcs5ZidUMUwshKktC4hrks0nAoKB71xQSASJeB+1BF0Pj08wUSpVp4Bp1ebsoGUVGMiO1bsoraeEQ9u98G237NkqDgxQ+/xwjna51s15N1Z7u0DJxq+I7tCiqyv2//Z958PP/hcG//teM37pqTZMQYpnmFp4cPXqUc+fOMTAwwNTUFAcPHuR73/vess95GamdKER9kpC4BsVMAlPRMHQnHsfaVydvb/Jg01XuPknXxbCzQ3fS5mnjyfQIhXL9L4CY2wva8a23wYTC519QGhysi6/lSy1a0PKo4gta9vz93+fA8f+YPb/1j8A0Gb91hWx8tGLvLza23t5eent7Fx0LBoNcuHCBcDhMLBZb1jkLxePxRR/LDixC1CcJiWtQziYw0DBs7jX3JALomsquNh+pmSKPp+pj2Lnd04GqqDzKPKx1U5ZN9fuxv/MttG1bKfYPULh8GSObrXWzXq1KC1r2/Nbv89qJ/4T9v3OK5r2HcAabmRq8ydTQbRl+FsuyVLkaYNHCkuWcA9DT00MikahY24QQ1SMhcQ1KM0kMVGxOX8X2Yg647XSG3AxNTJPN137YWVd12r0dTMxMkC81SKkZZnsVd+7EcfgQFIoUfvFLSsPDtW7Wq1V5QYuq6YS79xHctoeZqVEZfhZLejbELbVyee6cuQUmyzkHrCHpWCw2/1xfX99zQTKRSDzX2yiEWH9STHsNjJk0ZVRcnpfvsrJS25o8xDN57j5J88bWYM0LIre4W3gyPcKjzCOiwZ6atmWl1FAI+7ffoXT3LsVbtzEmJ7Ht349iq49dbpZkc0HX34EnX83u0LIHIjsreglPUzt2j5/4wHXGb10huHU37khbRa8hGksikeDkyZP09vaSSCQ4ceIEp0+fJhwO89lnnxGLxYjFYvPhLZFIzPceLuecOQcPHuTChQucOXOGQ4cOAVY9RbCGnc+cObOoDdFodH5XFyHE+lLMup+wVV9SqRSBQICJqQlG/6f/jGxiHN9v/mN276xsuYZktsC14QTbmjx0RTwVfe/VeDI9wnB6mNeaXq+L7fpWozw2RvHGDdB07K8dQA2Fat2klzNNmLxr3Xyd0PaGtdClggyjTPLBXbKTT3CFWwlu3YWqNd7fjvm8NYfT4XDUuCUwcugIxpMnqG1ttF+9XOvmCCHq0Ep/Zs1lj2Qyid/vr2bTFpHh5lXKF2Ywc1nKmhOvu/KhaW7Y+cFkluk6GHZudregqzYeZRq33p7W0oLjnXdQXE7yV65SitX5ohZFgabd0H4IMk9g+OdQrOyQv6pqhLbvJdS9j1xiQopvCyGEmCchcZUKMxlryzPdjdfjrso1tjV5cNpUBkZrX8pFUzQ6PB3Ec5PMlOpjUc1qKC4X9sOH0bu7Kfb3U7jaV//b+vk7YOvfgVIOHvwNzCQqfgl3uJWWfUfQbHbG7/SRHhmq7wAthBCi6iQkrlJhJoNRLIDDjcdVneFXVVXobvaSnCkyka5cSZTVanY3o6s2HjdwbyJYO8rYdvRgP3wIczpD/vPPKdf7JHlnALb+Kugua0FLsvKrzXWHi6Zdb+Fr307q8SATd76kVKjzAF2HvKdO4vvD/wveUydr3RQhhFiTxpt8VCeKM9OYpTKa041ms1ftOmGvg5DHzuB4hrDHjqrWbhGLqqh0ejsZSg7S7unAbatOD+p60cJh1HfeoXDtOoWrfejRKHq0u+YLhV7I5rRK5Ix+A0++hHwKmvdaw9IVoqgq/o5uHP4QU7GbjN34gkDXTjxN7RW7xkbnO32q1k0QQoiKkJ7EVSrkpqFcwu7xoajV/TJGm73ki2Ue18GWfU2uJhy6o6HnJi6kOBzYDx1E7+6mNDBAse9LzEIdFw5XVWh/E5r3w1QMHn4O5WLFL+PwBmnZdwRXqJnE/dtM9l+jXKzjr4sQQoiKk5C4SoWZDIYCDmf1Vx67HTrtQRcPJqcplGq7J7GqqHR4OpjKxckWG6BA9TLMDz8fOoiRTpH/5S8xpqZq3ayXC0dhy7cgl7AKb+czFb+EqtsIbd9LuOcAhUySsRtfMJMYr/h1hBBC1CcJiauUSydB1XB6fOtyva0RDwpwf6L2hY8jriYcurOhdmFZDi0SmV397LJWPz94UOsmvZyn2Sq8rajWgpYKF96e4wo207L/bezeAPGB68QHb2KUKt97KYQQor5ISFylYjqNoqxfSLTpKlsjHkaTM0znalsS52lv4tSG6U2cozid2A8dQt/aRfH2HQrXrmGW63jrOrsHtv4KuJutwtvjt636ihWm2exEdrxGaPtecslJRm98wczUWMWvI4QQon5ISFyl8kweTdewu9evqGV70IXTphEbr/zQ4kpFXBFsmp0n0yO1bkrFKaqKbfdu7K+/hjE2TuHzL+p772dNh87D0LQX4vfg0RdVmacI4I600brvbWu3ltgNJgeuy1xFIYTYoCQkrpKRzWKz62iO9dsNRVUVulu8JLIFJjO1LYmjKirtnjYmc5Pky7Uvz1MNWlsb9m+9DYZB4ZefUx6r856zyA7Y8g7MTM3OU6xOfU3N7rB6Fbv3U8gkGb3xOdMTtf1jQVEUDKO283WFEGK5DMOo30oaC0hIXCUzX8Bms6Pa13eLuojXQchtZ3Asg2HUtthxs6sFVVEZna7OXLh6oHq92N/5Fmo4ROGrryn2D9R3kem5eYqqZgXFVPVWobvDLbTsP4IzECFx/zYTd7+kmKvNnFld1ymVar8zkRBCLEepVELX678KoYTE1SqX0J0uVN227pfubvYyUywznq5toWNN1Wh1tzI+M07J2Li/oBVdx/7mm9h27KAUi1H86ivMeg4kc/MUvW0w0gej16FKvWyabifcvY/Ijtcp5XOM3bxM6vEg5jr36qmqKj2JQoiGYRgGapXL51VC/bewTulmCZvThaZXr5D2i3icOk1eBw8mszXvTWx1t2GYBmPZ0Zq2Yz3o0W7sb72JMZWo/3mKqgYdB6HlNUjer8q+zws5AxFa9r+Nr3Ur6Sf3Gb35BbnU+u5i43Q6Sadrv4WlEEK8TDqdxulc31HI1ar/vs46pVNCd3pq0pMIVkmcvvtxxlI52oKumrQBwKbZaHI1M5odpdXThqZoNWvLetCam1G+9TbFr76m8MvPsb3+GlpTU62b9WKh7daWfo+vwP2/gvZD4KlOe1VVw98ZxRVpJXH/LpP3vsYVasG/pQd9HaZlOBwOAJLJJKqqout6Q/ylLoTY+AzDoFQqYRgGTqdz/udVvZOfoKvkUA1wuFCruCXfy8z1Jg7Ha9+b2O5pp2SUmJyZrGk71ovq8WD/1tuooSCFvi8pxQZr3aSXc4Vg26+DIwAPfwmT96pSJmeOzemhefdbhLbvJZ9OcPW//r9z5b/6v3Hr//tfVe2acxwOB4FAAI/Hg6Zt7D9YhBCNQ9M0PB4PgUCgYQIiSE/iqjlsCopmR9Vq9yWsl95Ep+4k5AwxMj1Cs6u5IVZsrZWi69jefBNlIEaxvx9jOoNt3z6Ueg0mut3aoWXyLkzchuwktL8FevV+WLkjbTiDTfzivzhDPjmJ3ReklJsGRcHm8rLnt36/atdWVRW7vTZ/wAkhxEYhPYmrZLMpoDurvm/zy3icOs2++pib2OZuJ1/KEc+t7zy0Wprfzm+unuLly5i52i4meilFgabdVpmcfAru/7UVFqtI1fT5P6QUReXun/9Lrl/8L7j1P1a/V1EIIcTaSEhcJYdNR7WvX43EF+mKeMiXyoymahtOvHYvfoefkQ1YXPtVtLY27EcOYxaK5D//HCORqHWTXs7TDNt+DexeGP5F1Yef9/z93+fA8f+Yff/g5PwcXqNUJB67QSk/U7XrCiGEWBsZbl4lh92GYnfXuhl4HFZv4vBklla/E1Wt3VBvu6eDO/HbpPJJ/I5AzdpRC6rfj+Nbb1P45hsKV65i27cPraO91s16MZvT6lFcOPzc9qZ1vMKeHVYuZNNgGOTTCUZvfI63uRNf+/aaLQITQgixNAmJq6TqGso67rbyMl0RD31DcUZTOdprODcx4AjgtrkZzY5uupAIoDgc2A8donTrNoXr19EzafSdO+t3jubc8LMrAk++tFY/t70J3taqXXJhYDSMMpknD8iMDjM9+QRf21Y8LVtQ1Tqd1ymEEJuMDDevkqkqqHUSEq3eRCfDdTA3scXdSiKfIF+q47l5VaSoKrb9+7Dt3kXp/oP6L7wNVkmcbb8OzpC17/PoNTDKVb+sqmr4O7ppPfAO7nArqceDjF77JZnR4XUvxi2EEOJ5EhJXyVQUNEfth5vnbI2462JuYsQVQVVUxmbqfJ/jKtO3bXtaePuLyxgzdT73TrfDlrdni28/sLb0y6XW5dKazU5w6y5a938LZyBC8mE/T67/gszYIwmLQghRQxISV8nU9LqaQ+WenZv4KJ6t6d7CmqLR7GpmPDtO2ax+b1Q905qasH/rbTAMCp9/Tjk+VesmvVpou7WoRVHgwd9APFbVRS0L6Q4Xoe17aN3/LRy+EMnhu4xe/6WERSGEqBEJiatk6o66CokAHSE3M8Uy8Uyhpu1ocbdSMkrEZzZPOZwXUT0e7G8fQfH6KPb1UXr0qNZNejWHD7b+KgS3w/gNqwB3cf16QnWnm3D3Plr2v43dGyA5fJcn135B+skDjHKdD90LIcQGIgtXVsm01Wbf5pfxu2z4XTYeTmWJ+GpX0d2pOwk4AoxmR2l2N9esHfVCsduxH3yL0u07FG/cxMxk0Hftqt8FLQCqCi37wdMKT76Cob+ElgMQ6Fq3JticHsLR/ZRy3aSfPCD1OEb6yX3GbnyBous4PIGqFuQWQojNTkLiamn2uutJBNgScnPzcZLUTBG/q3bta/W0cTd+h0whjdfuq1k76oWiqtj27UXxeijeuYs5ncX2+msoep3/J+hpgu2/DmM3rLCYeQKtr1d1p5Zn6U43oe178HVsZ3p0mIGfXaSQSeDwh9lx9F10Z/3MDRZCiI1EhptXS7XVbN/mlwl77bhsGo/i2Zq2I2AP4NCdjGZHa9qOeqNv3Yr94FsYiQZZ0AKg2aD9Teg4AjNTVq9iav2HzXW7k0DXzvn/7kzDYPTG50z2XyOfSax7e4QQYqOr826M+mWqWl32JCqKQmfYTf9omplCCZe9Nv/EiqLQ6m5hOD3M1nIRm1Z/X6ta0SIRlG+9TfHLryh8/jn2N95ADYVq3axX87WBKwxj12GkzwqKra9XpQD3y+z9rX9EcSaD7nAT3LabzJNhJu58ic3tw9uyBVeopabbZQohxEahmBVaCjs0NERvby8DAwN89NFHACSTSa5evcp3v/vdSlyiLqRSKQKBAPcu/b/ZcfTfr3VzllQ2TC7HJmnyOdjRWruh3pJR4quxL2n3dtDp7axZO+qVWShQ+PobzGQSfe9e9M6OWjdp+dJPYOwaGCVo3geBrdaK6BowTZN8Kk5mdJh8egpVs+Fp7sDd3IFuX98AK4QQ1TCXPZLJJH6/f92uW5E/tz/44AOi0Sjvv/8+58+fnz8eCAQIhUL8yZ/8SSUuU1cUW+12NnkVTVXoCLoYS+YolmpXOkRXdSKuJsayYximlDB5lmK3Yz90EK29neKNGxTv3q1p+aIV8bXB9r8Lvg4Y/cbaA7owXZOmKIqCMxChadebtOx/G1e4hczYQ0av/YLJ/mvkUvHG+boKIUQdWfNY5I9+9CNisRgDAwN0d3fzk5/8ZNHzb731FtFolB//+Md8//vfX+vl6oZaB/s2v0x70MXDeJaRxAxbm2q3M0yru5Xx7BhTuSkirkjN2lGv5nZoUbweinfvYU5PY3utARa0gDVXse0N8HXC6NfWXMXwTgjvsFZH14DN6SG4dRf+zigz8VGmxx8xee9rHvX9Faqm4w63su+3T9akbUII0WjW/JtoYGCATz/9dP7jpcp6zPUobiSKsz625HsRm67SEnDyODHDlrAbVa3NUKDb5sZn9zGaHZWQ+BL6tm0oHg/Fb65R+OIytrfeRHXVb2/1Ip4mq1dx8h7E70HqIbS+Bp7alT9SNR1Pcyee5k4KmSSf/4v/K/lUHLs3SMv+t/E0teMMNMncRSGEeIk1/4Ts6elZ9PGLhnUGBwfXeqm6otjquycRoDPkolg2GKvxVn0t7hYyhTS5Tbqf83LN7dBilssUPv8cY6oBdmiZo2rQvMfarUV3WgW4H/dBKV/rlmH3BuYXmSmajlkuEY/dYOSbvyXx4C6FTLLGLRRCiPq05pC43ILAAwMDa71UXdFs9b9a12XXafI6eDRV2636Qs4wmqoxMTNeszY0CtXjwfGtt1E8XgpXG2SHloUcPtj6HWh7E7LjMPgziA9AnWyrp2oazXsO0bLvbTxNHcwkxhm/08eT678k9XiQUr4BShIJIcQ6WfNw89TUFD/96U/53d/9XWDp0PiDH/yAQ4cOrfVSdaUey98spSPk4pvhBIlskZCnNnUdVUUl4mxifGaCTu+W+t5ppA7MLWhpqB1anhXoAm8bTNyG8VuQfADN+8Hb8vScn/8p5FPg8MN3/rCqzdnz93+f4kwGm8sLgM3lIbClB39nlEI6QTb+hMzoA9IjQ9jcPtzhVlyhFjR77XYuEkKIWqtICZzDhw+jKArvvfceAwMDnDhxgkQiweXLlzl//jxHjx7lk08+qUR7a25uGfrIYD9t23te/YI60DcUx2nT2NcZqFkbMoUMNydvsCu8m6AjWLN2NJrSgwcU79xFi0SwvXYApQF6sJ+TS1l7QGcnwNNibe9n98A/2wLpR9bCl3/8sNatxDDK5BKTzMRHyaUmwTSxewK4wq24gk0SGIUQNVOrEjgVq5P48ccf88EHH1hvqiiYpkkwGOTs2bOcPLlxVhPO/UONDQ/RvGVbrZuzLCOJGQZG0xyJRnDYtJq149r4N7hsLnYEd9asDY2oPDlJ8ZtvUOwOa0GLu/7nwy4pPQLjN6E4A6Fu+O9+AzKP6yYkLmSUiuSSk2Tjo+RTcQDsngDOYBOuUDO6o0EWFQkhNoSGD4lzBgcHicViRKNRuru7K/nWdWHuH2riyWMire21bs6ylMoGX8Qm6Qy52VbDcjgj0yM8TA/zVstBdLUBSrzUESObpfjlV5iFPLbXX0eLNOhKcaMMUzGI98NP/0OYmQRvB/xR/c69nAuMM1Pj8z2MNpcHZ7AZZ7AJu1v2JhdCVFetQmLFf1N3d3dvyHD4LLWBhv10TaXF7+RJYoauGpbDaXI28TA9zOTMJK2e1pq0oVGpbjf2t49QvHadQt+X2HbtRN/WGD3Zi6gaRHZCYBsos73a5Twkh8G/pWa7tryMqttwR9pwR9owyiXyqTgzU+Nkxh6SHhlCtdlxBppwBiI4/CFUtXa99UIIUUnr1p3zgx/8gH/+z//5el2u6jS9NotAVqs96GIkMcNkJk+zvzZbldk0G0FHkPGZcQmJq6DYbNjeehOlv5/inbsY6TS2ffsas9afbodfOWPNUyzNwJOvrFXQTbutBS91GBbBqr/oCrXgCrVgGgaFTJJccoKZ5CTZiccoqordG8QZiOD0h9GdDTo1QAghWMFw81dffbXqi8TjcU6cOMHk5OSq36Ne1KrLtxK+eWDV3Xt9a+0Km0/l4tybuseBptdwN0CtyXpVHhmheOMmis+L/Y03UJwNvkfxTMJaCZ0dt1Y7N+0Bb2P9IVHMTZNLTJJPTZLPJME00RwunP4wDn8Yhy+Iqsk0CyHEytX9nMRwOEwymVyy3t5caY6XPQdQLpdX28660cghcTyV4/ZIioPbw3gctfllZZgGX419SZOria3+BhwurSNGKkVh9o83++uvowaDNW1PRWTjVlicmQRnECK7F5fNaRBGuUQ+nSCfipNLTlIuWIXk7Z4ADn8Ihz+E3e1vzF5gIcS6q/uQuGPHDi5cuEA0GiUQeFpKJZlMcubMGU6fPk00Gn3udZcvX+bixYt8/PHHDReqltLIIdEwTM7+TzdRFeiKePj3vrO9Ju14kLrPxMwEb7a8harIL8m1MPN5Ct98g5lMoe/dg97ZWesmVcb0hBUWc1NWWAzvtHoW63QY+lVKuSz59BS51BT59BRmuWQNTXsCOHxBHL4QNrdPQqMQYkl1v3Dl+PHjvPXWW88dv3DhwksD4NGjRzly5Aiffvop3//+91ffUrFmqqrgc+qMpXIYRu12YGlyNfNk+gnJfIKQM1yzdmwEisOB/dAhSndmC28nU+h7djd+2PA0gedXYHrc2hP68WWw+6xFL76OhguLutON7nTjae7ENE2KMxnyqTj5dIL0kwekHg8uCo12bwCbxy+LYIQQNbXskPjDH/5wyeOmab4y1QYCgZpuCyee+v7f28GV2CQ9rbUr2+G2uXHbPEzMTEhIrABFVbHt3Yvi91O6dRsjk8b++uuNP08RwNNs3bJxiN+DkT6rhzG8w1oN3YAhSlEU7G4fdrcPX9s2KzRm09bwdHqK9JMHmEYZ5s7zBqybJ4Bma6wFc0KIxrbmiWnL3SqsobYU28CcNo2w18FIYob2YO0KAje7mnmQvk+xXMSmNU45oXqmd3aier0Uvv6a/OefW/MUQ7VbpFRR7jC4vwW5pNWzOPqNFRaD262b3ri7oSiKgt3jx+7x42vbimmalGamyWeSFDIJsvFRMqPDAGh2p3WuN4Dd48fm8jZ+r7EQom6tOST29/cv67yBgYG1XkpUSHvQxfWHCZLZAgF3bXomIq4IDzPD3Jm6za7Qbuya9JBUghoI4PjWtyh88w35K1ex7d6FvnVrrZtVOc4AdB6GwjRMDVplc+L9Vq9iKAqOxi9srSgKNrcXm9sLLdYc01IhRyGTojidpDCdIvlwHExzvrfRNtfj6PahOVzyR7kQoiLWvOPKZ599xsWLF19aA/HDDz/kyJEj/O7v/u5aLlUXGnnhyhzTNLk6GMfnsrG7vXafQ7aY5e7UHUxgV2gXHlvtdoPZaEzDoHTvHqX7D9Da26zhaH0Dll8pFyFxHxKDUMqBuwmC3Q29yGU5TMOgkE1RnE5TmE5RyKYp52cAUDQdu9uLbTY82txedIdbgqMQDazuF668yPe+9z3+4i/+gkgkMr9IJRgMkkgkGBgY4NNPP+Xdd9/dEAFxo1AUhZaAk4eTWUotBrpWm+Eqt83N/sgB7k7d4dbkTXqCPTJHsUIUVcW2ezdqIEDx5i0K6cvY3ngd1bPBgrhmg8gOqxcxM2L1Lj6+DLoLgtsgsLWhh6JfRFFVHN4gDm9w/phRKlLIpilOpyhkM9auMLPD1IqqYXNZvZM2lweby4vu9srCGCHES1Vs7+be3l4++OAD+vr65o9Fo1HOnj3L7/3e71XiEnVhI/QkAuSLZS7HJtnR6qOthnMTAcpmmVgixlQuzlb/Vto8jbEndqMwpqcpfvU1Zj6Pbd9etLa2WjepunJJSAxB6hFggrfdCovuyIbuXVyKUSpSzGasXsdshuLMNKXc9PzzusON7vLMBkcrPMpwtRD1p+7rJK7E4ODght2/eaOERIDrDxOUDZM3argDyxzTNHmYGWYkM0LIGabJFcHvCKAp0tNRCWapRPHmTcpPRtG3bUXfuXPjL3goF609oZP3oZABm9sKi/4usG2Ald+rZBhlSjPTFGcyFLMZSrlpitlpjHIRsHopdYdVssfm8lgh0umR8ChEDW2okLiRbaSQOLcDy6HtYdw12oHlWRMzE4xkHjNTmkFRFPx2PyFnmKAjKItbKqD04AGlu/dQfF5sr7+O6qptL/K6ycYh+QDSj8E0wNMCgS7wtMJGD8vLVC4WKM5kZgPkNKVclmJuGrNcsk5QlKfhcbbuo+5woTvdqLpUKBCimjZ8SHzvvff45JNPVvSaRCLBRx99RCKRIBaLEY/H+fDDDzl+/Piazl2LjRQSDcPk84EJ2oIuupu9tW7OIrlSjkR+iqncFJliBtM0afO0scXXJbu0rJGRTFL45hqUitj27UNrbaw9ktekXLSCYvIB5BKg2qzi3P4tVpkd8ZxyIT8fGEu57Pxjo1iYP0fVbOhOF5rDNR8c9dnHEiCFWLuGXbjys5/97JXnxGIxent7V/S+iUSCM2fOcPbsWYKze9L29fVx6NAhjh8/zoULF1Z1rnhKVRVa/E7Gkjm2N3nqaijJqTtp09tp87RTMkqMZcd4lHlIqpCiJ7gDl75JesCqQA0EcLzzLYo3b1L4+hv0rV3ou3Zt/OFnsBa6BLdZt3wGUg+tW/I+2Dzg7wRfJzjq64+mWtLsDjS7A4d/8bQUo1yilJ+hlJuhnM9SzGUp52fIJ+PzQ9dgrbbW7U50pxvN7kR3ONHsTjSHE93u2hzfd0I0qDX3JIbDYRKJxIsvoCjzu60YhrHs9z1z5gwffvjhfOib8/HHH3PmzBkuXbrE0aNHV3zuWm2knkSATK7Il/en2NcZIOKt71Wg2WKW/sQ98uU8W33baPVsoh6wKikND1O6cxfF67GGn93uWjdp/ZkmzMSt+YuZETBK4PBbPYy+DrBvsBXh68AoFSkVclaALMw8vc/nKBdyi85VdfvT4Gh3PH+vyzQTIRp2uHnHjh1cunRpyYUqyWSSS5cuMTU1xcmTJ1f0vj09PcDzRbjneghPnTrFuXPnVnzuWm20kAjQNxTHZdPY2xmodVNeqWyWGU4NM5YdJeQMsT3QjU2V4ay1MFIpit9cw8zn0ffsRu/srHWTascwYHrMGpLOPAGzDP1/DopqrZL+tQ9r3cKGZ5om5WKecv5paJy7lQp5K0Qu/LWkKFZgtFk9mk/v7Wg2B+pskJQeSbGRNexw8+nTp1+4kjkQCHD8+HGSySQ//vGP+f73v7/s941Go1y5cuW543O9hfF4fFXniue1BZwMjk9TLBnY9Pr+QaspGtsD2wk4/AwmB7k2/g3b/NuIuJpq3bSGpfr92N/5FqU7dyneuIkxOWkV37ZtwvCtquBrs25GGTKj8K//ALLj4ApD17fB124FRufG+CNxvSmKYg0/251LbpBjmiZGqWgFx2Ke8mxwLBfylPM5Cpkk5WJ+cZDEGtbWbA40mw3VZoVIVbdb9zY7mm7dq7qtrqbWCFHP1hwS/8k/+SevPCcQCLDSDstLly4teXyuDuORI0dWda54XrPPColjqRyd4cYYbgw5w3htPu6nhxhIDDA+M0G3fzsOffOWNlkLRdex7d+HGglTvHWL/C9+if21Axtn7+fVUDXwd8DcqnrVbm0LODUIk3etkjqeVmt3F1dEVklXiKIos72ELx9mtoJk3gqSxQLG3H2paIXJ6RRGsYBplJ+/hqZboVG3WTebHW3usW6bfX72sa6jahIsxea0bnVPkslkRd7n3LlzBINBTp06VdFzXyafz5PP5wGry3ejsekqYY+d0QYKiQA2zcaO4E6mnFPcTw1xbeIaHd5O2jxtsgJ6lbS2NtRAgML1G+SvXEXv3o4ejcpQHlihsf0ta0h6ZtIajs48sbYEVHVwN4O3DTzNG3KXl3ozF+hsvHyRkWGUMYoFjGKBcskKkUapOPux9bg4nSJfLmGUikuGSrCCparp1nVn7xVVm22HjqLqqLp1jqJpVrDUNFRVQ9F0CZmiIa1LSEylUnzxxRdrfp/e3l56e3u5cOHCc4tU1nLuq3z00Uf88R//8Zreo961BpzceJQkkyvidTbWMGPIGcJv9/Mo85BHmYfEc5N0B6KyF/QqKS4X9sOHKA8OUYrFMMbHsR04gOpbYmxwM1JVKwh6mqH1NcilYHrUGpp+8qV1jiNg1WL0tIArtOl2eqknqqqhOlzgWF5FBNMwngbJ8tx9afbY0yBpfZy1Qmip+LSe5AsoqjYbGudCpGYdmw2Riz5W555XZz9WF5yroSjq02PyvSWqaM0LV141lDtXt/DSpUt897vfXculCIVCnD17dlk9gys591We7Uns6uraUAtXwJoH9EVskojXwY7Wxg0D08VpBpMxZkozdHg7afe0S6/iGhjpNMVr1zGzWfSeKNr27Zvvl9I/2wLpR1ZpnH/88OXnlvIwPW4tfsmOQ7lg1WJ0N1k3T7Oslt7ADKOMWSpZobJcxCyXMcplTKM0+7iEWS7NHrNuRrmMWS5Z4dSw7l8VOBdRlNnQOBsmFRVFVeY/RlHnQyWKMh8s5861Xm8dn3vM3PkoKKoye6+ioICqWPez10VRQGH+GIqy4PESx+fbvLF/jpimOT9v1sScO7j4mMkzHz89b/752Y+TqSTN7Vsab+HK1NQUBw8efGFYDAaDvPvuuwQCa1s5e+LECT788MNlhb6VnLscDocDh2NjDx8pilUz8UlihmizF1VtzP+APTYP+yL7eZx5zOPMI6ZyU0QDUdy2xhlGryeqz2ctahkYoHivn/LEBLb9+zdnqZzl0B0Q2GLdTNPaR3ouMI5dB0zQXVZYdDdZcxk38RaBG42qamDX0Fj77wvTMGaDpHVvGOVFH1v3xoKPy1a4MMqYc/dzz5sGZqmEaSx43jSs800DDOt+7ti6mw2MT4Ojsui4deT5Y1Ti19Szn+6Cz9989sn552aD3YJQV23pzPSrT6qCNfckHj58eMmVxZV05swZIpEI77//fkXPXY2NWAJnznS+RN9QnH0dASK+xg/FmUKGwWSMXDlHp3cL7Z72Df/XazUZU1MUb9y0SuXs6EHbunVzfD1//qeQT1m1E7/zh6t/n3LJmsuYnbB6Gwtp67jNM9vTGJHQKOrCcyHS5LmP53u7DGO+F2xhaLIezx1nwePZ47PnPb0zFwe0ZwLZs0GskpvFLflzbNGxZ0PswueVpx8+85r58xec+/Shsvi8uaBsdb8+87FCKp2uSU/imkNiMplccy/hy5w/f55EIvFc6Dt//vxzPYUrOXe1NnJIBKtmotuusaej/msmLodhGjzKPGQkM4LX7qMn2INDa/wAXCtmqUSpv5/Sg2HUYADbvn2oXtmdZFVKechOzgbHycWh0RV+epPdX4TY9Bq2TmI1A2Jvb++SoQ94bpeXlZwrXqzZ52B4MkupbKBrjT+XT1VUunxbCTpCDCT6uT5xjWggSsgp+/SuhqLr2PbsQWttpXjzFoVffo4ejaJt3yYroFdKd1gldvwd1sfzoTFu3VIPAdMqweMKgzM0ex+wVloLIUSVrbkncWhoiO3btz93/Cc/+QmxWAxFUTh48OCKF63EYjGOHTvG0aNH51cnz4W9WCzG6dOnOX78+IrPXauN3pOYK5a5HJtkd5uflsDGGvYqGSUGkzGmclO0elrp8m2VRS1rYBoGpViM0uAQqtdj9SpW8Y/GTadchFziaXDMJawdYFCsoW9XyAqOzqC1GGYzDP0LsUk17LZ87733Hp988slLz/nyyy+5evXqinZc6enpIRaLvfD5q1evcvDgwRWfu1YbPSQCfP1gCk1VOLAlWOumVMWT6ScMpx/g0t3sCO7AKQW418RIpSjeuIGRzqB3bUHfsWNz7tZSbaYJ+TTkpmZ7G6egODuZXbVZPYxzodEZANvySr4IIepfww43LydjvvXWW1y9enVF7/vsPsyVOle8WovfycBomkLJwF7n2/StRpunDZ/dR3+inxuT19nu3y7b+q2Bta3fO5QfPKA0EKM8NoZt1y609vZaN21jURRrK0CnH4LbrGNzvY1zt9QwxO9Zz2l2q16jM/D03uaWHkchxLKtOSQuZ3VjKpVacU+iqJ0mr4OB0TQT6TwdoY3ZG+GxeTgQOcBQapCBxACpQppt/m0y/LxKiqKgb9tmzVW8c4fCteuojx5j27dXyuVUk2Z7Wth7TjEH+aRVfieftOY2lvqt51Qd7D5ruNoZsO7tPtDWbfOtjalSK+CFqDMr+snw5Zdfcu7cufm5hgBXrlzhN3/zN5c8Px6PzxfTvnDhwtpbK9aFTVcJeeyMp3MbNiQCaKpGT3AHPruf+6khpovTMvy8RorTif2NNyhPTFC8dZvCz3+Btm0renc3ii5BZF3YnNbN2/r0WClvhZh8ygqPM3FIPmC+vIjN/TQ8OnzWY7vX2l1GAtCr/eJPnxZcl6+R2EBWPSfx4sWLnDp1ikgkwve+970XntfT08PRo0d56623Vt3IerIZ5iQCjKVy3BlJcaQ7gtO+8VdSThen6U/0UzKKdAeihGX185qZ5TKlofuUBwfBpqPv2InWIbUq64ZhWGV38umnAbKQhlJu9gTFCo8X/6FVDNzTBv/JLSs8Ss/jYivZlUeIVWi4OYnHjx/n4MGDfPDBB/yLf/EvKtkmUQfCHjuaojCWzrE1svG3EfPYPOyP7GcwOUj/1D3aPG1s8XXJ8PMaKJqGrSeK3tFO8d49ijduUH44jG33btQ17qcuKkBVZxe7PLMivVy0gmMh87R2I4BRgAd/Yz3WnVZYtHmse7vHemxzW+8rhNgQ1vTnYDQa5dixY5Vqi6gjuqYS9joYT+U3RUgE0FWdnaGd86ufM7PDz3bNXuumNTTF5cL++uuUt0xRunOb/BeX0drb0HfsQHVt3OkMDUuzgTts3cBaAAOgOWHbr1nhMZ+2VlbnpqzFMqYx+2LFWlVt8ywOjja39bHUdxSioax5zODkyZPLOu+rr77izTffXOvlxDpq9jsYf5RjOlfC49w8w0ttnjY8Ns/86ucdwZ347L5aN6vhaeEQ6jvvUH78mNLAAIW//Tla1xZrvqJdgnjdU5Slex5N0xqiLkxbwXHufiYOyeHZ2o6zNMdsiHQvvuku67iESCHqyrr95j9z5gx//ud/vl6XExUQctuxaSpj6Rzdzs21NZjP7uNA5AD9iXvcjt+iy9dFm0dKuqyVoijonZ1obW2U7z+gNDRE+fFj9O3brb2gNQkJDUeZ6z10AUuUkirlZ4NjdvEtl4DiDPOLZ8DqtZwLjPP3zqePNYcMZwuxjlYUEn/wgx8QCoX4p//0n84fO3LkyCtfN7fCWTQWVVWIzA45b2/ybLoFBzbNxu7wHh6mh3mQesB0cZrtgW40RYLMWimahh7tRtvSSSkWo9g/QGl4GL27G62zU7b420h0h3VjicVghmH1QpZmrMBYmpkNkTOQnbA+NkqLX6M5ZoPjwtvssbnnNLuESSEqYEUh8S/+4i9QVXVRSJyamuLo0aP09PS88HWmaXL27NnVt1LUTIvfwZPkDKmZIgH35hsSVBWVrf5teGxeBpMxZkoz7AzuxCFlcipCsdutvaC3bqU0MEDx1m1Kg0PoPVG09nYJixudqoLdbd1epFyywmIpNxskc09vuYR1X84//zrNPhsaHYvvNfvsY/vTc2S1thBLWtF/GUvtbBKNRpe1uvnKlSsruZSoE36XDbuuMpHOb8qQOCfiiuDSXdxL3OXG5A12BHfgd8g+xZWiut3YX3sNo7vb6lm8cZPy4BB6NIra3rbperHFApoOms+q3/gipmkNa5fz1v1ccJw/NltgvJQHo/j86xV1QWi0W9sczn9se+aYbfaxbfPMoZRamU9tsq/FmvduXq5kMkkg0Pi/VDdLncSFBkbTTGTyvB2NbPpf1iWjRH+in3QhJfMUq8hIpykNDFAeG0d1u9G6t0vPYq1stF+KhmEFxVIeyoXZW94q/TP/8YKbUXx+yHuOolqB8Sf/PsxMgLsZ/oP/1drZRrXN3uuzgVJf+qZo9T80LnUgn6rR16Lh6iSu1EYIiJtVk8/B48QM6VwJv8tW6+bUlK7q7A7tZnh+nmKW7YHtMk+xwlSfD/ubb2KkUpQGByneuElpYMDa+m/LFlngsp42QjBcSFVBnZsnuUxzwXIuSBql2Y9nH7/1j6wgrbusc4pZ694oWbeFK7yXoqgLAqNu9VDOf6w9c6/Pnv/M8flj6oLj6oJbA4RRUXfWLSS+9957fPLJJ+t1OVFBfpcNu6Yymc5v+pAI1grdrf6t/M2jv+avH/4VXruP/+OBP5B6ilWg+v3Y33gDY3qa0uAQxbv3KA0Oom/dhralU0rniPXxqmD5Gx+//PWm+TQwzoVHs/z0mFECo/zM8dn7csF6bJafv5+vT7kCi4Lj7A1lduhceea52Y/LBeu15QI8+Xr2POXpaxVlZffwisc8/XjRsZcdX4PnBlTNFzxnWl97sO5Tj5Y4z3z+45c+Xnj+M48XnptaUNh+HVUsJH711VfE4/Eln0skEvT19VXqUmKdKYpCxOdgPJ2nu2VzlcJ5md77f8FkbhKfzce32t9hZ3AnXrt8fapB9XiwH9iP0ROlPDREKRajNDiI1t6Otn0bqvslCx+EqDVFmZ3baJstFVQhpmkFxYWh0Zi9f+5jY/Fx03j6+mdvmIvPWXi9fHrBOeYz98YSx14QhBqVUXh6P7KcXLNE+H1VWF7q2Ey2Ao1fuTWHxMHBQQ4dOkQikXjpeZt9Lluja/I6GEnMkJ4p4pPeRAD+wY7fIVvM4tAcODQ7t+I36Q5EaXItUStOVITqcqHu3Yve00P54SNKww8oPXyI1tyEtm07WjhU6yYKsX4U5enQcjXNjZLoDtj2K2t/P9PkuR66ZwPkC3vmWHx8qfdeykszyBLPLdV7qSjWzkNg3e/4t5/v+VSeDYUVkkpV9v2Wac0h8cyZM/zoRz/i6NGjL513+Bu/8RtrvZSooYDbhk1TmcjkJSTO+u0dvzP/2DANhpKDxBIDzJSybPF2yR9GVaTY7Vadxe3bKI88oXz/PoUrV1B9XrSuLrS2NhRdypoIUZcUpfIhar0sDIGboHTSmj/DI0eO8Hu/93uvPE/2eG5siqIQ8dqZTOfpbpYh1Wepiko02IPb5mY4PcxMaYaewA60zVIio0YUVUXv7EDv7KA8OUl5eNiqtXj3HlpHO1pXF6pnc+w9LoQQlbbmkBgKLW9455/8k3+y1kuJGot4HTxJbr69nFeizdOOS3fTn7jHzckb7AztwimFt9eFFomgRSKYMzOUHj2i/PAhpQfDqOEwemcHamurlNARQogVWPNPTNM0SS1jrPynP/3pWi8laizotqOrChOZJXY3EPMCjgD7IvsxMLkxeZ1UoTZzSTYrxeXCtmMHjl/7NewHDoBpULh2nfxf/RXFO3cwMplaN1EIIRrCmkPiyZMn+fTTT/nqq69eep6Uv2l8qqoQ9jqYSOdq3ZS659Jd7I/sx2PzcCd+m/HsWK2btOkoqorW0Y7jyBEcf+c7aJ2dlEdGyP/8F+S/+ILS8DBmcYndN4QQQgAVGG7+zd/8TcBawJJIJIhGowSDwUXnJBIJYrHYWi8l6kCTz8FYKsd0voTHIUPOL6OrOrtCu3mQesBgcpBsKctW3zZZ0FIDqseDumsX+o4dGOPjlB89pnj7DqU7d1GbImgdHahNTTIcLYQQC6z5t/zly5c5evQoJ0+eJBKJLHnOxMQEP/7xj9d6KVEHQm47mqowmc5LSFwGVVHZHtiOS3fxIH2fmVKOHcEd6Kp87WpBUVW01la01lbMfJ7ykyeUR0YofPU1ik1Ha21DbWtDDQUlzAsx59t/+HRrxs1uk30t1vybKhqN8umnn77yvMHBwbVeStQBVVUIexxMZPJsbZJVo8vV6mnFpTu5N7ugZVdotyxoqTHF4UDftg192zaMTIbyyAjlkSeUHj5EcTrmw6T6zMiIEJvORtuacS022ddCMc0XVZ5cnqGhIbZv3/7K85LJ5IbYv7lWm2zXk4l0nluPkxzuDuOyS4/YSuRKOe5M3aFkFNkZ3Inf0fj/TWw0RiJBeXSU8ugoZi6P4nSitbagtbSgBKWHUQix/mqVPdY8AefUqVPLGkoOBAL89Kc/5fDhw/zmb/7mKxe6iPoV8tjRFIWJtKxyXimn7ny6oGXqjixoqUNqMIht924cv/qr2A8fRmtuovzkCfnLV8j/9V9TvHmL8uQkprGKfXOFEKKBrLkb6MSJE3z/+99/5Xlffvkl3//+9/nRj35EIBDgzJkzXLhwYdP2xjUyTVUIeexMZgp0RWTIeaWeXdAyU5qhy7dVeqjqjKIoaOEQWjiEvmcPZiJBeXyc8uiYNSRt01EjEbTmZmvRi012IhJCbCwVGyv86U9/yieffEI0GuUf/sN/yBtvvLHo+Q8++ICzZ8/O787S09PD+fPn+aM/+qNKNUGso4jPwZ2RFPliGYdNdhVZqecXtMzQIwta6paiKCihEGoohG3XLoxUivL4BMb4GIVr10FRUINBtOYm1KYmVK/sSiSEaHxrnpMIEA6HSSQSHDx4kEQiweDgIOfPn+cP/uAP5s+JRCJcvXp10fzF//Q//U8bbicWmZNoKZUNftk/QU+rj/agq9bNaWjJfJL+xD1sqk0WtDQgM5ejPDGBMT6OEZ/CLJeteYyRCGpTBDUcll5GIcSaNOycxB//+MecPXsWwzC4cuUK/f39TE5O8ud//ueLdmJJJBLPLXB5tp6iaBy6phJw25mU3VfWbG6HFhNkh5YGpDidzPzP/wu5n/3/KNy9i/3QQbTWFoxEgsLX35D7y78if/kyxYEYxtSUzGUUQjSMNY9t9ff388Mf/nDRsWAwyNmzZ+nt7eV3f/d3SSaTS75W5mA1trDHztDENKWyga5JEeK1mNuhpT9xjzvx22z3b6fZ3VLrZollypz/EcaTJ6htbfh/8H9Ci0RgNxgzMxiTkxiTk5Qf3Kc0MICi66ihIGo4bPUyer3ys1AIUZfWHBJfVEBbURReNZKdSCTWenlRQ2Gvndh4hmS2SMTnqHVzGp4saNl4VJcLdcsW2LIF0zQxUymMyTjleJzSvX5Mw7AWwMzOd5TQKISoJ2sOicFgkD/5kz9ZtABlaGiI06dPc+HCBQBisdiSgXFiYmKtlxc15LLruO06k5m8hMQKeXZBS7Y0Izu0bBCKoqAEAqiBAHq0G9MwMBJJjKkpjGdDYzA4f1MCAdkuUAhRE2v+zXPy5EneffddNE0jGo0Sj8dJJBKcPHmS3t5eJicnOXv2LGfPnuXHP/7xfLmcn/70pxw5cmTNn4CorbDXzmgyh2ma0vtRQa2eVpy6k37ZoaUheE+dxMxkUFawqllR1fkSO/RE50OjmUhgTE1RGhzCLJVQVBUl4EedDZhqIIDilO8FIUT1VWR1M0Bvby+9vb309PRw9OhRuru7SSaTXLlyhcOHDxMIBPjggw8YHBwkHA4Ti8X48z//80pcel3J6ubFUjNFvn4wxetdQQJue62bs+HkSjnuTt2haBTZEdxJQHZo2TRM08RMpzESidlbEjOXA6zFMmpwNjT6/Sh+P4ompaiE2KhqlT0qFhKX8tlnnxGLxThy5AhvvvkmAD/5yU+Ix+OcPHmyWpetKgmJi5mmyecDk7QGnHQ3S224aigZJQYS/aQKKbp8W2nztNW6SaJGzFwOI5mcv5mpNGa5bNVp9HpQ/QGr19HnQ/H5ZJhaiA1iQ4bEOV9++SVXr15d1s4s9U5C4vPujqRI54oc6l56EZNYO9M0GU4/4Mn0E5rdzWzzb0dVJABsdqZpYmYyGKkUZiKJkUpiZKbBNBcHR78P1eu1gqMu81uFaDS1yh7r8tMiGo3ywQcfbIiQKJ4X8ToYTeWYKZRw2eUXUDUoisJW/zZcupuh1CAzpRw7gzuxaVKkeTNTFAXF50P1+aCzEwCzXLaCYzqNmUxhpFOYIyPz9RlVtxvF50Xx+lB9s8HR6ZQ5xUKI51TkN/rPfvYzzpw5QywWe+65uTI3Z8+ercSlRB0KeuyoikI8U6AzLCGxmprdzbh0J/cS/VyfvM6u0C48Ntk/WzylaNr8Kmq2WMdMw8CcnsZIpTEzaYx0BuPBfUrFkvUaXUfxelA9XhSvB8XrRfV4ZIGMEJvcmoebv/zyS773ve9x6tQpenp6uHr1KocOHSIcDhOPx7l69SrHjh2b37O50clw89KuP0xgGCavbw3VuimbQqFc4N7UXbKlLNFAlIirqdZNEg3IzOUwMhmr53H23sxMz/c6KjYdxbMgPLrdKB4PisslPY9CrKOGHW4+f/48g4ODBALWqsvDhw8TCoXmt+A7efIkg4OD/OxnP+O73/3uWi8n6lTE62BgNC27r6wTu2Znb2Qfg8lBBhIDTBezdPm65Be3WBHF6URzOqHp6R8ZpmlizszMB0YjO20NWT95Yi2SwSrfYwVGt3XvdlvD2G43ikNqpgqxUaw5JB48eHA+IIJVXLu3t3fR/MPu7m4+++wzCYkbWNhjpx+ITxdo8csQ1XpQFZWeYA8em5vh9DDZUpaeYA82VeYpitVTFAXF7Qa3GxbsDGmaJuRyGNks5vS0NXw9ncUYeTJfmgdmh67dLqu30e1efC9zH4VoKGsOic/+B9/d3b1hVjKL5XPYNHxOG/FMXkLiOmvztOPW3fQn+rk5eYOdwV24be5aN0tsMIqigMuF5nLBM9uxmuWy1fuYzVq36ay1b/WTUStAzs1qUhQrKLpcqC4XOJ3z94rLieJwbKqyPelz5+eLsPtOn6p1c4R4zppD4tTUFGBtxReLxfjud79LIBDgv/wv/0v+4A/+YP68S5cuSXDc4MIeOw+nshiGiapKb8F68jsC7G86wL0pa4eW7kC3zFMU60bRNGu3mSV2nDENAzOXs8JjLmeFyZkZjEwac3yMUqG44I0UFIf9aa+jw/HcPQ7HhumNzJz/EcaTJ6htbRISRV1ac0g8deoUH3zwARcvXmRqaorJyUk++OADotEo586d43vf+x69vb0cPny4Eu0VdSzidXB/cprkTJGQR3ZfWW8OzcHeyF6GkkOz8xSn2eLrknqKoqbm5i/iXrp32yyXrfA4d5uZmb3PYSZT1uPZhTTWG84GSYdjwc0JDjuKffa43b6hwqQQtVKxYtpffvklwWCQ7u5uAGKxGCdOnODLL7/k2LFjDbkF31JkdfPLfT4wQbPPSbRFdl+ppSfTTxhOP8Bj87IjuAO7JqFdNC6zUHgaIgsFzFwe8nnMfB4zn7PuF/ZIzlLsNhS7A+x2K1ja7Sg2uxUobTYrTNpnH9vWfy7vyKEj8z2J7Vcvr/v1ReNo2NXNc956661FH0ejUa5evVqptxcNIuyxMzWdByQk1lKbpw2PzUN/op8bk9fpCe7Ab5c/akRjUuxWwOMlvxxNw4BCYTYwFjALBcgXMAv5+WBpptLWx7P1IRdfRLFCpc0ONpv1WLeB3WbtUmO3W/c2u1UaSNdB12sSLoVYL1L5WFRUyOPgSTJHrljGadNq3ZxNzWf3cSBygIFkP3fit9ni66Ld0w7Av+r/M7LFLG6bm9/e8Ts1bqkQa6eoqrUAZhkFwE3TtAJloYBZLFqPi0WrN7I4d6yIkctBsYRZLCwdLMEKl7oGuu3pvW02QM4FSU17+vHcY02bH0Y3TRPTNGV4XNSdZYfEDz/8kI8++mjVF/qTP/kT/uiP/mjVrxeNIei2oQBT0wXag65aN2fTs2k2dof28DAzzHDqAZlChu5AN/+6/8+YzE0ScUYkJIpNR1EUa87iCmo6mqYJpZIVIItFzFLp+ftS6en9zAxmqWSFy3IJs1R+usp7lvPv/T3IzYDTRe5SrxUmNRU0/ZnHqhUuVW3xY1WdPUdDUTVQldnn1NnnrMeoqvU5LzymKBJKxSstOyT29fWt6UKXL8t8i81A11R8LpuExDqiKApdvq14bV5iyRg3Jm9gmMarXyiEmKcoijUMvYbhZbNctgLk7L398GEwFhwrl6FUxjSWeFwsYRp5K4SWDeu4YViPZ4ucr/hzmg2LqLMBUlGXOKaAsiBUqsrs83PPKcAzzz17s6729D0Unh5feL/UscUtnj2sLHV45Z5ZkbF4icZzT778fu7x/PFn3tM0Fx5c/PpnXzf72Jx7jWlSyGRW+MlVxrJD4qVLl9i5cyfBYHDFF0kkEkvu6yw2ppDHzsO4lMKpNyFnmP2z9RRLxguGzoQQVWP1DmqrzjQvYxqGFSYNw3q86GMTzKcfY5iY5uw5szfrnNnzDMN6bJizQ+Kzzy04TrlkhZjZIGPOPWeaT4PO7DA6C8KO1dglQtGiYMRzva51QVkipC4Kuc+c+1zYXRCUl3rts4F5wevMfKHCn8zyLDskBgKB+ZXLK2UY0muxmYQ9du5PTJPOFQm4ZVVtPXHqTvZF9nF02zEmZiYIOUOUjBK6KtOThWhkyuywMqy+Y61emc/21i190urefKkh9wXH6mVI3pFK1eS6y/7NMDQ0xKeffgrA0aNHVxwY33333ZW1TDQsj0PHrqnEpwsSEuuQqqj8owN/QDwXZzAZ48bkdXYEd+KxeWrdNCGEeM6SPXdiXayoJ/HkyZOAVRPxxz/+MeFwmKNHjy6rZs+RI0dW30rRUBRFIeixMzVdoLu51q0RLxJ2hvEs2M6vy9dFq7utbv5yFkIIUVtrLqb92WefEYvF6Onp4bvf/W6l2lW3pJj28oylctwZSfF2NIJDSuHUNcM0eJh+yJPpEfyOANFAVIpvCyFEHalV9qjYjivJZJLe3l6mpqY4fPgwb775ZiXetu5ISFyeYsnglwMT7Grz0RqQVc6NIJlPEksOYJom3YEoIWeo1k0SQgjBBgiJCw0ODtLb24uiKBw9epTt27dX+hI1IyFx+b66P4XTprKnI1DrpohlKhpFBpMxErkELe5WuvxdaIr0BAshxIukz53HzGRQvF58p09V5RobKiQu9N5775FIJDh9+jS/+7u/W81LrQsJict3f2Kax1NZ3tnRJPPcGszo9CgP0vdxaA6igR68dtlmUQghlrIee3DXKnuo1XjTn/3sZ/zgBz9A0zQuXLjA5ORkNS4j6lzIY6dkmKRmirVuilihVk8r+yMHUBWNW/GbPEw/lALcQgixyVSsONrQ0BDnzp3j/PnzTE1NEY1G+eEPf8jx48dXXV9RNDafU0dXFRJZqZfYiNw2N/si+3iceczI9GMS+QQ9wR5cuswxFUKIzWBNITGVSvHpp59y7tw5+vr6CAQCvPvuu5w+fZq33nqrUm0UDUpRFEKzpXC2NUkNvkakKipbfFsIOoLEkgNcn7gmpXKEEGKTWFVI/OlPf8onn3zCxYsXMU2To0eP8umnn/J7v/d7lW6faHAhj527T9IUSwY2vSqzG8Q68Nq97G86wKP0Qx6kHhDPxdnu78Ztc9e6aUIIIapk2SHxq6++4ty5c3z66adMTU1x8OBBfvjDH3Lq1CkCgVevXv3www/56KOP1tRY0XiCs8PMU9kCLX5njVsj1kJTNLb6txFyhud3aun0dtLmaUdV5A8AIcTm5D11cn5180az7NXNqqoSCoU4deoUp06dWtE8w2QySTQa3RALWGR188r1DcXxOHR2t8vXa6MwTIPHmUeMTI/g0l10B6KyrZ8QQlRJrbLHioabT548STgc5ic/+cn8sbmM+aL5SRMTE5w/f55kMrmGZopGFvLYGUvmat0MUUHWXMWu+V7Fm5M3aHW30entRFOlrqIQQmwEy+5J/I3f+A3+4i/+YtUX2rFjB/39/at+fb2QnsSVm5oucP1hgoPbw3gcFVtQL+qEYRo8mX7C48wjdFVn2+yQtBBCiMqo+zqJx44dW9OFTp8+vabXi8bld9lQFYVktlDrpogqUBWVDm8HrzW9hkt3cW/qHnen7pAv52vdNCGEEGtQ9R1XNhrpSVydbx5MoWsq+zpli76NLp6Lcz91n7JRosPbSZunTRa2CCHEGjTEnEQhVivgtvN4KotpmlJfb4MLO8ME7AEeZR7yKPOQ8Zlxtvq2EnKGat00IYQQKyB/3ot1EXTbKBkm0/lSrZsi1oGmWuVy9kcO4NAc3Ju6y+34bWZKM7VumhBCiGWq65CYSCQ4c+YMp0+f5tixYxw6dIiLFy8u+/WHDh2qYuvESvicNjTF2qJPbB5um5s94T3sDO0kX85xfeIaD1L3KRnyx4IQQtS7uh1unguIZ8+eJRgMAtDX18ehQ4c4fvw4Fy5ceOnrT58+TV9f3zq0VCyHqir4XTYS2QJbwrJLx2YTcoYJOILzq6DHZ8bp8HbS6m6V+YpCCFGn6van80cffbQoIAIcPHiQs2fPcvHiRXp7e1/42r6+Pq5cubIOrRQr4XfbSM0UMQxZK7UZza2CfqP5TcLOCA/Tw3wz/jUTMxPI+jkhhKg/dRsSL168uORw8dGjRwFe2pP4ySef8N5771WtbWJ1gm47ZcMkI/MSNzWbZqM70M2Bptfw2DzEEgPcmLxBMi8F94UQop7UbUiMRqPE4/Hnjs/1LC71HMDHH3/Mhx9+WM2miVXyOXU0VSEh9RIF4NJd7AztYm9kH6qicid+m5uTN0kVUrVumhBCCOp4TuKlS5eWPD43z/DIkSNLPheNRhcNUYv6oSgKAZeNZLYIkVq3RtQLn93Hvsg+pnJTPMo84vbkLfwOP1u8W/DafbVunhBCbFp125P4IufOnSMYDHLq1Knnnvvkk084fvx4xa+Zz+dJpVLzN7F6Qbdd5iWKJYWcIQ40HWBnaCfFcpGbkze5E79NppCuddOEEGJTaqiQ2NvbS29vLz/60Y+e6y2s5jDzRx99RCAQIBAI0NXVVZVrbBZBtx3DNEnnpBSOWFrIGeZA02v0BHeQL+e5OXmT2/FbpGTOohBCrKuGCoknTpzg3Llzz/UWVnuY+cMPPySZTJJMJhkeHq7KNTYLt0NDV6Veong5RVGIuCK81vQ6O0I7KRllbsdvc2PiOlO5uKyGFkKIdVC3cxKfdeLECT788MMXDjOfPXu2atd2OBw4HI6qvf9moigKQbedRLbANjy1bo6oc4qiEHaGCTvDJPNJHmcecW/qHi7dRZunjYirSeosCiFElTRESDxz5gxHjhzh/ffff+65ixcv0tfXx+nTpxcdn6uTOHf82ZqLonYCbjuxsTRlw0RTZR9nsTwBR4CAI0C6kObJ9AiDyUEeZh7S4m6lxd2CTbXVuolCCLGhKGadj9ucP3+eRCLxXEA8f/78kr2Kc06fPs358+crPiyVSqUIBAIkk0n8fn9F33uzyOZLXB2Kc2BLkJDHXuvmiAaVK+V4Mv2EiZlxTEyaXM20ultx22RHHyHExlKr7FHXPYm9vb1LBkSwtu0Tjcnt0LFrKolsQUKiWDWn7mR7YDudvk7Gs2OMZkcZz47hs/tocbcScoY2zFD0v+r/M7LFLG6bm9/e8Tu1bo4QYpOo25AYi8U4ffo0R48e5cyZM8DTYDj33Mu8qNi2qA8Bt92qlyjEGtlUGx3eTto87fwPt/8l49kJFAV+dcuv0+JuodnVjF1r7D9G/nX/nzGZmyTijEhIFEKsm7oNiceOHSMWi3H+/Pkln3/RQpXz589z6dIlLl68CMChQ4c4fPgw586dq1pbxcoF3Tb6R3OUyga6tjF6e0RtqYrKpft/wWRukpAjxL/T8w8YyTzmceYRQUeQJlczQUcQRZF5sEIIsRx1GxIHBgZW9bpTp069dK6iqA9Btx0TGE3m6AzLHDJRWaqi0h3opsvXRTwXZyw7xr2pu9g0O82uJppczTh1Z62bKYQQda1uQ6LY2Jx2jc6Qm9h4Bpuu0uKXX9ii8nRVp8XdQou7heniNOPZcUazozzOPMZr99HsaiLkDKOr8qNQCCGeJT8ZRc10N3solg3ujqTQVYWwV2pRiurx2Dx4Ah66/F0kclOMz0wwlBpiKDVEyBki4mwi4AhsmMUuQgixVhISRc0oisLOVh9lw+T24xT7twQIuF++wOBf/nyI6XwJj0Pn3/vO9vVpqNhQNEUj4moi4mqiUC4wmZtkcmaCe1N30VWdkDNMxBXBZ/PJ/EUhxKYmIVHUlKoq7Gn3c/1hgpuPkrzeFcLjfPG35X//iyHGU3ma/Q4JiWLN7Jqddk877Z52ssUsk7kJJmcmGc+OYdNss7u9RPDavBIYhRCbjoREUXOqqrCvM8C14QTXHyZ4rSuI2yHfmmJ9uW1u3LatdPm2kimkiefiTObijE6PYtPshBwhws4QPrtfAuM6kzqRQtSG/CYWdUHXVA5sCfLN8BRXh+LYNRWnXcNl03DaNdx2nbAU3hav8A92/M58mFgLr92H1+6zAmMxw1QuPrtKehRN1Qg5QoScIfyOAJqiVaj14kWkTqQQtSEhUdQNm67yxtYQ8ekCuUKZ7OxtMpOnZJj4nLaKb7MoNpZKBwhFUfDZffjsPrb6tzFdnJ4PjBMzEyiKQsARIOgIEXQEG75otxBCLCQhUdQVXVu6HE5qpsitx0mKZQmJonY8Ng8em4ctvi5ypRxT+SkSuSnup4YYMk3cNg9BR4CAIyjzGIUQDU9ComgIfpeNt7aFmfuVWzaWFxZlNbSoFqfupF23Fr2UjBLJfIJEPsFYdozHmcdoqobfHiDoCOB3BHBoqy/xVKlhdCGEWAkJiaJh2HUVXbNiYqlscPtxkh2tvpdu6yerocV60FV9vqyOaZpMF6f5yb2LJHJTKIrKO+3v4NSd+O0BAg4/Prt/RQW8ZR6eEKIWJCSKhjI3fKdrKvHpAl/dn6Kn1UdIFrWIOqEoCl67l78c/hmTuUnCjjD/wb7/kGQ+SbKQZCw7iqIouHU3/tnA6LP50FRZACOEqC8SEkVD0lSFt7aFuPckzfWHCZp9TqItXuy67JYh6ouiKLP1FsMA5Es5UoUUqUKKiZkJRjIjz4VGr80rWwUKIWpOfgqJhuWy67y+NcRYMkdsPMPVwUm2NXlpDzplwYCoWw7dSbPupNndAkBuNjSmCykmZiYZyYwAVt1Gr82Hz+7Fa/etaU6jEEKshoRE0VD+3W9vn1+IMqcl4CTstTM4nmFgLM1YKsfOVt9Ld24Rol44dSdO3UnLgtCYLqTJFNOkZoenAWyaHa/Ni9fuxWvz4rF5ZJ9pIURVyW9R0VBetPhE11R2tvlp8TvpH83w1YMpuiJuqasoGs5caGymGYCiUSRTSJMpZsgUMjxKP8Qwjfkhas9sYPTYPLh0l/SiCyEqRkKi2FACbjtvbQvxYHKaBxPTUldRNDybaiPkDBOandNomAYzpRkyhQyZYmZRb6OmarPB0QqNbpsHpybTL4QQqyMhUWw4qqqwvdlLxPt0DlepbGCapvyyFA1PVdT5ENhKKwBlo8x0MUOmOG3tCpOf4sn0k/nz3TbPbHh049LduGwu2U5QCPFKEhLFhuVz2bDN1lUsGyZfP0iwp92P0y6/HMXGoqka/tmi3XNKRolscZrpUpZscZpUIcn4zNj8H0sOzTE/RO3SXbh1N3bNIX9ICSHmSUgUG9rcLzybrlIsG3x5P86uNj8Rn6wUFRubrurPBceyWWamOMNMKct0MUu2lCWRT1A2yoAVNudCo3P23qW7sKt2CY9CbEISEsWmoCpWXcW7T9LcfJykI+iiu9mLqsovPrF5aIpmrY62e2eXxVgK5QLZUpaZ0gzZohUgJ2cmMUwDsIasreDoxKW7cGjO+Y/XY4W1bEsoRG1ISBSbhq6p7OsM8Hgqy+D4NKmZIns7AjL8LKqikYKNXbNj1+wEHcH5Y6ZpUijnmSnnmCnNMFPMkivnF/U8zr3Wqbtwak6cugOHZq3OdmiOigVI2ZZQiNqQkCg2tKXqKnaE3PhdNm49TtE3O/zcJMPPosIaPdgoioJDd+LQnYvCI1hleXKlHLnSDLlyjlwpT6aYZmJmfL73EawA6dAc86Hx6c2JTbOt82ckhFgpxZRCciuSSqUIBAIkk0n8fn+tmyPWoFQ26B9NM57O0xlys73JI8PPQqyBaZqzAXKGfDlPvpwnV86RL1n3C3sgVUWdD5F2zYFjwWO7Zpd5kEIsUKvsIT2JYtPSNZU9HQH8C4afZfWzEKunKMr80PVSSkaJQrlAvpybD5GFcoHpYoZ4Lr8oRCqKgq7acGh2bKodx+z72lQbds2BTbVh02xSykeIKpKQKDa9jpAbn9PG7ZGUtfq53b+oxqIQojJ0VUdX9RfO0ywbZSs4GgUK5blbnoJRJJFPUCgXFg1ng7Ui267asWk2Kziqduzzj23os2FSV3TpmRRihSQkCoFVU3F+9fOjJFtCbrY3e+SXihDrSFM13KobNy9e7FM2yhSMAsVygYJRpFDOUzSKFI0i+XKBTHGa4hJh0uqZ1GeD49N7fTZMaoqGbfbjuTAre2OLzU5CohCz5lY/P4xnGRrPkM4V2dMRwK5X9hfF/6v3Lr//q1Hcjs31n181P+9n33st13rRazfrv9vL/Hc3/1tO7HoXl+5at2tqqoZLdb3ymmWjPB8eS0aRolGyPi4XKZnW8WwxS8ksUTJKS+7zrirqfGDUVR1N0dFV7bl7XbWhKSqaqqMpGpqqyTC42BDkp50Qz9gSduNz6tx+bA0/7+0I4HdVbiXm39we43cPd226sFHNz/vZ917LtV702s367/Yynz/+Jf/29n97XUPicmmqFdacOF95rmmalM0yJaM0HyhLRonybIAsGWXKZomiUaJQzFMyy5SNEmWzvGS4BKvncmFgVBVtNlhqqIo6f2zueU1RURXVOvbMvTr7nIxsiPUmP+2EWELAbefNbSFuj6S4Npygu9lDR6j+690JIVZOURR0xeotZBmhcqGyUaZkligb5fmgaZjl+WBZNsuUDYOyWcYwrd7NuXPLpoFhljFM44Vh89l2qguCo6aoKCjzIXI+TKLMnqPMhsu5Y9bxucfWcVAUFRVlwXnW6MlcMJ3/n/L8vdjYJCQK8QIOm8ZrW4IMTUwzMJYhNVNiZ5sPTcrkCCFmaaqGhgZrHF0uz4ZFw1gcHq3H1sdPHxuYPH288Db/PvPnmU8fmyYm5rIC6XLNBcW5YPk0PCqzAdR6zOz/L3x+/th82FQW/D/PhNDFP3eX+im80tD67Nfh+a/K0yNz5y71isXPLX6NyfPPWcefvsfC8xa/39PzM+np5X1SFSYhUYiXUFWFaIsXr1On/0marx+U2NcZwGmT+UZCiMqxhpw1UKtfZNw0TQxmQ+PCx7OBZS5cmgsCjGkamLDg3mQuxliLhMznjj99bGLlnqfhxzRZEIyWDk9PH809v+izWPPX4dlw+vT4gscsHVSfD7LKonMVZeFrlSWPzb1WWfDxwiss7K3NKOmVfGoVIyFRiGVo8Tvx2HVuPk7y1f0p9nb4CbiXrgUnhBD1TFEUq/dTBkUahrO0smkQlSLr+4VYJo9T582tITwOnWvDCR5PZWvdJCGEEKJqJCQKsQI2XeXAlgAdITcDYxnuPUlhGLKzpRBCiI1HhpuFWCFFmZ2n6NC5N5pmOl9mb4cfh8xTFEIIsYFIT6IQq9QScPJ6V5B8qcxXD6ZIzxSXPM80TXKF8pLPlcoGxZKx5HONrpqf97PvbZom+WJ5/r0XTnB/1bVe1E7TNMnkiku+diP/u72MaZrkS7klnyvNFqsWQmwcEhKFWAOfy8abW0M4dI1vhhOMpZ7/BTo0nuG9//zf8Mc//YbEdAEAwzQ5/7N7vPef/Rs+/fz+ejd7XVTz8372vYtl4//f3t3EtnHu+x3/US+k3jWk/JI4tpMMT85J7nlJDmUDLdCuQvYCRTcFSHtxF8XBram2FwdFUUCEF+05d+VSxUW76YJMUXTRokioVdF2QwZFi+KigC2enPd7k3Ic27Hj2BY5epcsWexCIY+oISVSosS37wcwEg3Jmed5Zubhn888L/qz/3hPf/Lv/o/C//Z/a2lju+5j1Urnv/7vv9Mf/6v/qX/47/9vT523w3y18pX+8af/SP9m4S+0tLUkSXpV3NV//v1/0j/JzOi/5f5ri1MIoJlcxWZOmNQDlpeXNTk5qaWlJU1MTLQ6OWgTu7tFffHNip4tb+qKb0Rvnqtc93ltc0f/4X/l9JdfPNf2q11dMka0s7urf/H3f6jXjfZbraJZTjPfB/d9cWJIj/Lr5WUUGzlWrXT+87/7nv7HZ0967rwdZn17XR//9X/Rvad3tVPc0cWRi9rZfaV/Ov3PdHHkYquTB3SlVsUetCQCTdDX59L3Xp/Q2+fH9Ci/rt8/WdbOqz88jhwdGtBP//h7+hvfOadny1tae7mjv/iTQNcHGqeZ74P73tzZ1cc//Vv629+70PCxaqXTf3G8J8/bYUYGR/STH/ypAhen9WLjhdZ3NvQv/+bPCRCBLkSQCDTRZd+Ivv/GpOz1l/rlQ9vR1+3PQt/Vxckh/fTvfFfD7t4ZN3aa+T6475Mcq9Zne/W8HeYffP8nOj98Xj/5wZ9qaKA1c7gBOF0EiUCT+cY8ev+qV7vFoj57WNDS+svyawP9ffp7H7yhH13xtjCFZ+80831w3yc5Vq3P9up5O8xA34CCb4b0nu+9VicFwCmhT2KD6JOIem3v7Or3T5a0srmj71wc08XJ3n1ECQA4PvokAl1mb+JtQxcmhvT50xVZz1YdC8oDANCu6FwDnKK+PpfeeW1co55+Wc9Wtf5yR+++PqGBfn6fAQDaG99UwBm45B3R9y8bWtnYrjqgBQCAdkOQCJwR76hb71/1qlhlQAsAAO2GIBE4QyOeAb1/1asRd79+89WSntobrU4SAABVESQCZ+wPA1o8+uKbFd1nQAsAoA0xcAVogb0BLRMacQ/o/vO9AS3fY0ALAKCN8I0EtNAbvhH90RuTWmJACwCgzRAkAi3mG/PoAwa0AADaDEEi0AYODmj5ZokBLQCA1iJIBNrE/gEtnz9lQAsAoLUYuAK0kYMDWtZYoQUA0CJ88wBt6A3f/hVaCtp4udPqJAEAegxBItCm/rBCi/TZg4LsNQa0AADODkEi0MZGPAP64E2vxocG9ZuvbD0pMKAFAHA2CBKBNjfQ36fvX57U68awcs9W9MXTFe3uMqAFAHC6GLgCdACXyyX/xXGNegaUe7a3Qst7lyblHuB3HgDgdPANA3SQ14xh/eDypDZfvtJnDwpa3dxudZIAAF2KIBHoMJMjbr3/pleD/S796qGt58ubrU4SAKALESQCHWhosF8/uuqVb8yjv/p6WfefM/E2AKC56JMIdKj+PpfevTShsfyAvny+qrXNHb17iYm3AQDNwbcJ0OEulybe3tzWLx4UtLbJxNsAgJMjSAS6gHfUrR+/6VN/n0u/fFignyIA4MQIEoEuMeTu1/v7+yk+o58iAOD46JMIdJFSP8Xx/IDuP1/Vyua23mU+RQDAMfDNAXShN3wj+uEVQxsvX+kXD/JaWmfdZwBAYwgSgS41OeLWB296NTTYr18/svU4v97qJAEAOkhbP262bVt37tyRbduyLEv5fF63b99WOByu+ZlkMqlcLlf+2+/3KxqNnkVygbbjGezXDy8bevBiTda3j5+/c3GcaXIAAEdq2yDRtm3FYjHF43EZhiFJymazmp6eVjgcViqVcrw/EokoEokoHo+X3x+JRHTt2jUFAoGzzgLQFvr6XHr7wpjGhwf1+dNlffagoPcuTWp0qG1vfwBAG2jb5oQ7d+5UBIiSFAgEFI/HNT8/r0wmU/H+W7duyTTNilbDfD6vfD5fsQ+gV50b9+jHb3rV1+fSZw8L+treaHWSAABtzFVs0zky/H6/JFU8Opb+0JoYjUaVSCQqtuVyOZmmearpWl5e1uTkpJaWljQxMXGqxwJOw+5uUdbzVX1tb+j8uIfHzwDQ5loVe7TtN4Npmsrn847tpVbB/a/duXNHhmGceoAIdIO+Ppe+c3Fc774+ofzaS332oKDVze1WJwsA0GbatlNSOp2uuj2bzUqSrl+/Xt5mWZZM05Rt20omk5KkxcXFpg1a2dra0tbWlqS9aB7oBucnhjQ2NKC/erKsXz609da5Ub3hG2l1sgAAbaJtHzfXEgqFdO/ePd2/f7/cquhyuRQMBhUKhTQ7O1t+byQSkc/nKz+WPq6f//zn+vM///OKbTxuRrfY3S3qyxdrelxYl3fUre++NsHk2wDQRlr1uLmjgsRMJqNQKKRUKlWeBse2bXm9Xkly9Em0LEt+v1/pdFrBYPDYxz3YknjlyhWCRHSd/OqWvni6Ikl657Vx+cY8LU4RAECiT2JdIpGIEolExTyJpdZE0zQdfRJLf5+0JdHj8WhiYqL8D+hGvjGPfvyWT6NDA/rt4yVZz1a1u9sxvyEBAE3WMUFiJBLR7du3a/YxrDXNjWEYsizrFFMGdA/3QJ9+cNmQ/8KYvrY39AsGtQBAz+qIIDEWi+n69esV/Q33Kw1aqabWdgC1XfKO6IOrXrlc0mcPCnq4uKYO6pkCAGiCtg8Sk8mkpqamHAFiaRSzJAWDwarT5ZRcu3bt1NIHdKvRoQF9cNWry74RPXyxpl8+tLW+tdPqZAEAzkhbB4mZTEa2bVdtQdzfQhiLxcrrO+9Xmi5nZmbmVNMJdKu+PpfeOj+mH131aufVrj57UNCTwjqtigDQA9p2dLNlWQqFQgoGg+X+hqXA0LIszczMVAxgicViymazFfMrTk9PKxgMltdybgZWXEGverVb1JfPV/XE3tDE8KDeuTiuEU/bTrUKAF2DKXAO8Pv9hw44WVhYUCAQqNiWTCaVTqfl8/mUz+d18+bNikCyGQgS0euW1l/qi6cr2trZ1dWpEV32jcjlcrU6WQDQtQgSOwRBIrDXqvjw2wm4Rz0D+u5rExodolURAE4D8yQC6Bj9fS69fWFM71/1qijpFw/yuv98Va+YVxEAugZBIoBjGx8e1AdXvbo6NaonhQ0t3M9rcXWr1ckCADQBQSKAE+nrc+nquVEF3vJqxN2v3z1e0u8eL2lz+1WrkwYAOAGCRABNMewe0A+uGHr39QmtbG4rez+vr/LrLO0HAB2KnuYAmur8xJC8o249eLGmL5+v6pulDb19fky+MU+rkwYAaAAtiQCabqC/T/6L4/rxmz4N9vfpt4+X9JuvbG28ZMUWAOgUBIkATs3o0IB+dNWr9y5Nan3rlbJfFnT/2ap2Xu22OmkAgCPwuBnAqTs37pF31K3HhXV9tbiub5Y3dcU3oteNYfX1MRE3ALQjgkQAZ6K/z6WrU6O6ODGkh4truv/tEn9vnhvV+XEPq7YAQJvhcTOAM+UZ7Nc7r03ox2/5NOLu119/vaxfPrRlr71sddIAAPsQJAJoiVHPgL5/2dCPrhiSpF9/ZevXj2wtrRMsAkA7IEgE0FKTI2598KZXf3RpUtuvdvWrR7Z+88jWysZ2q5MGAD2NPokA2sLUuEe+MbcWV1/qwYs1ffawIN+oW1enRjU+PNjq5AFAzyFIBNA2XC6Xzo17NDXm1vOVLT38Nlj0jrh1ZWpEkyPuVicRAHoGQSKAtuNyuXRhYkjnxz16sbKlR/l1/eqRrYnhQV32jWiK1VsA4NQRJAJoWy6XS+cnhnR+Ykj51b1g8XePlzTqGdAb3mGdHx9inkUAOCUEiQA6gm/MI9+YR0vrL/Uov67Pn67oyxdret0Y1iVjWAP9jMMDgGYiSATQUSZH3JoccWt9a0ePCxt6tLi3isuFySFdMoY14qFaA4BmoDYF0JFGPAN657VxvXVuVF/bG+V/xohbrxvDmhpzs4oLAJwAQSKAjjY40Ker50Z12TeiF6tbelLY0O+fLMk90KfXjWG9Njks9wCPogGgUQSJALpCX9/eiOgLE0Na3dzW1/aGvlpc18MXa5oa8+ji5JC8o7QuAkC9CBIBdJ2xoUG989qg3j6/q2fLm3q6tKnfPt5rXXxtclgXJ4Y05O5vdTIBoK0RJALoWgP9fbrkHdEl74hWNrb1dGlTjwvreri4JmPErQsTHk2NeRgZDQBVECQC6Anjw4MaHx6UeWFML1Y29Wx5S58/XVG/a1W+MY8uTHh4HA0A+xAkAugp/X0uXZwc1sXJYW1uv9Lz5U09W97Ubx9vyt3fp6lxj86PezQxPEjACKCnESQC6FlDg/26MjWqK1OjWt3c1vPlLT1b2dTX9obcA306Pz6kc98GjADQawgSAUB7g13Ghgb11vlRLW9s68XKlp4v7/Vh9Az0a2rMXQ4YaWEE0AsIEgFgH5fLVV7VxbwwpqX1bS2ubu3NwWhvaLC/T1Njbk2NeWSMuFk7GkDXIkgEgBpcLpeMUbeM0b2AcWVzR4urW1pc2dLTpU31f/u6b9Qt76hbnkGm1QHQPQgSAaAOLpdLE8ODmhge1Nvnx7S2taP86pbyay/1/75ZUVHS+NCgvN8GjONDAzyWBtDRCBIB4BhGPQMa9QzoytSotnd2VVh7qfzalp58Ow/jQJ+rHDAaI7QyAug8BIkAcEKDA326MDmkC5NDKhaLWt7Ylr2+rcLaS33+dEWSNOIekHfUrcmRQU0ODzKBN4C2R5DYoGKxKElaXl5ucUoAtCuXJK9b8rr7tb3jkr3+Ussb6/ryybJevnoll6SxoYFvH1+7NTE8qH4GwACooRRzlGKQs0KQ2KCVlb1WgStXrrQ4JQAAoJcsLi5qcnLyzI7nKp51WNrhdnd39eTJE42Pj/d0p/Tl5WVduXJFjx490sTERKuT03KURyXKoxLl4USZVKI8KlEelZaWlnT16lUVCgUZhnFmx6UlsUF9fX26fPlyq5PRNiYmJriB96E8KlEelSgPJ8qkEuVRifKo1Nd3tn2Z6TkNAAAAB4JEHIvH49HPfvYzeTyeVielLVAelSiPSpSHE2VSifKoRHlUalV50CcRAAAADrQkAgAAwIEgEQAAAA4EiQAAAHAgSAQAAIADQSIAAAAcCBIBAADgQJAIAAAAB4JEAAAAOBAkAgAAwIEgEQAAAA4EiQAAAHAgSAQAAIADQSIAAAAcBlqdALQH27Z1584d2bYty7KUz+d1+/ZthcPhuvcRCoUUCAR08+ZNBQIBWZalRCIh27aVSCROMfXN16y8zM3NaXFxUVNTU8rlcgqFQg2VaTvIZrO6c+eO4vG4TNM89n46+fqwbVu3bt3SzZs3Dz1/zTjfnXDN1FMevVSn1FMevVanHFUmvVSvNHovtFU9UkTPKxQKxWg0WiwUCuVtCwsLRUnFcDhc934CgUBRUsW/YDB4Cik+fc3ISzQaLcbj8YptwWCwmEgkmpnUU5dKpRxlcfBfPXnqxOsjHA4Xo9FoMRqNFiUVU6lUzfc243y3+zVTb3n0Sp3SyPXRK3VKvWXSK/VKo/dCu9UjrmKxWGw8tEQ3icViun37tgzDqNg+NzenWCymdDqtYDB45H4ikYiuX7+uu3fvyjRNhUKhuj7Xjk6al2w2q+npaR28vWptb2dzc3NKp9NVf+3n83lZlqWFhYUj99PJ14dlWfL7/UqlUjVbRU56vjvpmjmqPHqtTjmqPKTeq1OOKpNeqVcauRfash5pOKxE1zFNs2iapmN76ddONBqtaz/1vq8TnDQv0Wi0GAgEqr6mI1oc2s3s7GzN1w7+Qj5MJ18fuVzu0PPWjPPdSdfMUeXRa3XKUeVRLPZenXJUmfRKvdLIvdCO9QgDVyDTNJXP5x3bS798qr2Gw2UymZr9bAzDUDqdPuMUHV8oFKq6fX5+XqFQyPELuRc143x30zVDndJ83XR9SL1TrzRyL7RjPUKQCKXTaRUKBcf2bDYrSbp+/XpD+8tms0omk+XPd7Lj5sWyLPl8vqqv+Xw+3bt3rxnJOxPVHt1YlqW7d+8eqyN0N10fJc043910zVCn1EadsqdX6pVG7oV2rEcIElFTIpGQYRiKRqN1vT+fzysWiymfzysajSqfz2t6elqWZZ1ySpvvJHmxbfvQ1w3DOPI97S4Wiykejzf0mW66PvZrxvnuhWtGok6hTjlcL9UrB++Fdq1HmAIHVWUyGWUyGaVSqbqb/UOhUEXlHwwGdfPmTYVCIeVyuVNK6enoprw02/z8fMMtQRJl2uuoU7onL6ehl+qV49wLrUJLIqqKRCJKJBINNftXax0IBoOyLEvz8/PNTN6pO0lejrrpO/0XfywWO9bjoG66PvZrxvnu9mtGok6hTjlcL9Ur1e6Fdq1HCBLhEIlEdPv27bofCR2m1IG20zpVV9OsvOTz+bb/9VhLNpuVZVknmvx2v266Pmppxvnu5GtGok6phTplTy/VK8e9F1pVjxAkokIsFtP169c1Ozvb0OcikYimp6drvt5JoxmbkRfDMGq+z7ZtXbt27djpa6VEInGsirybro9qmnG+u/WaoU6hTjlKr9QrR90L7ViPECSiLJlMampqynEBJ5PJIz+bzWarjqgqXazH6WvSKs3Iy40bNw7tOF1r+od2l8lkjvVrtpuuj2qacb678ZqhTtlDnXK4XqhX6rkX2rEeIUiEpL2b1Lbtqr9w6unHEA6Hqzbtl/qENOMx01lpRl4ikYiy2ayj7DKZjKTq0z90gsOmVzhMN10f1TTjfHfbNUOd8gfUKYfr9nql3nuhLeuRhqbeRlfK5XJF0zSL0Wi0ODs7W5ydnS2vvRkMBh0ztFebQb60PuV+CwsLRcMw2m4lgKM0mpdaM+qHw2HHqgLVyrNTFAqFutZF7cbro7Q6wmFrnzZyvjv9mjmqPHqtTjmqPHqxTqnnnikWu79eafReaLd6hLWbIb/ff2jz9MLCggKBQPnv6elp+Xw+xy8427YVi8VkGEZ5f7dv3674bKdoJC+1ykPaW59zcXFRU1NTyuVyCoVCxxrB1y68Xq+i0eihc5l10/URi8VkWVa5Y71hGAoGg/L5fEokEo7313u+O/Waqbc8eqVOaeT66JU6pdF7RurueqXRe0Fqr3qEIBEAAAAO9EkEAACAA0EiAAAAHAgSAQAA4ECQCAAAAAeCRAAAADgQJAIAAMCBIBEAAAAOBIkAAABwIEgEAACAA0EiAOBUzc/Py+/3y+/3KxQKtTo5DZuZmSmnf25urtXJAc7MQKsTAADofjMzM5qdnW11Mo6ltObw/Py87t692+LUAGeHIBHoAH6/X/l8Xj6fT6ZpyjAMSVImk5EkmaYp0zQl7S14f+/ePUlSPB5XNBptSZrhZFmWIpGILMvStWvXlE6nW52kthOLxTQ/P698Pq9PP/1UgUCg1UkCehZBItABLMtSIpFwBHyxWExzc3NKpVLlIFHaCxQjkYhyuVzF+/1+v4LBYLllpJN1Yl5M09TCwkJHPnLNZrNnErDF43HdvHlT09PTp34sAIejTyLQ5mzbVjgcPrRF0OfzVfxtGIYSiYRs23a8t9QK2Q06NS+d2Dp2lsH4/h88AFqHlkSgzeXzed28ebPhz+1/LF1ysGWxk3VTXjqBZVlndqxODf6BbkNLItDmbNs+dsvK1NRUk1ODXjQ/P9/qJABoAYJEoM0ZhnHsILETH2uivViWpVu3brU6GQBagMfNQJs7Sf+sYDBY/v9QKCTLsmQYhhYWFiRVH22bTCYl7T3OzWazisfjCgQCymaz5VHTCwsL8vv9Nac0mZubKz8yzOVy8vv9dY+ytixL8/PzMk1T+Xy+3JJ69+5dxePxY+XFtu3y52uVZ2l6k6mpKS0uLtZM80nydpBt2/rkk0/Kfy8sLCgWizkGIX344Yeybbsiv6W0JBIJ5fN5pVKp8vneXxY3btxQPB4vl0U6nS6Pqj6qrJPJpNLptHw+n+7du6dIJFI+9mFl2YhMJqN0Oi2/31/eduPGjZrvt21byWSyfA4WFhY0MzNT/kFUGh1duj4++ugjhcPhijKLxWIyDEPRaFTxeLyuaw7oSUUAHWt2drYoqVgoFOp6fzQaLQYCAcf2cDhcDAaDxUQiUbGvVCpVNAyjuLCwUEyn0xWfMU2zmEgkHPsKBALFVCrl2P/s7OyR6SsUCsVwOOzYnkgkHNuPm5dqZmdnHfsvFAqONJ8kbwePd1j5HSzr0nGq5TeXyxUlVf1MMBgshsPhYjweLxaLxWI6nS5KKuZyuYbLOhgM1p2/g1KpVDkN+83Ozhaj0WjFtkKhUIxGo0VJxYWFhaqfOfh+wzAq8l8oFIqSqh6zWNwrl1wuV35vveWQSqUaPtdAJ+NxM9BD9rfW7Geapu7du6dr165VDBoIBoOybVuJRKKiVbL0WiqVqthWWo1if8uNtNfqNDc3d+Tgh0wm4xipLUnRaNSx/bh5yWazjmPOzc3po48+qtheal0qpfmkeaumWgtkPB6vaLEruX79etV9HNaaV2oBLp27YDCoQqEg0zQbKuvTkMlklEwmHaOmDcPQzMxM1c9ks9mKc1J6fzQaVSwWc2z7+OOPq+4nFAqVy63V5QC0M4JEAJqampJt244+jKUgq1pAZhiG8vl8xbZYLFZ1JHZppHVp8u9aTNPUJ5984gjkJFUNnKo5Ki/V0hwOhx0javP5fMV7T5q3g2qN4C0Fs81a/i2fz1eURem4zSjrk4jFYo4fHiW1At/SNXcwIPf7/Y5tMzMzymazVX8U7A/0W10OQDujTyIASUe3SB2l9CWdy+VqjoY9atqaQCCgYDCo6elpmaapYDCoUCikcDhcM6CoppG+cvtb2vYrtbpJzclbI0r94Zq1r2qaVdbHlc1mG57ayTTN8jmR9s6LbdsV/TRLAoGAAoGAEolERWvlwfPd6nIA2hlBIgBJJ5+brhRIRSKRql+uBx/T1pJKpZTJZMr/TSaTMk1TqVSq7tHa9eallOajpgpqVt7qZRhG0+YlPOyRaTPK+jhKeTvONWfbtu7cuSPbtjU9PV0O8PYPACqZmZnRzMxMOUisNZ1Uq8oBaHc8bgbQFPvXjj6uUvBQWm4vl8upUCgoEAjoww8/bEYyK5TSvLi4WNf7TpK3Rpxkbsx6nbSsTzJ34nHL07Isvf322/L7/eVlKk3TrBkIl/p8lkZ2J5NJR0B/1tcc0EkIEgE0Ralv3mGPSY8KCjKZjKNvn2EYSqVS8vl8p7Lqh2mah+63FLCdNG+NsCyr5kCVZh230bI+2JfzpOfCNM2GH9FHIhGZpukY8HMwbfvzFY1GD11SsBXXHNApCBIBNM3++fgO2j/P4mEOjpguOa3HfvF4XPPz81WDLcuyygFEM/K2X63gbn5+XoZhOOagNAyj6mcaHTCzX71l7ff7m96KGo/Hqz4illSzLGv1H11YWKhI3/7ArjSAZW5urma3gLO+5oBOQZAIdLDSF+PBlpR6PrPf4uLioUFAvQFCNBpVMBh0TGFi27YymUxdAwGqjTQtHf/gI9hm5CUcDisajVYdyTo/P18OLJqRt4PHPfjIttTf7tNPP3W8/8aNG+WBGvvfXzr31Vq89r9eTb1lHQwGK46dzWZPHECFw2HduHGjYuqa0vFLQdvBtAeDQUdQbFmWQqFQ+bPZbFbXrl0rv14awPLxxx/XfITfyDUH9BJXsVgstjoRAOqXzWbLHffv3btXfhwaCATk8/mqPlqzbVu3bt1SJpORbdsKBoPlFTP2bw+Hw7p586bC4bDm5uaUTqeVyWTK+0+lUrIsS7FYrOIzMzMzFUFSMpksr8pSWsWinlVJSqtelAKeUpBg23a5Za2ZeTl47NLKH4el+bh5O7iP0qoiB1dcicfjNQd0ZDIZJRIJXb9+vfyeaDQql8slwzDKc1dWO0emaVasHlJPWVcrn+np6fJx61Wa27DWfu/evetYccXr9cowDN24caPimp6ZmVE+ny8HhqURyaVzHIlEHGkrtQBXS3Mj5VBKK6uwoFcQJAIATtVhQeJZmJuba8qxCRLRa3jcDAAAAAeCRABA1ygt3VdycIUVAPUjSAQAdI1EIqE7d+6U/7Ysi8EnwDGx4goAoGuUBuwkk8mag3AA1IeBKwCAUzU/P69bt27J5/PJNE2l0+lWJ6khMzMzymQy5ZHsDFxBryBIBAAAgAN9EgEAAOBAkAgAAAAHgkQAAAA4ECQCAADAgSARAAAADgSJAAAAcCBIBAAAgANBIgAAABz+P/EVfMHRPADBAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Make a dictionary for colors on the plot\n",
    "band_colors = {'lsstg':'#4daf4a', 'lsstu':'#377eb8', 'lsstr':'#e41a1c', \n",
    "               'lsstz':'#a65628', 'lssti':'#ff7f00'}\n",
    "ax = kn_object.plot_data(show=False, band_colors=band_colors)\n",
    "ax.set_ylim(27, 18)\n",
    "\n",
    "upper_limits = kn_sim.observations[kn_sim.observations['detected'] != 1.0]\n",
    "\n",
    "for band in band_colors.keys():\n",
    "    up = upper_limits[upper_limits['band'] == band]\n",
    "    plt.scatter(up['time (days)'], up['limiting_magnitude'], s=100, marker=r'$\\downarrow$', color=band_colors[band])\n",
    "\n",
    "    \n",
    "# We can also plot the true data \n",
    "tt = np.linspace(0.1, 20, 100)\n",
    "# specify output_format \n",
    "parameters['output_format'] = 'magnitude'\n",
    "for band in band_colors.keys():\n",
    "    parameters['bands'] = band\n",
    "    out = redback.transient_models.kilonova_models.one_component_kilonova_model(tt, **parameters)\n",
    "    plt.plot(tt, out, color=band_colors[band], alpha=0.3)\n",
    "\n",
    "plt.xlim(0.1, 20)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6fa2a373",
   "metadata": {},
   "source": [
    "Now we plot the light curve in the Rubin ToO format for 350 Mpc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d3e1765f",
   "metadata": {},
   "outputs": [],
   "source": [
    "kn_object = redback.transient.Kilonova.from_simulated_optical_data(name='my_kilonova', data_mode='magnitude')\n",
    "\n",
    "# Dictionary for filter colors based on the ToO simulations code\n",
    "band_colors = {'lsstg':'g', 'lsstu':'b', 'lsstr':'r',\n",
    "               'lsstz':'k', 'lssti':'yellow', 'lssty':'orange'}\n",
    "\n",
    "# Use the plot_data from the Redback Kilonova Object documentation\n",
    "ax = kn_object.plot_data(show=False, band_colors=band_colors)\n",
    "ax.set_ylim(27, 18)\n",
    "ax.set_xlim(0, 7)\n",
    "\n",
    "# Plot detected data points (detections) as markers\n",
    "for band in band_colors.keys():\n",
    "    detected_times = kn_sim.observations[kn_sim.observations['band'] == band]['time (days)']\n",
    "    detected_magnitudes = kn_sim.observations[kn_sim.observations['band'] == band]['magnitude']\n",
    "    plt.scatter(detected_times, detected_magnitudes, color=band_colors[band], marker='o', s=100, edgecolor = 'k', linewidth = 1.5,label=\"Det.\" if band == \"lsstz\" else  \"_nolegend_\")\n",
    "\n",
    "# Save the upper limits (non-detections) as a variable\n",
    "upper_limits = kn_sim.observations[kn_sim.observations['detected'] != 1.0]\n",
    "\n",
    "# Plot the upper limits (non-detections) for all bands using a loop\n",
    "for band in band_colors.keys():\n",
    "    up = upper_limits[upper_limits['band'] == band]\n",
    "    # Plot non-detections as downward-pointing triangles\n",
    "    plt.scatter(up['time (days)'], up['limiting_magnitude'], s=100, marker='v', edgecolor = 'blue', linewidth = 1.5,color=band_colors[band], label=\"UL\" if band == \"lsstz\" else  \"_nolegend_\")\n",
    "\n",
    "# Plot the true data with a line to match the ToO simulations notebook formatting\n",
    "tt = np.linspace(0.1, 7, 100)  # Adjust the time range for better visibility\n",
    "parameters['output_format'] = 'magnitude'\n",
    "\n",
    "# Plot the light curve for the Kilonova object\n",
    "for band in band_colors.keys():\n",
    "    parameters['bands'] = band\n",
    "    out = redback.transient_models.kilonova_models.metzger_kilonova_model(tt, **parameters)\n",
    "    plt.plot(tt, out, color=band_colors[band], alpha=0.8, linewidth=2,label=f'{band[-1]}')\n",
    "\n",
    "# Add horizontal dashed lines for each filter with their fiveSigmaDepth magntidue\n",
    "plt.axhline(y=24.5, color='r', linestyle='dotted', linewidth=2)\n",
    "plt.axhline(y=23, color='yellow', linestyle='dotted', linewidth=2)\n",
    "plt.axhline(y=25, color='g', linestyle='dotted', linewidth=2)\n",
    "\n",
    "# Get current handles and labels\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# Filter out unwanted labels (like 'lssti', 'lsstg', 'lsstr'), these were appearing in the final legend when calling the plt.legend function\n",
    "filtered_handles = []\n",
    "filtered_labels = []\n",
    "for handle, label in zip(handles, labels):\n",
    "    if label not in ['lssti', 'lsstg', 'lsstr']:  # Exclude these labels\n",
    "        filtered_handles.append(handle)\n",
    "        filtered_labels.append(label)\n",
    "\n",
    "# Rebuild the legend with only the desired items\n",
    "plt.legend(filtered_handles, filtered_labels, loc='center left', bbox_to_anchor=(1.05, 0.5), fontsize=10, frameon=True)\n",
    "\n",
    "\n",
    "# Set labels and title of the plot\n",
    "plt.xlabel(\"Time from start of Rubin observing window (days)\", fontsize=13)\n",
    "plt.ylabel(\"AB Magnitude\", fontsize=13)\n",
    "plt.title(\"Gold Event: Metzger 2017 (350 Mpc)\", fontsize=15)\n",
    "\n",
    "# Add the annotation \"Faint KN\" at the top-right corner\n",
    "plt.text(6.5, 19, \"Faint KN\", fontsize=20, weight='bold', color='black', \n",
    "         horizontalalignment='right', verticalalignment='top')\n",
    "\n",
    "# Show the plot\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0ce7b770",
   "metadata": {},
   "source": [
    "Now we plot the light curve in the Rubin ToO format for 700 Mpc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e76e98a4",
   "metadata": {},
   "outputs": [],
   "source": [
    "kn_object = redback.transient.Kilonova.from_simulated_optical_data(name='my_kilonova', data_mode='magnitude')\n",
    "\n",
    "# Dictionary for filter colors based on the ToO simulations code\n",
    "band_colors = {'lsstg':'g', 'lsstu':'b', 'lsstr':'r',\n",
    "               'lsstz':'k', 'lssti':'yellow', 'lssty':'orange'}\n",
    "\n",
    "# Use the plot_data from the Redback Kilonova Object documentation\n",
    "ax = kn_object.plot_data(show=False, band_colors=band_colors)\n",
    "ax.set_ylim(27, 18)\n",
    "ax.set_xlim(0, 7)\n",
    "\n",
    "# Plot detected data points (detections) as markers\n",
    "for band in band_colors.keys():\n",
    "    detected_times = kn_sim.observations[kn_sim.observations['band'] == band]['time (days)']\n",
    "    detected_magnitudes = kn_sim.observations[kn_sim.observations['band'] == band]['magnitude']\n",
    "    plt.scatter(detected_times, detected_magnitudes, color=band_colors[band], marker='o', s=100, edgecolor = 'k', linewidth = 1.5,label=\"Det.\" if band == \"lsstz\" else  \"_nolegend_\")\n",
    "\n",
    "# Save the upper limits (non-detections) as a variable\n",
    "upper_limits = kn_sim.observations[kn_sim.observations['detected'] != 1.0]\n",
    "\n",
    "# Plot the upper limits (non-detections) for all bands using a loop\n",
    "for band in band_colors.keys():\n",
    "    up = upper_limits[upper_limits['band'] == band]\n",
    "    # Plot non-detections as downward-pointing triangles\n",
    "    plt.scatter(up['time (days)'], up['limiting_magnitude'], s=100, marker='v', edgecolor = 'blue', linewidth = 1.5,color=band_colors[band], label=\"UL\" if band == \"lsstz\" else  \"_nolegend_\")\n",
    "\n",
    "# Plot the true data with a line to match the ToO simulations notebook formatting\n",
    "tt = np.linspace(0.1, 7, 100)  # Adjust the time range for better visibility\n",
    "parameters['output_format'] = 'magnitude'\n",
    "\n",
    "# Plot the light curve for the Kilonova object\n",
    "for band in band_colors.keys():\n",
    "    parameters['bands'] = band\n",
    "    out = redback.transient_models.kilonova_models.metzger_kilonova_model(tt, **parameters)\n",
    "    plt.plot(tt, out, color=band_colors[band], alpha=0.8, linewidth=2,label=f'{band[-1]}')\n",
    "\n",
    "# Add horizontal dashed lines for each filter with their fiveSigmaDepth magntidue\n",
    "plt.axhline(y=24.5, color='r', linestyle='dotted', linewidth=2)\n",
    "plt.axhline(y=23, color='yellow', linestyle='dotted', linewidth=2)\n",
    "plt.axhline(y=25, color='g', linestyle='dotted', linewidth=2)\n",
    "\n",
    "# Get current handles and labels\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# Filter out unwanted labels (like 'lssti', 'lsstg', 'lsstr'), these were appearing in the final legend when calling the plt.legend function\n",
    "filtered_handles = []\n",
    "filtered_labels = []\n",
    "for handle, label in zip(handles, labels):\n",
    "    if label not in ['lssti', 'lsstg', 'lsstr']:  # Exclude these labels\n",
    "        filtered_handles.append(handle)\n",
    "        filtered_labels.append(label)\n",
    "\n",
    "# Rebuild the legend with only the desired items\n",
    "plt.legend(filtered_handles, filtered_labels, loc='center left', bbox_to_anchor=(1.05, 0.5), fontsize=10, frameon=True)\n",
    "\n",
    "\n",
    "# Set labels and title of the plot\n",
    "plt.xlabel(\"Time from start of Rubin observing window (days)\", fontsize=13)\n",
    "plt.ylabel(\"AB Magnitude\", fontsize=13)\n",
    "plt.title(\"Gold Event: Metzger 2017 (700 Mpc)\", fontsize=15)\n",
    "\n",
    "# Add the annotation \"Faint KN\" at the top-right corner\n",
    "plt.text(6.5, 19, \"Faint KN\", fontsize=20, weight='bold', color='black', \n",
    "         horizontalalignment='right', verticalalignment='top')\n",
    "\n",
    "# Show the plot\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f2582882",
   "metadata": {},
   "source": [
    "You can now use the simulated object and do parameter estimation. There are multiple examples available at \n",
    "https://github.com/nikhil-sarin/redback/tree/master/examples. Alongside other examples to simulate full survey or single lightcurves for Rubin or ZTF \n",
    "https://github.com/nikhil-sarin/redback/blob/master/examples/simulate_survey.py\n",
    "https://github.com/nikhil-sarin/redback/blob/master/examples/simulate_single_transient_in_rubin.py\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}