Climate/02_one_layer_atmosphere_model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b> One-Layer Atmosphere Model </b><br>\n",
    "Reference: Walter A. Robinson, Modeling Dynamic Climate Systems"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.style.use(\"seaborn-dark\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Step size\n",
    "dt = 0.01\n",
    "\n",
    "# Set up a 10 years simulation\n",
    "tmin = 0\n",
    "tmax = 10\n",
    "t = np.arange(tmin, tmax + dt, dt)\n",
    "n = len(t)\n",
    "\n",
    "# Seconds per year\n",
    "seconds_per_year = 365 * 24 * 60 * 60"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Albedo\n",
    "albedo = 0.3\n",
    "\n",
    "# Emissivity\n",
    "e = 0.85\n",
    "\n",
    "# Atmospheric Absorption\n",
    "atmospheric_abs = 0.1\n",
    "\n",
    "# Stefan-Boltzmann constant\n",
    "sigma = 5.6696e-8 # W/m^2*K^4\n",
    "\n",
    "# Solar constant \n",
    "solar_const = 1367 # W/m^2\n",
    "\n",
    "# Density of water\n",
    "water_density = 1000 # kg/m^3\n",
    "\n",
    "# Depth of the mixed layer\n",
    "depth_mixed_layer = 50 # m\n",
    "\n",
    "# Specific heat capacity of water\n",
    "spec_heat_water = 4218 # J/kg*K\n",
    "\n",
    "# Specific heat capacity of atmosphere\n",
    "spec_heat_atm = 1004 # J/kg*K\n",
    "\n",
    "# Gravity accelleration\n",
    "g = 9.81 # m/s^2\n",
    "\n",
    "# Atmospheric pressure\n",
    "atm_press = 1e5; # Pa (kg/m*s^2)\n",
    "\n",
    "# Atmospheric mass\n",
    "mass_atm = atm_press/g;  # kg"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Heat capacity water\n",
    "heat_capacity_water = water_density * depth_mixed_layer * spec_heat_water # J/K\n",
    "\n",
    "# Heat capacity atmosphere\n",
    "heat_capacity_atm = mass_atm * spec_heat_atm # J/K"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initialize temperature\n",
    "ts = np.zeros((n,)) # Surface temperature\n",
    "ts[0] = 273.15\n",
    "\n",
    "ta = np.zeros((n,)) # Atmospheric temperature\n",
    "ta[0] = 273.15"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Inflows (Solar to earth, IR down) - Surface\n",
    "\n",
    "# Absorbed solar energy\n",
    "solar = solar_const/4 * (1 - albedo) * seconds_per_year\n",
    "\n",
    "# Solar to earth\n",
    "solar_to_earth = solar * (1 - atmospheric_abs)\n",
    "\n",
    "# IR Down\n",
    "ir_down = sigma * e * ta[0] ** 4 * seconds_per_year\n",
    "\n",
    "# Outflows (IR) - Surface\n",
    "ir = sigma * ts[0] ** 4 * seconds_per_year"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Inflows (solar to atmosphere, IR) - Atmosphere\n",
    "\n",
    "# Solar to atmosphere\n",
    "solar_to_atm = solar * atmospheric_abs\n",
    "\n",
    "# Outflows (IR down, IR to space) - Atmosphere\n",
    "\n",
    "# IR to space\n",
    "ir_to_space = ir * (1 - e) + ir_down"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Flows of energy\n",
    "\n",
    "# IR norm\n",
    "ir_norm = np.zeros((n, ))\n",
    "ir_norm[0] = ir/seconds_per_year\n",
    "\n",
    "# IR down norm\n",
    "ir_down_norm = np.zeros((n, ))\n",
    "ir_down_norm[0] = ir_down/seconds_per_year\n",
    "\n",
    "# Solar to earth norm\n",
    "solar_to_earth_norm = np.zeros((n, ))\n",
    "solar_to_earth_norm[0] = solar_to_earth/seconds_per_year"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "for k in range(1, n):\n",
    "     \n",
    "    # Inflows (Solar to earth, IR down) - Surface\n",
    "    solar_to_earth = (solar_const/4 * (1 - albedo) * seconds_per_year) * (1 - atmospheric_abs)\n",
    "    ir_down = sigma * e * ta[k-1] ** 4 * seconds_per_year\n",
    "    \n",
    "    # Outflows (IR) - Surface\n",
    "    ir = sigma * ts[k-1] ** 4 * seconds_per_year\n",
    "    \n",
    "    # Calculate the temperature - Surface\n",
    "    energy_ts = ts[k-1] * heat_capacity_water + (solar_to_earth + ir_down - ir) * dt\n",
    "    ts[k] = energy_ts/heat_capacity_water\n",
    "    \n",
    "    # Inflows (solar to atmosphere, IR) - Atmosphere\n",
    "    solar_to_atm = (solar_const/4 * (1 - albedo) * seconds_per_year) * atmospheric_abs\n",
    "    \n",
    "    # Outflows (IR down, IR to space) - Atmosphere\n",
    "    ir_to_space = ir * (1 - e) + ir_down\n",
    "    \n",
    "    # Calculate the temperature - Atmosphere\n",
    "    energy_ta = ta[k-1] * heat_capacity_atm + (ir + solar_to_atm - ir_down - ir_to_space) * dt\n",
    "    ta[k] = energy_ta/heat_capacity_atm\n",
    "    \n",
    "    # Calculate IR norm, IR down, Solar to earth norm\n",
    "    ir_norm[k] = ir/seconds_per_year\n",
    "    ir_down_norm[k] = ir_down/seconds_per_year\n",
    "    solar_to_earth_norm[k] = solar_to_earth/seconds_per_year\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Convert to °C\n",
    "ts = ts - 273.15\n",
    "ta = ta - 273.15"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAGoCAYAAABL+58oAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAACGlUlEQVR4nOzdd3gUVdvA4d9sT7KpJCDdhF6kSw0gNkBBUUEQxMKHHRQsFOkCChYsqIAVX31FpNj1VcGCFBERFJSO9NASUjZty5zvj91sEkJJSNkkPPflXjN7ZubMM2cX98mZMzOaUkohhBBCCFGJGAIdgBBCCCFESZMERwghhBCVjiQ4QgghhKh0JMERQgghRKUjCY4QQgghKh1JcIQQQghR6UiCU8lt3ryZoUOH0rdvX/r06cPw4cPZtWtXketxOBwMGjSI66+/nm+//bYUIi05W7Zs4corryzSNgcPHmTkyJGlFFHhXHnllWzZsqVYdSxatIg33nij0OuXh+MuirS0NO64444y3++xY8cYNGhQkbbJ27aHDh2idevWpRFaiVm+fDn33XdfoMPwa9SoEUlJSaxcuZIZM2aUWL0TJ05k69atJVbfhUhKSqJRo0bnXW/cuHG8/fbbZRBR5WQKdACi9DidTu677z7eeecdmjVrBsBnn33GPffcw8qVKzEajYWua9u2bSQmJvL999+XVrgBdeTIEf79999Ah1Fst912W5HWr2jHnZKSUuwk8EJUq1aNjz76qEjbVLS2La+uuuoqrrrqqhKrb+3atQwcOLDE6hPllyQ4lVhmZiZpaWlkZGT4y2644Qbsdjsej4fff/+d6dOn8+WXXwKwfv16//u5c+eyefNmjh8/ToMGDdiyZQvHjh3jxhtvZPHixSxcuJAVK1aQnZ1NZmYmY8eO5ZprrsHtdvPcc8/x008/YTQaad26NVOmTMFisTBv3jy+++47dF2nZs2aTJkyhWrVquWLOSMjg6lTp7Jv3z5SUlIICQnh+eefJy4ujqFDh9KqVSv++OMPEhISaNu2LbNnz8ZgMPDhhx/y3nvvYbfbadiw4VnbZP78+QXivvLKK5k4cSLHjh3j//7v/5g2bRp33nknHTt2ZPPmzbjdbsaMGcPixYvZu3cvzZs3Z86cORgMBlasWMGrr76Kx+PBbrczfvx4WrRowZ49e5gwYQJOpxOlFP3792fIkCHMnTuXXbt2cfLkSRITE2ncuDEzZ87EbrcDsHjxYqZMmUJSUhI33ngjo0ePBuCHH35g3rx5uFwubDYbY8eOpXXr1vk+p0aNGlG3bl1OnTrF5MmT+ffff5k8eTJJSUkYDAYeeOABrrvuOn9beDyefMf99ttv88cff/D888+TmZmJpmmMHDmSHj16sHz5cr777juysrI4fPgw1atXZ8iQIXzwwQfs27ePu+++m2HDhrF8+XK++uordF3n2LFjVKtWjVmzZlGtWjXS0tKYOXMmO3fuxOVy0alTJ8aMGYPJZKJ58+ZcddVVbN++neeff54dO3awePFiXC4XKSkp3HPPPQwePJjx48eTlZXFjTfeyPLly2natCnr1q0jKioK8P7Vv27dOnbt2sXMmTMJDg4mIyODpUuXsnr16jO24dk+q7wOHTpE37592bRpE3PnzuXw4cOcOHGCw4cPExUVxYsvvpjvu3x6206bNg2Px8PkyZPZsmULqampjBkzhp49ewIU6t8GwGuvvcZXX32F0WgkNjaWSZMmERMTc85/G2f7TM9l8+bNPPfcczidTk6cOEHnzp15+umnmTdvHrt37+aFF14AYOPGjUyfPp1PP/30nN+dpUuXkpmZid1u5/3338+3rz///JMZM2aQmZmJ2WxmzJgxdOrUyb98+fLlfPvttyxYsIChQ4fSrFkzfv31VxITE7njjjtITEzkt99+IzMzk5deeolGjRqdNf4XX3yR48eP8/jjj/Pss88SFxd31u9kXuPGjcNqtbJlyxZOnjxJ7969iYqK4scff+TEiRPMmDGDTp06kZaWxrRp09i+fTuaptG1a1ceffRRTCYT3333HS+++CJBQUE0b948X/1Llixh0aJF6LpOREQEkyZNol69euf8jEQhKFGpvfPOO6pFixbqyiuvVI8//rhasmSJysjIUEop9euvv6rrr7/ev27e96+88orq2bOncrlcBZYdOnRIDR06VGVmZiqllPryyy9Vnz59lFJKvffee2rIkCEqMzNTeTwe9cgjj6hPPvlEffLJJ2rUqFH++j766CM1fPjwAvF+8803avr06f73kyZNUk899ZRSSqnbb79dPfzww8rj8ai0tDQVHx+v1q1bp/755x/VqVMndfz4cf82PXr0KFD3ueLOe3wHDx5UDRs2VCtWrFBKKTV58mTVo0cPlZaWprKyslSXLl3Uxo0b1e7du1Xnzp3VgQMHlFJKrV27VnXp0kWlpaWp8ePHqwULFiillDp+/LgaNWqU8ng86pVXXlHdunVTJ06cUB6PRz366KNq1qxZSimlevTo4T/W48ePq+bNm6sjR46of//9V/Xp00clJSUppZTauXOn6tKli0pPTy/wOb3yyitq2rRpSiml+vXrpz744AOllFJHjhxRV111lUpLS8vXJnmPOzk5WV177bXq4MGDSimljh49qrp166YOHz6sli1bptq2bauOHDmiPB6Puu6669TIkSOVx+NR27ZtU5dddpnyeDxq2bJlqlWrVmrv3r1KKaWee+45NXLkSKWUUuPGjVP/+c9/lFJKud1u9fjjj6s33nhDKaVUw4YN1SeffKKUUsrhcKhbb73Vf7ybNm1SrVq18n82OfM52yUmJhZ4/+uvv6rGjRurQ4cOKaXUOdvwbJ9VXnn3+8orr+Rry/vuu0+9/PLL6nRn+k7973//U0op9d1336mrrrpKKaUK/W9j6dKlauDAgSo9Pd0fx7Bhw5RSZ/+3ca7P9HTLli1T9957r1JKqdGjR6tff/3V/3l06NBBbdmyRZ08eVK1adNGnTp1Siml1BNPPKEWLVp03u/O5ZdfXuC7p5RSTqdTdenSRf34449KKaW2bNmi+vTpozwej/+zzBvX7bffrkaMGKGUUmrz5s2qYcOGauXKlUoppWbOnKkmTpx4zviV8v47++uvv5RS5/5O5jV27Fg1YMAA5XQ61fHjx1XDhg392y1cuFDdfffdSimlxowZo6ZPn650XVfZ2dlq2LBhasGCBerEiROqbdu2ateuXUoppebPn68aNmyolFJq/fr1avDgwf7/L//yyy+qd+/e/v2+9dZbBeIRhSM9OJXc3XffzYABA9iwYQMbNmzgzTff5M0332Tp0qXn3bZVq1YF/pIBqFmzJrNnz+aLL75g//79/Pnnn6SnpwPe7t8bb7wRm80GwEsvvQTAI488wpYtW7jlllsA0HWdzMzMAnX36tWL2rVr8/7777N//35+++23fGMXevTogcFgwG63U7duXVJSUvjnn3/o0qULMTExAAwcOJDVq1cXKe7Tmc1m/zieOnXq0Lp1a38vS9WqVUlJSWHbtm107NiR2rVrA9CpUyeioqLYunUr11xzDWPHjuWvv/6iU6dOTJw4EYPB4D/G6OhoAPr378/TTz/N2LFjAejTpw8AMTExREdHk5iYyJ9//snx48e56667/PFpmsaBAwfO+jklJyezfft2BgwYAED16tVZsWLFGY81x+bNmzlx4gQPPfRQvv3s2LEDgMsuu4zq1asDUKtWLeLj4zEYDNSuXdvfIwbQpUsXYmNjAbj11lu58cYbAfjpp5/YsmWL/7uXlZWVb//t2rUDICQkhPnz5/Pzzz+zb98+tm/fnq8XsrCqV69OzZo1AVizZs1Z2/Bcn9XZtG/f3v99aNq0KSkpKeeNx2w2+3tsGjduTGJiIgA//vhjof5trFq1iptvvpng4GAA7rjjDubPn4/T6QTO/G/jXJ9pjRo1zhrrrFmzWLVqFfPnz2fv3r1kZWWRkZFBlSpVuOKKK/jss8/o168fq1evZsqUKfz+++/n/O40atTI31557dy5E4PBwBVXXAFA8+bN+eKLL87Zjtdccw2A/99d165dAe+/099+++2c8Z/ufN/JvHr06IHZbCYmJobg4OB8+01OTga8n9GiRYvQNA2LxcKgQYN47733qFu3Lg0bNqR+/fqA9/9Rc+bM8cewf//+fGO8UlJS/HWKCycJTiW2ceNGNm3axPDhw+nRowc9evTg0UcfpW/fvqxZs4aoqChUnkeRuVyufNvn/I/0dH///TcPPvggd911F126dOHyyy9n2rRpAAV+aE+ePImu6+i6zvDhwxk8eDDgHR90ph+FDz/8kI8//pghQ4bQt29fIiIiOHTokH95TuIE3v+BKqX80xxnG1t0rrhPZzab0TQt3/vTqTM8xk0phdvtpkePHnz77besXbuWdevW8dprr/nHcOSNT9f1fD+medsv57h0XadTp07+ZBEgISGBqlWr8v3335/xc8qpJ+8x7N27lxo1auRrw7w8Hg/16tVjyZIl/rJjx44RFRXFF198gcViOeM+Tnf68eW813Wdl19+2d/1npqami++nOM4evQoAwcO5NZbb6Vt27b06tWLH3/88Yz7yivnh/70+nL2fbY2bNy48Rk/qzp16px1X2f6Hp5P3u9Q3uMu7L+N0/eh6zput/ucMZ3rM50wYYJ/sO2gQYOwWq3+dYYMGULjxo3p2rUrvXv35s8///Tvf8iQIUydOhWTycS1115LSEjIeb87Z/t/idFozNcW4E164uLizrg+UOB7eKZ/m+eKP6/zfSfPtd8zff91XS/w3u12F/iO5N1W13VuvPFGnnjiCf/748ePEx4efsY4ROHJVVSVWFRUFPPmzeP333/3l504cYLMzEwaNmxIVFQUR44cITExEaXUef/Cz7FhwwaaN2/O3XffTfv27Vm5ciUejwfw9mJ8+eWXOJ1OdF1n6tSpfPXVV8THx7N06VIcDgcAL7/8MmPGjClQ9+rVq7npppsYMGAAsbGx/PDDD/66z6Zz586sWbOGo0ePAvDJJ58UOW6j0VggwTufjh07smbNGg4ePAjAunXrSEhIoGXLljz22GN8/fXXXH/99UyZMgW73U5CQgIAK1euJC0tDV3X+fjjj887HiJnP3v27AHg559/5oYbbiA7O/us29jtdpo1a8ann34KeH/Mb7vtNtLS0vKtl/e4W7Vqxf79+9mwYQPgHVjes2dPjh8/XqR2+fXXXzl27BgAH330kf/44uPjWbhwIUopnE4nDzzwAB988EGB7bdu3UpUVBQPPvggXbt29Sc3Ho8Hk8mEx+Px/1hERUX5Bx2fawD8udrwXJ9VcRT2O1XYfxvx8fEsX77c3xPx/vvvc/nllxf44c3rXJ/pzJkz+eyzz/jss8/yDU5PSUlh69atPP7441x77bUcO3aMAwcO+H+827Rpg8Fg4O233/Zvd6Hfnbi4ODRNY82aNYD3j5A777yzQKJQFOeL32g0+hPDwn4nCys+Pp7//ve//vo+/vhjOnfuTLt27di9ezfbt28HvOOKcnTp0oWvvvrK31aLFi3izjvvvOAYRC7pwanEYmNjee2113jxxRc5evQoVquV0NBQnnrqKf9fSIMGDeKWW24hJibG3018Pn369OG7777juuuuw2w206lTJ1JSUvyXkh8+fJibb74ZpRTt27dn6NChGAwGjh07xq233oqmaVSvXp1Zs2YVqHvYsGFMnjyZ5cuXYzQaadasGTt37jxnPI0aNeKJJ57gzjvvJCQkhBYtWhQ57gYNGmA0Gunfvz8vvvhiodqhfv36TJkyhREjRuDxeLDZbMyfP5/Q0FAefPBBJkyYwOLFizEajVx99dW0b9+e3377jejoaO655x5OnTrF5Zdfzv3333/O/TRo0ICnnnqKRx99FKUUJpOJefPmnfWv4hwvvPAC06ZN4/3330fTNGbOnOk/jZe37pzjXrJkCa+88grPPvss2dnZKKV49tln/ad5CqtatWo88cQTnDhxgvr16/PUU08BMGHCBGbOnEnfvn1xuVx07tyZ4cOHF9i+S5cuLF26lF69ehEUFESLFi2Iiopi//791K1bl6ZNm9K7d28WLVrExIkTeeqppwgLC6Nz584Fjq8wbXi2z6q4CvudGjBgQKH+bfTv35+EhAQGDBiAruvUrVuX559//pwxREVFFfkzDQ8P59577+Wmm24iIiKCyMhI2rRpw/79+/2Df2+++Wa+/vpr/6XOF7If8PaKzJ07l6effppnn30Ws9nM3Llzz5m0nc/54r/66qsZPXo0M2bMKPR3srAmTpzIjBkz/PV17dqV+++/H4vFwvPPP8/jjz+O2Wzm8ssv92/TtWtX7rnnHoYNG4amadjtdl599dWz9iSJwtNUYfpWhRAlYu7cuf6rnCqjvFe8iMrJ7XYzYsQIbrjhhnxX5QlR3sgpKiGEEIWye/duOnXqREhICL169Qp0OEKck/TgCCGEEKLSkR4cIYQQQlQ6kuAIIYQQotKp0FdRnTiRdv6VisFut+JwnP1SXFHypM0DQ9q97EmbB4a0e9kr7TaPiQk9Y7n04JyDyVT4h1GKkiFtHhjS7mVP2jwwpN3LXqDaXBIcIYQQQlQ6kuAIIYQQotKRBEcIIYQQlY4kOEIIIYSodCr0VVRCCCGEKF+UUqAU+B5wqorx8NTikARHCCFEQCmlwOPxv5R/3u39ccxX7s7zXgfdA7oC3ZNbT86Pq657t/e/PKQFmclKy/Rv458q5a1f173be06vz7uu0j156jvDPnK2zVemg8c79ZcrvOU5yYBS3mU55TlJgZ67jjdx8L3XFeCb+spVnvm89eauq+cmHznr6apgDKfVmxtTThyn139arKdJb9MW+8vzyuS7lJckOEIIUUEojwfcbpTHDS6Xd97tzp26XN6k4FzLfPP5lrnd4HblX+bx+Ovx7tOXdOh5E5CcRMOXdOg6yn2G8pxleRKUfHWU4V/4JXL3NE0DgyHfSzMYQDOA0TfNKTOcoUwzgEHz1uN/GUDT0Ay58xg0QMuzrgHNYMzdxuAbZeKr07uK4Qz1+tYhp34tNwZ89Wjkj8FfflqcBoP3SeenxZ03RjTQcso1CG91Ga6SaPcikgRHCCHOQSkFTicqOwuVnQ3Z2aisLO+PvzPbN3V6py4nKtuJcjm9iYIzZ5rtf5+zXRY6zvTM3O2cTpTTN59Tp9Ppe+9NQEo9ETAawWhCM5nAbPLOm81gNHrLDIbcdYxG37wRzWQGqw2MRgx5y41GMJz23j9v8v7w51tmyrOdIc96Ofvxbe+rU8uJJ+fH3mj0/vgajN4fXN9Uy/M+LCKYVIczN/nIk5j4tzXmLTMUTGY0rXQ/h0omJCKY5OSMMt+vJDhCiArJm3hkozIyUVmZqMxMVGaGd5qVlZuIOLO9iUmWN0FR2b5lvoRFZWWDM2eZdxuceeedZ+x2LxKDASwWNLMFzGY0iwXdZkM3GtEsvjJbEFpoGJrF4k0qzBbvMovFl3CYvYmEf96XDJi9Zf5lpjxJiX+ZCXzLvctM3mThtG21nB6BSswSEYwpAD+2ouxJgiOEKBPK6USlp6PSHaj0dPSMdFRGBmRm5ElSfAlK5ukJS6Z3ncwMyMxNaIrco2E0olltYLWi2axoljzzwSEYIqPAakOzWn0v33KrFc3mLcfie5+TfFgsvmQkNynRzBawmHMTGlPB/9VGBOivWiEuFpLgCCHOSem6NzFJS/UnKHp6er5kRaWnozLSUQ6Hd5qejp5nXqU7vGNACsNsRgsKQgsK9vZqBAehBQVhqFoVbEG+ZUG56wQF5VkvGM1mA5sNLSdRseUkKbYzJhpCiMpJ/rULcRFQSnl7Q9LSUGmp6KmpufNpaXicGWSeSEL3leWUq7RUlMNx/p4SoxEtJAQtxO6bhmCIicEQcilacAia3Z5/GuKbDwryvrfZcpMVSUKEECVA/k8iRAWksrPQk5NRycnoyafQU1JQKd557zQ5d3lKMiot1Xu1yllkGI1o9lC00FAMoWFo4RGYa9X2jgkJDUULDcMQ6l3uTVDsGHyJjBZi9/aQyMBLIUQ5IgmOEOWA8nhQqSnoiYnoSb6Xbz4nifEnLinJkJV15ooMBrTwCAzh4WjhERgvjcUcEY4WGu5LXkLRwsJzE5nQULSwMCKrR5OSklmmxyyEEKVJEhwhSolSyjsWJfEkKimpQOKS971KPnXm00A2G4bISAwRkWiRUZhj43wJTASGiAjvfEQkWoTvvT30gq6Ekd4XIURlIwmOEBdA6bq3R+X4cTwnT6AfP4Z+4jj6iRPe6ckT6MePe6/6OZ3JhCEqCkNUFQwxVTE1aoKhShUMkb6yKlW806gqaMHBZX9wQghRCUiCI8QZ6A4H+tEEPEcT0I8moB87iufECfQTvkTm5MmCVwUZjRiqRGOIicEYG4e5fUcMMVUxREf7ExZDlSrecS3SYyKEEKVKEhxx0VFKece7HE3Ac/RovkTGOz2Kcpx2Q3eLFUPVqhhjqmJu3hJD1are5CUmBkNMNYwxVdEiI713YBVCCBFwkuCISkl5PN5elyOH8Rw6iH74EJ6cV8IRyMw/oFYLCsZQvTqGS6pjbtEKwyXVMV5yCYZLamC85BK0iEjpdRFCiAqkXCU4f/75J88//zzvv/8++/fvZ9y4cWiaRoMGDZgyZQqGi+A24qLwlMuFJ+FIbvJy6CCew4fRDx/0JjFud+7KFgvGGjUx1qyFue3lGKtdgqF6dYyX1MBwySVy2kgIISqZcpPgvPnmm3z++ecEBQUB8MwzzzBq1Cg6dOjA5MmTWblyJddcc02AoxSBoKel4tm/L/d1YD/u/fvQE47ku7eLFhSMoVYtjPXqY+nWA2PNmhhq1cZYsxaG6JiL4jk7QgghvMpNglOnTh3mzp3LmDFjAPj7779p3749AN26dWPNmjWS4FRiStfRjx0lfetRMv/e7k9iPAf2o04l5a5oNmOsXQdT/QYYr7waY+06GGvWxlirlpxGEkII4VduEpyePXty6NAh/3ullP/HKiQkhLS0tALb2O1WTKbSG9RpNBqIiJDLdEua59Qpsnftwul7Ze/ahXP3Lu+DF30M4eFY4uII7nEF5thYLLGxmGPjMNesKQN5S4F818uetHlgSLuXvUC1eblJcE6Xd7xNeno6YWFhBdZxOLJLNQZ52m/xqOwsPPv+xb1nD+69u/Hs3YN7z25UUqJ/HS0sDGNcfay9rscYV4+I5k3IjLoEQ2RkvrpcvhdppfuZX6zku172pM0DQ9q97JV2m8fEhJ6xvNwmOE2bNmX9+vV06NCBVatW0bFjx0CHJM5BT3fg2bUT984duHftwL1zB579+3LHyFismGJjsbTviLFefUxx9TDF1UerUiXfaaWgiGCy5X8+QgghiqncJjhjx45l0qRJzJkzh7i4OHr27BnokISPnpLsTWT8r+3ohw76lxuqRGNs2AhL1ysw1W+AKa4ehpq15NSSEEKIMqMppVSgg7hQJ04UHJdTkqQrE1RWFu6d23H98zfuf7bi3vYP+tEE/3JD9eqYGjbG1LCRd9qgIYYq0Re8P2nzwJB2L3vS5oEh7V725BSVCDil696rl3zJjOufv/Hs3e0/zWSoXh1T02aYbh6AqUFDTA0bYQgLD3DUQgghREGS4FzE9Ix03Fu34PrrT9x/b8G9/R+UwwGAFhKCqUkzLEPuwNS0OeYmTTFEVQlwxEIIIUThSIJzEdGTEnH99Seuvzbj+nMznt07QdfBYMBYrz7Wq67B1KQ5pmbNMdapKzfGE0IIUWFJglNJKaXQDx/C9edmb0Lz1+bcgcBWK+ZmzbEMvRtzy1aYmjXHEBwS2ICFEEKIEiQJTiXiOX4M1x+/49r4O64/fkc/fgzw3mvG3KIVphtuwtyiJaaGjdHM5gBHK4QQQpQeSXAqMP3UKVybN+La+DvOjRv8PTRaeDjm1m0x334n5lZtMNa9VE43CSGEuKhIglOBqOxsXH9txrl+Ha6NG/Ds3gWAFhyCqVVrgvrdgrltO4xx9SWhEUIIcVGTBKec8xw6iHP9Opy/rsO16XfIzgaLBfNlLbHe8wDmtu0wNWqCZpKPUgghhMghv4rljMrMxLVpozepWb8O/bD3AaSGWrWx9bkRS8fOmFu1QbPZAhypEEIIUX5JglMOeI4dw7lmFc7Vq3Bt/gNcLrDZMLdpR9Ctt2Hp0AljzVqBDlMIIYSoMCTBCQClFJ5dO3GuWUX2L6vw7NoBgLF2HWw3D/D20rRohWaxBDhSIYQQomKSBKeMKJcL1+Y/cK729tTox4+BpmFqdhnB94/A0rU7pjp1Ax2mEEIIUSlIglOKlNOJc8N6nD+uwLnmF+9jEKxWLJd3wDLsHiyd4zFERgU6TCGEEKLSkQSnhPmTmh9W4FyzCpWejmYPxdL1CizdrsDSrr0MEBZCCCFKmSQ4JUC53bg2rCd7xXf5k5ruV2LtcRXmtpfLnYOFEEKIMlSsBGfnzp389ttvJCcnExUVRadOnYiNjS2p2Mo9966dZP3vK7JXfItKSkILDZOkRgghhCgHLijB2bNnD7Nnz8Zms9GwYUOqVq1KSkoKc+bMwe128+ijj9KgQYOSjrVcUJmZZH37NVmfLsOzZzeYTFg6x2PteR2WTl0kqRFCCCHKgQtKcL7++mteeOEFQkNDCyxLSUlh4cKFPPLII8UOrjzRU5LJXPQBWZ99gnKkYWzYmJBHx2C98moM4RGBDk8IIYQQeWhKKVXcSjweD0ajsSTiKZITJ9JKtf6IiGBOJaaR+fEiMv/zDiojA0v3HgT1H4TpshZomlaq+78YRUQEk5ycEegwLjrS7mVP2jwwpN3LXmm3eUxMwc4WgAt6IuPRo0cZNGgQKSkpgLdH59Zbb+XYsWMXHmE55D55kpRRD5Hx+iuYLmtJxML/EvbUM5hbtJTkRgghhCjHLugU1ZQpUxg+fDjh4eEA9O3bF5PJxJQpU5g/f36JBKbrOlOnTmXHjh1YLBZmzJhB3bpldyM8/VQShx9+APfRBOwTpmDteZ0kNUIIIUQFcUE9OOnp6Vx99dX5ynr37u3v0SkJK1aswOl0snjxYh577DFmzZpVYnWfj9J1/nzmBdwJRwh//mVsva6X5EYIIYSoQC6oB+dsw3ZKYDiP38aNG+natSsArVq1YuvWrQXWsdutmEwlP/Zn69a93FutJ/95OJ563buUeP3i7IxGAxERwYEO46Ij7V72pM0DQ9q97AWqzS8owWnRogX/+c9/uOOOO/xl77//Po0aNSqxwBwOB3a73f/eaDTidrsxmUx51skusf3llWz0fhDZHbvJYLQyJgMAA0PavexJmweGtHvZC9Qg4wtKcEaPHs3MmTOJj4+natWqpKamEh8fz/jx44sVZF52u5309HT/e13X8yU3pcls8J65c7n1MtmfEEIIIUrWBWUMFouFadOmMXHiRJKTk4mMjCzx5KNNmzb8+OOPXHfddWzevJmGDRuWaP3nYjZ5x9u4PJLgCCGEEBXRBWUla9euZdWqVXTv3p1OnTqVdEwAXHPNNaxZs4ZBgwahlOLpp58ulf2cSU4PjlMSHCGEEKJCuqAE59VXX+WDDz5g6NChpZbgGAwGnnrqqVKp+3zMxpwenJIbNC2EEEKIsnNBl4n36NGDWbNmccUVV5RwOOWD2ejrwZExOEIIIUSFdEE9OPfcc09Jx1GuWIxyikoIIYSoyC6oB2fKlCns3LnzjMu2bdvG5MmTixVUoOWcopIeHCGEEKJiuuDLxF966SW2bt1KbGws0dHRpKamsn37di677DJGjRpVwmGWLZOvB0euohJCCCEqpgtKcCIiIpg6dSoOh4M///yTU6dOUaVKFSZMmEBwcMW/Q6TJoGHQ5BSVEEIIUVEV6+Y1drudLl0q56MMzEaDXEUlhBBCVFAXNAbnYmA2anKKSgghhKigJME5C7PBIIOMhRBCiAqqWKeojh07xnPPPUdSUhK9evWiUaNGtGzZsqRiCyjpwRFCCCEqrmL14EyaNIlbbrkFl8tFu3btmDlzZknFFXBmo/TgCCGEEBVVsRKcrKwsOnXqhKZpxMXFYbVaSyqugLPIIGMhhBCiwipWgmO1Wvnll1/QdZ3NmzdjsVhKKq6AMxk1uUxcCCGEqKCKleBMnz6d5cuXc+rUKd555x2mTp1aQmEFnvcycUlwhBBCiIqoWIOMFy5cyIsvvlhSsZQrFqMmY3CEEEKUC3PnvsiOHdtISkokKyuLGjVqEhERyYwZs4tdd2pqKqNGPUBYWDgvvfR6CURbPhQrwdm9ezepqamEhYWVVDzlhkl6cIQQQpQTI0eOBuDrr79g//59PPDAyBKre+/e3VSvXoOZM58rsTrLg2IlOHv27KFDhw5ERUWhad4HVK5evbpEAgs0i1HD4ZQERwghRH5Z//uKrK++KNE6bdf3xdbr+mLV8fPPP/DBB+9hMpmIjo5h2rSneffdN6lSpQr9+vVn//59PPfc07z66hsMHXortWvXBeDIkcOcPHmCt99eQI8eVzF37ovouk5ycjKPPz6Oyy5ryZdffsonnyxD1z3Ex3fn//7vPn74YQWLF/8Xg8FAixatSjTpKgnFSnB+/PHHkoqj3PHe6M8d6DCEEEKIQvn++28ZPHgoPXpczTfffEl6evpZ183MzOSuu/6Phg0b88cfv/PZZ8v4v/+7j5Urv2PEiNHUq1ef7777H19//QW1atXmgw/e4733FmGxWJk//1WOHj3KO+8s4K233sdmszF9+iQ2bPiVyy/vWIZHfG7FSnDGjx9foOyZZ54pTpXlhtlokKuohBBCFGDrdX2xe1tKw8iRo3n//YUsW/YxdeteSrduV+RbrlT+W5/UqXNpgTqio6uycOFbWK1WMjIyCAkJ4fDhw8TG1sNqtQHwwAMj+eefrSQnn+Lxxx8GICMjg8OHD3H55aVyaBekWAnOddddB3gb7Z9//uH48eMlElR5YJbLxIUQQlQgn3/+Cf/3f/cSGRnFs8/OZNWqn7BYrCQmJgKwc+f2fOvnDC3J6+WXn2Py5Blcemksb7+9gISEI9SsWYsDB/bhdDqxWCxMnDiGESNGU7VqNV566XVMJhNff/0FDRo0LJPjLKxiJThdu3b1z3fr1o1hw4YVO6DyQh7VIIQQoiJp0qQZY8aMIjg4hKCgIDp3jic9PZ3Jk8ezadNGGjVqct46rr22N5MmjSU0NIyYmKqkpCQTGRnJkCF3MmLEvWiaRpcuXbnkkuoMHDiEESPuxePxUL16Da688poyOMrC09TpfVZFkHdA8YkTJ3j33Xf5/PPPL6iu77//nv/973+88MILAGzevJmZM2diNBqJj49nxIgRBbY5cSLtwgIvhFkrdvHDrpN890CnUtuHKCgiIpjk5IxAh3HRkXYve9LmgSHtXvZKu81jYkLPWF6sHpyvvvrKP2+xWC54/M2MGTNYvXo1TZrkZpdTpkxh7ty51K5dm3vvvZd//vmHpk2bFifcIjHLoxqEEEKICqtYCU6bNm0YMGCA//1//vMfmjVrdkH1XH311SxevBgAh8OB0+mkTp06AMTHx7N27doCCY7dbsVkMhbjCM4uNNiC06MTERFcKvWLMzMaDdLmASDtXvakzQND2r3sBarNLyjB+fLLL/nhhx9Yv349v/76KwC6rrNz507uuOOOs263ZMkS3nvvvXxlTz/9NNdddx3r16/3lzkcDux2u/99SEgIBw8eLFCfw5F9IeEXisftweXROXUq/YwDsUTpkO7jwJB2L3vS5oEh7V72KtQpqq5duxITE0NycjIDBw4EwGAwULt27XNuN2DAgHw9Pmdjt9vzXb+fnp5e5ndLthg1lAKPApPkN0IIIUSFckEP2wwPD6dDhw688847XHrppdSqVYvq1atz5MiREgnKbrdjNps5cOAASilWr15Nu3btSqTuwjIbvE0jV1IJIYQQFU+xxuA8+eSTbN68mczMTDIzM6lTpw4ff/xxiQQ2bdo0Hn/8cTweD/Hx8bRs2bJE6i0ssyk3wQkyl844HyGEEEKUjmIlONu3b+err75i8uTJjB49mkceeeSC6+rQoQMdOnTwv2/VqlWJJUsXwmzwnpdyypVUQgghAqw4TxP/+ecfadasOdHRMaUe59tvL/A/+yrQipXgREZGomkaGRkZREVFlVRM5YLZ6E1w3HKKSgghRIAV52niS5Ys4tJLnyyTBKc8KVaC06xZM95++22qVq3K6NGjycrKKqm4As5s9J6ikh4cIYQQeX319zE+33q0ROu8ofklXN+sWrHq2Lt3d4EngaelpbF7905mzJjMpEnTmT59MtWqVSMhIYGrrrqWf//dw86dO+jcOZ777nuInTu38+KLz2E0GrFYLIwZM5HIyEgmTx5Heno6WVlZ3Hvvg7Rv35EBA26kadNmHDlyiNjYeowbNwmAX35ZxY8/riQlJYXhw++nT59eZ3zy+NtvL2Dr1r/IzMxk3LhJ/P77er7//ls0TeOqq65lwIBBxWqPYiU4/fr1o2rVqthsNlatWkWLFi2KFUx5kpPgyCBjIYQQFcG//+4t8CTwsWMnUr9+Q5544knMZjMJCYd58cXXyM7OYsCAG/n006+xWm3079+X++57iNmzZzJu3EQaNGjEL7/8xKuvzmHYsPtISUnhhRde4dSpUxw8uB+AEyeOcc89r1KrVm0mTRrHL7/8BEBMTAzjxk3ijz9+58MP/0PXrh3P+ORxgLp1Yxk16nH+/XcvK1d+z+uvvwXA6NEP0aFDxzM+ELSwipXgTJgwgUWLFgFw5ZVXFqeqcsfiO0UlCY4QQoi8rm9Wrdi9LaXhTE8CP1316jX9VypHRUURFhYO5D548+TJEzRo0AiAli3bMH/+q8TF1ePGG29m6tQJuN1u+vf39qxUq3YJtWp5bw9z2WUtOHDAm/jkPPOqSpVosrKyOHDgwBmfPA5Qp05dAPbu3cOxY0d55JEHAEhLS+PgwYOBS3CCg4N5+umniY2NxeC7rDrnvjgVnc13h+QslyQ4Qgghyr8zPQkcvPep03Xvb9n5blwbHR3D7t27qF+/AZs3/0Ht2nXYs2c3GRnpPPfcy5w8eZIHHhhGly5dOXHiBImJJ6lSJZq//vqTXr2uY8eO7Zy+i5o1a53xyeOrVv2EwXdBT506dbn00jheeOEVNE1j8eL/Uq9eg2K1R7ESnNatWwP4H8VemdjM3oQt0+0JcCRCCCHE+Z3pSeAAzZu3YMaMKYwZM+G8dYwdO4EXX3wWpRRGo5Fx4yYRHR3Du+++wQ8/rEDXdf7v/+4DwGIx8+KLz3Ls2DGaNbuMLl26sWPH9gJ1RkVFnffJ4w0aNKRdu8t58MH/w+l00aRJM2JiijcoulhPEwdYu3YtBw8epGXLlsTGxmK1WosVUFGU5tPEdxx3cPv7f/DsDU3p0SC61PYj8pPbqAeGtHvZkzYPDGn3knPDDT35/PNvz7tehXpUQ445c+Zw9OhR9uzZg8Vi4Y033mDOnDnFqbLcyLm5X6ZLenCEEEKIiuaCHtWQY+PGjTz77LMEBwdz0003cejQoZKKK+CCfKeosiTBEUIIIQooTO9NIBUrwfF4PGRnZ6NpGh6Pxz/QuDLwDzJ2yyBjIYQQoqIp1imqO++8k5tvvpmkpCQGDBjAXXfdVUJhBV5OD46cohJCCCEqnmIlOL1796Zz584cOHCAWrVqERkZWVJxBZzJaMBs1MiUy8SFEEKICqdYCc4ff/zBtGnTSExMpGrVqsycOZMmTZqUVGwBZzMbZQyOEEIIUQEVK8GZMWMGL7zwAvXr12fnzp1MnjyZjz76qKRiC7hgs1Fu9CeEEEJUQMUaFRwaGkr9+vUBaNiwITabrUSCKi+CLEYZgyOEEEJUQMXqwalSpQoTJkygY8eO/P333+i6zuLFi4HK8cgGm9koV1EJIYQQFVCxEpy4uDgA9u/fj91up3379pw4caJEAisPgqUHRwghhKiQipXgPPDAA+zatQun0+kva9GiRbGDKi+CzEZSM5znX1EIIYQQ5UqxEpx7770Xp9NJWFgY4H1K6auvvloigZUHQWYjx+UUlRBCCFHhFCvByc7O5oMPPiipWModGWQshBBCVEzFSnDatWvHL7/8Qr169fxlNWrUKFIdaWlpPPHEEzgcDlwuF+PGjaN169Zs3ryZmTNnYjQaiY+PZ8SIEcUJ9YIEW+QycSGEEKIiKlaCk5iYyNNPP53vFFVR74Pz7rvv0rFjR+666y727t3LY489xieffMKUKVOYO3cutWvX5t577+Wff/6hadOmxQm3yGxm6cERQgghKqJiJTh79+7lm2++KVYAd911FxaLBfA+vNNqteJwOHA6ndSpUweA+Ph41q5dWyDBsdutmHwPxSwNIVYTWW6d8PAgNE0rtf2IXEajgYiI4ECHcdGRdi970uaBIe1e9gLV5sVKcBo1asTmzZvzJR45ycqZLFmyhPfeey9f2dNPP02LFi04ceIETzzxBE8++SQOhwO73e5fJyQkhIMHDxaoz+HILk7452UzGfDoipNJ6ZiNledJ6eVZREQwyckZgQ7joiPtXvakzQND2r3slXabx8SEnrG8WAnOhg0b+Omnn/zvNU1j5cqVZ11/wIABDBgwoED5jh07ePTRRxkzZgzt27fH4XCQnp7uX56enu4/DVaWbGZv71CmyyMJjhBCCFGBFCvB+eKLLwA4deoUERERF3QaZ/fu3TzyyCO89NJLNG7cGAC73Y7ZbObAgQPUrl2b1atXB2yQMUCmSyescj2FQgghhKjUit2DM23aNDweD7169aJGjRpn7KE5lxdeeAGn08nMmTMBb3Izb948pk2bxuOPP47H4yE+Pp6WLVsWJ9QLYrd6myfDKQONhRBCiIqkWAnOSy+9xAcffMDIkSO5//77ue2224qc4MybN++M5a1ateLjjz8uTnjFFmbzNk9qliugcQghhBCiaIo1sMRgMPhPTVmtVkJCQkoqrnIhLMgMQFq2O8CRCCGEEKIoipXg1KlThxdeeIHk5GTeeOONIt/kr7wLs0mCI4QQQlREF5TgjBo1CoBp06ZRo0YN2rZtS1BQENOnTy/J2ALO34OTJQmOEEIIUZFc0BicpKQk78YmE7fddluJBlSe5I7BkQRHCCGEqEguKME5ePAgc+bMOeOyRx99tFgBlSdmo4Egs0FOUQkhhBAVzAUlODabjdjY2JKOpVwKtZrkFJUQQghRwVxQghMdHc1NN91U0rGUS2E2s/TgCCGEEBXMBQ0ybt68eUnHUW6F2kyS4AghhBAVzAUlOGPHji3pOMqtUKtJBhkLIYQQFYw8QfI8Qm0yBkcIIYSoaIr1qIaLQZhVTlEJIYQoSCmFAnQFKIWuQPeVqZx5BQrvMuVbR5Fn3l+H8m3jK1Og+7YjT106vjrPUJeuFJB3v7n157zHF0/ucm/h2dZXvjr9733rkvcYz1SvbwUFtLq0CvXCraX9cRQgCc55hFpNpDs9eHSF0VD0p6ULIUSgKaXw6Aq3rvD45nNe7pzy05a5T5vmLHN78sznqcdb5v3B03N+cH3lulLoeZZ58izL+ZH06N5tvGW56+i+ehW56+Ssr3zre+v3bqNOqz9vLLpSKE3D49HzJRy5iUHexCE3Nt33g54zn5OAqAB+phVJ8xrHefe2VmW+X0lwziPvAzcjgy0BjkYIEUhuXeF0696XR8el67g83h99l677pzllbt98bpmOW1e4PAqTxURqerY3wfD41stZnqfMu75eoMyt6wUSFY8Ct0cvkMR4ytEvsUEDg6ZhNGj+ee8LX1mecoOGMe86Bt+2vmX+ujQwGDRMebYtUJemYTSAzWrG7fagAUZNQ9NA862j4ZtqoPm21dB8y/OW5Zk/fbs88/nq8m2Try7Ity+Dpvna6PS6c+fzreury0D+48BXb0455N3eGyOnvdfy1KHhLcxdrvm2yT12/7I8bZTTBXD6vmpXCyMzPbt0v1hnIAnOeVQJ8SY1iemS4AgRaEp5k4Mst4csl06WWyfL5SHLrZPtK3N6dLLdOi6PTrZbeRMRX1lOYuL05Mx7E5Zs3zpOd862imzfe/+2vr/6S5oGWEwGTAYNk0HDbDRgNnqnRoOGOU+ZyaBhMxv965oM3h/ynKnR4P3x95YZfFPyLTP56jVq+ZflrJ+7zJsQ5Cv3JREF9mnITTJyE5A8iUWeZYEWERFMcnJGoMO4qFjNRjIDsF9JcM4jxu5Nak6kZ1M/pnI9LV2IkqaUIsutk+H0kOnykOH0vVy57zNdvmTErZPlyk1Mstwesn1lOQlMtlvPk8x4p8XJMawmAxZfspAzb8kzDbYYiQgyYzUZMPvKrEbvvNWk5a5rNGA2GbAYNcwGAyajN0HwJia5Zf6p0eBLXnLLoqNCSE/LklPfQpQSSXDOI6cH56TDGeBIhCgd2W6dk45sjpzKxJHt9r6cHtKz3XmSEg/p/qRF907zLMs7LUoCYjZq2ExGbGYDNpMBa575aLvFW2Y2+pYZsPnmbWaj932+Mu/2Vn/SkpuQWH09JFo56EHIYTMbyZLkRohSIwnOeUTnJDjpkuCI8setKxxZblKyXKRlu0nLdpOe7fEnKacnLP5pthtHtgeH042rEAM0jBoEWYwEm40EW4wE+aYxdou/zF9uNuZbN+/7IIuBILMRm8mboEjvhRCitEiCcx42sxG71UiiJDiiFGW7dVKzXKRmuX0vFyl55lOz3KRkuknL9s37yh3ZnnPWqwHBFiN2qwm71YjdYiIq2EKdyCDsVhMhvmUxEUEYPHq+9UKsRkLMJoIsRizG8tX7IYQQ5yMJTiHEhFg5IaeoRCEppUh3ejiV4eJUpotTGS6SM50kZbhI9r0/lekiOWea6SLbrZ+1PoPmfSZamM1EuM2boNSNCibcZiLMZvIvC7OZsFtMuUmK1USwxViogZ0y8FIIUdlIglMIVewWOUV1kdOVIjnTxUmHk8QMp3ea7iQxw8WpDGdu0uJ7ne20j81kIDLYTESQmSohFurFhBBhMxMelD9ZCc8zX9gkRQghRK6AJzgZGRk89thjpKamYjabmT17NtWqVWPz5s3MnDkTo9FIfHw8I0aMCFiM0SEW/jqSGrD9i9LjdOucSM8mMd3FyXRv0pIzTUx3+hOapHTnGe8lEmIxEhlsJjLIzCWhVppUsxMRZPGXRfimOe9tZmPZH6QQQlyEAp7gfPzxxzRr1owRI0awfPly3nzzTSZOnMiUKVOYO3cutWvX5t577+Wff/6hadOmAYkxOsTCSUc2SikZh1CBZLk8nHA4OZaWzXFHtnea5ps6nBxPy+ZUpqvAdgYNIoMtRIdYqBJipmHVEN+8xT/NmZeERQghyqeAJzh33XUXHo93oOSRI0cICwvD4XDgdDqpU6cOAPHx8axdu7ZAgmO3WzGZSu8Hxmg0EBERTL1LQnF6FC6TkaqhtlLbn8ht8/PRdcUJRzaHTmVyKDmTw6cyOZqaxdHULBJSsjiWmsWpjILJS0SQmUvCbFwSEUTrupFcEmajaqiVqqFWYkKtxNitRIVYLrqrewrb7qLkSJsHhrR72QtUm5dpgrNkyRLee++9fGVPP/00LVq04I477mDnzp28++67OBwO7Ha7f52QkBAOHjxYoD6Ho3Rv/Zwz8DLG6m2mLfuSaFs7olT3ebHLaXNdKZLSnRxJzSYhJYsjqVkcSckiITWLhNRsElKzCoxziQgyU9VuoWqolWbV7FQLtVLVbqVqqIWqdivVQq3n73HxeEhLDcQ9NwNLBhmXPWnzwJB2L3ul3eYxMaFnLC/TBGfAgAEMGDDgjMv+85//sGfPHu677z4+/fRT0tPT/cvS09MJCwsrqzALqBMVBMCBU5mS4JQgt0fnSGo2B5MzOejriUlwONl/Mp2jadkFriyKDDJTI9xGwxg7V9SvQvUwG9XDbdQMs3FJWCGSFyGEEBeNgJ+iWrBgAdWqVaNfv36EhIRgNBqx2+2YzWYOHDhA7dq1Wb16dUAHGVcLtWIxahw8dfH9ZV9cLo/OkZQsbxKTnMXBU5kcTPYlMylZ+QbuBpuNXBodQr3oELrW8yYwNcNtVA+3Uj3MRpAkMEIIIQop4AnOLbfcwtixY1m2bBkej4enn34agGnTpvH444/j8XiIj4+nZcuWAYvRoGnUigjigCQ4Z5WW5ebfpAz2JWZ4p77XkZSsfA8oDLEYqR0RRJNqoVzbKIbakUHUjgiiVkQQUcFmIiNDpPtYCCFEsQU8wYmOjubtt98uUN6qVSs+/vjjAER0ZnUig9ifdHEnOEopTjic/JuUwf6kDP5N9CYx/yZl5rvTs8WoUScymMZVQ7m2cVXqRARRK8JGncggIoLMciWaEEKIUhfwBKeiiK0SzC97k8hyeS6KsR6ObDd7Tqaz+2Q6u06k++fzPhrAbjUSGxVM50sjuTQqmEurBBMbFUyNcNtFdxWSEEKI8kUSnEJqdkkoHl2x47iDljXDAx1OiXHrigOnMth9wpvA5EwTUnOvULNbjTSIDqFX46rEVgkhtkoQsVHBVAmxSG+MEEKIckkSnEJqVt17FdfWhLQKm+Bku3V2n0xn+7E0th9zsP2Yg72J6Th9I32NBo26kUG0qBHGTS1CaBATQv3oEKqFWiWREUIIUaFIglNI0SEWaoRZ2ZpQMR7ZkOXysPNEnmTmuIO9iRl4fCN+w20mGlW1c2vrmv5E5tKoYCwmQ4AjF0IIIYpPEpwiuKxGGL8fTEFXqlw9/NDt8fbMbE1IY2tCKtuOOdiXlOG/eikyyEzjanbi46JoXC2UxlXtVA+TXhkhhBCVlyQ4RdAlLopvt59ga0IaLWoE7saDx9Ky+TshlS0JafydkMo/xxz+m+JFBZtpekkoVzaIpnE1O42rhVLVLmNlhBBCXFwkwSmC+NgqmAwaP+w8WWYJjtOts+1YGn8dSfX30Bx3eC/JNhs1Gle1c3OL6jSvHkrz6mHSMyOEEEIgCU6RhNpMdLw0km+3H+fB+EtLZbxKutPNliOpbDqcyuZDKfx9NM3fO1Mz3EbrWuE0rx7GZdVDaRBjlzEzQgghxBlIglNEg9rUZMTSLXz591Fublmj2PUlZ7jYfDiFTYdT2HQohZ3HHXgUGDVoWNXOLS2r07pmOC1qhhEVbCmBIxBCCCEqP0lwiqh9nQguqx7K/DX76VavCtF2a5G2P5qaxabDKWw+lMqmQyn8m+R9LIHVZKDZJaHc2aEObWqG07xGKCEW+XiEEEKICyG/oEWkaRoTezbkjg828einfzPnpuZEh5y5Z0Upxb6kTH/vzOZDKRxN895AL8RipFXNcK5rWpXWtcJpUi1UTjcJIYQQJUQSnAsQVyWEZ/o0YfyX2+j/zgaua1qN+tHB2MxGUrPcHErOZJfvjsCpWW7Ae3VT61rh3N6uFq1qhVM/OkQeZyCEEEKUEklwLlDXelV4//Y2vLF2P19sPUqWbyAwQJDZQP3oEK5qGE3zS8JoVSuc2hE2ubpJCCGEKCOS4BRDbJVgnunbBLeuSEp3ku3WCbWZCLOZytWNAIUQQoiLjSQ4JcBk0KgaWrTBxkIIIYQoPTKqVQghhBCVjiQ4QgghhKh0JMERQgghRKUjCY4QQgghKh1JcIQQQghR6WhKKRXoIIQQQgghSpL04AghhBCi0pEERwghhBCVjiQ4QgghhKh0JME5ja7rTJ48mYEDBzJ06FD2798f6JAuCi6XiyeeeILBgwfTv39/Vq5cGeiQLhqJiYl0796dPXv2BDqUi8aCBQsYOHAgN998M0uWLAl0OJWey+XiscceY9CgQQwePFi+62Xgzz//ZOjQoQDs37+f2267jcGDBzNlyhR0XT/P1iVDEpzTrFixAqfTyeLFi3nssceYNWtWoEO6KHz++edERETw4Ycf8tZbbzF9+vRAh3RRcLlcTJ48GZvNFuhQLhrr169n06ZNLFq0iPfff5+jR48GOqRK7+eff8btdvPRRx/x0EMP8dJLLwU6pErtzTffZOLEiWRnZwPwzDPPMGrUKD788EOUUmX2B6wkOKfZuHEjXbt2BaBVq1Zs3bo1wBFdHHr16sUjjzwCgFIKo9EY4IguDrNnz2bQoEFUrVo10KFcNFavXk3Dhg156KGHuP/++7niiisCHVKlFxsbi8fjQdd1HA4HJpM8hrE01alTh7lz5/rf//3337Rv3x6Abt26sXbt2jKJQz7l0zgcDux2u/+90WjE7XbLP4hSFhISAnjb/+GHH2bUqFGBDegisHz5cqKioujatStvvPFGoMO5aJw6dYojR44wf/58Dh06xAMPPMD//vc/NE0LdGiVVnBwMIcPH6Z3796cOnWK+fPnBzqkSq1nz54cOnTI/14p5f9+h4SEkJaWViZxSA/Oaex2O+np6f73uq5LclNGEhISuOOOO7jxxhvp27dvoMOp9JYtW8batWsZOnQo27ZtY+zYsZw4cSLQYVV6ERERxMfHY7FYiIuLw2q1kpSUFOiwKrWFCxcSHx/Pt99+y2effca4ceP8p09E6TMYclON9PR0wsLCyma/ZbKXCqRNmzasWrUKgM2bN9OwYcMAR3RxOHnyJMOGDeOJJ56gf//+gQ7novDf//6XDz74gPfff58mTZowe/ZsYmJiAh1Wpde2bVt++eUXlFIcO3aMzMxMIiIiAh1WpRYWFkZoaCgA4eHhuN1uPB5PgKO6eDRt2pT169cDsGrVKtq1a1cm+5WuidNcc801rFmzhkGDBqGU4umnnw50SBeF+fPnk5qayuuvv87rr78OeAeqyeBXUdn06NGDDRs20L9/f5RSTJ48WcaclbK77rqLJ598ksGDB+NyuRg9ejTBwcGBDuuiMXbsWCZNmsScOXOIi4ujZ8+eZbJfeVSDEEIIISodOUUlhBBCiEpHEhwhhBBCVDqS4AghhBCi0pEERwghhBCVjiQ4QgghhKh0JMERQpQrDz/8MAsWLPC/dzgc9OzZk+3btwcwKiFERSOXiQshypWkpCRuueUW3nzzTerXr8/kyZO59NJLGTZsWKBDE0JUIHKjPyFEuRIVFcWkSZOYOHEio0eP9j+vafjw4WRnZ2O1Wpk+fTrVq1fnhRdeYOvWrSQnJ9O4cWOeeeYZ5s6dy6ZNm8jIyGDmzJk899xzOBwOMjMzGT16NPHx8YE+RCFEGZAERwhR7lx55ZV8//33jB8/nkWLFvHMM88wdOhQunfvzrp163j++eeZNm0aYWFhvPvuu+i6zvXXX8+xY8cAiIuLY+LEiezatYvk5GTeeustEhMT2bdvX2APTAhRZiTBEUKUS/369SMrK4tq1aqxc+dOFixYwFtvvYVSCpPJ5H9I5aOPPkpwcDAZGRm4XC4AYmNjAWjQoAEDBw7k0Ucfxe12M3To0EAekhCiDEmCI4Qo9+Li4hg2bBht2rRhz549bNiwgVWrVpGQkMBLL71EUlIS33//PTlDCnOeXrxjxw7S09N54403OH78OIMGDaJHjx6BPBQhRBmRBEcIUe6NHTuWqVOnkp2dTVZWFhMmTKBWrVq8/vrrDBkyBE3TqF27NsePH8+33aWXXsprr73GN998g67rPPzwwwE6AiFEWZOrqIQQQghR6ch9cIQQQghR6UiCI4QQQohKRxIcIYQQQlQ6kuAIIYQQotKRBEcIIYQQlY4kOEIIIYSodCTBEUIIIUSlIwmOEEIIISodSXCEEEIIUelIgiOEEEKISkcSHCGEEEJUOpLgCCGEEKLSkQRHFEqjRo3o27cvN954o/81YcIE/7KkpKQARwjz5s3jiiuuYPz48YEOpciGDRsWkDa87777WL58eYHyhIQE+vTpww033MCmTZtKPY68x3/llVeyZcuWUt/nuXg8Hh544AF69uzJBx98UGr72bJlS5GfcP7XX38xefJkANavX0+fPn1KI7QSM3fuXJ566qlAh1FsTz31FHPnzj3nOocOHaJ169ZlFJE4H1OgAxAVx3vvvUdUVFSgwzirpUuX8vzzz9OuXbtAh1Jka9asCXQI+axfv57o6GgWLlxYJvsrb8d/7NgxVq9ezebNmzEajaW2n8suu4xXXnmlSNvs3r2bY8eOlVJEQlQekuCIEvXaa6/x1VdfYTQaiY2NZdKkSfz555+8/fbbLFq0CIBevXrRu3dvHnnkEY4ePUr//v354YcfmDlzJn/88Qdms5latWrxzDPPEBISkq/+o0ePMnXqVA4fPoxSin79+jF8+HBGjRrFsWPHmDBhAo888gjXXXddvu2WLFnCokWL0HWdiIgIJk2aRL169Rg3bhx2u50dO3Zw9OhR4uLimDNnDiEhIezZs4eZM2eSnJyMx+Nh6NCh9O/fn/Xr1zNz5kyCg4PJyMhg6dKlLFy4kKVLlxISEkK7du1YuXIlX331Fd26dePjjz8mNjYWgLvvvpshQ4Zw9dVX+2PL6XG68847eeONN3A4HDz11FMkJyejaRrDhg2jX79+Bdr6xx9/ZMGCBTidTpKSkujXrx+jRo1i/fr1vPjii9SuXZtdu3bhdDqZPHkyHTt25NixY4wbN47jx49To0YNEhMTC9T766+/8tJLL5GWlsbQoUMZMWJEgeP95JNPeP/99zEYDERHRzNp0iRiY2MZN24cVquVLVu2cPLkSXr37k1UVBQ//vgjJ06cYMaMGXTq1Cnf/k4/foDFixczZcoUkpKSuPHGGxk9ejQAP/zwA/PmzcPlcmGz2Rg7dmyBv5jdbjfTp08v8F06deoUffv29fdIHTp0yP9++fLlLF26lMzMTOx2O4mJibjdbm6++Wbmzp3Lb7/9xuLFi3G5XKSkpHDPPfcwePBgABYsWMAnn3yCyWSibt26zJo1i9DQ0LN+5/Jav34906dP58svvzzndzFHQkICr7zyCmlpaYwfP55+/fqRkZHB6NGj2bt3L9nZ2cyYMYN27drhdDp5/vnn2bBhAx6Ph6ZNmzJx4kTsdnu+GFwuF7NmzWLdunUYjUZatGjB+PHjsdvtXHnlldx0002sW7eOhIQEevfuzZgxYwr9WRT2Oztx4kSioqJ49NFHAfj888/59ttvee211866n7lz57J582aOHz9Oo0aNeP755/37OXToEHfeeScdO3Zk8+bNuN1uxowZw+LFi9m7dy/Nmzdnzpw5GAwGVqxYwauvvorH48FutzN+/HhatGiBw+FgwoQJbN++napVq2I0Gmnbti3gTYCfeuopEhIScLlcXH/99dx///3nPHYRAEqIQmjYsKHq06ePuuGGG/yvkydP+pclJiaqpUuXqoEDB6r09HSllFKvvPKKGjZsmMrMzFRt2rRRKSkp6uDBg6pLly5q4MCBSimlPvjgAzVlyhS1YcMG1atXL6XrulJKqWeffVZt3LixQBxDhgxR77zzjlJKqdTUVNW3b1/15ZdfKqWU6tGjh/rrr78KbLN+/Xo1ePBglZGRoZRS6pdfflG9e/dWSik1duxYNXDgQJWdna2cTqfq16+fWrp0qXK5XOq6665TW7du9e+rd+/eatOmTerXX39VjRs3VocOHVJKKbVq1SrVs2dPlZKSonRdV+PHj1c9evRQSik1Y8YMNXv2bKWUUvv371fdu3dXbrf7jO2bmJioXC6Xuuqqq9S3336rlFLq6NGjqmvXruqPP/7It76u6+r2229X//77r3+9Jk2aqMTERPXrr7+qJk2aqH/++UcppdTbb7+thgwZopRS6sEHH1QvvviiUkqpffv2qVatWqlly5YViGfZsmXq3nvvVUqpAse7du1adfXVV6vExET/ur1791a6rquxY8eqAQMGKKfTqY4fP64aNmyo/vOf/yillFq4cKG6++67C+wr7/HnfI5PPfWUUkqp48ePq+bNm6sjR46of//9V/Xp00clJSUppZTauXOn6tKli//7luNs36WDBw+qVq1a+dfL+37ZsmXq8ssvV2lpaQWWORwOdeutt/r3u2nTJv+yFStWqGuvvVYlJycrpZR6+umn1euvv37O71xev/76q7r++uuVUmf/Lp7vs2nSpInavHmzUkqpd999V91xxx1KKaXmzp2rZs2a5W+HF154QU2ZMqVAfS+//LIaMWKEcjqdyuPxqHHjxqlJkyb5P4tZs2Yppbzfscsuu0wdOHCg0J+FUt7/D0ybNu2c39l//vlHdenSRblcLqWUUoMHD1arVq06535eeeUV1bNnT/82eR08eFA1bNhQrVixQiml1OTJk1WPHj1UWlqaysrKUl26dFEbN25Uu3fvVp07d1YHDhxQSnm/2126dFFpaWlq5syZasyYMUrXdZWYmKi6deumXnnlFaWUUkOHDlUrV65USimVlZWlhg4dqr766qsC3zERWNKDIwrtfKeoVq1axc0330xwcDAAd9xxB/Pnz8dgMNC5c2fWrFlDcnIyAwcOZPHixaSlpfHDDz8wfPhwGjZsiNFoZMCAAcTHx9OzZ09atGiRr/6MjAz++OMP3nnnHQBCQ0O5+eabWbVqFddff/1Z4/rpp5/Yv38/gwYN8pelpKSQnJwMQNeuXbFYLAA0bNiQlJQU9u3bx4EDB3jyySf922RlZfHPP/9Qr149qlevTs2aNQH4+eef6dWrF2FhYQAMGTKEX3/9FYDBgwdz++23M3r0aBYvXkz//v3Pecpj3759ZGdnc+211wJQrVo1rr32Wn755Zd8fx1rmsb8+fP56aef+PLLL9mzZw9KKTIzMwGoUaMGTZo0AaBp06Z88sknAKxdu5axY8cCULduXTp06HDWWPLKe7y//PIL1113nf+7cPPNNzNz5kwOHToEQI8ePTCbzcTExBAcHEzXrl0BqFOnjr/NzydnXElMTAzR0dEkJiby559/cvz4ce6666587XDgwAEaN27sLzvbdyknvrNp1KhRgd4NgJCQEObPn8/PP//Mvn372L59OxkZGQCsW7eOXr16ER4eDuT2Rj377LNn/c5FREScNYYzfRfPp3bt2rRs2RKAxo0bs2zZMsD7vU9LS2Pt2rWAt6emSpUqBbZftWoVo0ePxmw2AzB06FAeeugh//KrrroK8H4Xq1SpQkpKSqE/i7zO9Z1t0qQJtWrV4qeffiI2Npbjx48THx/Phx9+eNb9ALRq1QqT6cw/Y2azmSuvvBLwfvdat27t/3yrVq1KSkoK27Zto2PHjtSuXRuATp06ERUVxdatW1m3bh1PPvkkmqYRFRXFNddcA3j/P7RhwwZSUlJ4+eWX/WXbt28v8P8sEViS4IgSo5TK917XddxuNwDXXHMNq1atIjU1leHDh7N3715WrFjBzp07ufzyyzGZTHz22Wf88ccf/Prrr4waNYo77rgj3//YdF0/5z7ORtd1brzxRp544gn/++PHj/t/lGw2m39dTdNQSuHxeAgLC+Ozzz7zLzt58iShoaFs3rzZn8QBmEymfHHlTWBiY2Np1KgRK1eu5IsvvmDJkiXnjfV0SqkCx5iRkcFNN93E1VdfTbt27bjllltYsWKFP44zHdPp8zmxF0be4z39Mzg9xpwf6KLu42zb5MSs6zqdOnXipZde8i9LSEigatWq+bbN+dxO/y5dc801+WJ3uVxnPca8jh49ysCBA7n11ltp27YtvXr14scffwS8n7Wmaf51U1NTSU1NPe937mzO9rmdS05icvo2uq7z5JNP0r17dwDS09PJzs4usP3p3zld1/O1jdVqLVD/uT6Le+65h+PHjwPkG0B9vu/skCFDWLZsGZdeeim33normqadcz/ff//9WT+znHbJ+9nkbacc5/sun+nfdc7/hz766COCgoIASEpKwmq1curUqbPGI8qeXEUlSkx8fDzLly/3/3X7/vvvc/nll2OxWLjiiitYt24d27Zto0WLFnTp0oWXX36Zbt26YTKZ+PHHH7nrrrto3bo1I0eOpF+/fmzfvj1f/Xa7nZYtW/Lf//4XgLS0ND799FM6d+58zri6dOnCV1995f+f7qJFi7jzzjvPuU1sbCxWq9Wf4ORcVbR169YC63bv3p3vvvuOtLQ0wDvYOa/Bgwfz7LPP0rJlS6pVq3bG/RmNRtxuN7GxsZjNZr777jvAe67/22+/LXCM+/fvx+FwMGrUKK688kp+++03nE7nGROkvLp27crixYsBOHLkCOvXrz/n+mcSHx/P119/7b/qadmyZURERFC3bt0i15Uj5/jPpWPHjqxZs4Y9e/YA3p6zG264ocCP9tm+S2FhYbhcLnbv3g3A999/X6jYtm7dSlRUFA8++CBdu3b1Jzcej4fOnTvz/fff43A4AO8VQwsXLryg71xhFaatwPs5/fe///V/LyZNmsScOXMKrNe1a1c++ugjXC4Xuq7z3//+ly5dupyz7nN9Fm+++SafffYZn332mb/3B87/ne3Zsyfbtm3ju+++45ZbbjnvfkpCTv0HDx4E8I81atmyJV27dmXp0qXouk5KSgorV64EvP8fatWqFe+++y7gTWpvu+02/3JRfkgPjigx/fv3JyEhgQEDBqDrOnXr1vUP/AsLC6NevXoEBQVhNBqJj49nwoQJ/lMx3bp1Y9WqVfTp04fg4GDCw8OZPn16gX08//zzPPXUUyxfvhyn00nfvn25+eabzxlX165dueeeexg2bBiapmG323n11Vfz/XV3OovFwuuvv87MmTN56623cLvdPPLII7Rt27ZAUtCpUyduvfVWBg4ciM1mo0GDBv6/7MB7ymbixIn5Tlec7pprrmHw4MG8/vrrvP7668yYMYO5c+fi8Xh46KGH6NixY771GzVqxBVXXEHv3r0JCwujTp061K9fn/379xfoQclrypQpjB8/nt69e3PJJZec9XTCuXTp0oW77rqLO++8E13XiYqKYsGCBRgMF/73Ut7jP5sGDRrw1FNP8eijj6KUwmQyMW/evAJ/xZ/tuxQaGsoTTzzBPffcQ1RUFL169Sr08S5dupRevXoRFBREixYtiIqKYv/+/XTv3p3du3dz2223AVC/fn2mT5+O3W4v8neusFq3bs1LL73EQw89xB133HHW9R588EFmz57NTTfdhMfjoUmTJowbN67Aeg888ACzZ8+mX79+uN1uWrRowaRJk84ZQ2E/i7zO9Z2tU6cOFouFnj17cvLkSf/pzwvZT1HUr1+fKVOmMGLECDweDzabjfnz5xMaGsrIkSOZMmWKf6B8w4YN/ds9//zzTJ8+nb59++J0Ov23VDjfaVBRtjRVmD5QIcRZbdmyhU2bNvl/bN59913+/PNPf7f6H3/8waRJk/jyyy9L5AdOiMooIyOD22+/nSlTpvjHFAlRHNKDI0QxxcbG8uabb/Lxxx+jaRrVq1f39z6NHTuW3377jdmzZ0tyI8RZ/PLLLzz22GPccsstktyIEiM9OEIIIYSodGSQsRBCCCEqnYAkOImJiXTv3p09e/awf/9+brvtNgYPHsyUKVP8I+pfffVV+vfvz6BBg/jrr78CEaYQQgghKqgyH4PjcrmYPHmy/34PzzzzDKNGjaJDhw5MnjyZlStXUqNGDX777TeWLFlCQkICI0eO9N+8Kq8TJ9JKNVa73YrDUTKXI4rCkTYPDGn3sidtHhjS7mWvtNs8Jib0jOVl3oMze/ZsBg0a5L85199//0379u0B7+Wda9euZePGjcTHx6NpGjVq1MDj8QTkScsmU+k9ZE+cmbR5YEi7lz1p88CQdi97gWrzMu3BWb58OVFRUXTt2tX/UD2llP/qkpCQENLS0nA4HPluZ55TfvpjAux2a6k2nNFoICKiZO63IApH2jwwpN3LnrR5YEi7l71AtXmZJjjLli1D0zT/HW3Hjh2br2cmPT2dsLAw7HY76enp+cpDQwt2QZV2N2NERDDJyRmlug+Rn7R5YEi7lz1p8/NTukIp7x/C+ed9j1FQOev4lulnnpJn+5AQK2lpWd5tc+pQuY9l8K+vOG2d3OV5t0MpFLnr5tvuTOtzprpz93n2dfLGl2dfnLkef3lObOSv72zL/fvOV55n3/7yvPXgn8mtJ/dYatSPpGHXM9/FvSSc7RRVmSY4ObfYB+8D3aZOncpzzz3H+vXr6dChA6tWraJjx47UqVOH5557jv/7v//j6NGj/rulCiFEZZXzw63rCt2jozwKPeelK3S3d6p8r/zzvh9zz5nLC25XhHLP2dbHP58vwdDJ997/Y3yGBCU3CcmbrOSuI0DTAM37HDC8/4Gm5SvPucWWfx3vSnnWP225v2L/xFeu5Vk3N4D8deSpu5Db6p7AfJgBv9Hf2LFj/c9IiYuLo2fPnhiNRtq1a8fAgQPRdZ3JkycHOkwhRAWnlC9ZcCs8bh3do+NxK3R3nqnntPdu5VtP923nSz70PMmH/6X75/Mv1/Otl5s05MThS2A8uX8RlyXNAJpBw2DQ0Hwvg6+sMOUGk4amGdA037KcH1GDlqfM90N8Whma5t1/3mmeenLm/T/kp6+br27Osk3+fYaG2nCkZ3vfQ/7kIWcecuPLmZ5tXY1zrgM5MeatLzchyFued33yxFLRBaq3skLf6K+0r6KSLuSyJ20eGGXd7jk/7m6Xjuf0l1vH7fROvWUKt8uDx6XwuDx43Mq/rtuln5ag6Ogu5Z26cxOInASlRGhgMHp/4A1GDYPR4Jv6fvCN3h997/LcZfmWGw3Ygs243Z7TluVfP/dl8CYUedctQhJy9vLK8yNaWPL/mLJX2m1eLk5RCSEqDqW8iYTbqeN2enBn67hdHu/7bN80Z1neac4y37oFkhfXhScbBqOG0WTAaPa9TBpGswGDyTtvspoxmjTfewMGk4bReNp7f3me90ZfPcac5QbvdmZDnu1zE5CSID+0QpQuSXCEqESU8p6CcWV7vK+s3Kk724PTN81b7sr2gEeRleHOn6y49CKdMjFZDJgsxtyp1YDZZiQo1OxNSiyGfMmJyXza/HnWyUlEhBCiMCTBEaKc8bh0nFlunBkenFluXJkenJlu38s7709QcpKUnMQl2+Mdy3EeRpOG2WbEZDVithmxBZsJDrd4kxOrsWCyYjGcMYHJeW80GbzjHoQQopyQBEeIUuJ26WSnu3BmuMlOd5OdkZukuLLyJyzOTA8uXxLjOdfpGw3MViOWIBNmqwGzzURQmJmwGBtmmxGz1fsy2YxYrLkJjDnP1GQ1YjTlv8ennC4RomJISDjClClPUrfupezcuZ3Q0DAAUlKSGTTodq6//oYAR1h+SIIjRCHousKVmZuoZGe4yU53+aZuXxKT+z47w43HdfZLI00WA5YgE5YgI+YgE6FVrFiCQ7DYcsssQUbfOrnzZqtRekqEEAA88MDDdOzYGYDU1BSGDr2V667re9ENHD8bSXDERUvpCmemmyyHy/fyzae58pR5X84MN2e73tBkMWANMWMNNmGzmwmrGoQ12OQvs4aYsAabsAR7kxWzrWAPihCiYtm36ST//nGyROuMbRPNpa2jL2jbxMRELBarJDd5SIIjKh3do8hyuMhMdZKZ6iQj1ZukZOdNWtJcZKW7vTciO43RpGGzm7GFmrFHWalSx44txOxPVKwhJizBvvdBJoxmSVaEEGVv3rxX+M9/3uHo0QQuvTSO6dNnBTqkckUSHFGhuJ0eMtN8yUuKN3nJTPPOZ6a5yExxkuVwFeht0QxgC/EmLdYQM+GXBHuTGLuZIHtOuQlbqNl7Gkj+ChJCnMOlrS+8t6Wk5JyiWrduNfPmzaVGjVoBjae8kQRHlBtKV2SkZJN40EH6qWzSk51kpGSTkewkI8XbG+PM9BTYLudS5KAwC2ENwv3zweEW/7w12CRjV4QQlVKnTvFs3bqFZ5+dyYwZswMdTrkhCY4oM7pHkZnqJD3Zm7SkJzvJSM72TzNSnAUucbYEGQmOsBISaSW6bijBYd6ExfvyzputpfdEeSGEqAjuums4d989hLVrV9O5c3ygwykX5FEN5yCXzhad2+nBkZRNelI2aUlZpCdl4/C9MlKyCzxAz2Y3ExxhISTCQnCElejqdjSr5n8vyUvZkO962ZM2Dwxp97Inj2oQFUZ2hhtHUhaORG/ikn4qy5vEJGaT5XDlW9cSZMQeZSOqVgh1WkQREmH1JjSRVoLDLAUG6Mr/fIQQQpQESXDEGekeHUdSNmkns0hLzCLtRBZpJ7NIPZmFM8Odb92gMDP2KBuXNAzHHmnFXsWKPcqGPcqKJUi+YkIIIcqe/Ppc5LIz3KSdzPQmL74kJu2kt0cm7yXUNruJ0OggajWNJDTam7zYo2yERFkxyWXSQgghyhlJcC4SrmwPqcczSTmeScqxTO/8scx8p5QMRg17lJWwGBs1m0YSFm0j1PeSnhghhBAVifxqVTJul07aicx8yUzK8Uwykp3+dYxmA2FVbVzSIJywqjbCYoIIjbYREmGVpzULIYSoFCTBqcCyHC6SEzI4lZBB8tEMkhMycCRm+W9yZzBqhEbbiK5tJ6xdEOFVgwivFkRIhFXuCSOEEKJSkwSnAlC6wnEqm+QjGZzyJTLJRzPISss9vRQcYSGyejC1m0cRXs2bzNirWDEYZXyMEEJUFu+/v5Dff/8Nj8eNpmk89NAoGjducsZ1c548/sYbC4u8n88+W87119+AyVR6aUJ2djbfffcNffv24+23F1ClShX69etfYvVLglPO6Loi7WQWpw6nk3Q4nVMJGaQczcDt9N5ARjNohMXYqFYvjIjqwURWDybikmAZIyOEEJXcv//uZc2aVcyb9zaaprFr1w5mzJjKe+8tKvF9vf/+u/TqdX2pJjhJSYl88cWn9O3br1Tql1/FAFJKkZ7s5NQhbzKTdDidU0fS/cmMyWIgonowl7aO9iczYTFB8nBHIYQIMOv2pdi2fVSidWY1GUR247P3YNjtdo4dO8pXX31Ghw6dadCgEW+++R4AO3du58UXn8NoNGKxWBgzZmK+bX/8cQXLly/B7fb2/Dz99PPs3bubefPmYjabueGGm+jV63oAvvzyU5KSEpk69UmeeeYF5s59kb/+2gzANdf04tZbb8tX96ZNG3njjdcxGo3UqFGTMWMmkJ2dxaxZM3A40jh1KpEbb+zPTTf1Z8SIe4mMjCI1NZUaNWqwb9+/vPvumwD88ssqfvxxJSkpKQwffj/x8d2K1Z6S4JShLIeLxIMObyJzOJ2kIxn+e8oYjJo/mYmqGUJkzRBCo20YZKyMEEIIICamKrNmzWHZssW8886b2Gw27r33Qa644ipmz57JuHETadCgEb/88hOvvjqHhx4a5d/24MEDPPfcy9hsNp59dia//baO6OgYnE6nP0nK0adPPxYufJupU59mzZpfSEg4whtvLMTj8fDAA/9H27aXU69efcD7h/rs2TOZN+8tIiOjePPNeXz99Rc0atSEq6++lu7dr8TpdHDHHUO56SZv8nb11T3p3r0HCQlH2LNnN3fffQ9vv72AmJgYxo2bxB9//M6HH/5HEpzyStcVqccyOXnQQeIBBycPeB8gCaBpEFY1iJqNI4isGUJUzRDCqwVhNEnPjBBCVATZjfufs7elNBw6dJCQkBCefHIKANu3/8Pjjz9MmzbtOHnyBA0aNAKgZcs2zJ//ar5tIyOjmDFjCsHBwezfv4/mzVsAUKdO3XPuc//+f2nZshWapmEymWjW7DL27dvrT3CSk0+RmHiSSZPGAd5xNZdf3oFOnbrw8ccf8vPPPxIZGY7bnXuD2LPts1Ej71iiKlWiycrKKmrzFCAJTglxZrpJPJRO4gEHiQe9r5xTTTa7iSq17dRrH0OV2nYiqwdjssgzloQQQhTenj27+OyzT5g9ew5ms5natetgt4diMBiJjo5h9+5d1K/fgM2b/6B27Tr+7RwOB2+/vYBly74EYPToh8h5DOXZzhJomgGlFHXrxvL1158zcOAQ3G43W7f+Re/effzrhYdHULWqt2fJbrezevXPBAUF89FHH9C8eQtuuqk/O3du4aeffvJvYzAY8uxDz7PPEmsqQBKcC+bMcnNyv4Pje1M5/m8ayUczQHk/oPBLvKeaqtS2U6WOnZAIC1pJf3JCCCEuKt27X8m+ff8yfPgdBAcHoeuKBx98BLvdztixE3jxxWdRSmE0Ghk3bpJ/u5CQEC67rCX33383RqOJ0NBQTp48QfXqNc66r5YtW/H44w8zd+4CNm3ayH333Y3L5eLKK6+mUaPG/vUMBgOPPPI4TzzxCEopgoNDmDRpGpqm8eKLz7Jy5XdERoZjNBpxOp359hEZGYnL5eb111/BarWWeHsV6Wniv//+O++99x4bN27EbDZjNBpp3bo1Q4YMoU2bNiUe3PmU9dPE0xKzOLztFIe3JZN00IFSYDBpVKltp2psKNF1Q4mqGSJPwC4GedhmYEi7lz1p88CQdi975f5p4tOnT8dutzNy5Ejq16/v72LasWMHn3/+OZ9//jlTp04tkWDLE49b5+CWJHb/dpykQ+kARNYIpnG36lSLCyOqtl2exSSEEEKUM4XuwUlMTKRKlSpnXX7y5Emio6NLLLDCKO0enMyTTn7+cAepJ7IIjbYRd3kMtZpGEhJR8l1pwkv+ugoMafeyJ20eGNLuZa/c9uCsWbOGb775hiFDhlClShUWL17MwIEDC6xX1slNafK4df7+4TA71hwjONxCl8H1qdE4QsbRCCGEEBXEeROcZcuWMXXqVObNm0dycjLbtm0ri7gCypXlYd/mRGo3iaRtv0sx22RMjRBCCFGRnHfwSEhICGFhYYwdO5Y1a9awZcuWsogroGx2MzeMaUWv+y+T5EYIIYSogM6b4HTv3t0///jjj9OvX7/SjEcIIYQQotjOe4rq6quvzvd+6NChpRaMEEIIIc4s79PBZ86cys6d2wkNDQMgJSWZQYNu5/rrbzjr9vPmzaVu3Uu57rq+ZRVyQBX5Rn8ulwuz2XxBO/N4PEycOJF///0XTdOYNm0aVquVcePGoWkaDRo0YMqUKRgMBl599VV++uknTCYTTz75JC1atLigfQohhBCV0QMPPEzHjp0BSE1NYejQW7nuur5yQYxPkRKciRMn8tNPP2E0GqlatSqNGjWiUaNGhe7V+fHHHwH46KOPWL9+PS+++CJKKUaNGkWHDh2YPHkyK1eupEaNGvz2228sWbKEhIQERo4cybJly4p+dEIIIUQp+O7QN3xz6MsSrbN3rT5cW6v3BW2bmJiIxWItkNz89NNK3nvvbSIiInG5XNSteylAgSeE9+zZm0ceeZCFCz9k69YtPP74w3z99UpOnjzBrFnTufrqnqxbt4bs7CwOHz7EkCF3lvueoCIlOBs2bODnn3/GaDRy7Ngxtm/fzo4dOwq9/dVXX80VV1wBwJEjRwgLC2Pt2rW0b98egG7durFmzRpiY2OJj49H0zRq1KiBx+MhKSmJqKiofPXZ7VZMptIbBGw0GoiICC61+kVB0uaBIe1e9qTNA6Ok2j04yYKphB+QHBxsOWds6elBmEze+C0WEwsWvMqHHy7kyJEj1KtXnxdffCnf9i6Xi9dee4klS5YRHh7Ogw/eT3Cwhc2bf+PkyWN8/PHHuN1u7rjjdq64oitVqkSSmZnC5s2/UaNGDQ4f3svff/9Nz57XYrVayM7O5I033mT//n2MGPEQgwcXvGXMmQTqu16kBKdly5YkJydTpUoVqlWrRrVq1fINQi7UDk0mxo4dy/fff88rr7zCmjVr/BlnSEgIaWlpOBwOIiIi/NvklJ+e4Dgc2UXad1HJDaHKnrR5YEi7lz1p88AoqXaPj7qK+KirSiCi/M4VW2pqJm63TnJyBk6nm/vuG0HHjp1Zt2418+bNJSwsOt/2J0+ewG4PBSykpGTSqFEzMjKcHDmynaZNLyMlJROARo2asWXLP3Tu3I1vv13B779v5LbbhvLDDz+zdetfjB8/mV9/Xcull9YjOTkDmy2czMysQrdjoG70V6T0c+DAgdx+++28/fbb/P7776SlXdidhGfPns23337LpEmTyM7OTVLS09MJCwvDbreTnp6erzw09MwHIIQQQlzMOnWKp2vXK3j22Zn5yiMiInE4HJw6dQqA7dv/AaBu3Vj/6amcJ4TXqlWHrl2vYMWKbwkJCaFDh8788svPOJ0uoqK8TzGoaGN7ipTgPPHEE/Tr1w+Px8OHH37IoEGDClxldS6ffvopCxYsACAoKAhN02jevDnr168HYNWqVbRr1442bdqwevVqdF3nyJEj6LpeoPdGCCGEEF533TWcffv+Ze3a1f4yk8nE6NFjeOyxETzyyIO4XG4AunTpSvXqNbnvvru59967uOKKK2nUqDFVq1YjOzubtm3bExYWhtFopHPnLoE6pGIr0tPEBw8ezIcffpivzOl0YrFYCrV9RkYG48eP5+TJk7jdbu655x7q1avHpEmTcLlcxMXFMWPGDIxGI3PnzmXVqlXous748eNp165dgfrK+mniovRJmweGtHvZkzYPDGn3sheoU1RFSnCmT59OnTp1uPPOO0sssOKQBKfykTYPDGn3sidtHhjS7mWv3D5sM6+TJ0+ydu1a3nzzTZo1a0bjxo1p3LgxvXtf2GVtQgghhBCloVAJzqZNm2jVqhUvv/wy4D0ttWvXLnbu3Mlff/0lCY4QQgghypVCJTiffvopTz31FJdeeildu3ala9euNGvWjGbNmpV2fEIIIYQQRVaoBGfatGkA7Nmzh19++YVx48bhcDjo0KEDXbt2pU2bNhiN8tRtIYQQQpQPhRpknJGRQXBw/rsQZmVlsX79elatWsWmTZtYvnx5qQV5NjLIuPKRNg8MafeyJ20eGNLuZa9cDzLu1asXV155JYMGDaJx48YA2Gw2unfvXuQ7GQshhBBClLZC3ejvf//7H02bNmXixIkMGjSITz/9FKfTWdqxBYzDlcbE38ewbPfSQIcihBBCiAtQqAQnODiYW2+9laVLlzJlyhQ2bdpE7969mTVrFv/++29pxxgAGmuPr+Z/+74JdCBCCCGEuACFOkW1f/9+HA4H6enpOBwO2rZtS506dViyZAnvvfce27ZtK+04y5TdbKd9TCfS3amBDkUIIYQQF6BQCU7Pnj2pVq0a11xzDWFhYQQHBxMSEsKDDz5ISEhIaccYEDajlZPOzECHIYQQQogLUKgE55NPPuGjjz5i/fr1XHfddfTr14/o6OjSji2gbMYgMt2S4AghhBAVUaESnCZNmjBt2jQyMjL48ssveeihh7jkkksYNGgQnTp1Ku0YA8JmtJHlzgp0GEIIIYS4AEV6FpXBYOCqq66iQ4cO/Prrr0ybNg1N0/jmm8o3GNdqtEkPjhBCCFFBFSrBufzyy9E0jaCgIOx2O3a7nZCQEBo2bFhpx+AEGYPI8mShKx2DVqiLzYQQQghRThQqwfntt9/QNK20YylXrEYrAE7dic1oC3A0QgghhCiKQnVNaJrGlClT2LVr1xmXb9u2jcmTJ5doYIFmMwYBkCWnqYQQQogKp9BjcEaPHs1LL73E1q1biY2NJTo6mtTUVLZt20aLFi0YNWpUKYZZ9nJ6bbI8MtBYCCGEqGgKneBEREQwdepUHA4Hf/75J6dOnaJKlSpMmDChwIM4KwNJcIQQQoiKq0hXUQHY7Xa6dOlSGrGUK0Em7ymqTI88dVYIIYSoaOTyoLMIMdkBcLgcAY5ECCGEEEUlCc5Z5CQ46e70AEcihBBCiKIqcoJz3333sWLFCjweT2nEU26EmL3398mQBEcIIYSocIqc4IwZM4Y//viDm2++meeee459+/aVQliB5+/BkVNUQgghRIVT5ASnXr16jBkzhnfffZejR4/Sp08f7r77bjZt2lQa8QVMsCkYDU1OUQkhhBAVUJGvovr555/55JNP2LNnDzfccANPPvkkbrebe+65h88//7w0YgwIg2YgxByCwy09OEIIIURFU+QE54svvmDw4MG0b98+X/nIkSNLLKjywm62yykqIYQQogIqcoKzfv16vvnmGyIjI0lOTsZisRAdHV3pHtUAEGoJI82VGugwhBBCCFFERU5wLr/8ckaMGEFcXBwHDhzg1Vdf5aGHHuKJJ54gPj6+NGIMCOv2pUSe2k+6+pfwT/oHOpyLhtFkJNxdua/QK4+k3cuetHlgSLuXPa3tHVDnhjLfb5EHGR89epS4uDgA6tSpQ0JCAnXr1sVoNJZ4cIEWoQwkowc6DCGEEEIUUZF7cGJiYnj++edp3bo1mzZtIjo6mjVr1mA2m0sjvoDJbtyfCONhkvZ9TUqfpYEO56IRERFMSrI8HqOsSbuXPWnzwJB2L3sREcEQgDYvcg/Os88+S9WqVVm1ahXVq1dn1qxZBAcHM2fOnNKIL6AirBGkudLw6O5AhyKEEEKIIihyD86mTZu444478pW1bt26UNu6XC6efPJJDh8+jNPp5IEHHqB+/fqMGzcOTdNo0KABU6ZMwWAw8Oqrr/LTTz9hMpl48sknadGiRVFDLbYIayQAKa5UoqxRZb5/IYQQQlyYIic4c+fOpWPHjhe0s88//5yIiAiee+45kpOT6devH40bN2bUqFF06NCByZMns3LlSmrUqMFvv/3GkiVLSEhIYOTIkSxbtuyC9lkcMUHRACRln5QERwghhKhAipzgaJrGQw89RGxsLAaD9wzXo48+Wqhte/XqRc+ePQFQSmE0Gvn777/999Tp1q0ba9asITY2lvj4eDRNo0aNGng8HpKSkoiKKtsk45KQSwA4lnmU+mENy3TfQgghhLhwRU5wbrnlFjRN49SpU4SFhRXp6qmQEO8DLB0OBw8//DCjRo1i9uzZaJrmX56WlobD4SAiIiLfdmlpaQUSHLvdislUildvOWsCkEqSd5CUKHVGo0HaOgCk3cuetHlgSLuXvUC1eZETnJo1a/LUU0/h8Xjo1asXNWrUKNL2CQkJPPTQQwwePJi+ffvy3HPP+Zelp6cTFhaG3W4nPT09X3loaGiBuhyO7KKGXyTh4RFYDBb2JR0kWUbdl4mIiGBp6wCQdi970uaBIe1e9kq7zWNiCuYHcAFXUb388st88MEHREdHc//997No0aJCb3vy5EmGDRvGE088Qf/+3pvnNW3alPXr1wOwatUq2rVrR5s2bVi9ejW6rnPkyBF0XS/z01PgPR1XNegSjmUeLfN9CyGEEOLCFbkHx2AwEBERgaZpWK1W/2mnwpg/fz6pqam8/vrrvP766wBMmDCBGTNmMGfOHOLi4ujZsydGo5F27doxcOBAdF0P6GMgqtmqcTzzWMD2L4QQQoii05RSqigbTJgwgaioKH7++Wf69OnDnj17mD17dmnFd04nTqSVav0REcFMXDWJdcfXsOzqL0t1X8JLuo8DQ9q97EmbB4a0e9mrMKeopk2bRo0aNWjbti1BQUFMnz692MGVZ1WDqnHKmYTTU7rjfYQQQghRcop8ispkMnHbbbeVRizlUh37pQDsc+yjYXijwAYjhBBCiEIpcoLz6aefsmDBApxOJ0opNE1j5cqVpRFbudDAd/+bXak7JMERQlQISikUCqV0PEpHodCVjkJHVwqFjlL4pr51UXm2y1N2pvJzbIOvXFe+eaXQOW3+TNsXqMu3Vb54c4/LvyffPP5jVOBbV0fPnVfeBycHJZpJz8j2txP+rfO3ne8dOYM4cvZNnrW96+Yuz1n5tBpyKvbP58R4er2565KvLO9IktPXPf0YlDqtrrxzBY7ttGMnT13nODb/NvmOt+Cx5WhVvQVXx1xPWStygvPmm28yf/58qlevXhrxlDs1gmsSYrKzK2UH1A50NEJUHkopdOXBrTy4dTdu5cajPOjK45vqZ3h/prLc97nL9UKt41HuPMv1s6yTf3+6P2Hw+JMF3Xcsuv8HV+FRHt8Puu47Vt037z0Og1HD5Xb72yFnW91Xd95tPcrj/zHXz7JvfxLj205UXhpavimafw4NzX9vOfzr5dlG0/KUa+SseloNOdXmL9e0PLXmLMutN3dd8sQAutFVMRKc2rVrU7du3dKIpVzSNI0GYQ3Zlboz0KEIcUa60nHrLly6G5fuxKV8U92FW3fh1F24dCdu3Y0zT7lLd+H0lZts4EjPxK3cuHU3HuX2Jx4e5cHjS0DyLde963h867j9Ze5863t0jz95yb9t+XuIrVEzYtCMGDUjRs3ge2/AqJkwaIbcFwY0zYBRM6ChYdAMaOQs13zzGgbN6F9uNBixYEHTvO8tZhMeN3nWN/jq0U7bj4bBt1zTDBjQvNMC8759++fPsq1/m5yfPe3M874fLC3v/JnKfPX5f0Q1zReXt8R7DHjnfevl7ss7DDQnTnw/qN5578+ldx3NX79/X752y/mR1nxtlVN/Tlviq8MbG4SHBZOWmnXW5CAnnsIkB77NC5Uc5Gyv5an37OU5ceVNRHIThoomUAO7i5zg2Gw2hg8fTpMmTfwNXthHNVRUDcIb8tn+5Tg92ViM1kCHI8oRpRQu3Um2no3Tkzt16tlke7Jx6k7/1OnJ9pfnXS9nnWxPNi5fMuLSXbiUNwlxeVy4Vf5ExaU7cfoSldJIFAyaEZNmxKiZMBmMmDQTRoMJk2bKP2/wrmPUjJgNZoKMQXmWmfLU4V3HO8273ITRV3/OOnmTjJwf+pwEI3/iYcyTgORPTE4vO9M6Oe9zlxX5motikat5AiMiJJggl7T7xaDICU737t0v+FENFVWbKpez5N+P2Jy0ifYxF/agUVH2dKWT7cki05NFtieLLE8mWZ4sMt3eqbcsi0xPpn9ZljuTbE+2r8y7TU7yka17E5TTE5jiMBssWAwWrEYrVoMVs8GM2WDxTc1YDVbsJjsmgxmLb5nJYPJtZ8akedez5CnP2dZsMBfYzuJbfvp20ZFhOFJd/iSjrH/shRCipJX5oxoqotZV2mAzBrHm2C+S4JSSnJ6QU1lZJGQkkunOIMOdQbo73Tvv8b7PKc/0LcvweOfzJi/ZerY3UdGLdmm/hobVaCPIaMNmDMJmtGE12rCZbESaIrEYrFiNFixGq2/e6k1ODFZfmcVfZvElLBajxT+1+Naz+rY1G8zlJpEItQTjMclftUKIyqPICU7OoxpGjhzJ/fffz2233caAAQNKI7Zyw2K0cnl0B9YcW8XIpqMxGYrcbJWarnQy3Ok4XA7S3Q4cLgcOd5p36koj3Z2Ow5WGw+0gw51ORp4kxZu4eMs8ylOo/VkNVoJNIQSbggkyBRNsCibcEkk1o42gvImJ0YbNFORPWHLKgk5LXrzvg7AYLBX6PLcQQohcZfqohoqsd+3r+eXYT/yQ8D3X1uwd6HBKnEt3kepMIdWVSqorhTRnKmm+pCT9jAlLblmGOz3PpYdnFmQMJsQcQogphCCjNymJslbJTVJ8ZVVCI9FcpnxlOclMsCnYP8ZDCCGEOJci/1LUqVOHF154geTkZN54442L4hQVQIeYzsTa4/hg93t0v+RKrOV0sLFHeXC40khxppDmSiXV6U1YUl2ppOaU+eZzy9LI9Jz79ESIyY7dbMduCiXEHEK1oGrUM9fPLTOFYDeHEmK2YzfZvfO+MrsppNBJiQy8FEIIURKK/Cwqt9vNkiVL2LlzJ3FxcQwcOBCLxVJa8Z1TWTyLKu+P7YYT6xm7YTR9at/I6OZjyuR0hq500lxpJDtPeV/Z3ukp3/T0+TTX2dvEgIFQSxih5jDCzGGEWcK9U3MYoZYwwszhhFvC/ctDzWHYzXaCTMEYtbIZTC4JTmBIu5c9afPAkHYve4F6FpU8qqEILo/pwKC4IXy0979keTJ5sMkjRFgji1SHUoosTyanfMnKKV9ikpKd7J/PW57sTEY/y9iUMHM4kdZIIiyRxIXWJ8ISQbglgjBfshKWJ1kJs4QRYrKXm0GtQgghRGmSwQxFdE+jB7EZg/jP7nf5KeEH2kZfTsPwxkRaorCZbL6bqWWT5c7yJygpzmT/NMWZfNare7yDZSOItERySdAlNI5oQqTFm8BEWCOJtET55iMIN4fLWBQhhBDiLAr9C5mWlkZo6Jm7gS4mmqZxR4NhdK9+JV8f/Jxfj69lw4n1Z7w1us1oI9wSQYQlgkhrFLGhcf4EJsLX85J3vryO6xFCCCEqmkKPwbnttttYtGgRU6ZMYdq0aaUdV6GU9Rics3HpLhyuNLI8WZg0Exaj1XcZsiQsRSXnxwND2r3sSZsHhrR72Sv3Y3BMJhO33HIL+/fvZ8eOHQD+p4l/9NFHJRNlBWU2mIm0RgU6DCGEEEL4FDrBWbhwIceOHWPq1KlMmTKlNGMSQgghhCiWQl9SYzQaqVGjBm+88QYhISEkJiYSFBREzZo1SzM+IYQQQogiK/I1w9988w0DBw5k/vz5DBw4kM8++6w04hJCCCGEuGBFvs544cKFLF++nJCQEBwOB3feeSc33nhjacQmhBBCCHFBityDo2ma//lTdrsdq1WuFBJCCCFE+VLkHpzatWsza9Ys2rVrx++//06dOnVKIy4hhBBCiAtW5B6cZ555htq1a7N27Vpq167N9OnTSyMuIYQQQogLdkHPohoyZEhpxCKEEEIIUSLkyYtCCCGEqHSKnOC4XK7SiEMIIYQQosQU+RTVwIEDiY2N5dprr6V79+7YbLbSiEsIIYQQ4oIVOcFZvnw5e/bsYeXKldx1111UqVKF1157rTRiE0IIIYS4IEVOcLZt28batWtZv349APXq1SvxoIQQQgghiqPICc7tt99O7dq1GT16NN27dy+NmIQQQgghiqXICc769evZuHEjq1ev5p133qFKlSrMmTOnNGITQgghhLggRb6KKjU1laNHj3LkyBEyMzOpUaNGkXf6559/MnToUAD279/PbbfdxuDBg5kyZQq6rgPw6quv0r9/fwYNGsRff/1V5H0IIYQQ4uJV5B6c4cOHc/XVV3PbbbcRFxdHVFRUkbZ/8803+fzzzwkKCgK8d0YeNWoUHTp0YPLkyaxcuZIaNWrw22+/sWTJEhISEhg5ciTLli0raqhCCCGEuEgVuQdn+PDhfPbZZ7zzzjsMHDiQzz77rEjb16lTh7lz5/rf//3337Rv3x6Abt26sXbtWjZu3Eh8fDyaplGjRg08Hg9JSUlFDVUIIYQQF6ki9+C89957LF++nJCQEBwOB3feeSc33nhjobfv2bMnhw4d8r9XSqFpGgAhISGkpaXhcDiIiIjwr5NTfnpvkd1uxWQyFvUQCs1oNBAREVxq9YuCpM0DQ9q97EmbB4a0e9kLVJsXOcHRNI2QkBAA7HY7Vqu1WAEYDLmdSOnp6YSFhWG320lPT89XHhoaWmBbhyO7WPs+n4iIYJKTM0p1HyI/afPAkHYve9LmgSHtXvZKu81jYgrmB3ABp6hq167NrFmzWLFiBbNmzaJOnTrFCqxp06b+e+qsWrWKdu3a0aZNG1avXo2u6xw5cgRd14s81kcIIYQQF68i9+A888wzLF68mHXr1hEXF8fjjz9erADGjh3LpEmTmDNnDnFxcfTs2ROj0Ui7du0YOHAguq4zefLkYu1DCCGEEBcXTSmlCrPiCy+84B8rA96xM+A9ZfXoo4+WTnTnceJEWqnWL12ZZU/aPDCk3cuetHlgSLuXvUCdoip0D05cXFyJBSOEEEIIUZoKneAkJCTw4IMPAnD8+HGqVq1aakEJIYQQQhRHoQcZ//rrr/754o67EUIIIYQoTYVOcPIO1SnksB0hhBBCiIAodIKTd4Bx3nkhhBBCiPKm0GNw/v77bwYNGoRSit27d/vnNU3jo48+Ks0YhRBCCCGKpNAJzueff16acQghhBBClJhCJzg1a9YszTiEEEIIIUpMkR/VIIQQQghR3kmCI4QQQohKRxIcIYQQQlQ6kuAIIYQQotKRBEcIIYQQlY4kOEIIIYSodCTBEUIIIUSlIwmOEEIIISodSXCEEEIIUelIgiOEEEKISkcSHCGEEEJUOpLgCCGEEKLSkQRHCCGEEJWOJDhCCCGEqHQkwRFCCCFEpSMJjhBCCCEqHUlwhBBCCFHpSIIjhBBCiEpHEhwhhBBCVDqS4AghhBCi0pEERwghhBCVjiQ4QgghhKh0JMERQgghRKUjCY4QQgghKh1ToAM4G13XmTp1Kjt27MBisTBjxgzq1q0b6LCEEEIIUQGU2x6cFStW4HQ6Wbx4MY899hizZs0KdEhCCCGEqCDKbQ/Oxo0b6dq1KwCtWrVi69atBdaJiQkt9TjKYh8iP2nzwJB2L3vS5oEh7V72AtHm5bYHx+FwYLfb/e+NRiNutzuAEQkhhBCioii3CY7dbic9Pd3/Xtd1TKZy2+EkhBBCiHKk3CY4bdq0YdWqVQBs3ryZhg0bBjgiIYQQQlQUmlJKBTqIM8m5imrnzp0opXj66aepV69eoMMSQgghRAVQbhOcQJHL0wPD5XLx5JNPcvjwYZxOJw888ABXXXVVoMO6KCQmJnLzzTfzzjvvyB8RZWTBggX88MMPuFwubrvtNgYMGBDokCo1l8vFuHHjOHz4MAaDgenTp8t3vZT9+eefPP/887z//vvs37+fcePGoWkaDRo0YMqUKRgMpX8CqdyeogoUuTw9MD7//HMiIiL48MMPeeutt5g+fXqgQ7oouFwuJk+ejM1mC3QoF43169ezadMmFi1axPvvv8/Ro0cDHVKl9/PPP+N2u/noo4946KGHeOmllwIdUqX25ptvMnHiRLKzswF45plnGDVqFB9++CFKKVauXFkmcUiCc5rCXJ4uSl6vXr145JFHAFBKYTQaAxzRxWH27NkMGjSIqlWrBjqUi8bq1atp2LAhDz30EPfffz9XXHFFoEOq9GJjY/F4POi6jsPhkAtWSlmdOnWYO3eu//3ff/9N+/btAejWrRtr164tkzjkUz7N2S5Pl38QpSskJATwtv/DDz/MqFGjAhvQRWD58uVERUXRtWtX3njjjUCHc9E4deoUR44cYf78+Rw6dIgHHniA//3vf2iaFujQKq3g4GAOHz5M7969OXXqFPPnzw90SJVaz549OXTokP+9Usr//Q4JCSEtLa1M4pAenNPI5emBk5CQwB133MGNN95I3759Ax1Opbds2TLWrl3L0KFD2bZtG2PHjuXEiROBDqvSi4iIID4+HovFQlxcHFarlaSkpECHVaktXLiQ+Ph4vv32Wz777DPGjRvnP30iSl/e8Tbp6emEhYWVzX7LZC8ViFyeHhgnT55k2LBhPPHEE/Tv3z/Q4VwU/vvf//LBBx/w/vvv06RJE2bPnk1MTEygw6r02rZtyy+//IJSimPHjpGZmUlERESgw6rUwsLCCA313kk3PDwct9uNx+MJcFQXj6ZNm7J+/XoAVq1aRbt27cpkv9I1cZprrrmGNWvWMGjQIP/l6aL0zZ8/n9TUVF5//XVef/11wDtQTQa/isqmR48ebNiwgf79+6OUYvLkyTLmrJTdddddPPnkkwwePBiXy8Xo0aMJDg4OdFgXjbFjxzJp0iTmzJlDXFwcPXv2LJP9ymXiQgghhKh05BSVEEIIISodSXCEEEIIUelIgiOEEEKISkcSHCGEEEJUOpLgCCGEEKLSkQRHCFGuPPzwwyxYsMD/3uFw0LNnT7Zv3x7AqIQQFY1cJi6EKFeSkpK45ZZbePPNN6lfvz6TJ0/m0ksvZdiwYYEOTQhRgciN/oQQ5UpUVBSTJk1i4sSJjB492v+8puHDh5OdnY3VamX69OlUr16dF154ga1bt5KcnEzjxo155plnmDt3Lps2bSIjI4OZM2fy3HPP4XA4yMzMZPTo0cTHxwf6EIUQZUASHCFEuXPllVfy/fffM378eBYtWsQzzzzD0KFD6d69O+vWreP5559n2rRphIWF8e6776LrOtdffz3Hjh0DIC4ujokTJ7Jr1y6Sk5N56623SExMZN++fYE9MCFEmZEERwhRLvXr14+srCyqVavGzp07WbBgAW+99RZKKUwmk/8hlY8++ijBwcFkZGTgcrkAiI2NBaBBgwYMHDiQRx99FLfbzdChQwN5SEKIMiQJjhCi3IuLi2PYsGG0adOGPXv2sGHDBlatWkVCQgIvvfQSSUlJfP/99+QMKcx5evGOHTv4/3bu0AZCMAbAaNkDElZhCxJ2wCIwLAGCRdjp9yRgyOlbgLs07/mKus+053nGvu9RSom+76Prul+uArxE4AB/b5qmWJYl7vuO67pinueo6zq2bYthGKKqqmiaJkopX3Nt28a6rnEcRzzPE+M4/mgD4G2uqACAdPzBAQDSETgAQDoCBwBIR+AAAOkIHAAgHYEDAKQjcACAdD5AuDG75gZyewAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 576x432 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot Ts and Ta\n",
    "fig, (ax1, ax2) = plt.subplots(2, 1, figsize = (8, 6))\n",
    "ax1.set_xlabel(\"Years\")\n",
    "ax1.set_ylabel(\"Temperature (°C)\")\n",
    "ax1.set_title(\"Surface and atmospheric temperatures in the one-layer climate model\")\n",
    "ax1.plot(t, ts, c=\"tab:red\")\n",
    "ax1.plot(t, ta, c=\"tab:blue\")\n",
    "ax1.legend((\"T surface\",\"T atmosphere\"))\n",
    "ax1.grid()\n",
    "\n",
    "# Plot the Flows\n",
    "ax2.set_xlabel(\"Years\")\n",
    "ax2.set_ylabel(\"Flow of energy ($W/m^{2}$)\")\n",
    "ax2.set_title(\"Flows of energy to and from the surface in the one-layer model\")\n",
    "ax2.set_ylim(0, ir_norm.max() + 10)\n",
    "ax2.plot(t, ir_norm, c=\"tab:purple\")\n",
    "ax2.plot(t, solar_to_earth_norm, c=\"tab:orange\")\n",
    "ax2.plot(t, ir_down_norm, c=\"tab:green\")\n",
    "ax2.grid()\n",
    "ax2.legend((\"IR\", \"Solar to earth\", \"IR down\"))\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Final Temperature, T=16°C\n"
     ]
    }
   ],
   "source": [
    "print(f\"Final Temperature, T={round(ts[-1])}°C\")"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "7b11ec4e63c99d8162a346e44061eab97ff71e663e3c205d951a4ad2f12c3cfa"
  },
  "kernelspec": {
   "display_name": "Python 3.7.10 64-bit (system)",
   "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.7.10"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}