Tutorial4_OnLine.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tutorial 4: Continuous degradation modelling, Prognostics and Remaining Useful Life\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "## IMPORT\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Description of the model and the assumptions\n",
    "\n",
    "- We consider an item which is inspected every 15 months (perInsp$=15$) and the total number of inspections is 35 (nInsp$=35$).\n",
    "\n",
    "\n",
    "- We consider the degradation model $X(t)$ which is deterministic:\n",
    "$$\n",
    "X(t) = 0.01\\,t+0.01\\,t^2\n",
    "$$\n",
    "\n",
    "\n",
    "- The observation $Y(t)$ of $X(t)$ is defined as $Y(t) = X(t)+\\epsilon(t)$ where $\\epsilon(t)$ is a Gaussian noise (ie. a noise having a normal pdf with mean and standard deviation respectively equal to 0 and $\\sigma$).\n",
    "\n",
    "\n",
    "- The item is considered to be failed when $Y(t)$ reaches a level L.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first define all the parameters we need"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "## PARAMETERS\n",
    "\n",
    "# `mdl`is a dictionary to define the model Y(t)\n",
    "mdl = {}\n",
    "# Parameters for X(t) function\n",
    "mdl['para'] = np.array([0.01, 0.01])\n",
    "# Mean of the noise\n",
    "mdl['muEps'] = 0\n",
    "# Standard deviation of the noise\n",
    "mdl['sigEps'] = 100\n",
    "\n",
    "# Period of inspection\n",
    "perInsp = 15\n",
    "# Number of inspections\n",
    "nInsp = 35\n",
    "# Vector time of inspections\n",
    "vTime = np.arange(0, (nInsp+1)*perInsp, perInsp, dtype=float)\n",
    "\n",
    "# Failure level\n",
    "L = 5000\n",
    "# Time to end the simulation of one path of Y(t) - Should be large enough to be sure that Y(t) overpasses L\n",
    "tEnd = 1000\n",
    "# Prediction horizon for the RUL: from the last inspection date `vTime[-1]+1` to `tEnd` with time step `1`\n",
    "vTimeProj = np.arange(vTime[-1]+1, tEnd, 1, dtype=float)\n",
    "\n",
    "# Number of trajectories for projection simulation\n",
    "nRnd = 1000\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We creat a function GetState that simulates a given number of paths of $Y(t)$ \n",
    "\n",
    "- First we calculate the values of $X(t)$ at any times in vTime. This is put in vStateExact\n",
    "\n",
    "\n",
    "- Second we duplicate $X(t)$ to prepare for the number of simulated paths we want (nbHist)\n",
    "\n",
    "\n",
    "- Third we add the random noise at any times to all duplicated $X(t)$, in order to get nbHist path of $Y(t)$. This is put is mState. This is a matrix with nbHist columns. One path is in a column. The number of lignes is the number of inspections.\n",
    "\n",
    "\n",
    "- As output, we keep  vStateExact and mState.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# GetState is the function that generates `nbHist` paths of Y(t)\n",
    "def GetState(mdl, vTime, nbHist):\n",
    "    # Duplicate vTime\n",
    "    mTime = np.tile(np.asmatrix(vTime).T, (1, mdl['para'].shape[0]))\n",
    "    # Calculate the power of t \n",
    "    mTime = np.power(mTime, range(1, mdl['para'].shape[0]+1))\n",
    "    # Calculate X(t) by applying a product of matrixes\n",
    "    vStateExact = mTime@(np.asmatrix(mdl['para']).T)\n",
    "    # Duplicate X(t)\n",
    "    mStateExact = np.tile(vStateExact, (1, nbHist))\n",
    "    # Simulate Y(t)\n",
    "    mState = mStateExact+np.random.normal(loc=mdl['muEps'], scale=mdl['sigEps'], size=mStateExact.shape)\n",
    "    # mState[mState<0] = 0\n",
    "    return vStateExact, mState"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use GetState to simulate one path of $Y(t)$\n",
    "\n",
    "- $Y(t)$ is in vObs for any time $t$ in vTime, that is to say from $0$ to the last inspection date, at every inspection date.\n",
    "\n",
    "- $X(t)$ is in vState for any time $t$ in vTime"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "## History simulation: simulate the learning data set - `vState` is X(t) and `vObs` is Y(t)\n",
    "vState, vObs = GetState(mdl, vTime, 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We plot the obtained path of $Y(t)$ in black and the curve for $X(t)$ in red from $0$ to the last inspection date\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Degradation')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4EAAAGpCAYAAAAzyjO+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeZyPVf/H8deZMcZe3BhLfiFuERLinnRnSHFnS4SaQosxqewZa9FYY6TFvmTPIPsSGo11kFTWQoXIFoURs57fH98vTRnM9p3vLO/n4/F9mDnXua7rM93XPQ9v51znGGstIiIiIiIikj14uLsAERERERERST8KgSIiIiIiItmIQqCIiIiIiEg2ohAoIiIiIiKSjSgEioiIiIiIZCM53F2AKxQuXNiWLl3a3WXc5MqVK+TNm9fdZUgmpmdIUkvPkKSWniFJLT1Dklp6hpLm66+//s1aWySxY1kyBJYuXZpdu3a5u4ybhIeH4+fn5+4yJBPTMySppWdIUkvPkKSWniFJLT1DSWOMOXarY5oOKiIiIiIiko0oBIqIiIiIiGQjCoEiIiIiIiLZiEKgiIiIiIhINqIQKCIiIiIiko0oBIqIiIiIiGQjCoEiIiIiIiLZiEKgiIiIiIhINqIQKCIiIiIiko0oBIqIiIiIiGQjCoEiIiIiIiLZiEKgiIiIiIhINqIQKCIiIiIiko0oBIqIiIiIiGQjCoEiIiIiIiIpsG3bNjp16kRERIS7S0kWhUAREREREZFkWrJkCU88+ijxkydTz88vUwVBhUAREREREZEkstYyZcoUurRpQ5i1TAGqx8YSHh7u7tKSLIe7CxAREREREckMjh8/TseOHTm5bh07cubkbqCVhwffensT4ufn7vKSTCOBIiIiIiIit3F99K9y5crk2riR3blyUbxQIY5Mm0aNIUMICwvD19fX3WUmmUYCRUREREREbuH66N+6desYU7483X76CfPAA7ByJVVLlaKquwtMAY0EioiIiIiI/EPC0b+ILVv4tl49uh8+jGnUCLZsgVKl3F1iiikEioiIiIiIJHD8+HEaNWpEQEAAdR58kNP/+Q8PfvkldO0Ky5ZB/vzuLjFVFAJFRERERET4++jf1q1bmTFkCKsjI8kTHg4ffwxjx4Knp7vLTDW9EygiIiIiItlewnf//Pz8mNO9OyUDAyEyElatgkaN3F1imlEIFBERERGRbMtay9SpU+nZsyfx8fGMGzeOwGLF8HjuOShcGLZtg8qV3V1mmtJ0UBERERERyZYSvvtXo0YN9u7ZQ+c//8SjVStH8NuxI8sFQFAIFBERERGRbOaf7/6NGzeOsM8/p8zw4fDWW9CqFYSHQ7Fi7i7VJTQdVEREREREso1/vvs3ffp0ytx9NzRuDGFh0L8/vPsueGTd8bKs+5OJiIiIiIg4bdu2jWeeeYaKFSv+NfoXFkaZ+Hh45BHYtAlmzIAhQ7J0AASNBIqIiIiISBb35Zdf0qBBA+Lj4/Hw8GDBggW0bNkStm6Fp5+G+HhYvx7q1nV3qekia0dcERERERHJ1uLj4+nSpQvx8fEAGGM4dOgQzJ0L9etDwYKwfXu2CYCgECgiIiIiIllYcHAw+/btw8vLC09PT3J6efHCkSPwwgvg6+sIgOXLu7vMdKUQKCIiIiIiWdLChQsZNGgQ7du3Jzw8nGFvv83x//6XUtOnQ4cOsG4dFCrk7jLTnd4JFBERERGRLGf37t20b98eX19fJk2ahPelSzzy1luOzd+HD4egIDDG3WW6hUKgiIiIiIhkKadPn6Z58+YULlyYJUuW4H3kCDRtCqdOwcKFjn0AszGFQBERERERyTKuXbtGixYtuHDhAlu3bsVn5054/nnIlw82boRatdxdotvpnUAREREREckSrLUEBASwfft2Zs+aRbW1a6F5c6hQAb76SgHQSSOBIiIiIiKSJYwaNYrZs2cz7O23eWbpUpgzB9q0genTIU8ed5eXYSgEioiIiIhIprdixQr69OlDYPPm9Fm7FnbsgOBg6N8/2y4AcysKgSIiIiIikqnt27eP559/Hv/772f8V19hLl6ExYuhRQt3l5YhKQSKiIiIiEimde7cOZo2bYp/jhyM//lnTNGisHUrPPigu0vLsLQwjIiIiIiIZErR0dE827IlHX/5hYl//IFHjRqOBWAUAG9LIVBERERERDIday09AgLosnkz/eLi4KWXICwMihZ1d2kZnqaDioiIiIhIpjNj8GA6zpxJFWMgJAS6ddMCMEnkspFAY0wpY8yXxpiDxpj9xpiuzvZBxpiTxphvnZ+nEpzT1xhzxBjzgzGmYYL2Rs62I8aYPq6qWUREREREMr6d779P48GDKZcjB6xcCd27KwAmgytHAmOBntba3caY/MDXxpj1zmPvW2tHJ+xsjKkEtAUeAEoAXxhj/u08PA54AjgBfGWMWW6tPeDC2kVEREREJAPKOWcO1aZN41TOnOTdts3xHqAki8tCoLX2FHDK+fVlY8xBoORtTmkOzLfWRgE/G2OOALWcx45Ya38CMMbMd/ZVCBQRERERyS5iY7nWpQuPTJvGRi8vyu7cSV4tAJMi6bIwjDGmNPAQsMPZ9IYxZo8xZroxpqCzrSTwS4LTTjjbbtUuIiIiIiLZwR9/EN+4MbkmTOBjY8ixfj2lFABTzOULwxhj8gGfAd2stZeMMROAYMA6/wwBXgYSm8RrSTyo2kTuEwAEAPj4+BAeHp4m9aelyMjIDFmXZB56hiS19AxJaukZktTSMyTJlfvECar060fOEycIBPJ06UJla/UcpYJLQ6AxxgtHAJxrrV0MYK09k+D4FGCl89sTQKkEp98D/Or8+lbtN1hrJwOTAWrWrGn9/PzS5odIQ+Hh4WTEuiTz0DMkqaVnSFJLz5Cklp4hSZb16+HNN7kaE8Pj1lKzRw+aNm2qZyiVXLk6qAGmAQettWMStBdP0K0FsM/59XKgrTHG2xhTBigP7AS+AsobY8oYY3LiWDxmuavqFhERERERN7MWPvoI/vc/IgsWpOq1a+Rp1Ij33nvP3ZVlCa4cCawDvAjsNcZ862zrBzxnjKmGY0rnUaATgLV2vzFmAY4FX2KB1621cQDGmDeAtYAnMN1au9+FdYuIiIiIiLtER8Mbb8CUKVx5/HEe2L2bPOXLM3/+fDw9Pd1dXZbgytVBt5D4e36rb3POUGBoIu2rb3eeiIiIiIhkAWfPwrPPwqZNRPXoQe01a7gMhK1YwV133eXu6rKMdFkdVERERERE5LZ27YIaNWDnTuJnzaLlDz/w/aFDLFq0iHLlyrm7uixFIVBERERERNxr1ix49FHw8GDPhAn4TZnCqlWr+Oijj6hfv767q8tyXL5FhIiIiIiISEIxMTFcvnyZS+fPk2/QIArPm8dvVavyScOG9OvYkdjYWHLkyEG1atXcXWqWpBAoIiIiIiLJtnnzZmbPnk358uUpUaIEly5dSvRz8eLFm9quXr1KEWAB4Ae8D7y1Zw9xe/bcuL517gXo6+vrnh8wC1MIFBERERGRZNm2bRv16tUjLi7upmOenp4UKFCAu+66iwIFClCgQAF8fHwoX778je//HRlJ2/nzyX3lCl8FBlKlWTO2FSjAjz/+yMsvv0xMTAw5c+bUfoAuohAoIiIiIiLJMn78+BsB0NPTk65du/LWW29RoEABcufOjWPL8FuYPRsCAqBIEQgL4+Hq1W8cqlWrFqVLlyY8PBw/Pz+NArqIQqCIiIiIiCTZxYsXWbt2LcYYPDw8yJkzJ61ataJYsWK3PzEmBt56Cz74APz8YMECRxD8B19fX4U/F1MIFBERERGRJOvfvz/nz59nypQpnD17NmkjdufOQevWEB4OXbvCqFHg5ZUu9crNFAJFRERERCRJtm/fzvjx43nzzTd55ZVXknbS7t3QogWcOQMzZ0K7dq4tUu5I+wSKiIiIiMgdxcTE0LFjR0qWLMmQIUOSdtLcuVCnDlgLW7cqAGYQGgkUEREREZE7CgkJYd++fSxdupT8+fPfvnNsLPTuDe+/D3XrOt7/K1o0fQqVO1IIFBERERGR2/rxxx8ZPHgwLVq0oHnz5rfv/Ntv0KYNbNgAb74JISF6/y+DUQgUEREREZFbstby2muv4eXlxUcffXT7zt9843j/7/RpmDED2rdPlxoleRQCRURERETklubOncv69ev5+OOPKVmy5K07zpsHr74K//oXbN4MDz+cfkVKsmhhGBERERERSdT58+fp3r07tWvXJjAwMPFOsbHQsyf4+zuC365dCoAZnEYCRUREREQkUW+99RZ//PEHkydPxtPT8+YOv/0GbdtCWBi88QaMGaP3/zIBhUAREREREbnJl19+ySeffEKfPn2oWrXqzR2+/dbx/t+pU/DJJ9ChQ7rXKCmj6aAiIiIiIvI3165do1OnTpQtW5a333775g4zZoCvL8TEON7/UwDMVDQSKCIiIiIifzNs2DAOHz7MunXryJ07918Hrl2Drl1h8mSoXx8+/VT7/2VCGgkUEREREZEbDhw4wIgRI/D39+eJJ57468CxY/Df/zoCYJ8+sHatAmAmpZFAEREREREBID4+noCAAPLnz8+YMWP+OrB2LTz/vGMl0KVL4U4bxkuGppFAEREREREBYOrUqWzdupXRo0dTtGhRiI+Hd9+F//0PSpaEr79WAMwCNBIoIiIiIiKcPn2a3r174+fnR4cOHeDCBXjhBVizBl58ESZOhDx53F2mpAGFQBERERERoVu3bly9epWJEydivvkGWraEkydh/HgIDARj3F2ipBFNBxURERERyeZWr15NaGgoAwYMoMKWLfDII473/zZvhtdeUwDMYhQCRURERESysStXrtC5c2cerFCBfj//DK++6lgFdPduqF3b3eWJCygEioiIiIhkY++88w7m2DG2WIvnJ59A//7w+edQpIi7SxMX0TuBIiIiIiLZ1O7duzk4Zgz7vL3Je+YMLF8OTZu6uyxxMYVAEREREZFsKDYqil1NmrDCWuLLl3fs/3fffe4uS9KBpoOKiIiIiGQ3589zolo1Ak6d4ljduuTYuVMBMBtRCBQRERERyU527SK2WjWKf/89Hz7wAKU3bIDcud1dlaQjhUARERERkezAWpg8GVunDhfOn6eBtzfNVq7EeCgSZDf6X1xEREREJKu7ehVefhk6deJMxYpUvHqVp4cOpXTp0u6uTNxAIVBEREREJCv78UfH5u8zZnCtd29qnjnDvQ89RNeuXd1dmbiJQqCIiIiISFa1aBFUrw7HjsGqVfS4fJlTZ88yefJkcuTQRgHZlUKgiIiIiEhWExUFb74Jzz4LFSvCN9+w7e67mTBhAl26dKFmzZrurlDcSPFfRERERCQr+eknaN0avv4aevSA4cOJBgIaN6ZUqVIEBwe7u0JxM40EioiIiIhkFZ99Bg895HgPcOlSCAkh4uuvadq0Kfv372fcuHHky5fP3VWKm2kkUEREREQks4uKgt694cMPoVYtCA2F0qWJiIigXr16REVF4enpSeHChd1dqWQAGgkUEREREcnMfv4Z/vtfRwDs1g02b4bSpfntt994/fXXiYqKutE1PDzcfXVKhqEQKCIiIiKSWS1Z4pj+eegQLF4M77+P9fJi5syZ3H///ezZs4ccOXLg6elJzpw58fPzc3fFkgFoOqiIiIiISGYTHQ1BQTB2LNSsCQsWQJkyHD58mMDAQDZs2ICvry+TJ0/m8uXLhIeH4+fnh6+vr7srlwxAIVBEREREJDM5ehTatIGdO6FLF3jvPaKN4b0hQxgyZAje3t5MmDCBgIAAPDwcE/8U/iQhhUARERERkcxi2TLo0AGsdawE+swzbN26lYCAAA4cOMCzzz7LBx98QPHixd1dqWRgeidQRERERCSji46Gnj3h6afhvvtg927+qF+fwMBAHn30USIjI1mxYgULFixQAJQ7UggUEREREcnIjh2Dxx6DMWPgjTewW7awYNcuKlasyJQpU+jevTv79++nSZMm7q5UMglNBxURERERyahWrID27SE2FhYs4FitWnRu2ZLVq1dTvXp1Vq5cSY0aNdxdpWQyGgkUEREREcloYmLgrbegWTMoXZrYnTsZ88svVKpUiY0bNzJmzBh27NihACgpopFAEREREZGM5PhxaNsWIiKgc2d2+/vT0d+f3bt306RJEz7++GPuvfded1cpmZhGAkVEREREMopVqxybv+/bx7UZM+ieMycP//e/nDp1ioULF7J8+XIFQEk1jQSKiIiIiLhbVBT07Qvvvw/VqrEhMJAOAwdy4sQJAgMDGT58OHfddZe7q5QswmUjgcaYUsaYL40xB40x+40xXZ3thYwx640xh51/FnS2G2PMh8aYI8aYPcaY6gmu1d7Z/7Axpr2rahYRERERSXeHDsEjj8D77/NT48ZUj4ri8cBA7rrrLrZu3cr48eMVACVNuXI6aCzQ01pbEfgP8LoxphLQBwiz1pYHwpzfA/wPKO/8BAATwBEagXeA2kAt4J3rwVFEREREJNOyFmbOxFavTsyRIwTXqMF9q1bxzcGDeHl5MW7cOHx9fd1dpWRBLguB1tpT1trdzq8vAweBkkBzYKaz20zgaefXzYFZ1mE7cLcxpjjQEFhvrb1grf0dWA80clXdIiIiIiIud+kSV1u2hA4d2BkfT5lLlxiydy/GGADi4+PZunWrm4uUrCpd3gk0xpQGHgJ2AD7W2lPgCIrGmKLObiWBXxKcdsLZdqv2f94jAMcIIj4+PoSHh6fpz5AWIiMjM2RdknnoGZLU0jMkqaVnSFIruz9DsbGxHF+0iIYzZuATFcUAYOV99/FC48YULVqUfv36ERMTQ44cOShQoEC2/m91K9n9GUoLLg+Bxph8wGdAN2vtpev/upFY10Ta7G3a/95g7WRgMkDNmjWtn59fiup1pfDwcDJiXZJ56BmS1NIzJKmlZ0hSK7s+Q99//z2fTJtGvokT6RMZyRkPD6b4+9Pu7bcZ8u9/3+jn6+t747+RpoImLrs+Q2nJpSHQGOOFIwDOtdYudjafMcYUd44CFgfOOttPAKUSnH4P8Kuz3e8f7eGurFtEREREJLUiIyNZuHAh06ZN48etW5kNNAB+feQRii1dymtFitx0jq+vr8KfuJwrVwc1wDTgoLV2TIJDy4HrK3y2B5YlaG/nXCX0P8BF57TRtcCTxpiCzgVhnnS2iYiIiIhkKNZaIiIiePXVVylevDgvv/wyFY8e5UjevDyeKxdMnkyJLVvIkUgAFEkvrhwJrAO8COw1xnzrbOsHjAAWGGNeAY4DzzqPrQaeAo4AfwIvAVhrLxhjgoGvnP3etdZecGHdIiIiIiLJcvbsWWbNmsX06dM5ePAgefPm5fmWLXknKoqSoaFQpQrMnw+VKrm7VBHXhUBr7RYSf58P4PFE+lvg9VtcazowPe2qExERERFJuYiICDZs2IC3tzfbtm1jxYoVxMbG4uvry9SpU2nz0EPk69gRdu+GN96AUaMgVy53ly0CpNPqoCIiIiIiWcW2bdvw8/MjJiYGgIIFC9KtWzdeeuklKlWsCLNnw2OPgbc3LF0KzZu7uWKRv1MIFBERERFJoqtXr/Laa6/dCIAeHh706NGDAQMGwKVL8OKLMHcu1K0Lc+bAPfe4uWKRm7lsYRgRERERkazk1KlT1K1blz179uDl5YWnpyfe3t48/vjjsHMnPPSQ472/4GAIC1MAlAxLI4EiIiIiInfwzTff0LRpU/744w+WLFmCj4+PY7+6xx7Dd/Nm6N8fSpSAjRuhTh13lytyWwqBIiIiIiK3sXjxYl588UX+9a9/sWXLFqpVqwaAb5ky0K4drF8PrVrB5MlQsKCbqxW5M00HFRERERFJhLWWoUOH0rJlS6pUqcLOnTtvBEDWrIGqVWHLFkf4W7BAAVAyDYVAEREREZF/uHbtGi+88AIDBgzA39+f8PBwihUrBn/+CW++CU89BcWKwa5d0LEjmFvtjCaS8SgEioiIiIgkcPr0afz8/Jg3bx5Dhw5l9uzZ5MqVy7HnX82a8PHH0L27YzEYbf4umZDeCRQRERERcfr2229p1qwZ58+f57PPPuOZZ56BuDgYPRoGDoQiRRzvADZo4O5SRVJMI4EiIiIiIsDSpUupU6cO1lq2bNniCIBHj0K9etCnj2PT9717FQAl01MIFBEREZFszVrL8OHDadGiBZUrV2bnzp08VK2aY7P3Bx+Eb7+FmTMdi78UKuTuckVSTSFQRERERLKta9eu0a5dO/r160fbtm0JDw+neK5c8Nxz8OKLjhVA9+xxbAWhxV8ki1AIFBEREZFs6cyZM9SvX585c+YQHBzMvHnzyB0R4Qh+n30Gw4ZBeDiULu3uUkXSlBaGEREREZFs57vvvqNZs2acO3eOhQsX0qppU3jrLQgJgQoVYPt2qFHD3WWKuIRCoIiIiIhkK8uWLcPf35+7776bzZs3UyNnTnj4YceiL507w6hRkCePu8sUcRlNBxURERGRbMFay8iRI2nRogUVK1Zk5/bt1Ni82REAz5yBVatg3DgFQMnyNBIoIiIiIlleVFQUAQEBzJo1i9atWzNjyBByd+gAYWHQrBlMmQJFi7q7TJF0oRAoIiIiIllWREQEK1euZMWKFezdu5dBgwbxdsWKmNq1ISoKJk+GV1/Vyp+SrSgEioiIiEiWERMTw+nTpzl58iQbNmxg0KBBxMTEADCyf396//QTDBoEtWo59gEsX969BYu4gUKgiIiIiGR41louXbrEyZMnb/s5c+YM1tqbzq/r4UHA+PFw8SK88w707w9eXm74SUTcTyFQRERERNwuIiKC+fPn4+PjQ758+fj1119vCnhXrly56bxChQpRsmRJSpYsSbVq1W58XbJkSc6dO0eXwED6RkfTOz6emNy5HYu/+Pq64ScUyTgUAkVERETErSIiIqhbt+6NaZsAXl5elChRgpIlS/Lggw/y1FNP/S3glSxZkhIlSpA7d+5bX3jvXlrdey/5Dh3iTNOm+MydC/nzp8NPJJKxKQSKiIiIiFstXrz4RgD09PQkKCiI4OBgPDxSuJtZbCyMHg1vv02+ggVhyRJ8nn46DSsWydwUAkVERETEbeLi4li/fj3gCIA5c+akSZMmKQ+AP/wA7dvDjh3QqhWMHw9FiqRhxSKZn0KgiIiIiLjN2LFj+e677xgwYAB58uTBz88P35S8sxcfDx9+CH37OjZ7//RTaNNGWz+IJEIhUERERETcYv/+/fTv35/mzZvz7rvvYlIa2H7+GV56CTZuhCZNHHv/FS+etsWKZCEpHGcXEREREUm5mJgY2rVrR/78+Zk0aVLKAqC1MGkSVKkC33wD06fD8uUKgCJ3oJFAEREREUl3w4YNY/fu3SxatAgfH5/kX+DECXjlFVi3Dho0gGnT4P/+L+0LFcmCNBIoIiIiIunq66+/ZsiQIfj7+9OyZcvknWwtzJoFlSvDli2OhV/WrVMAFEkGhUARERERSTfR0dG0a9cOHx8fPvroo+SdfOYMtGjhWP2zalXYswdee02Lv4gkk6aDioiIiLhZREQE4eHhKV8ZMxOZNm0aBw4c4PPPP6dgwYJJP3HhQkfgi4yEkBDo2hU8PV1XqEgWphAoIiIi4kYRERHUr1+fqKgocuXKRVhYWJYNgps2bWLhwoV06tSJhg0bJu2k8+fh9dchNBQefhhmzoSKFV1bqEgWp+mgIiIiIm4UHh7OtWvXsNZy9epVVq1a5e6SXOLy5ct06NCB4sWLM3r06KSdtGKF492/xYth6FDYtk0BUCQNKASKiIiIuFGtWrX+9v306dM5ePCgm6pxnV69enH06FGCgoLIly/f7TtfvOjY969ZMyhaFL76Cvr1gxyaxCaSFhQCRURERNzo3LlzAHTs2JGpU6cSHx+Pr68vYWFhbq4s7Xz++edMnjyZnj17UrVq1dt3Xr/eMfo3ezYMGOAIgA8+mD6FimQTCoEiIiIibhQaGkqJEiWYOHEir7zyCjt27OCee+6hUaNGTJs2zd3lpdrvv//OK6+8wgMPPEBwcPCtO0ZGQufO8OSTkC8fRERAcDDkzJl+xYpkEwqBIiIiIm5y8eJFVq9ezbPPPouHh+OvZffeey9bt26lfv36vPrqq/Tp04f4+Hg3V5pyb7zxBmfPnmXmzJnkypUr8U7r1sEDD8DEidCzJ+ze7VgERkRcQiFQRERExE2WLVtGdHQ0bdu2/Vv7XXfdxapVqwgMDGTkyJG0bt2aP//8001VptyiRYuYN28eAwYMoEaNGjd3+P13ePllaNgQcud2bP4+erTjaxFxGYVAERERETcJDQ3l3nvvpXbt2jcdy5EjB+PHj2fMmDEsXryYevXqcfr0aTdUmTKnT58mMDCQGjVq0K9fv5s7LFvmGP2bNQv69oVvv4VHHkn/QkWyIYVAERERETe4cOEC69ato3Xr1hhjEu1jjKF79+4sWbKEffv2Ubt2bfbt25fOlSaftZaAgAAiIyOZNWsWXl5eN455/fEHtG0LTz/tWPlz504YNgxuNVVURNKcQqCIiIiIGyxevJjY2NibpoImpnnz5mzatImYmBgeeeQR1q5dmw4VptyMGTNYsWIFw4YNo1KlSo5Ga2HePGq1bw9LljgWffnqK6he3b3FimRDCoEiIiIibhAaGkq5cuV46KGHktS/Ro0a7Ny5k7Jly9K4cWMmTJjg4gpT5tixY3Tt2pXHHnuMbt26ORpPnoTmzcHfn6slS8I33zi2f0gwQigi6UchUERERCSdnT17lg0bNtCmTZtbTgVNzD333MPmzZtp1KgRnTt3pkePHsTFxbmw0uSJj4/n5ZdfJj4+nk8++QQPY2DqVMe7f198ASEh7P7oI7g+OigibqEQKCIiIpLOFi1aRHx8PG3atEn2ufnz52fZsmV06dKF999/n2eeeYbIyEgXVJl848aNY8OGDYwZM4ayxjj2/OvYER56CPbsgR49wNPT3WWKZHsKgSIiIiLpLDQ0lEqVKlG5cuUUne/p6ckHH3zARx99xMqVK3nsscc4efJkGleZPIcOHSIoKIinGjWi49WrULky7NgBEyZAWBiUK+fW+kTkLwqBIiIiIuno119/ZXs6WmUAACAASURBVPPmzcmeCpqYN954g+XLl3P48GFq167Nt99+m0ZVJk9sbCzt2rWjipcXS377DdOtG9StC/v3Q2AgeOivnCIZif4fKSIiIpKOFi5ciLU2RVNBE9O4cWO2bNmCMYZHH32UVatWpcl1k2PU8OH47djBtqtXyfnjj469/1atglKl0r0WEbkzhUARERGRdDR//nyqVatGhQoV0uyaDz74IDt37uT++++nWbNmfPjhh2l27Tv5YcECGr79NiMAj2bN4MABePFFSOUop4i4TpJDoDGmpDHmEWPMY9c/rixMREREJKs5duwY27dvT7NRwISKFy/Oxo0badasGV27duXNN98kNjY2ze9zQ1QUsf36UbZNG0p5eHBp2jTMokVQrJjr7ikiaSJHUjoZY0YCbYADwPV1iC2wyUV1iYiIiGQ5CxYsAKB169YuuX7evHlZtGgRQUFBhISEsGvXLp588kkaNWqEr69v2t1oxw54+WVyHDjALKDorFk08vdPu+uLiEslKQQCTwMVrLVRSb2wMWY60AQ4a62t7GwbBHQEzjm79bPWrnYe6wu8giNkdrHWrnW2NwI+ADyBqdbaEUmtQURERCQjmT9/Pg8//DBly5Z12T08PT0ZPXo0np6evPfee2zfvp3g4GCqVq1K1apVuffeeyldujSlS5fm3nvvpVSpUnh7eyft4pcvQ//+8PHHRBUpQktjKNqhA9MVAEUylaSGwJ8ALyDJIRCYAXwMzPpH+/vW2tEJG4wxlYC2wANACeALY8y/nYfHAU8AJ4CvjDHLrbUHklGHiIiIiNsdOXKE3bt3ExISki73u/vuu/Hw8CA+Ph6AixcvsnHjRk6cOHGjDcAYQ/HixW+EwoQBsXTp0vzf//0fuXPnhqVL4Y034NdfOdGsGY9u2sS1okU5NHZsuvw8IpJ2khoC/wS+NcaEkSAIWmu73OoEa+0mY0zpJF6/OTDfOdL4szHmCFDLeeyItfYnAGPMfGdfhUARERHJVEJDQwF49tln0+V+fn5+eHt7Ex0dTc6cOZk3bx6+vr7ExMRw8uRJjh49ytGjRzl27NiNP7dv387ChQv/9i5hSWCytzdPRUVx7K67mP7kkwxduZK4uDi8vb3Zv39/2k41FRGXS2oIXO78pIU3jDHtgF1AT2vt7zh+v2xP0OeEsw3gl3+0106jOkRERETSTWhoKHXq1KFUOm2b4OvrS1hYGOHh4fj5+d0Ial5eXjdG+xITFxfHr7/+ytEffyTX9OlUXbAAYmOZWq4cIdZy+IsviItzLBERGxtLeHi4QqBIJpOkEGitnWmMyQlcn6L5g7U2JgX3mwAE41hUJhgIAV4GEltD2JL46qU2sQsbYwKAAAAfHx/Cw8NTUJ5rRUZGZsi6JPPQMySppWdIUkvPUMocPXqUvXv38uabb6b7fz9fX1+ioqKSdd+8R47wYEgIBb7/ngs1a3Koe3fKlSjBBGDv3r306tWL2NhYcuTIQYECBZJ1bT1Dklp6hlIvqauD+gEzgaM4AlspY0x7a22yVge11p5JcM0pwErntyeAhP8sdg/wq/PrW7X/89qTgckANWvWtH5+fskpLV1c/5c4kZTSMySppWdIUkvPUMq88847eHh40K9fP4pl5C0UrlyBwYNhzBj4179g7lwKPfcc/0mw55+fnx81a9a8aYQxqfQMSWrpGUq9pE4HDQGetNb+AOBctOVToEZybmaMKW6tPeX8tgWwz/n1cmCeMWYMjoVhygM7cQTO8saYMsBJHIvHPJ+ce4qIiIi4k7WW0NBQ6tatm7ED4Oefw2uvwdGj8OqrMHIkFCqUaFdfX19NARXJxJIaAr2uB0AAa+0hY4zX7U4wxnwK+AGFjTEngHcAP2NMNRxTOo8CnZzX22+MWYBjwZdY4HVrbZzzOm8Aa3FsETHdWrs/6T+eiIiIiHt99913/PDDD3Tv3t3dpSTuzBno1g3mz4f774eNG+Gxx9xdlYi4UFJD4C5jzDRgtvN7f+Dr251grX0ukeZpt+k/FBiaSPtqYHUS6xQRERHJUEJDQ/H09KRly5buLuXv4uNh2jTo3Rv+/NMxDTQoCJK6Z6CIZFpJDYGvAa8DXXBM0dwEjHdVUSIiIiJZwfWpoA0aNKBw4cLuLucvBw5Ap06wZQvUrQuTJkGFCu6uSkTSSVJXB40Cxjg/IiIiIpIEu3bt4ueff2bgwIHuLsXh2jUYOtTxvl/+/DB9OnToACaxhdpFJKu6bQg0xiyw1rY2xuwlka0ZrLVVXVaZiIiISCY3f/58vLy8ePrpp91dCnz5pWP07/BheOEFxwqgRYq4uyoRcYM7jQR2df7ZxNWFiIiIiGQl8fHxLFiwgEaNGlGwYEH3FXL+PPTqBTNmwH33wbp18MQT7qtHRNwusc3Yb0iwnUNna+2xhB+gs+vLExEREcmcIiIiOHHiBG3atHFPAdbCrFmOFT/nzIG+fWHvXgVAEbl9CEwgsd8W/0vLQkRERESykvnz55MrVy6aNWuW/jffu9exzUP79lC+POzeDcOGQe7c6V+LiGQ4tw2BxpjXnO8DVjDG7Enw+RnYkz4lioiIiGQucXFxLFq0iMaNG5M/f/70u/GlS9CjBzz0EBw8CFOnOlYArVIl/WoQkQzvTu8EzgPWAMOBPgnaL1trL7isKhEREZFMbNOmTZw+fTr9poJa69jsvWdPOH0aAgIcq4D+61/pc38RyVRuGwKttReBi8BzAMaYokAuIJ8xJp+19rjrSxQRERHJXEJDQ8mbNy+NGzd2/c0OHoTXX3es/lmjBixdCrVquf6+IpJpJemdQGNMU2PMYeBnYCNwFMcIoYiIiIgkEBMTw6JFi2jatCl58uRx3Y0iIyEoCKpWhW++gfHjYccOBUARuaOkLgwzBPgPcMhaWwZ4HNjqsqpEREREMqkNGzZw/vx52rZt65obWAuffQYVK8J778GLL8KhQ/Daa+Dp6Zp7ikiWktQQGGOtPQ94GGM8rLVfAtVcWJeIiIhIphQaGkqBAgVo1KhR2l/88GFo1AhatXK877dlC0yfrk3fRSRZ7rQwzHV/GGPyAZuAucaYs0Cs68oSERERyXyioqJYvHgxTz/9NN7e3ml34T//dGzxMGoU5MoFH3wAnTtDjqT+VU5E5C9J/c3RHLgGdAf8gbuAd11VlIiIiEhmtG7dOi5evJi2U0GXL4euXeHoUXjhBccU0OLF0+76IpLtJCkEWmuvJPh2potqEREREcnUQkNDKVSoEA0aNEj9xX76yRH+Vq6ESpUgPBzq1k39dUUk27ttCDTGXAbsrY5bawukeUUiIiIimdDVq1dZtmwZbdu2xcvLK+UXunbNMdo3fLhjoZdRoxxhMDXXFBFJ4E77BOYHMMa8C5wGZgMGx5TQ/C6vTkRERCSTWL16NZGRkanbIH7NGnjzTfjxR2jdGkJC4J570q5IERGSvjpoQ2vteGvtZWvtJWvtBKClKwsTERERyUxCQ0MpWrQofn5+yT/52DF45hl46inH6N/69RAaqgAoIi6R1BAYZ4zxN8Z4GmM8jDH+QJwrCxMRERHJLCIjI1m5ciWtWrUiR3JW7LxyBd55B+6/Hz7/3LEC6J49kBbvFIqI3EJSQ+DzQGvgjPPzrLNNREREJNtbsWIFV69eTfpUUGvh008d4e/dd+Hpp+H776FvX0jLrSVERBKR1NVBj+LYJkJERERE/iE0NJQSJUrw6KOP3rnz7t3QpQts3QoPPQTz5sF//+v6IkVEnJIUAo0xuYBXgAeAXNfbrbUvu6guERERkUzh4sWLrFmzhs6dO+PhcZtJVmfPQv/+MG0aFC4MU6bASy853gEUEUlHSZ0OOhsoBjQENgL3AJddVZSIiIhIZrFs2TKio6NvPRU0Otqxymf58jBjBnTvDocOwauvKgCKiFskNQSWs9YOBK5Ya2cCjYEqritLREREJHOYP38+9957L7Vr17754OrVUKUK9OoFderAvn2OQHj33elfqIiIU1JDYIzzzz+MMZWBu4DSLqlIREREJJM4f/4869evp02bNhhj/jrwww/QuLHjA7BqlSMQVqjgnkJFRBJIagicbIwpCAwAlgMHgJEuq0pEREQkE1iyZAmxsbF/TQW9eBF69oTKlWHLFhg9Gvbudez/JyKSQdxxYRhjjAdwyVr7O7AJKOvyqkRERCTDiYiIYO7cuXh7e+Pr6+vucjKE+fPnU65cOR6qWhWmToV+/eC33+CVV2DIEPDxcXeJIiI3uWMItNbGG2PeABakQz0iIiKSAUVERFC3bl1iYmKYM2cOGzZsyPZB8MyZM3z55ZdMfOEFTO3ajq0f6tSBNWugRg13lycicktJnQ663hjTyxhTyhhT6PrHpZWJiIhIhvH5558TE+NYIiAqKorw8HD3FpQBrJ06lTnx8XScNcux/cOnn8LmzQqAIpLhJWmfQOD6foCvJ2izaGqoiIhItnD06NEbX1trKVeunPuKcberV2HUKFoPHgzGwMCB0Ls35M3r7spERJIkSSOB1toyiXwUAEVERLKBX3/9lYULF/LEE0/g7++Pt7c3n332mbvLSn/x8Y7Rvvvvh3feYXl8PJO6dYPBgxUARSRTSdJIoDHmmUSaLwJ7rbVn07YkERERyUgGDx5MbGwsEydO5Pjx45QtW5bg4GB69epFzZo13V3eLUVERBAeHo6fn1/q31/ctMmx199XX0G1aixs2pQ248bxfadOaVOsiEg6Suo7ga8AUwF/52cK0APYaox50UW1iYiIiJt9//33TJs2jcDAQMqWdUwC6tWrF4ULFyYoKAhrrZsrTFxERAR+fn4MGDCAevXqpfwdxkOHoEULqFsXfv0VZs6Er79mzNdfU61aNSpo3z8RyYSSGgLjgYrW2pbW2pZAJSAKqA0Euao4ERERca/+/fuTO3duBgwYcKOtQIECDBw4kA0bNrB+/Xo3Vndrn3zyCdHR0cTHxxMVFUW9evUoV64czZo1o2/fvsyePZvdu3fz559/Jn6Bc+fgzTfhgQfgiy9g6FBHIGzXjqPHj7N9+/a/9gYUEclkkrowTGlr7ZkE358F/m2tvWCMiXFBXSIiIuJm27dvZ/HixQwePJiiRYv+7VinTp0YO3YsQUFBNGjQAA+PpP67sutFRUWxbt06jDF4eHjg6enJiy++yKVLlzhw4ABr1qwhNjYWAGMMZcqUoVKlSlSqVIkq5cvj9913lJw1C3PlCgQEwDvv/G2/vwULHLtmtW7d2i0/n4hIaiU1BG42xqwEFjq/bwVsMsbkBf5wSWUiIiLiNtZagoKCKFq0KD169LjpuLe3N0OGDMHf35/58+fz/PPPu6HKxI0cOZJjx44REhJCVFTUTe8ExsTEcOTIEfbv38+BAwfYv38/B/fvp9CaNXSOi+MeYDnwYYkS5Dp+nEohITzwwANUqlSJihUrMm3aNO655x7OnDlzY4qsiEhmktQQ+DrwDPAoYICZwGfW8SJAPRfVJiIiIm6yZs0aNm3axLhx48iXL1+ifdq2bcuoUaPo378/LVu2xNvbO52rvNmhQ4cYOnQobdq0STS8Anh5eVGxYkUqVqzoaNi40bHoS1wc1ypVYlPr1nzn4UERZ0Bcv3490dHRf7uGMYbHH3+csLCw1C86IyKSzpIUAq211hizC7horf3CGJMHyAdcdml1IiIiku7i4uLo06cP5cqVo2PHjrfs5+HhwciRI2nYsCGTJk2iS5cu6Vjlzay1BAYGkjt3bsaOHXvnE374wbG/3/LlUKoUzJ5Nruef5zEPDx5L0C02NpYff/yRAwcOMH78eMLCwrDWEh0dTXh4uEKgiGQ6Sd0ioiMQABQC7gNKAhOBx11XmoiIiLjD3Llz2bt3L6GhoXh5ed227xNPPMHjjz9OcHAwHTp0oECBAulU5c1mz57Nl19+ycSJEylWrNitO54759jbb+JEyJMHhg2Dbt0gd+5Eu+fIkYMKFSpQoUIFihUrxtatW4mOjiZnzpz4+fm55ocREXGhpL7F/TpQB7gEYK09DBS97RkiIiKS6Vy7do2BAwdSs2ZNWrVqdcf+xhhGjBjBb7/9xujRo9OhwsT99ttv9OjRA19f31uPXl69CiNGwH33OQJgp05w5Aj07XvLAPhPvr6+hIWFERwcrKmgIpJpJfWdwChrbbQxBgBjTA4gY24MJCIiIik2fvx4jh8/zieffJLkFT9r1qxJmzZtCAkJoXPnzrcfhXOR3r17c/HiRSZNmnRz3fHxMG8e9OsHv/wCzZrByJFw//0pupevr6/Cn4hkakkdCdxojOkH5DbGPIFjldAVritLRERE0tsff/zB0KFDadiwIfXr10/WuUOGDCE6Opp3333XRdXdWnh4OJ988gm9evWiSpUq/zwIDz8ML74IRYvCl1/CsmUpDoAiIllBUkNgH+AcsBfoBKwGBtz2DBEREclU3nvvPS5cuMCIESOSfW65cuXo1KkTkydP5tChQy6oLnFRUVF06tSJMmXKMHDgwL8O7NvnGPGrV8/xDuCcObBzJ+gdPhGRpIVAa208sBTobK1tZa2d4tweQkRERLKAkydPMnbsWPz9/alWrVqKrjFw4EBy5crFgAHp9+/EI0aM4NChQ0yYMIE8efLAzz9Du3ZQtSps2gTDhztWAfX3hwy0ob2IiDvd9rehcRhkjPkN+B74wRhzzhjzdvqUJyIiIulh8ODBxMbGEhwcnOJr+Pj40KtXLxYuXMjOnTvTsLrE/fDDDwwbNoznnnuOhtWqQZcuUKECLFwIb70FP/0EffokedEXEZHs4k7/JNYNx6qgD1tr/2WtLQTUBuoYY7q7vDoRERFxue+//55p06bRuXNnypQpk6pr9ezZkyJFihAUFIQrJw1d3xOwWO7cTPbxcaz4OX48vPSSY8XPkSOhUCGX3V9EJDO7UwhsBzxnrf35eoO19ifgBecxERERyeT69etH3rx56d+/f6qvlT9/ft5++23Cw8NZu3ZtGlSXuLlTp1I9PJwfYmPJN3YsNGkCBw7ApElQsqTL7isikhXcKQR6WWt/+2ejtfYccPvdY0VERCTDi4iIYMmSJfTu3ZsiRYqkyTUDAgIoW7YsQUFBxMfHp8k1b4iN5fL771MvMJAQwLtOHfj6a5g/H/7977S9l4hIFnWnEBidwmMiIiKSwVlrCQoKwsfHh+7d0+4tj5w5czJ06FD27NnDvHnz0uai1sKiRVC5Mvl79OCX+Hh+nj4ds3YtVK+eNvcQEckm7hQCHzTGXErkcxmocodzRUREJANbtWoVmzdvZtCgQeTNmzdNr926dWuqV6/OgAEDiIqKSt3FvvgCatWCZ5/lyrVrNAdW9O1LmZdeSpNaRUSym9uGQGutp7W2QCKf/NZaTQcVERHJpOLi4ujTpw/ly5fnlVdeSfPre3h4MHLkSI4dO8aECRNSdpGdO+Hxx+GJJ+DsWaKnTKFmzpzsK1uWAQn3BBQRkWRx2YY5xpjpxpizxph9CdoKGWPWG2MOO/8s6Gw3xpgPjTFHjDF7jDHVE5zT3tn/sDGmvavqFRERyU7mzJnD/v37GTZsGF5ervl33QYNGvDEE08wZMgQLl68mPQTDx6Eli2hdm3Yuxc++AAOHWLYiRN8f/gwEyZMILe2fRARSTFX7po6A2j0j7Y+QJi1tjwQ5vwe4H9AeecnAJgAjtAIvINjW4pawDvXg6OIiIikzLVr1xg4cCC1atWiZcuWLr3XiBEjOH/+PKNGjbpz5+PH4eWXoXJlWL8eBg+GH3+ELl34/uefGT58OM8//zxPPvmkS2sWEcnqXBYCrbWbgAv/aG4OzHR+PRN4OkH7LOuwHbjbGFMcaAist9ZesNb+Dqzn5mApIiIiyTBu3Dh++eUXRo4ciTHGpfeqXr06zz33HGPGjOHUqVOJdzp3Dnr0gPLlYe5c6NrVEf7efhvy57+xJ2CePHkYM2aMS+sVEckOcqTz/XystacArLWnjDFFne0lgV8S9DvhbLtV+02MMQE4RhHx8fEhPDw8bStPA5GRkRmyLsk89AxJaukZksjISAYPHkytWrUAkv08pOQZatKkCQsXLqRTp0706NHjRnuOixcptWAB9yxejEd0NKcbNuRo+/ZE+fjA/v03+q1Zs4aNGzfSq1cvDh48yMGDB5N1f8lY9HtIUkvPUOqldwi8lcT+GdLepv3mRmsnA5MBatasaf38/NKsuLQSHh5ORqxLMg89Q5Jaeoakb9++REZGMnnyZB588MFkn5/SZ2j79u2MHz+eUaNGUaFIERgzBj78ECIjoU0bePttilesSPF/nHfu3DlatmzJo48+ysiRI/HwcOWbLJIe9HtIUkvPUOql92/SM85pnjj/POtsPwGUStDvHuDX27SLiIhIMp08eZKxY8fi7++fogCYGgMGDKBYrlx816IFlCkDQ4dCo0aOhV8+/RQqVkz0vF69enH58mUmTZqkACgikkbS+7fpcuD6Cp/tgWUJ2ts5Vwn9D3DROW10LfCkMaagc0GYJ51tIiIikkyDBg0iPj6e4ODg9L3xxYsUnTiRI3FxtD54kPPVq8N338GCBfDAA7c8bcOGDcyaNYvevXtTqVKldCxYRCRrc9l0UGPMp4AfUNgYcwLHKp8jgAXGmFeA48Czzu6rgaeAI8CfwEsA1toLxphg4Ctnv3ettf9cbEZERETu4ODBg0yfPp0uXbpQunTp9Lnp5cuOKZ8hIfD77+Ro3Jj627YRFx9PeJUqib7zcd21a9cIDAzkvvvuo3///ulTr4hINuGyEGitfe4Whx5PpK8FXr/FdaYD09OwNBERkWynX79+5M2bN30CVWQkfPwxjBoFFy5AkyYwaBA5atSg1fjxvP7666xZs4annnrqlpcYNmwYhw8fZv369doTUEQkjWlyvYiISBa3bds2li5dSlBQEIULF3bdja5ccQS/MmWgb1/HZu87dsCKFVCjBgAdO3akXLly9OnTh7i4uEQvc/DgQUaMGMELL7xAgwYNXFeviEg2pRAoIiKShVlrCQoKolixYnTr1s01N/nzT8dqn2XLQu/eUL06RETA6tXg3IriOi8vL4YOHcrevXuZM2fOTZeKj4+nU6dO5MuXj5CQENfUKyKSzSkEioiIZGErV65ky5YtDBo0iLx586btxa9ehbFjHeGvZ0+oWhW2boW1a+E//7nlaa1ataJmzZoMHDiQa9eu/e3YjBkz2Lx5M6NGjaJo0aK3uIKIiKSGQqCIiEgWFRcXR58+ffj3v/+/vTuPs7H8/zj+us7stuxUlqwxJEo0lMZMZL6J7GXJVwutvpW0WPp+K5RKiyKVhFBkjUIaJtRR+Mkue2VPhoxl1uv3x31mjCWGGeecmXk/H4/7ce5zn/vc53PqYrznuu7rqs7999+fcxc+eRLeew+qVIGnnoLwcFi8GBYsgEaNLvh2l8vF0KFD+eOPPxgxYkTG8QMHDvDMM89w66230qNHj5yrV0RETqMQKCIikkeNHz+eDRs2MGTIEIKCgrJ/wcRErpo5E6pWhd69nceFC53t1lsv6lJRUVHccccdDB48mMOHDwPOmoAJCQlaE1BE5DLT37AiIiJ50IkTJ3jxxRdp0KABbdu2zd7Fjh1zhn1WqUL1d9+Fa66B776D77+Hpk0v+bKvvfYa8fHxDB06lO+++47PPvuM559/npr/sHC8iIjkDIVAERGRPKhv377s2rWL7t27Y8z5VuQ7jyNHYMgQJ/Q99RRUrcrqN96AJUsgOhou9boedevWpUuXLrz11lt06NCBq6++mn79+mXrmiIicmEKgSIiInnMxx9/nHGv3TPPPIPb7b64Cxw8CAMGQMWK0L8/3HSTE/zi4oivXz/b4S+zNm3akJSUxOHDhzl48CCrVq3KsWuLiMi5KQSKiIjkETt27OCee+6hZ8+eGceSkpKIi4vL2gX27HFm+axYEQYPdnr7Vqxwlnq45ZbLUvPmzZszeipTUlKyXquIiFwyhUAREZFc7vDhwzz77LPUqFGDr776ih49ehAWFkZAQADBwcFERkae/wI7d8IjjziLvL/7LrRtC+vXw7RpGYu8Xy6RkZGEhoZmvVYREcm2QF8XICIiIpcmOTmZDz/8kP/9738cOnSI++67j0GDBlGuXDkeeugh4uLiiIyMJCIi4twX2LQJXn0VJk4Elwt69IDnnnPW/fOSiIgIYmNjL1yriIjkGIVAERGRXMZay+zZs3n22Wf59ddfadq0KcOGDaNevXoZ50RERPxzoPrlF2fCl6lTITQUHn8cnnkGypXz0jc43XlrFRGRHKfhoCIiIrnIqlWriI6OpnXr1gB89dVXxMbGnhYA/5HbDS1bQr16MG8ePP+8MxT0nXd8FgBFRMT7FAJFRERygV27dvHvf/+bG2+8kbVr1/L++++zdu1a7rrrrvMvAWGts5h7VBQ0agTLlsErr8Dvvzu9gaVLe+9LiIiIX9BwUBERET+WkJDA0KFDGTZsGKmpqfTt25d+/fpxxRVXnP+N1sLXXzuzfC5bBmXLwrBh0LMnFCrkneJFRMQvKQSKiIj4odTUVD799FMGDBjA/v376dSpE6+++iqVKlU6/xuTk+HLL+H112H1ame5h5EjnUlfQkO9U7yIiPg1hUARERE/8+233/LMM8+wdu1aGjVqxMyZM7n55pvP/6ajR2H0aOf+vt9/hxo1YOxY6NwZgoK8UreIiOQOuidQRETET6xfv56YmBjuuOMOjh07xpdffsnSpUvPHwD37oUXXoAKFeDpp+Gaa+Crr5x1/rp3VwAUEZGzqCdQRETEh9xuN3PmzGH9+vXMnj2bIkWK8Oabb/L4448TEhLyz2/cuBHefBMmTICUFGeB92eegYYNvVe8iIjkSgqBIiIiXmCt5ciRI+zfv599+/axf/9+3G4377//PikpKQB06NCBDz74gBIlSvzTRWDxYif8zZkDYWHw4INOD2CVKl784j6klQAAIABJREFUNiIikpspBIqIiGSB2+0mLi6OyMjI0xY2T0hIyAh1+/btO23/zMfExMR/vH5AQAD16tU7dwBMTYXp0+GNN2D5cihZEl56CR591NkXERG5CAqBIiIiF+B2u2natClJSUm4XC7Cw8M5duwY+/bt4/jx42ed73K5KFWqFGXLlqVMmTLUqFGDMmXKULZs2YxjZcuW5ffff6ddu3YkJSURHBxMZGTk6Rc6fhw+/RTeegu2b4eqVeGDD5x7/cLCvPPlRUQkz1EIFBEROY/U1FT69++f0YuXmprK8ePHufnmm08LdJn3S5YsSUBAwAWvXatWLWJjY8/uYTxwAN5/31na4a+/4OabnV7A1q0hC9cVERE5H4VAERGRf3Do0CG6du3KokWLMkJdcHAwn3322WlDQrMjIiLi1LW2bHEWdB83DhIToVUr6NsXGjfOkc8SEREBhUAREZFzWr16NW3atGHXrl188MEH1KlTh++///6sewJzhNvt9PTNnAnBwXDffdCnD1x7bc5+joiICAqBIiIiZ5kwYQI9e/akWLFiLF68OGOdvkaNGuXchyQnw4wZ8O678OOPUKwY9O8Pjz8OZcrk3OeIiIicQSFQRETEIzk5mT59+vDee+/RpEkTpkyZQpmcDmQHD8LHHzv3++3a5SztMHw43H8/FCyYs58lIiJyDgqBIiIiwN69e+nYsSNLly7lqaeeYujQoQQFBeXcB6xZ44S9iRPh5Em4/XZnps+YGE32IiIiXqUQKCIi+d6PP/5I+/btOXLkCJMmTeLee+/NmQunpsLs2c6Qz7g4Z1mH7t2hd28ID8+ZzxAREblILl8XICIi4ivWWkaMGMFtt91GgQIFWLZsWc4EwPh4Z5bPqlWhTRtnjb/XX3eGf44apQAoIiI+pZ5AERHJl06cOMHDDz/M+PHjadmyJZ999hlFixbN3kU3boT33nOWeDh+HJo0ccJgq1YQqB+5IiLiH/QTSURE8p0dO3bQtm1bVq9ezUsvvcSAAQNwuS5xcExaGsyd69zv9+23EBICnTs7Qz7r1s3ZwkVERHKAQqCIiOQr8+fP595778Vay5w5c/jXv/51aRf6+28YO9bp+du6Fa66CgYNgp49oVSpHK1ZREQkJ+meQBERyRfS0tIYPHgwMTExlC9fnhUrVlxaANy6Ff7zHyhXznksVQo+/xx27nTW+VMAFBERP6eeQBERyfOOHDlC9+7dmTVrFp07d+bjjz+mQIECWb9AWhosWADvvw9ff+3c39exoxMCb7rp8hUuIiJyGSgEiohInrZhwwbatGnD9u3beffdd3niiScwxmTtzX/+CZ9+Ch9+6MzwWbo0DBwIDz8MV155eQsXERG5TBQCRUQkz/ryyy/p0aMHhQoVYuHChdx6660XfpO1sHSps5D7tGmQlOTM8jl4sLPcQ0jI5S9cRETkMlIIFBGRPMXtdrNw4UI2bNjApEmTiIiIYOrUqVx11VXnf+Phw/DZZ846fhs2wBVXOD1+vXppXT8REclTFAJFRCTPcLvdREVFcfLkSQDatm3L559/TnBw8D+/acUKJ/h9/rmztl+DBjBmDHTqBBdz36CIiEguoRAoIiJ5xvTp0zMCoMvlon79+ucOgMeOOaFv1ChYudIJe126OL1+N97o5apFRES8SyFQRETyhJUrV/LJJ58ATgAMCQkhMjLy9JPWrXMmeRk/3lnnr3ZtZ8bPrl2d4Z8iIiL5gEKgiIjkevPnz6ddu3aULFmSESNGsHPnTiIjI4mIiIDERJg61en1W7oUgoOd5R0efhgaNYKszhQqIiKSRygEiohIrjZu3DgefPBBatWqxTfffHNqApitW+HZZ50lHg4ehCpV4I034N//hpIlfVqziIiILykEiohIrmSt5dVXX6V///5ER0czffp0ioSEwJdfwujR8O23EBAArVvDI49AVBS4XL4uW0RExOcUAkVEJNdJTU2ld+/ejBw5ki5duvDpk08SNHAgTJgAhw5B+fLw8svwwANwoaUhRERE8hmFQBERyVVOnDhB586d+X7mTGbcfjutN27E3HSTc69fmzZw//0QHe30AoqIiMhZFAJFRCTX+OvAAV5p2pSOGzbwZWAggd99B/XqwXvvQefOULy4r0sUERHxewqBIiLi/3bs4PA775A4ciTvpKSQVLAggfffDz16OCFQREREskwhUERE/NPx4zB9OowZA4sWUQRYERjI4f/+l/Dnn4fQUF9XKCIikispBIqIiP+wFpYvd4Lf55/D339z/MorGRYczOzixfn0u+8Ir1XL11WKiIjkaj4JgcaYncBRIBVIsdbWN8YUByYD1wA7gY7W2nhjjAHeBf4FHAf+ba39P1/ULSIil8mBA87MnmPGwPr1EBYGHTqwoHx5Wg4dSvUaNZg7dy7lypXzdaUiIiK5ni8XTGpqra1rra3vef48EGutrQbEep4DxADVPFtP4AOvVyoiIjkvORlmz4a2beHqq6FPHyhcGD76CLt3L2/Urk3zwYNpdMstLFmyRAFQREQkh/jTcNDWQKRnfxwQBzznOT7eWmuBZcaYosaYK621e31SpYiIXDprYdkymDgRJk+GgwehdGl48klnkpfwcFJTU3n66acZPnw4HTt2ZPz48YSEhPi6chERkTzDONnKyx9qzA4gHrDAh9baj4wxh621RTOdE2+tLWaMmQO8Zq1d6jkeCzxnrV1xxjV74vQUUqZMmRu/+OILb32dLEtISKBQoUK+LkNyMbUhyS5ftaGwP/6gzHffUea77wjbs4fU4GD+atyY/bffzqEGDbCBzu8kk5KSGDJkCN9//z3t27fnkUceweXy5aAVOZP+HpLsUhuS7FIbypqmTZuuzDTq8jS+6glsbK3dY4wpDSwwxmw6z7nmHMfOSq7W2o+AjwDq169vIyMjc6TQnBQXF4c/1iW5h9qQZJdX29D+/U5v34QJzmQvxjiLuA8ZQkCbNpQuUoTSmU6Pj4/n7rvvZvHixQwbNoynn37aO3XKRdHfQ5JdakOSXWpD2eeTEGit3eN5PGCMmQE0APanD/M0xlwJHPCcvgson+nt5YA9Xi1YRESyJiEBZs50hnsuWACpqc46fsOGwT33wFVXnfNtf/zxBy1atGDLli18/vnn3HPPPV4uXEREJP/wegg0xhQEXNbao5795sDLwFdAd+A1z+Msz1u+Ah43xnwBNASO6H5AERE/kpLiBL6JE2HGDGd9v4oV4bnnoEsXCA8/79vXrl1LTEwMR48eZd68eURFRXmpcBERkfzJFz2BZYAZzsoPBAKTrLXzjDHLgSnGmAeA34EOnvO/wVkeYivOEhE9vF+yiIicJn09v4kT4YsvnCUeihWDbt2c4Ne4MWThXr4RI0bQp08fChcuzJIlS6hTp44XihcREcnfvB4CrbXbgevPcfwvIPocxy3wmBdKExGRC9m2zQl+EybAli0QEgItW0LXrhAT4zw/h/j4eDZs2MDGjRvZsGEDGzZsYNWqVRw44Iz8N8Zw7Ngxb34TERGRfMuflogQERF/tHcvTJ0KkyY5yzsYA7fd5gz3bNcOijoTO1tr+fPAgbPC3oYNG9i3b1/G5cLCwqhZsyZXXnklf/75J9ZakpOTiYuLIyIiwlffUkREJN9QCBQRkbPt3QvTpsGUKbB0KVjLgbJlOfHoo1R47jn2BAQ4AW/cuNMC319//ZVxicKFCxMeHk5MTAw1a9YkPDyc8PBwKlasiMvlwu12Ex0dTVJSEsHBwZrpTURExEsUAkVExHGO4EetWvzWowetPvuMtfv3Y0aNInTsWI4fP57xtmLFilGrVi3atWtHeHh4RuC7+uqr8dz/fU4RERHExsZmTPWtXkARERHvUAgUEcnP/iH48d//QocOuI8coV27duxNTgacIZ916tShW7duGYGvdOnS5w175xMREaHwJyIi4mUKgSIi+c2ePU7w+/LLcwY/wsM5cOAAzz//PJ9++iklS5YkODiY1NRUgoODeeuttxTcREREcjGFQBGR/MAT/OqOHg1r154Kfv/7nxP8atYEICUlhVHvv8/AgQNJSEigb9++DBw4kHXr1mnYpoiISB6hECgikledo8cv6Jprzgp+6X744Qcee+wxVq9eTXR0NO+99x41Pedo2KaIiEjeoRAoInIBbrc79/SC7d4N06efPdTTE/yW799/1iyc+/fv57nnnmPcuHGUK1eOKVOm0L59+0u+z09ERET8m0KgiEgmaWlp7N69m82bN7Nlyxbi4uKYOnUqaWlpBAcHs3DhQho1auTrMk+xFtatg1mznG3FCuf4OYZ6ArB/f8ZuSkoKH3zwAQMHDuT48eM899xzDBgwgEKFCnn3O4iIiIhXKQSKSL5jreXPP//MCHqbN2/O2N+yZQsnT57MODcwMJDU1FQAEhMTiYmJ4d5776Vly5ZERUVRoEAB73+BlBSnly89+O3Y4Rxv2BCGDIG77z5rqOeZli5dymOPPcaaNWto1qwZ7733Htdee60XihcRERFfUwgUkTzjzGGbhw8fzgh2mYPe5s2b+fvvvzPeFxQUROXKlalevTrNmjWjevXqVK9enWrVqvHbb7/RrFkzkpKSCAgI4IYbbmDixIl8+OGHhIaGEh0dTcuWLWnZsiXlypW7fF8uIQHmz3dC39dfw6FDEBIC0dHw/PNw111w5ZUXvMyhQ4fo3r0748ePp3z58kydOpW2bdtq6KeIiEg+ohAoInmC2+0mMjKS5ORkjDEUKVKEw4cPZ7xujKFixYpUr16dbt26ZYS86tWrU7FiRQIDz/3XYbly5c5a0DwxMZHFixczZ84cZs+ezddff80jjzxC3bp1MwLhTTfdhMvlyt6X2rsXZs92gl9sLCQmQvHicOed0Lo13HEHZHHoZkpKCiNHjqRfv34kJSXxwgsv0L9/fwoWLJi9GkVERCTXUQgUkTxh2LBhJCUlAc5wz6pVq9KxY8eMoFe5cmVCQ0Mv6dpnzowZEhJCs2bNaNasGe+88w6bNm3KCIRDhgxh0KBBlC5dmjvvvJOWLVvSrFkzChcufOEPshY2bjw1zPOnn5zjlSrBI484we+WW+AfAus/WbJkCY8//jhr1qyhfv36TJw4kerVq1/UNURERCTvUAgUkVxv1qxZzJgxA5fLhTGG4OBghg8f7pWZPI0x1KxZk5o1a9K3b18OHTrEvHnzmDNnDjNmzODTTz8lODiYyMjIjF7CSpUqnbpAair8+OOp4Ld1q3O8fn145RUn+NWuDZcwXHPfvn307duXCRMmUKFCBaZNm0axYsUUAEVERPI5hUARydW+//57OnXqRP369RkyZAg///yzT5dyKF68OJ07d6Zz584kJyfz448/MmfOHObMmUPv3r3p3bs3Da+9lntLl+bGvXu5cd8+whISsEFBJDZuTNrDDxPQpg3BlSpd0n16brebhQsX8ueff/Lpp59y8uRJ+vfvT79+/ShQoABxcXE5/6VFREQkV1EIFJFca9WqVbRq1YrKlSvzzTffUKJECaKjo31dVoagoCBuu+02bmvShDc6d+avCRM4MXMmZX/9lcBffyUemA7MAuYlJ3M0Lg7i4uCZZzDGEBoaSlhYWJa20NBQ4uPjmTJlCikpKQA0bNiQzz77jGrVqvnuP4KIiIj4HYVAEcmVtmzZQosWLbjiiiuYP38+JUqU8HVJp4uPh2+/hXnznG3fPkoA3HADP0RG0m/xYn5IS8O6XHTs2JH3WrTgxIkTGdvJkydPe37mdujQobOO/f333xkB0OVy0apVKwVAEREROYtCoIjkOnv27KF58+akpaXx7bffUr58eV+XBGlp8H//5wS+uXNh2TLnWLFi0Lw5xMQ4s3mWLYvL7WZ5dDQkJRESHEzv3r1zZPiq2+0mOjqapKQkgoODadq0aQ58MREREclrFAJFJFeJj4/njjvu4ODBgyxatIgaNWr4rpiDB0/19s2fDwcOOBO41K8P/fs7wa9BAwgIOO1tERERZy07kRMu13VFREQkb1EIFJFc4/jx47Rs2ZLNmzfzzTffUL9+fe8WkJoKK1ac6u37+WdnWYcSJZxevpgYp9evdOkLXurMZSdyyuW6roiIiOQdCoEikiskJyfTvn17li1bxpQpU7w3AcyePc5C7em9fX/95fT2NWgA//2vE/xuvPGs3j4RERERf6UQKCJ+Ly0tjR49ejB37lw++ugj2rVrd/k+7NAhZ4bO2FhYuBA2bXKOlyoF//rXqd4+f5uIRkRERCSLFAJFxK9Za3nqqaeYOHEigwcP5qGHHsrZD0hIgCVLnMC3cCGsWuUM8SxYEJo0gQcfhKZNoW5dcLly9rNFREREfEAhUET82uDBgxk+fDhPPfUUL7zwQvYvmJjozNyZ3tP300+QkgLBwdCoEbz0EkRFOcM9g4Ky/3kiIiIifkYhUET81qhRoxg4cCDdunXjzTffxBhz8RdJTYWVK0/19C1dCidOOL169evDM89AdLQTAAsUyPkvISIiIuJnFAJFxC9NmTKFRx99lDvvvJNPPvkEV1aHYloL69c7gS82Fr7/Ho4ccV6rXRt69nR6+m67Da644vJ9ARERERE/pRAoIn5nwYIFdO3alcaNGzNlyhSCzjcsMzUV1q51eviWLHEmdTlwwHmtShXo2NHp6YuMhDJlvFG+iIiIiF9TCBQRv/Lzzz/Tpk0batasyezZsylw5hDN48ed9fmWLnW2H3+Eo0ed18qVg2bNnNAXFQUVK3r/C4iIiIj4OYVAEfEbGzduJCYmhjJlyjBv3jyKFi0KBw/CDz+cCn0rV0JysvOG2rWhSxe49Va45RaoUMG3X0BEREQkF1AIFBG/8Pvvv9O8WTOqGsPshx+m9IsvOqEvfZ2+4GBnxs4+fZzA16gRFCvm26JFREREciGFQBHxnZQUWLOGhHnzWP/aa/yckMCV1sKzz0LRotC4MXTv7vT03XgjhIb6umIRERGRXE8hUES858ABWL7cuafP7Xa2hAQKAeFAQLNm0Lat09MXHq7F2UVEREQuA4VAEbk8EhLg//7PCXzp22+/Oa+5XFC7NqlduvDq0qWM3riR92bO5K677vJtzSIiIiL5gEKgiGRfcjKsW3d64NuwAdLSnNevuQYaNoQnnnDu66tXj9SwMDp16sS09esZN26cAqCIiIiIlygEisg/crvdxMXFERkZSUREhHPQWti27VTYW77c6fE7edJ5vUQJJ+i1a+c83nQTlCp12nV//PFHnnzySZYvX86wYcO47777vPzNRERERPIvhUAROSe3201UVBRXnDzJTy4XbcqV47qTJ6kWH09hzxINJ10uthQuzMaSJdlQqBDrCxZkV2AgaX/9ReqcOaTOmkVaWhqpqakZj8eOHWPPnj1YawkMDDwVLkVERETEKxQCRcSRmgpbtsAvv8Dq1ZSaNo1tJ09yFUBaGim//87W0FAWFCnC+kKF2FS4ML8VLAiBgQQEBOByuQgICKCw5zH9eeZ9l8vFpk2b2L17NwDWWuLi4hQERURERLxIIVAkPzp6FNasgdWrM0Ifa9fCiRPO60FBBBYsyLfAGmP4JSiIV+fNo2HTptQA2mbjo91uN9HR0SQlJREcHExkZGT2v4+IiIiIZJlCoEheZi388cepoJf+uG3bqXOKF4e6deHhh53H66/n6+3badm2Le3atePGG29kcGQkDXOoty4iIoLY2Niz7zUUEREREa9QCPSSJUuWMHr0aEJCQvSPXrk8EhOdGTkzh73VqyE+/tQ51apBvXrQowdcf70T+q6+GozJOGXXrl3cFxVF3bp1mTBhAqGXYYH2iIgI/TkQERER8RGFQC9wu900a9aMxMREJk2aRPfu3XnooYdo2LAhAQEBvi5PcpuEBNi0CTZuPLVt2OD07qWmOueEhUGdOtChgxP06taF666DQoXOe+mUlBQ6d+5MYmIikydPviwBUERERER8SyHQC+Li4khJSQGciTDGjRvH2LFjKVasGM2bNycmJoYWLVpQpkwZH1cqfuWvv5xwlyns3fzLL7B//6lzAgOhalWoVQvat3eCX926zrFL+AXDyy+/zJIlSxg/fjzVq1fPwS8jIiIiIv5CIdALIiMjCQ4OJjExkZCQEGbMmMGRI0eYO3cu8+bNY/LkyQDceOONxMTEEBMTo17C/MJa2L37VG9e5t69P/88dV5YGNSowZHatQl94gmoWdPZqlaFoKAcKWXhwoUMGjSI7t27061btxy5poiIiIj4H4VAL0ifCGPMmDHcf//9GfdCdezYkbS0NFavXs3cuXOZO3cuQ4YMYdCgQeolzGvi453hmunbr786QW/TJmemznTFijnhrlUr5zE83HmsUAFcLjbGxVHmMsymeeDAAbp06UL16tV5//33c/z6IiIiIuI/FAK9JCIigsTExLMmw3C5XNSrV4969erRr18/4uPjWbBgwVm9hDfccMNpvYSBgfpf51fS0pwevW3bYPv20wPftm2nT84CcNVVTrjr3v30sFe69GmTtHin9DTuu+8+4uPjmT9/PoUucN+giIiIiORuShJ+plixYnTs2PGcvYSvvfYagwcPplixYjRr1oyYmBhKlizJ2rVrNdW+NyQmwo4dZwe8bduc44mJp84NCIBrroEqVeCmm6ByZWe/ShVn34+C1ptvvsn8+fMZOXIkderU8XU5IiIiInKZKQT6sTN7CQ8fPnxaL+GUKVMyzg0JCWHhwoU0atTIhxXncseOOWvq7dp16nHnzlNBb/du5x6+dIUKOaGuZk1o2fJUyKtSxRm+mQt6a5ctW0b//v1p164dDz/8sK/LEREREREv8P9/pUqGokWL0qFDBzp06IC1lt69ezNixAistSQmJtK1a1cmTpyoHsFzOX787ICX+fGPP+Dw4bPfV6aME+qaNj095FWpAqVKeX3oZk6Kj4/nnnvuoVy5cowePRqTi7+LiIiIiGSdQmAuZYyhc+fOfPLJJyQlJeFyuTh8+DCNGjXirrvuYvDgwVx33XW+LvPys9aZWGXfPifQZQ51mYPeoUNnv7dUKShfHipVgiZNoFw553n58s7+1VdDSIj3v5MXWGt58MEH2b17N0uXLqVo0aK+LklEREREvEQhMBdLn3U0Li6OyMhI6tSpw/Dhwxk6dCjXX389nTt35uWXX6Zy5cq+LvXiJCfDwYNw4ICz7d9/7v3055nvxUtXsqQT5CpWhMaNTw935cs7AS8fL4T+wQcfMH36dF5//XUaNmzo63JERERExIsUAnO5iIiI04Z/vvDCC/Tq1YvXX3+d4cOHM3nyZB566CEGDhzIlVde6f0CrXWGYh4+7Gzx8U7AO1+w++uvc18rONiZPTN9q1Xr1H6ZMqcHvLAw737PXOSXX37h6aefJiYmhj59+vi6HBERERHxslwTAo0xLYB3gQBgtLX2NR+X5LeKFy/Oa6+9xn/+8x9eeeUVPv74Y8aOHUvv3r157rnnKFasWNYvZq0zYUp6iMvqduTIqf2UlH++frFip4JceLhz713moFemzKn9K67I1ffg+YOEhAQ6depE8eLFGTduHC6Xy9cliYiIiIiX5YoQaIwJAEYAzYBdwHJjzFfW2g2+rcwH0tIgNdUJVqmpzpaUBCdOOD1u6Y/Hj3PliROMbNKE/1WowDfTprFh6FBGvfsukQ0aUD88nKDk5NPOP2v/6FEnxKWmnr+mAgWgaNFTW5kyUL366ccyb6VKOaGuVCmnd0+85rHHHmPLli3ExsZSqlQpX5cjIiIiIj6QK0Ig0ADYaq3dDmCM+QJoDeSeELhiBTc8+qgTmM4Mcefa/6fXMy9RkEWlgX+nPzl5kpOLF3NsyRICixShYIkSmIIFnbrCwpwAl75fqJDTU/dPYa5oUad3TkHuorjd7oz7OL05k+v48eMZP348L774Ik2bNvXa54qIiIiIf8ktIfBq4I9Mz3cBp81mYYzpCfQEKFOmDHFxcV4rLisKbt1KxbAwkgMDsSEh2IAArMuFDQgAlytjP+Mx/Xj6sUzHz3pPQACpISGkhYY6j2fsp7+WFhJCanAw6zZu5OOPP2bNmjVcWaAA/+7UiejoaAICArL2ZY4dc7bduy/vf7Q85OjRo6xbt45Fixbx3XffYa0lODiYt956i1q1amX5OgkJCZfUtn///Xd69epFnTp1aNKkid/9+RDvudQ2JJJObUiyS21IskttKPuMvYSeJW8zxnQA7rDWPuh53g1oYK194lzn169f365YscKbJWZJeu+PP7DWMn/+fPr168eqVauoXbs2gwYNolWrVlovLpustWzfvp0ffvghY1u/fj0ALpeLtLS0jHNr167NhAkTuP7667N07UtpQydPnqRhw4bs3r2b1atXc/XVV1/U+yVv8ae/hyR3UhuS7FIbkuxSG8oaY8xKa239c72WW2aF2AWUz/S8HLDHR7XkCcYYWrRowYoVK5g8eTJJSUncfffdNGrUiEWLFuF2u3n11Vdxu92+LtXvJScn8/PPP/P222/Tvn17rrrqKqpWrUr37t2ZPHkyFSpUYNCgQRm9gGFhYQQEBBAYGMj27dupW7curVu35nL94qJPnz6sWbOGcePGKQCKiIiISK4ZDrocqGaMqQTsBu4BOvu2pLzB5XLRsWNH2rZty9ixY3nppZeIiorKmDUyJCSE2NhYr9675u8OHz6M2+3O6OX76aefOHHiBACVKlXi9ttvp3HjxjRu3JhatWqdNQNn5rUda9asyfDhw3nnnXe46aabaNGiBQMHDqRRo0Y5Uuu0adMYOXIkffr04c4778yRa4qIiIhI7pYrQqC1NsUY8zgwH2eJiDHW2vU+LitPCQwM5MEHH6Rr1660a9eOb775BnCGEk6fPj1fhkC3282iRYu49tprOXHiBEuXLs0Y2mmtJSAggHr16tGzZ8+M0HfVVVdd8Lpnru344osv8uSTTzJy5EiGDRtG48aNiYp7xiUiAAAMc0lEQVSK4sUXX+S222675Pp37tzJAw88QIMGDRgyZMglX0dERERE8pZcEQIBrLXfAN/4uo68LjQ0lAEDBrBo0SJOnjyJtZa3336bpKQkBgwYkG+WFfjqq69o164dKZnWOCxSpAgRERF07NiRxo0b07BhQwoWLJgjn1ekSBGef/55nnjiCUaNGsUbb7xBZGQkt956KwMHDuT222+/qHs1k5OTuffee7HW8vnnnxOsGVxFRERExCO33BMoXhQREUFsbCyDBw9m1qxZPPDAA4wYMYIqVaowaNAgjh075usSL5vffvuNRx99lLZt22YEQJfLRe/evTl06BDz5s1j4MCBREVF5VgAzKxgwYL06dOHHTt2MHz4cLZv307z5s2JiIjg66+/JqsTOQ0YMIBly5YxevRoKleunON1ioiIiEjupRAo5xQREcELL7xAq1at+PDDD1m3bh3R0dEMHDiQatWq8dFHH53WS5bbbd26lQceeICqVasyevRo7rzzTkJDQwkICCAkJIR77rkn60to5ICwsDCeeOIJtm3bxqhRo9i3bx8tW7akV69ezJgx47QZRs80b948Xn/9dXr16kWHDh28VrOIiIiI5A4KgZIlNWrUYMaMGSxdupRKlSrRq1cvateuzYwZM7LcO+WPNmzYQNeuXbn22muZNGkSjzzyCNu2bWPWrFksXLiQV155xacT44SEhNCrVy+2bNnCmDFjOH78OG3btqVu3bpMmTKF1NTU087fs2cP9913H9dddx1vv/22T2oWEREREf+mECgXpXHjxixdupSZM2dijKFt27bccsst/PDDD74u7aL88ssvtG/fntq1azNz5kyefvrpjCGY5cs7q5Gk94b6w6Q4QUFB9OjRg3HjxjFhwgSSk5Pp1KlTxjqDKSkppKam0rVrV44dO8bkyZMJCwvzddkiIiIi4ocUAuWiGWNo3bo1a9eu5aOPPmLHjh3ccssttG7dmo0bN/q6vPP66aefuOuuu6hXrx4LFiygX79+7Ny5kzfeeIOyZcv6urwLCggIoEuXLqxbt47JkycTFBREt27dqFGjBnfccQeLFi3iySefpGbNmr4uVURERET8lEKgXLLAwEAeeughtmzZwuDBg4mLi6N27do89NBD7Nmzx9flnWbx4sU0b96cm2++mR9//JGXX36Z3377jUGDBlGyZElfl3fRAgIC6NixI7/88gszZswgICCA2NhYAN5++23cbrePKxQRERERf6UQKNlWsGBB+vXrx7Zt23jiiScYN24cVatWpX///hw5csRndVlrWbBgAU2aNOG2225j9erVvP766+zcuZOBAwdStGhRn9WWU1wuF3fffTfdu3fPWJQ+KSmJuLg43xYmIiIiIn5LIVByTMmSJXnnnXfYtGkTbdq0YciQIVSpUoV33nmHxMREr9VhrWXOnDlERETQvHlztm/fzrvvvsuOHTvo27cvhQsX9lot3tK0aVNCQkIICAggODiYyMhIX5ckIiIiIn5KIVByXOXKlZk4cSIrV66kXr16PPXUU9SoUYNJkybxww8/8Oqrr16W4YppaWlMnTqVevXqcdddd7F//35GjRrFtm3b6N27NwUKFMjxz/QX6Ws7+no2UxERERHxf4G+LkDyrhtuuIEFCxawYMECnn32Wbp06YIxBnDuaevQoQMVKlQgMDAwW9vmzZtZvHgxK1as4LfffqN69eqMHTuWzp07ExQU5OP/Ct4TERGh8CciIiIiF6QQKJdds2bNWLlyJffeey9TpkwBICUlhalTp+JyuTKWN8guYwwvvfQS/fv39+rC7iIiIiIiuYlCoHiFy+XiySefZPbs2SQlJREcHHzasEVrLampqaSkpFxwS05OztgfM2YMH330EWlpabhcLoKCghQARURERETOQyFQvCb9vrW4uDgiIyNPG7pojMkY3nkxUlJSGDduXEaw1IQoIiIiIiLnpxAoXpXT962dL1iKiIiIiMjZFAIl19OEKCIiIiIiWaclIkRERERERPIRhUAREREREZF8RCFQREREREQkH1EIFBERERERyUcUAkVERERERPIRhUAREREREZF8RCFQREREREQkH1EIFBERERERyUcUAkVERERERPIRhUAREREREZF8RCFQREREREQkH1EIFBERERERyUcUAkVERERERPIRY631dQ05zhjzJ/Cbr+s4h5LAQV8XIbma2pBkl9qQZJfakGSX2pBkl9pQ1lS01pY61wt5MgT6K2PMCmttfV/XIbmX2pBkl9qQZJfakGSX2pBkl9pQ9mk4qIiIiIiISD6iECgiIiIiIpKPKAR610e+LkByPbUhyS61IckutSHJLrUhyS61oWzSPYEiIiIiIiL5iHoCRURERERE8hGFQBERERERkXxEIdBLjDEtjDG/GmO2GmOe93U94p+MMWOMMQeMMesyHStujFlgjNnieSzmOW6MMcM9bWqNMeYG31Uu/sIYU94Ys8gYs9EYs94Y8x/PcbUjyRJjTKgx5mdjzGpPG3rJc7ySMeYnTxuabIwJ9hwP8Tzf6nn9Gl/WL/7BGBNgjFlljJnjea72IxfFGLPTGLPWGPOLMWaF55h+luUQhUAvMMYEACOAGCAcuNcYE+7bqsRPjQVanHHseSDWWlsNiPU8B6c9VfNsPYEPvFSj+LcUoI+1tiZwM/CY5+8btSPJqkQgylp7PVAXaGGMuRkYCrztaUPxwAOe8x8A4q21VYG3PeeJ/AfYmOm52o9ciqbW2rqZ1gTUz7IcohDoHQ2Ardba7dbaJOALoLWPaxI/ZK1dDBw643BrYJxnfxxwd6bj461jGVDUGHOldyoVf2Wt3Wut/T/P/lGcf4RdjdqRZJGnLSR4ngZ5NgtEAVM9x89sQ+ltayoQbYwxXipX/JAxphxwJzDa89yg9iM5Qz/LcohCoHdcDfyR6fkuzzGRrChjrd0Lzj/wgdKe42pXcl6eYVX1gJ9QO5KL4BnK9wtwAFgAbAMOW2tTPKdkbicZbcjz+hGghHcrFj/zDvAskOZ5XgK1H7l4FvjWGLPSGNPTc0w/y3JIoK8LyCfO9Rstrc0h2aV2Jf/IGFMImAY8aa39+zy/WFc7krNYa1OBusaYosAMoOa5TvM8qg1JBmNMS+CAtXalMSYy/fA5TlX7kQtpbK3dY4wpDSwwxmw6z7lqRxdJPYHesQson+l5OWCPj2qR3Gd/+pAGz+MBz3G1KzknY0wQTgCcaK2d7jmsdiQXzVp7GIjDub+0qDEm/ZfHmdtJRhvyvH4FZw9rl/yjMdDKGLMT5/aXKJyeQbUfuSjW2j2exwM4v4xqgH6W5RiFQO9YDlTzzIwVDNwDfOXjmiT3+Aro7tnvDszKdPw+z4xYNwNH0odISP7luZfmE2CjtfatTC+pHUmWGGNKeXoAMcaEAbfj3Fu6CGjvOe3MNpTettoDC621+g18PmWtfcFaW85aew3Ov3cWWmu7oPYjF8EYU9AYUzh9H2gOrEM/y3KM0Z8z7zDG/AvnN2EBwBhr7WAflyR+yBjzORAJlAT2A/8FZgJTgArA70AHa+0hzz/238eZTfQ40MNau8IXdYv/MMbcAiwB1nLqfpx+OPcFqh3JBRlj6uBMuBCA88viKdbal40xlXF6dooDq4Cu1tpEY0wo8BnO/aeHgHustdt9U734E89w0GestS3VfuRieNrLDM/TQGCStXawMaYE+lmWIxQCRURERERE8hENBxUREREREclHFAJFRERERETyEYVAERERERGRfEQhUEREREREJB9RCBQREREREclHAi98ioiIiHimJo/1PC0LpAJ/ep4ft9Y28klhIiIiF0lLRIiIiFwkY8z/gARr7Zu+rkVERORiaTioiIhINhljEjyPkcaY740xU4wxm40xrxljuhhjfjbGrDXGVPGcV8oYM80Ys9yzNfbtNxARkfxEIVBERCRnXQ/8B7gO6AZUt9Y2AEYDT3jOeRd421p7E9DO85qIiIhX6J5AERGRnLXcWrsXwBizDfjWc3wt0NSzfzsQboxJf08RY0xha+1Rr1YqIiL5kkKgiIhIzkrMtJ+W6Xkap37uuoAIa+0JbxYmIiICGg4qIiLiC98Cj6c/McbU9WEtIiKSzygEioiIeF9voL4xZo0xZgPwsK8LEhGR/ENLRIiIiIiIiOQj6gkUERERERHJRxQCRURERERE8hGFQBERERERkXxEIVBERERERCQfUQgUERERERHJRxQCRURERERE8hGFQBERERERkXzk/wHZaky2pjoj9QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "## History visualization: plot the data set, this is one path (X(t) is red, Y(t) with noise is black)\n",
    "fig = plt.figure(figsize=(15, 7))\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "ax.plot(vTime, vObs, '-k.')\n",
    "ax.plot(vTime, vState,'-r');\n",
    "ax.grid(True)\n",
    "ax.set_xlabel('Time')\n",
    "ax.set_ylabel('Degradation')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We consider now that one simulated path of $Y(t)$=vObs is a vector of observations, which is collected in the past, from time $0$ to the last inspection time. This is a learning data set.\n",
    "\n",
    "We creat a function to use it to estimate the parameters of the model.\n",
    "\n",
    "- We compute mean least square method to estimate parameters of $X(t)$. This is put in mdlEst['para'].\n",
    "\n",
    "\n",
    "- We estimate empirically the parameters of the noise $\\epsilon(t)$. This is put in mdlEst['muEps'] and  mdlEst['sigEps']\n",
    "\n",
    "\n",
    "- As output, we get mdlEst\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Estimation of the model parameters with mean least square method: vObs=mTime.mdlpara\n",
    "def GetMdlEst(vTime, vObs):\n",
    "    mdlEst = {}\n",
    "    mTime = np.tile(np.asmatrix(vTime).T, (1, mdl['para'].shape[0]))\n",
    "    mTime = np.power(mTime, range(1, mdl['para'].shape[0]+1))\n",
    "    # mdlpara is the solution of the mean least square estimation\n",
    "    mdlpara = (np.linalg.inv(mTime.T@mTime)@mTime.T*vObs).T\n",
    "    # mdlEst is for putting mdlpara in the same format as mdl (1 column vector)\n",
    "    mdlEst['para'] = np.asarray(mdlpara)[0]\n",
    "    vObsEst = mTime@(mdlpara.T)\n",
    "    # Estimation of the noise parameters\n",
    "    mdlEst['muEps'] = np.mean(vObs-vObsEst)\n",
    "    mdlEst['sigEps'] = np.std(vObs-vObsEst)\n",
    "    return mdlEst "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We apply the function to the simulated path $Y(t)$ which is in Vobs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Para: [0.14946854 0.00992559]\n",
      "muEps: -10.964145578980158\n",
      "sigEps: 113.19126095933196\n"
     ]
    }
   ],
   "source": [
    "## Model identification: estimate the parameters of X(t): i.e. 0.01 and 0.01 (noise variance is supposed to be known.\n",
    "mdlEst = GetMdlEst(vTime, vObs)\n",
    "print('Para: {}'.format(mdlEst['para']))\n",
    "print('muEps: {}'.format(mdlEst['muEps']))\n",
    "print('sigEps: {}'.format(mdlEst['sigEps']))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We consider now that the model parameters have been estimated on the learning data set with one path and we use the model to make predictions.\n",
    "\n",
    "- We simulate $Y(t)$ from the last inspection date to tEnd  by using the estimated parameters of $X(t)$ and of the noise as input for the function GetState\n",
    "\n",
    "\n",
    "- The number of simulated paths for $Y(t)$ is nRnd. We can simulated as much as path as we want on order to assess the uncertainty of the prediction.\n",
    "\n",
    "\n",
    "- The obtained $X(t)$ is in VstateEst and the obtained $Y(t)$ is in mEst\n",
    "\n",
    "\n",
    "- $X(t)$ and $Y(t)$ are calculated/simulated from the lastinspection date to tEnd with a given time step, that is to say for every times that are given in vTimeProj\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Model prediction\n",
    "\n",
    "vStateEst, mEst = GetState(mdlEst, vTimeProj, nRnd)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As a reference to make comparison, we make prediction with the true values of the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Theoretical prediction\n",
    "\n",
    "# Calculate the values of X(t) with the true parameters from the last inspection date to 750.\n",
    "vStateProj, vObsProj = GetState(mdl, vTimeProj, 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We plot the path from the learning data set and the paths for the prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4YAAAGbCAYAAACYpI5hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3xV9f3H8de5I3uThOywwoqIChXrxNGKddZRt1St1j2qdYuiYp1VHIA4WmsVcbSC2p+zUrVWbVVkCCSE7ITsdW/GXef3R0IKipV77gUy3s/Hg4e533vP53wuHnzw9nvO92uYpomIiIiIiIgMX7bd3YCIiIiIiIjsXgqGIiIiIiIiw5yCoYiIiIiIyDCnYCgiIiIiIjLMKRiKiIiIiIgMc47d3YBVSUlJ5rhx43Z3GyLf4Xa7iY2N3d1tiHyHrk0ZyHR9ykCla1MGsi+++KLRNM20cNQatMFw5MiR/Oc//9ndbYh8x4oVK5g5c+bubkPkO3RtykCm61MGKl2bMpAZhlEerlq6lVRERERERGSYUzAUEREREREZ5hQMRUREREREhjkFQxERERERkWFOwVBERERERGSYUzAUEREREREZ5hQMRUREREREhjkFQxERERERkWFOwVBERERERGSYUzAUEREREREZ5hQMRUREREREhjkFQxERERERkWFOwVBERERERGSYUzAUEREREREZ5hQMRUREREREBpHOzk68Xm9YayoYioiIiIiIDBKmabJu3ToaGhrCWlfBUEREREREZIDzeDy0trbi9/sBsNvtYa3vCGs1ERERERERCbu6ujrq6+v7X7e3t4e1voKhiIiIiIjIAGSaJhs2bKC7u/s77ykYioiIiIiIDAM1NTW43e5dci49YygiIiIiIjLA1NfXU1dXt8vOpxlDERERERGRAaaysnKXnk/BUEREREREZACoq6sjMTGRrq6uXX5uBUMREREREZHdzOPxUFVVRU1NDYFAYJefX8FQRERERERkN+no6MDlchEfHw+wW0IhKBiKiIiIiIjsNkVFRbu7BUCrkoqIiIiIiOwSfr8fj8dDIBBg1apVNDY2WqwT5sbYgWBoGMYzhmHUG4axZquxFMMw3jUMo7jvn8l944ZhGI8YhrHRMIxVhmHss9Uxs/s+X2wYxuytxqcZhrG675hHDMMwwv0lRUREREREdrfa2lo2bNhAQ0MDXq+X8vLyoGvU1ERw7rkTefvt5LD2tiMzhn8EZn1r7AbgfdM0C4D3+14DHAUU9P26EFgIvUESuA2YAewL3LYlTPZ95sKtjvv2uURERERERAa97u7u/kVmrPjkkwTOPnsSFRWRREeH91nEHwyGpml+CDR/a/h44Nm+n58FTthq/E9mr0+BJMMwMoEjgXdN02w2TbMFeBeY1fdegmma/zJN0wT+tFUtERERERGRIcPr9Vo6LhCAxYszufLKcYwdUc8btyzl4IPbwtqb1cVnRpqmWQtgmmatYRjpfePZwNY7MVb1jf2v8artjG+XYRgX0ju7SFpaGitWrLDYvsjO43K5dG3KgKRrUwYyXZ8yUOnalHCKi4sj2CfnWlvtzJkzmk8+SeToo5tYmHILGfe+xur9l4W1t3CvSrq9b2laGN8u0zQXA4sBJkyYYM6cOdNCiyI714oVK9C1KQORrk0ZyHR9ykCla1NCYZomhmHQ1dVFRUUFLpcrqOPXrYvhuuvG0NjoZO7Vq/jZGV4aes7HdcRMzJiYsPZqNRjWGYaR2TdbmAnU941XAblbfS4HqOkbn/mt8RV94znb+byIiIiIiMigY5omLpeLjo4O2tvb6erqCnpvQtOE115L5b77cslK6WDNAWeT9+pnrDv+TwTi4ugsLAx731a3q1gObFlZdDawbKvxc/pWJ90PaOu75fRt4KeGYST3LTrzU+Dtvvc6DMPYr2810nO2qiUiIiIiIjKo1NfXU1RURHt7O263O+hQ2N1tcMcd+cybl8/RhWtYlfhjxq94mdZDDyUQFbWTut6BGUPDMJbQO9uXahhGFb2ri94DvGQYxvlABXBK38f/BvwM2Ah0AucCmKbZbBjGncC/+z53h2maWxa0uZjelU+jgf/r+yUiIiIiIjLodHZ2AuB2u4M+tqoqguuuG0tRUQwPz3qeSz+5FMPvY+MDD9C2k29p/sFgaJrm6d/z1uHb+awJXPo9dZ4BntnO+H+APX6oDxERERERkYFqy8xgsDOEW3z4YSJz5ozCZoOHHyri/D/ej2dkOpvuu4+e3NwfLhCicC8+IyIiIiIiMuT5fD7WrFnD6NGjMQyD4uJiS3X8fli0KIs//CGT/QoqmHNbBekToymZfB/+2FjMnXj76NasPmMoIiIiIiIybPX09OD3+ykrK7McCpuaHFx+eQF/+EMm1x7yDivaZ/DjJ28BwDdixC4LhaBgKCIiIiIiEjS/3w/0zhxa8dVXcZx55iS+XhnL8qPv5t5PjgGbQe355+/Q8Xa73dJ5v49uJRUREREREdlBXq8Xv99PZWWlpeMDAXjuuZEsWJDN2MwWPp90LqPfXE7bj39M6Z134k9K2qE6UWGeTVQwFBERERER2QEej4fVq1dbPr6tzc7tt4/io4+SOOKIZuZe+hXZF31GzYUX9s4UBjELmJiYaLmP7VEwFBERERER+QGmaYYUCteujeGGG8bQ0OBk4S+WMP3qMRjOJNa+9BKB2Nig62VmZlruZXv0jKGIiIiIiMgPqK+vt3ScacJLL6Xxq19NwBHwsPLw87jopTNIf/UVAEuhcGfQjKGIiIiIiMj/UFpaSkdHR9DHud025s3L5513UjjxR2t4xn0miW+vov7UU2k88URLvUycOBGHI/wxTsFQRERERETkW0zTZMOGDbjdbkvHb9wYxXXXjaWqKpJHTljKRR/8GpvPR8k999B6xBFB1ysoKKC5uZnYnTTDqGAoIiIiIiKylUAgQFFRkeVQ+MYbKfzud/nExflZuLCIA2N9eIpzKL3jDnry83e4zujRo2lqaiI9PZ2EhAQSEhIs9bMjFAxFRERERGRYa2tro6WlhVGjRtHe3k59fb2lUNjdbXD//XksW5bKkXsW8dCBf6Bz2sl0MZH1f/wjGMYO1xo1ahQpKSmkpKQE3YcVWnxGRERERESGtY0bN9LU1ITb7aa4uJi2traga1RURHLuuRNZtiyVe496heVV+zPh6YeIqK3t/cAOhsIt+xNGR0cH3UMoFAxFRERERGTYMU2Tmpoauru7+8eKi4st1XrvvSTOPnsSDXUO/nn0b/jt26fiT4hn3bPP4glyW4nk5GSAnbLAzP+iW0lFRERERGRY6Onpoa6ujtzcXLxeL7W1tbS0tPS/7/f7g6rn9RrMn5/Niy+OZI89XLyddBJZb75D06xZVNx0E4GYmB2uZRgG48aNIz4+ntTUVCIiIoLqJVQKhiIiIiIiMiy0tbXR0NDAyJEjMU0T6A2LVtTURHDTTaNZsyaO00+v44orqun++yGUH7gPjSedFNTzhDExMUycOBGj75hdHQpBwVBERERERIaJLTOCzc3N/eFrS0AMxooVicydO4qAH947+mbGjrXR5DyBliOPDKqOzWZj9OjRJCUlBd1DuCkYioiIiIjIkNTS0kJ7ezv5fVtEeDweAGpqaizV83gMHn00myVLRrLv+Gr+mjKbrDffp9n7U5qOPz6oWcKUlBRGjx5tqY+dQYvPiIiIiIjIkNTW1kZTU1P/662fJwxWVVUE558/gSVLRnLTke/wD9cMMv+9gsqrr6Z03rygQiEwoEIhaMZQRERERESGKL/fj2maBAIBbDZb0IvLbPHee0nceecobDaTRbd8xAX3H4cvKYkNTz6Je889g643cuRIS33sTAqGIiIiIiIyJPl8PgDKysosbVjf02Pw8MM5vPxyOlMLW7nzd5VkZcVQ7riRtgMPxB/ks4FZWVlkBrl9xa6iW0lFRERERGRIqa2tZfXq1f3PFLa0tPT/vKO2bFj/8svp3DrrLf7VPpWC+s8BaD7mmKBD4aRJkwZsKATNGIqIiIiIyBDS3d1teXGZLd56K5m7787H6Qiw4oSbOeiN+/CmpmIGuel8bm4uLpeLpKQkYoLY03B3UDAUEREREZFBr6KigtbWVrxer+Ua3d0GDz6Yy1//msZBhRW8FDebjNdW0HrIIZTNmYM/MTGoevHx8aSnp1vuZ1dSMBQRERERkUGvoaEhpOPLyiK5/vqxlJREc+65tcxJepT0Rz+m4tpraTj11B1edTQuLg6XywWA0+kMqaddScFQREREREQGta6urpCOf+ONFO65J4+YKB9/mvMek49LpinwC9z7z6A7yG0lIiMjcblcpKSk4Ajy1tPdSYvPiIiIiIjIoFZaWmrpuK4uG3Pn5nP77aM5cHw568cexunzT8be2go2W9ChECAiIgKAhIQESz3tLoMnwoqIiIiIiHxLeXm5pRnDkpIobrxxDKWlUdx/9Ctc8dnFONrbqbzmGkvPEnZ0dACQkZGB0+kkJSUl6J52JwVDEREREREZdNra2vD5fDQ2NgZ1nGnCsmUjuP/+POJjvXz2syuZ/n+P05Oby/pHHqFr/Pigexk5cmR/MLTZbKSlpQVdY3dTMBQRERERkUEhEAjQ2NhIamoqJSUlmKYZ1PEul4158/J5990U9t23nTvmbmLsg+tonjWLihtuIBDklhKRkZH09PRgt9vZc889MXZwgZqBSMFQREREREQGtJqaGhwOBz09PdTX11NZWRl0jdWrY7n55tHU1UXw2HHPc9DZ0XjT8im7807MIFcPnTRpEg0NDcTExFBTU0NkZOSgWoF0exQMRURERERkQGtsbAx6dnCLQACefTaDRYuyyErv5MtZV7Dn8sU0eWZRdtddQYVCm81GIBAgOjqa/Px8gEF52+j2KBiKiIiIiMiAZJombrfb8qb1DQ1ObrttFJ9/nsAZB65iQetsEt9cScPPf07lNdcEVSs1NbU/DA5FCoYiIiIiIjIgVVRUBL24zBYff5zA7bePorvbxoJfvc4FS8/ECATYdNddtMyatcN1IiIiMAxj0G0/ESwFQxERERERGXC2LDQTLI/H4NFHs1myZCTjx3cyb14pYzJTaKs9mJoLLsCTk7NDdZxOJ16vlz322GNQLyqzoxQMRURERERkQGlpaaG8vDzo48rLI7nppjFs2BDDtUeu4EbXbZSn30kgKpayuXODqlVYWEhPT8+wCIWgYCgiIiIiIrtJRUUF0dHR/Qu4lJSU0NHRgd/vD6qOacKbb6Zw7715RDgDvHfSXA59/W788fFEVlXRNWFCUPXsdjt2u52YILevGMwUDEVEREREZLdobW3F6/USGxvLxo0bLS0y43LZuOeePN56awSHTS3jhejzGPnqB7Ttvz9lt9+OLyUl6JqRkZFBHzPYKRiKiIiIiMhu4fP58Pl8bNiwgUAgEPTxa9bEcPPNY9i8OYKLLqrmro0XkrziIyqvvpr6008Hmy2oeqNHj6anp4fU1NSgexnsFAxFRERERGSXM00T0zRxuVxBHxsIwHPPjWTBgmwyUjt55qFVFO5vo3rzldT9cjZdEycGXTMjI4MUC7OLQ0VwEVpERERERCRI69evZ/Xq1ZimSWtrK6Zp0tbWZqlWY6ODK64Yx6OP5nDKvqtZk3Ygx794FZgm3owMS6EwKSmJrKwsS/0MFZoxFBERERGRncrtdgPQ1tZGSUkJqamplrai+Oc/E5g7dxRut50Xfr6AU969HsM0Kb/xRrCwemhmZibt7e3k5OQMm9VHv4+CoYiIiIiI7BKlpaUAQYfCnh6DRx7JYenSdApHN7Fs74sY+9dXcBcWsmnevB3em3CL6OhoAoEAWVlZw36mcAsFQxERERER2Sk8Hg9VVVX9r60sMLNxYxQ33zyGkpJoTj+9jqvOWkPueR+xefZsqi++GBzBR5oxY8YQFRUV9HFDmYKhiIiIiIjsFCUlJXR2dlo61jRh6dI0Hnkkh/g4H/83+yHSLt4fHEmsfeklArGxQdUrKCjANE02bdqE0+m01NNQpmAoIiIiIiJh5/F46OrqsnRsU5ODuXNH8ckniRw3o4inzPNJe/ZjykbdRtOxxwYVCmNjY4mNjSUhIQHTNNlrr72G/fOE26NgKCIiIiIiITNNE4DNmzfj8/mor6+3VOfjj3sXmOnstPPMyX/krL9fg93louK3v6XpmGN2uI7T6WTKlCnbhEAFwu+nYCgiIiIiIiErLy+nqanJ8vHd3Qbz5+fw8svpFBR08vKPbmDqC7+ns6CAogUL6B43Lqh6KSkpCoJBUDAUEREREZGQtLW1hRQKi4qiueWW0WzaFM0ZZ9Rx2WXVODftQ51xJtWXXIIZGRl0TT1HGBwFQxERERERscTj8eDz+di4caOl4wMBWLIkncceyyYx3suHJ93EBF8JlRHX0TVxIlUWNqsHyMrKIjk52dKxw5WCoYiIiIiIWLJ69WrLxzY2Orj99lF8+mkiJ8xYz1O+8xjx6r9oPfhg8PksbUMRExNDXFwcmZmZlvsarhQMRURERETkB3m9XgzDwOFw0NXVFfQm9Vv7xz8SufPOfLq67Dx74pOc/v5vsXd3U37jjTSeeCIE+Wygw+HA5/MxduxYIiIiLPc1nCkYioiIiIjI9/J6vbS3t1NWVgZAfHw8LperfxXSYHR3Gzz0UC6vvprG+PGd3H/D5xx9+dX05OSwYd48ekaNCrrmyJEjSU9Pp7W1VaEwBAqGIiIiIiLyvZqamqiuru5/3dHRYanO+vXR3HrraEpLo/nt8f/k59dFEREZRdHChXQVFGAGsVhMSkoKUVFReL1e0tPTiYiIID093VJf0su2uxsQEREREZGBpaOjA7fbDYDP5wuplt8PzzyTwezZk+jsgH+e8FvueXMmmW8vB6Bz8uSgQiH0blqfmZlJXl6eZgnDRDOGIiIiIiKyjaKiIgDsdjuBQMBynaqqCObMGc2qVXGccdBqHms/j+TX/kPLYYfResghlmo6nU5SU1Mt9yTbp2AoIiIiIjJMeTwe/H4/0dHRANTV1bF58+b+9/1+v6W6pgnLl4/gwQdzsdlMlp65mBNf/y2G10vZnDk0HXts0AvMjBo1irKyMqKjo7HZdONjuIX0O2oYxtWGYaw1DGONYRhLDMOIMgxjtGEYnxmGUWwYxlLDMCL6PhvZ93pj3/ujtqpzY9/4BsMwjgztK4mIiIiIyI6orq5m06ZN/a+rqqpCvnW0udnBtdeO5c47RzF5spsXX/yGvffz0J2Xx7oXXqDpuOOCDoVAf3hNSkoKqT/ZPsszhoZhZANXAJNN0+wyDOMl4DTgZ8BDpmm+aBjGIuB8YGHfP1tM0xxnGMZpwL3AqYZhTO47rhDIAt4zDGO8aZrW/veEiIiIiIjsEK/X2x8Erawy+m0ffdS7DUVHh51HT3qJY0f/h8aMU2nP+DHtM2ZAkDN9ycnJtLS0AL3BcK+99tJs4U4S6u+qA4g2DMMBxAC1wGHAK33vPwuc0Pfz8X2v6Xv/cMMwjL7xF03T7DFNsxTYCOwbYl8iIiIiIrIdPp+PdevW0dXVhd/vx+/3093d3f9coRWdnTbmzcvj6qvHkZncwaqfnselfzmN9FdfwfB4ej8URKAbN24ce++9N2PGjOkfMwwDu92OYWG2UX6Y5RlD0zSrDcN4AKgAuoB3gC+AVtM0t8w/VwHZfT9nA5V9x/oMw2gDRvSNf7pV6a2P2YZhGBcCFwKkpaWxYsUKq+2L7DQul0vXpgxIujZlINP1KQPVULw2bTYbsbGxfPHFF0RGRmKz2Vi7dq3leqtXxzJnziiqqiK5+Zj3uembXxPzZgn1J59M9ZVXYu7gqqGBQKB/NvCrr77qH4+PjwcYcv8eBppQbiVNpne2bzTQCrwMHLWdj26Zk95etDf/x/h3B01zMbAYYMKECebMmTODa1pkF1ixYgW6NmUg0rUpA5muTxmohuK12dHRQVFRERMmTKCqqsryqqM+Hzz9dCbPPJNJWpqHPzzwGWfdciyBmBiK58+n/YADgqoXFxdHZ2cnwDa/51988cV3xiT8QlmV9Aig1DTNBgDDMP4C7A8kGYbh6Js1zAFq+j5fBeQCVX23niYCzVuNb7H1MSIiIiIiEkZbgmB9fb3lUFhWFsmcOaP55ptYTv3JRi6+uYO4OCdlc+fSsc8++INYICYjI4P6+nqys7MpLi7+zvv5+fl0d3db6lN2XCjPGFYA+xmGEdP3rODhwDfAB8DJfZ+ZDSzr+3l532v63v+72fuE63LgtL5VS0cDBcDnIfQlIiIiIiLb0dzcTGtrK4ClsGWa8PLLaZx55mSqqyJ48+R7ef6Tvcj+6h8AtB52WFChECA7O5u9996bhIQEEhIScDi2nbtKTU0lJycn6F4lOKE8Y/iZYRivAF8CPuArem/zfBN40TCMu/rGnu475GngOcMwNtI7U3haX521fSuaftNX51KtSCoiIiIiEj7l5eU0NjZiGIbl1UcbGx3ccccoPvkkkaOml/BMxIVkvPJ3OvbZh65x48LSZ0FBQVjqSPBC2uDeNM3bgNu+NbyJ7awqappmN3DK99SZB8wLpRcREREREdm+xsZGwPqWFO++m8w99+TR3W3j6ZOf5ey/X4O9o52qK6+k7owzwG4PuuaYMWOIjY211I+EX0jBUEREREREhq7WVjv33pvHu++mUFjoZu7cUqat3Ih3ZTLFjz9GVwgzfNHR0UTs4IqlsvMpGIqIiIiIDCF+vx+v10tUVBQej4fVq1dbqvPxxwnceeco2trs3Hvi6/x82je0jzqCpvzjaT76aEyn01JdwzAoLCwkMjLS0vGycygYioiIiIgMIWvXru0PhlYWmHG5bDz0UC7LlqUyaWwrKw66gcJlT9G9Mp9vjjgU7HZLoTAlJYXm5mYcDodC4QCkYCgiIiIiMgQEAgHKy8vxer2AtVVH//OfOObOHUVdXQS3Hv8+139zCbGvFdF4zDFUXntt0M8SRkZGYrfbGT9+PIFAgObmZpwWZxpl51IwFBEREREZpNrb2ykvL8dut+P1evH5fJbqdHcbPP54NkuWjCQvr5sl9/ydk28+Cn9CAhsffJC2Qw4JuqbdbmfChAn9QdBut5Ofn09CQoKlHmXnUjAUERERERmkysvL8Xg8IdVYsyaGOXNGU1ERxbk/L+a837iIjk6i4vrraZ05M+h9CceMGUNS3zG9253/V2pqaki9ys4Tygb3IiIiIiKymwQCgZBCoddr8PjjWZx33kQ83fDRyTfz5Dt7M6J0DQBNJ5wQVCg0DIOMjAySkpIwDOM7oVAGNs0YioiIiIgMIuXl5URGRlJdXW25RnFxNHPmjKK4OIZfHf4F9zddRNIr/6H1wAPxpKdbqhkZGUl2drblnmT3UjAUERERERkkOjo6+jert8Lng+eey+CJJzJJTPTz9qm/4/DX7wLDoGzOHJqOPRYszPRFR0eTl5dnuS/Z/RQMRUREREQGiaKiIsvHlpdHcvvto1i9Oo4jjmjmhhsqmPRyOe7CQsrnzMGTmRlUvfz8fHp6esjIyMAe5GqlMvAoGIqIiIiIDFCmaeLz+di0aRNpaWmWagQCsHRpOo89lk1khJ/XTnmIyQc5aE/an9pzz4XzzwdbcEuPTJgwgbi4OEv9yMCkYCgiIiIiMoCYpklLSwvd3d3U1tb2j7tcrqBrVVZGcscd+Xz1VTxH71vCU85fk/Hy+7Q0Hkr7/vuDw1ocUCgcehQMRUREREQGkKqqKurr60OqsfUsodMZ4MVTF3HiOzdhd7mouvxy6s46K0zdylChYCgiIiIiMkD4fL6QQ2FlZSRz5+azcmU8BxzQxoPHvsKMGy7GPXEiZQsX0j1unOXaW7akkKFHwVBEREREZDczTRPTNCkuLrZc49uzhA9e/QkHnxGJwVRK77iD5p/+1NKtoyNHjqSuro4xY8YQHR1NVFSU5R5l4FIwFBERERHZjerr66msrAypRkVF77OEK1fGc9SMUp6MuZyMJz7gm8OW4snMpPlnP7NcOzo6GrvdTmJiIrYgF6mRwUPBUERERERkN/D7/dTU1NDZ2Wm5RiAAL76YzuOP984SLjltMSe9cxOOtlZqzzsPb2pqUPWSkpJobW0FIDMzk4SEBGJjY0lKSlIoHOIUDEVEREREdoGenh4qKyvx+/1kZWXR3t4e0vOEW88SHnxAEy9EnEP2i3+jc/x4ih99hK4JE4KqZ7PZyM/P7w+GWVlZ/e9pn8KhT8FQRERERGQXWL9+PT6fDwhto/pvzxLefnspRx/djPO+WGp+/Ws2//KXmE5n0HVTU1NxOByMHz8ej8djuT8ZnBQMRURERER2gS2hMBRbzxIePaOEJ2KuoH3s6XQak6i87jowDEt1Y2JiSE9PByA+Pj7kPmXwUTAUEREREdmJXC4XGzZsCKmG3987S7hgQe8s4dJfLODn792Kvb0d7wH70DlpUtChMC4ujpiYGHJycjAsBkoZOhQMRURERER2gp6eHrxeb8ihsKwskjvvHMXXX8dxzIyNLI64lMyX3undl3DBAkv7EhqGQUFBgRaUkX4KhiIiIiIiYWKaJg0NDcTExIQcCH0++POfR7J4cRaRkb3PEp7fOJ+Ri/5O9cUXs3n2bEv7EkLvraMKhbI1BUMRERERkTDo6emhubmZmpqakGsVFUVzxx35rF8fy88PKOK6X3yO84BJ1PnOou2QQ+gePdpy7bi4OMaOHRtyjzK06H8TiIiIiIhY4Pf7MU0TgIaGBtatWxdyKPR4DBYuzOLssydRX+fktdMfYenaGex97/W9U4gOh6VQuPXWE+PHj8dhcaZRhi5dESIiIiIiQQoEAqxevZrc3FxGjBhBRUVFyDVXrYrlzjvzKS2N5uzD1/Bg1+WkLVmBe/Jkym67zdJto+PHjycmJga73U5tbS2maWqhGdkuBUMRERERkSD5/X78fj+dnZ10dXWFVKury8aCBVm8+GI66ele/njr+5z50AnYvF6qrrySutNPtxQKo6Ojt9l6YsqUKQQCgZB6laFLwVBEREREZAdVVVVRV1dHbGwsAPX19SHV++yzeObNy6emJpIzTizn11c2ExuTSEP1L2g69lh6csVD/VgAACAASURBVHODqhcREUFBQQFer5eYmJht3nNa2PRehg89YygiIiIi8j+YpkkgEGDz5s3U1dUB4Ha7Q6rZ0WHnzjvzufTS8UTYfXx6yvX88YO9SW6rAsOg5pJLgg6FAJMnTyYqKor4+HjsdntIPcrwohlDEREREZHvYbfb+fLLL4mMjKSnpycsNVesSOSee/JoaXFy83EfcH3x5cS/vJbWgw7CDHFWT2FQrFIwFBERERH5HhEREQBhCYXNzQ7uvz+Xd99NoWCcm/cPupopy5/Al5DAprvvpuUnP4EQFoYZP358yD3K8KVgKCIiIiLyLV6vl/b29rCs4Gma8NZbKTzwQC6dnTYuuqiaX/5yM5n31tL8059Sec01+JOSLNefOHEiUVFRmi2UkCgYioiIiIgALS0tGIaB0+lk/fr1QOi3Zm7e7OSee/L4+OMk9i2s47lR1+H/8eF0OgqpuOEGsFA/JSUFgObmZjIzM/sXwhEJhYKhiIiIiAiwadMmgLAELb8fXnopnQULsjBNWHjyEn75yW+JWltN1ZhEOgsLLYVCgJycHBoaGkLuUWRrCoYiIiIiMqx1d3dvs0F9qCuOFhdHc9dd+axdG8tP9y1ncdJvyH/lL3Tn5bFh8WJc++xjqa7dbsfv9+NwOIiPj6e2tpa4uLiQehXZQsFQRERERIaFQCBAU1MTqampGIZBa2srERERbNiwISwbv3d3Gzz1VCbPPZdBQoKPu+7axDkNj5Lz2DI2z55NzQUXYEZFWaptt9spLCwkEAhgGAbx8fHstddeeq5QwkbBUERERESGhfb2dioqKoiKiiIqKorS0lJiY2PDEgo//zyeu+/Oo6oqitlHrObyY7+AA6ZQ7zuNjhn70hXiiqFOp/M7G9QrFEo4KRiKiIiIyLDg9/sBqKuro6enh0AgQEdHR0g1W1vtzJ+fw+uvp5KX08mHp93G/q8/iHd9Mmv3exUcDkuhcMttoyK7ioKhiIiIiAwLW4KWx+MJeV9C04S3307mwQdzaW93cMsJH/DbjVeR8OIq2mfMoPzGGy0tLhMREYHH4yE3N5cRI0bg8/n4+uuvSUtLC6lfkR+iYCgiIiIiQ15PTw9tbW1A72IzpmlarlVTE8E99+TxySeJFBa6ef76N5h18yn44uMpveMOmo86yvJG9fHx8SQlJZGYmAiAw+FgH4uL1YgEQ8FQRERERIa8NWvW9P9sNRT6fPDii+ksWpSFYcDci/7NrHNt2G25VG2+nKZjjglpo3qbzYbD4SDpWzUMiyFTJBgKhiIiIiIyZHV0dNDc3BxynfXro5k3L59162I5er9NLIz7DVnPvsM3R7+EJzOT+rPOCrqmzWbrX/hm7NixOBwOIiMjQ+5VxAoFQxEREREZssrKyvB4PJaP7+42eOKJLF54YSRJiV6WnfowR71zB472djafcw7e5GRLdQ3DIC4ujvb2dgASEhKw2WyW+xQJlYKhiIiIiAwJHo+Hurq6/tsxq6qqQgqF//xnAvfdl0d1dSQnHlfL4s0nM2LpJ7gLCylesICuggLLtadOnUp9fT3t7e3k5+crFMpup2AoIiIiIoNOV1cXfr+fuLi4/rGqqipaWloA2Lx5s+X9CRsanDzwQC7vv59Mfl4XixZtYPp0F96Hx1Bx8P40nHKKpRVHt2a328nIyCA2NpaEhISQaomEg4KhiIiIiAw61dXVeDweJk+ejNfrpaKiArfb3f++lVDo98PLL6excGE2Pp/B7076G5etvZbKqOvpZA+qr7oqpJ7z8vJwu910dnYCvbeTKhTKQKFgKCIiIiKDSiAQoKenB5/Ph9/vp7y8vH8rCqvWrYvh7rvzWLculsOmV7I440bG/HUJ3pQU7C6X5bpbFphJSkoiLS1N+xHKgKVgKCIiIiKDSklJCd3d3QCsWrXK8i2jAC6XjUWLsnjppXSSk328cuZCjn3nNpxfNNJw0klUX3YZga1uV91RNpuNSZMm4ff7cTgcOJ1Oyz2K7AoKhiIiIiIyaHi93v6VPMHaLaMApgnvv5/EAw/k0tTk5OSTG7jkkhoKXlmLLzmZkvvvp3OPPSz3abPZiIqKsny8yK6mYCgiIiIiA55pmttsUh+K6uoI7r03j08+SaSwoJW3D7mdxH1zaI0/lM1nn83ms88Gh7W/JkdHR9PV1aX9CGXQUTAUERERkQHJ7/dj71v9s6mpKaStJwC8XoPnn0/nySezsNtNHjvlZc779/VEv1pKve0UWg891HIgnDBhApGRkfj9ftauXUtWVlZIvYrsagqGIiIiIjLgtLW1sXHjRtLT00lISKC8vDykeitXxnL33fls2hTNCQcU82jUteS8vJye7GyKH36Y9gMPtFQ3Li6OCRMm9L92Op1MmzYtpF5FdgcFQxEREREZUDweD1VVVQDU19fT0NBguVZrq51HH81h2bJUMjJ6eOihjRzX9QpZc/5G7bnnUnv++ZgWnwV0Op3k5eVZ7k1kIFEwFBEREZEBwev14vV6Wbdu3TbjpmkGXSsQgDffHMH8+dl0dDi46dgPOHvaZ7gP+gkt5k9wFxbiyc4Oqd/4+Hiio6NDqiEyUCgYioiIiMiAsHnzZurr60OuU1QUzb335vH113HM2KOOP4y+iYl/exbPf9JZe+RMTKfTUii02+34/f7+rScyMjJC7lVkoFAwFBEREZHdyufzUVxcTGdnZ0h1XC4bTzyRxdKl6SQk+PjzqU9x8j9uJXLNZhpOOIHqyy7DDGE/wdTUVJxOJ263G7/fr9lCGVIUDEVERERkl6urq6O1tZXY2FgMwwgpFJomvPVWCg8/nENzs4OTTmrgmp99zH7nXUDX2LGsf+op3HvtFVK/hmEQGRlJWloafr8/pFoiA1FIwdAwjCTgKWAPwATOAzYAS4FRQBnwC9M0WwzDMID5wM+ATuCXpml+2VdnNnBLX9m7TNN8NpS+RERERGTg8Xq92O12bDYbra2tuFwuXC5XSDVLSqK49948vvwynqkTm1n6y5dJOn06kE/xI4/Qvu++lragsNvtFBYW8q9//Yvp06f3B8Mt74kMNbYQj58PvGWa5kRgKrAOuAF43zTNAuD9vtcARwEFfb8uBBYCGIaRAtwGzAD2BW4zDCM5xL5EREREZAAwTRO3201rayurVq1i48aNeDyekAOh221j/vxszjhjMhs3RvPMaX/iU/feHPb7i4msrASgff/9Le9LmJeXh9PpJBAIEBMTQ3R0NDZbqH91Fhm4LM8YGoaRABwM/BLANE0P4DEM43hgZt/HngVWANcDxwN/MnuXlfrUMIwkwzAy+z77rmmazX113wVmAUus9iYiIiIiA0NtbS21tbX9rzs6Oli9erXleqYJ772XzEMP5VBfH8H5P1nJ3d3Xkv7i+3Tn5bHxkUfoyc21XN9ut5OTk0NysuYpZHgJ5VbSMUAD8AfDMKYCXwBXAiNN06wFME2z1jCM9L7PZwOVWx1f1Tf2fePfYRjGhfTONpKWlsaKFStCaF9k53C5XLo2ZUDStSkDma7PocnpdBIRERG2mbayskjuvz+Pzz5LYMKETh64Yy2n3XAY9q4uqi+5hLqzzsKMiLBU2+fzYRgGfr+fNWvW9I/r2pThIpRg6AD2AS43TfMzwzDm89/bRrfH2M6Y+T/GvztomouBxQATJkwwZ86cGVTDIrvCihUr0LUpA5GuTRnIdH0OHT6fj0AgQG1tLY2NjWGp2dVl45lnMnjuuZFERQWYf/Yb7HdpBnaHQcWNN9I5eTKezExLtePi4nC5XEydOpXS0lIMw2Dfffftf1/XpgwXoQTDKqDKNM3P+l6/Qm8wrDMMI7NvtjATqN/q81vP6+cANX3jM781viKEvkRERERkN6mursblctHd3R1yLdOEf/wjkQceyGXz5kjOPnQ19/muIeO5d9m49+9pO/hgWg8/PKRzTJgwof/n5ORketdLFBl+LM/rm6a5Gag0DGPLn6bDgW+A5cDsvrHZwLK+n5cD5xi99gPa+m45fRv4qWEYyX2Lzvy0b0xEREREBpHNmzfT2NgYllBYWRnJVVeN49prx5EU082nJ/6WP3w6g/TPP6T6ootonzEjDB1vKyMjg5EjR4a9rshgEOo+hpcDzxuGEQFsAs6lN2y+ZBjG+UAFcErfZ/9G71YVG+ndruJcANM0mw3DuBP4d9/n7tiyEI2IiIiIDCyBQICKigqysrKI6Huer6ioiPj4eGpqakKu39nZe9vo88+PxOk0ueqqSm7/+CQS//JvWg85hMrf/AZP9naXo9hhOTk5+Hw+LTAjspWQgqFpmiuB6dt56ztz+n2rkV76PXWeAZ4JpRcRERER2fm6u7tpamoiPj6eESNG4HK56OjoCHn7CdOEt99OZv78HBoaIjj78DWcdXU3IzIMGvJOp/6sM2k/8MCQzpGfn095eTmJiYlERUWFVEtkqAl1xlBEREREhhG/3w/0blbvdrvZtGkT0LtfoVVFRdHcf38uX30Vz5QJLbwx8xamvfk4NVMupP6ss2g7+OCQes7KyiIpKYno6GgSExNxOp0h1RMZihQMRUREROR/Mk2zf1GWLc8P1tTUhBQGAVpb7SxalMVf/pJGQoKPP572LL/4121Ev1xO60EH0RriaqAxMTHExcWRkZHR379Cocj2KRiKiIiIyPdqbm6mtLSUuLg44uLi2Lx5MxDaDKHfD3/9ayoLF2bT0WHnlFMauNe4nvwX/0h3Xh7FDz9s+bbR5ORkWlpaAEhJSdFiMiI7SMFQRERERIDesGea5jYb0jc1NQG9G72H+hwhwMqVsdx3Xx5FRTEcMLWWq6+sYNSeDrpWHkhVWjz1Z5yBGcKsXkZGBk6nk/r6ehwO/VVXZEdZ3q5CRERERIaW9vZ2Vq5cSWdnJz09PUDvhvXhUF/v5NZbR/GrX02krdXOslMf5u81Uzno/x4EwL3XXtTNnh10KIyOjt7mtdPpJCcnh6ysLBITE8PSu8hwoP+NIiIiIiJA7/ODpmmybt06AMaMGUNnZ2dINT0egxdeSOfppzPx+QxuO/5dri69jsSlK3FPmkTzUUeFVH/ChAnU1NTg8XhobW3F4XBgGAaZmZkh1RUZbhQMRURERAS3201VVdU2Y1tWHLXq448T+P3vc6moiOKQQ1p5uPAhpi66C19iImW33ELTcceBzdoNbFtud7Xb7eTm5hIIBPB6vf2LzIhIcBQMRURERIY50zRpbGwMW73y8kgefjiHjz5KYnSui8W/+5x9fmLHXrUP9a2nUXvBBfjj4y3VttvtpKSk4PF48Hg8/eM2m43IyMhwfQWRYUfBUERERGQYam1tpbS0lBEjRtDR0RHy1hMAHR12nnoqk6VL04iMNHnk5y9x/tc34X0jg40/mY8nJ4eq3/zGUu0twS8vL4+4uDg8Hg+BQCDknkWkl4KhiIiIyDCxJUj5/X5KS0sJBAI0NDSEXNfng9deS2XRoiza2hycd8TX3NV9HRl/fZeerCxqTzgBTBNCuM1z/PjxxMbG9r+OiIgIuW8R+S8FQxEREZEhzu/343a7KSkpCfss27//Hc+DD+awcWMM++zTwX0/WcqhD18GQPVFF1F31lmYUVEhn0dbT4jsXPoTJiIiIjLElZWV0draGtaalZWRzJ+fzYoVyWRldvPozZ+y3wlOHK5xNK+fRc0FF+DNyAjpHDExMf2roioYiuxc+hMmIiIiMoQ1NTXR0dERtnoul42nn85kyZJ0nE6T+09azq/X34LjxU6+OfZ5/PHxlN96a0jniIuLIzs7m5iYGCorK2lsbOxfhVREdg4FQxEREZEhyO/309zcTEVFRZjqweuvj2DBgmyam52ce8TX3OW/iaxX/4YnNZXqyy6zvPXEFiNGjGDUqFHbjKWkpPTvTSgiO4+CoYiIiMgQtG7dOnp6esJS64sv4njwwVyKimKYOtXFc1f8lVm/OweA2vPPZ/Ps2QRiYkI6R0xMDLm5ud8Zj4+PJ97i1hYisuMUDEVEREQGKZ/Ph8vlIikpqX/M7XazceNGfD5fyPWrqiKYPz+HDz5IJnNkFwuu+YAfnZaA4c+nfuMpNJx6Kp7MzJDOUVBQQHl5OdnZ2djt9pB7FhFrdLO2iIiIyCDV2NhISUkJfr+frq4u2tvbWb9+fcih0O228dhjWZxySiH/+lcC952wjA0pP+JXfzgZu9sFDgfVV10VUiiM6lupNCEhgSlTppCQkBBSzyISGs0YioiIiAwygUCAkpKS/pU6W1tbKSsrC7muzwfLl6fyxBNZNDU5mX3o19wduJGs1/4PT1oaVVdcEdItow6Hg4iICHw+H+PGjcPpdIbcs4iEh4KhiIiIyCDi8/no7Oykvb29fyzUUGia8MknCcyfn8OmTdFMnepi8fV/54Rbjweg5le/ou6cc0IOhXvuuSeAFpIRGYAUDEVEREQGuPb2duLj4zEMgw0bNoR164bi4mgefjiHzz5LIDe7kz9c+Tp7nJWFYaRSU3EBzbNmhbQfYWRkJMnJySQkJCgQigxgCoYiIiIiu5nH48Hv99PU1ERKSgoxW83MuVwuiouLcTgcTJw4ke7u7rCcs6HByaJFWSxfPoL4eD+Pn7yUc1bfTszCTaz5yV/wZmRQ98tfWq6/9957U11dzciRI4mIiAhLzyKy8ygYioiIiOxmlZWVeDweOjs7MQwDp9PZ//zdliDo8/lYs2ZNyOfq6rLx3HMj+dOfRuLzGVzzs4+5rvkW0l75kJ6MDMpvuQVvenrI57HZbNvdfkJEBiYFQxEREZHdzOPx9K8k6na7WbVqFZMmTQrbthPQu0H9m2+OYOHCLBoaIjj88BauOftrjrzwSEynk6rLLqP+tNMw+1YLtcpmszF27Niw9Cwiu46CoYiIiMhu5vP5CAQCQG9IBKirq8Pr9Yal/mefxTN/fg5FRTFMn1zPkpNfJOn8A4B4ym67jY5998WXnBzSOWw2G4FAgPHjxxMbGxuWvkVk11EwFBEREdmJqquriYmJIXk7wcvlclFaWtofBgF6enoAaG5uDvncmzZFMX9+Dv/8ZyLZmV0sP/kBjvzwASIW1rHm8FfoGTWKliOPDPk8AMnJyeTk5PRvoSEig8ug3eA+vqgIvr2y1bHH9o69/vp/xxYv7h278ML/jtXU9I5lZW17/LRpveNffPHfsdtv7x27/fb/jn3xRe/YtGnbHp+V1TteU/PfsQsv7B1bvPi/Y6+/3jt27LHbHm8Y+k76TvpO+k477TvNPPTQIfedhuK/J32nofedsnNytvudAvvsQ2dnZ38onDJrFtOmT8fZ0ND/0bx585g2fTqpf/lL/1jihx8ybfp0xl599bannz6dadOnA9DU5OB3v8sj+hc38PE/k3jjx9exPm5vjn3lt5j23r/+jXz++f5jnQ0NTJs+nSmzZm1Tc+JZZzFt+nRi1q3rH8t84gmmTZ9O9pNPkpCQwNSpU8mtr2fU6NE4ZszY9vdkEP17+r5rb9rWNYfIdxrMf570nb71ncJo0AZDERERkcHCNM3vjPV0d1NXVxf2cz31VAYnnrgHr72WSmZW7+zjrC8fwelqZ9Ndd7E5hJVGt5aRkUFBQQEOh4P0MCxWIyK7l7G9/1ANBhMmTDA3bNiwu9sQ+Y4VK1Ywc+bM3d2GyHfo2pSBbKhen6Zp8uWXX5KQkEB2djalpaVERkbS1tYW1vP4fLB8eSqLF2fS2BjBz/ffwB2jF9Jz1VlgGMSsWUNXQQFmZGRYzpeamkp+fn5Yag10Q/XalKHBMIwvTNOcHo5auglcREREZCfZEgC7u7tZ13c7Zrj2IQQwTVixIonHHsumvDyKGVPqeevAeez1zlPYPuth3VE/omviRDr32MPyOWJjY0lJSSEhIYGenh5qampISEgI23cQkYFBwVBEREQkTAKBAIZh4Pf7aWtro6ysDGCbxWXCZeXKWB55JIdVq+IYm9/BW6c8xKErHiJidQMthx5K9eWX05OXF9I5xowZs82iOVFRUSQmJobauogMQAqGIiIiImGybt06UlJSsNvtVFZW7pRzlJZG8dhj2fzjH0mkpnq4+eZyjjuyhqmnP44nI4NNv/sd7r32Cvk8ERER211JVUSGJgVDERERkTAwTZPu7m5cLhcxMTFhr9/Q4OSJJzJZvjyV6OgAD5zwV2a3PE7FUXMxo6LY8OSTeNPSvruSogVTpkwhIiIiDF2LyGChYCgiIiISgkAgQFdXFz6fD4DOzs6wbUwP4HLZePbZDF54YSR+P/xm1odc13Ibaa/9A09qKlEVFXSNH483TCuDjhgxQqFQZBhSMBQRERGxyDRNqqqqaGhoICoqCgCfz9cfEkPh8Ri88koaTz+dSVubg2OPqOZR83Ly3lpGIDqa6ksvpe700zH7zmtVYWEhAGvXrsVut5OZmRly7yIy+CgYioiIiATJNE38fj9ff/11/1i4VhsNBOCdd5JZuDCb6upIfvyjZi65oo5JE90kX1JF/emnU3vuufiTkkI+V0pKSn+g3WOPPXA6ndhs2uZaZDhSMBQREREJ0ldffbXdTetDYZrw6acJPP54FuvXxzJlXBMvnXAr0/71J9anPYvPSKX48cfBYnBzOBykpaURExNDSUkJhmEwatSo/vcjw7THoYgMTgqGIiIiIjvI7XZTUVER9lD49dexPP54Nl9+GU92RifLTriPWZ/8nojX6mjbf39sPT29HwxhNm/KlCnYbDYCgQBOp5O8vDyMMCxUIyJDg4KhiIiIyA/45ptvSEtLw+v10tnZGba6xcXRLFiQxUcfJTFihJcbryriluU/I+a1EtyFhZTeMRfX9OkhnSMiIoKJEyf23yJqs9nYc889w9G+iAwhCoYiIiIy7G0JfFtv3t7d3U1paWl/EKyoqAjb+aqqIli0KIu3304hLs7P3DM+5LCL44mODtDWdjC1v76Q1sMOC9vWEyIiP0TBUERERIa94uJiurq62GeffTAMg5aWFjZt2hT289TXO3nqqUyWLUvF4TCZc/Q7XF5/BylL/sW6Y1+gq6CAmksvDcu5IiIiyM3NDUstERn6FAxFRERk2AoEAng8Hrq6ugDYtGkTubm5YZ0dBGhttfPssxm89FI6fj9cfsQ/ud49l4w33sOXmEj1FVfQHaYQl5qaSmNjI5MmTcLh0F/1RGTH6L8WIiIiMiyZpslXX31FfHx8/1hra+s2m9WHqrPTxgsvpPPccxl0dto46qhmLjlnI7POnwWmSc0FF1B35pkE4uIsnyMuLo6kpCSqqqoAyM/PJz8/Pyz9i8jwoWAoIiIiw5Lf7wego6Njm/GeLSuAhsDjMXj11TSeeSaDlhYnP//xeuaMfRr/laeC4WDTvHl0FhbiS04O6Txjxowhua9GIBDYJuSKiARDwVBERESGBa/XyzfffENBQQExMTG0traG/Rw+H7z55giefDKTzZsjOWzvch7e/y4mv/9nbJ97+eaY6XSPG0f7gQeGfK7IyMj+UAiQmZkZck0RGb4UDEVERGTIM02TkpISfD4fjY2NdHZ24na7w1Y/EIB3301m8eIsysuj2HtiI6/teyPT/v4k9pVummfNovbCC+kJw3OEMTExdHZ2kpCQEIbORUR6KRiKiIjIkNXR0UFDQwMJCQn9QbChoSFs9U0TPvggiSeeyKKkJJqxY7u4//4SDt1vM1NO/DOuadOovvhiuseNC8v5xo79f/buOzyqKuHj+PfMTCa9J4TeQaSIIIoRXggoKkixYMG1rGXtddV3rauu666ujXV9LQhYFrEAioigQDCKGhQFaQLSayCQSupkZs77x4QICCowIRPy+zzPfWbunXPvORevV36ee89pR0JCAuvXr1cwFJGgUjAUERGRY47X68XpdFJQUEBBQUFQJ6WHQCCcNy+eV15pyqpVUbRtWcL0c//J/2yazJq+/4cNi+DHd97Bt9e8iEciJiaGtLS0mnkW27RpE5TjiojsoWAoIiIi9ZrH4yEsLAxjDB6PB5/Px6pVq0hOTqagoAAIzoAyEAiE8+fH8fLLTVm+PJqWzUp5//xnGfztaCKmbqa0SxdceXlUNW4clFDocrno2LEjbrcbp9MZhDMQETkwBUMRERGpt7xeL0uXLiUtLY3mzZuzceNGysvL8fl85ObmBrWu776L4aWXmrF4cQxNmlTy5G3zufXDS4h8fyNlHTuy5umnKerfH4wJWp2NGjUiMjIyaMcTETkYBUMRERGpt/ZMObFr1y7S0tKoqKigqqoqqHX88EM0L7/clO++iyMttYJn/pTJaVcnEuZ0Ur6qE9tuuonCAQPA4QhqvQBpaWlBP6aIyIEoGIqIiEi9tPfIoj6fjzVr1uDxeIJ2/GXLonj55abMnx9PcpKH8cPHctGP/yJi4jaWjvoIX1wc6x9/PGj1AbRv3578/Hzi4+MpLCzEBLH3UUTk1ygYioiISEgrLS2loqKC5ORkIDCRe15eHps2bdqnXLAGmFm5MpJXXmnKvHkJJMR7eGXY61y66klipq2komVLNt53H77o6COqIy4ujuLiYlJSUvB4PBQXF9OtWzfcbnfNADNJSUnBOB0Rkd9FwVBERERC2s6dOykqKqoJhrt27WLz5s1Br2fNmgheeaUpn32WSFycl5tu2srVJ33JSddcRWWzZqx/5BHyzz4bXEf+16eWLVuybNky3G43LVq0wOv14na7g3AWIiKHR8FQREREQprX68Xn87Fy5UoiIiKCOjE9wE8/RfLqq0347LNEoqO8PDvkHYa3/IbCqy8D2rP63/+muHfvoATCPdxuNy1btiQ+Ph6Hw6FQKCJ1TsFQREREQlZJSQnFxcVYayktLQ1qKFy5MpKxY5uQlZVIdLSPf50zhWs3/4PEGQupbNaMostGYiMiKO7T54jrcrlceL1ewsLCCA8PxxhDampqEM5CRCQ4FAxFREQkJPl8PvLz87HWBvW4K1ZEMWZME+bNSyAmxstf91nHiAAAIABJREFUR2Zxx5p7SPz4OzwpKWy65x52nXceNoi9eI0aNWLbtm3ExcXRunXroB1XRCRYjjgYGmOcwHfAVmvtUGNMG+AdIAlYCFxurfUYY8KBN4GTgDzgYmvthupj3AdcA/iA26y1nx5pu0RERKR+KC0txeFwUFxcTFxcHB6Ph9LSUtatWxfUepYvj+LVV5vw5ZcJxMVWccfVKzn3inISS0qJ+tNWNv/5z+w8/3xsRERQ6mvZsmXNADlhYWEAhIeHB+XYIiLBFowew9uBFUBc9fqTwHPW2neMMS8TCHwvVX8WWGvbG2MuqS53sTGmM3AJ0AVoCswxxnS01vqC0DYREREJYRUVFaxcuZLo6GhWr16Ny+XC5XJRUVERtDqWLYtizJimfP11PPFxVTx/zltcseFf8GMUa2JeoCqmMcumTgWnMyj1JSYm0qhRI6Kjo2uCYXJyMsYYjTQqIiHriIKhMaY5cA7wOPBnE5hsZyBwaXWRN4BHCATDEdXfASYDL1SXHwG8Y62tBNYbY9YApwDZR9I2ERERCW2bN28mNzd3n21erxev1xuU4y9ZEs2rrzYhOzsQCF88503+sO4p4j5ehictjZxhV4G1YMwRh8K4uDjS0tKIjIzE4XDg3O94xpiaUVVFRELRkfYYjgb+F4itXk8GCq21e+7oW4Bm1d+bAZsBrLVeY0xRdflmwPy9jrn3PvswxlwHXAeQmppKVlbWETZfJPhKSkp0bUpI0rUpoSQsLIyIID2yub8ffohmzJimfPttHAkJVdx66xZuDX+ZDk8/TmXTpmx84AHyhg7FVj/eebh8Ph9+vx+Xy0VOTg5bt279RZnY2MBfkfTvXv2le6c0FIcdDI0xQ4Fca+33xpiMPZsPUNT+xm+/ts++G60dA4wBOO6442xGRsaBionUqaysLHRtSijStSlHm8fjYenSpcTExNCuXTtKS0vJzc2luLi4VupbuDCGV19twoIFcSQnVjL2nFfpPbCKyv69KSkZxIYoJ3lDhgRt2om0tDRatWr1q2U2bdqE3+/npJNOCkqdcvTp3ikNxZHcGfsAw40xQ4AIAu8YjgYSjDGu6l7D5sC26vJbgBbAFmOMC4gH8vfavsfe+4iIiEg9Ya2lrKwMay1erxdXdQArKSlh8eLFtVQnfPNNLOPHN2HhwlhSkip57ZyXuXDls0R/vJrC3f1Y2783/pgY8oYPD0qd8fHxWGt/MxRCYAAaEZH64LCDobX2PuA+gOoew7uttX8wxkwCRhIYmfRK4MPqXaZVr2dX/z7XWmuNMdOAicaYZwkMPtMB+PZw2yUiIiJ1o6CggPXr19est2jR4ldKHxm/H774Ip7x45vw44/RNGrk4dXz32TUor8R/fFaKlq1Yv1jj5E/aFDQ6oyNjSU1NZWEhAQCwySIiBw7amMew78A7xhj/g4sAsZVbx8H/Ld6cJl8AiORYq1dbox5D/gR8AI3a0RSERGR+qeysnKf9c2bNwe9Dq8X5sxJ5LXXmrB2bSStm+7m4b+s4qwRpTSetR7HD5Z1jz9OwRlnBG2U0T1cLheJiYlBPaaISKgISjC01mYBWdXf1xEYVXT/MhXAhQfZ/3ECI5uKiIhIPeLxeLDWYq1l27baexOkqsrw8cdJvP56Y7ZsieD41vnMPOdJBnz3MrmVo8h1X0b+4MHkDxkCDkdQ6mzVqhUbN26sWd9/pFERkWNJbfQYioiISAPg9XpZunQpEJiuoTZUVBimTk3hv/9tzI4dbnp1zOGtIc+RPn8cYR/nU9K9O+XHHRcoHOTglpKSQmRkJJWVlWzatKlmknoRkWORgqGIiIgcsl27du3TmxbskUZLShxMnpzKW2+lUVAQRo8eu3nwwY1c9u41JMyYR1F6OtuvuoqSnj2DVmdMTAylpaVY+/Pg6NHR0URHRxMVFaVgKCLHNAVDEREROSCv14vD4cDr9ZKXl0fjxo0xxlBaWrpPKAymwkIn77zTiHffbcTu3S6GnbSSvyf/C+68iKrUVLYlXk/OdX+irHPnoNUZHR3NcccdR0lJCR6Phw0bNtC0adN9ytTWnIsiIqFCwVBEREQOaPHixcTExOD3+ykrK2PXrl00adKkVkLhzp1hvPVWI6ZMSaW83Mmlpy7i4ch/0v7LqRi/n/X92lFw9tmUd+oUlPoiIyNxOp00bdqUyMhIjDE1k9EnJycHpQ4RkfpEwVBERESAwEAy69ato1WrVpSWlgKBOQj3TM3g8XiCHgo3bAjnv/9tzIwZSfh8hrPP3MX/lV1Lqy+nY51O8oYPZ/vll+Np3jwo9aWkpBAdHU1iYiLGGBxBGqhGRKS+UzAUERERAMrKyigtLWXlypX4/f6a7Xu/cxcsy5ZF8cYbjcnKSsDtttx4+nwG3hBL8+Yeov/uIHfUKHZcdhlVqalBrTc8PJyUlJSgHlNE5FigYCgiItJA+Xw+CgoKcLlcOBwO1q5dC7BPKAwmayE7O4433mjM99/HEhdbxbODJnLF9udI+uR7lv/xHSpoz6YHHzziuiIiIqioqNhnW7NmzUgNctAUETlWKBiKiIg0ULt27WLLli21Xs+eSenffLMxP/0URdPUUiac/W9G/PQCMbPW4ElLY/Mdd+Bp0iRodXbp0oXly5fXhEOHw0Hjxo2DdnwRkWONgqGIiMgxrri4mNjYWIwxbN26FZfLRVpaGlVVVbVab0WF4cMPU3jrrTS2bQunTZtyHn54A+ekr6fnuffiadqU9Y88QsFZZ2GPcCoIl8uF1+sFoH379kDgfcItW7YQGRlJ5yCOYioicixSMBQRETlGVVZW4vV6Wb16NW3atCEpKYmCggLcbjcJCQnk5+fXSr1FRU4mTUrlnXcaUVgYRv/O65nc5Vk6F3/HmqH/ARPHyjffpKJ1awjS4C/Nmzdnw4YNAMTFxQHQqFEjkpOTcbn01x0Rkd+iO6WIiMgxat26dTUjiq5fv56CggIqKyvx+XwsW7Ys6PVt3x7GW2+lMXVqCuXlTi7stZhHYv9Fp6+n4PixksJ+/XCUleGPjqaibdug1dujRw+MMYSHh+8ziqoxRqFQROR30t1SRETkGFNaWkpBQcEvHhUtLCwEqHnkMlhWrYpkwoQ0Zs1KAuCss/K5u8dUMp64HmsM+YMHs+OKK6ho0+aI64qOjqZDhw788MMPNdv2TDkRExNDTEzMEdchItIQKRiKiIgcI0pKSoiIiGDVqlVYazHG1MpUExAYYfTrr+N46600vv02jqhIL4/0e5+z+mzFcW4fTOXx5Oy4il3nn09Vo0ZBq7dZs2Y4nc6a9T3vE4qIyJFRMBQREamnvF4vOTk5NKkezXPVqlX7/F4bodDjMcycmcRbb6Wxbl0kTVNLeePM/3De+peI/WwFJTu7surcPtjwcHJuuCGodTudTmJjYwHo3LkzFRUVxMfHB7UOEZGGSsFQRESkntq8eTP5+fnk5ubW+qTthYVOpkxJ5b33GpGXF0bHjmW8d9H/MTzrMcJn7aC8dWs2PvAAeUOGBKU+h8NBZGQkpaWlNdvC9hq5NDIyksjIyKDUJSIiCoYiIiIhz+v1snbtWtq0aVPzeKgxZp9RRXft2lUrdW/Z4uatt9L46KNkKiqcXNBzMef8pZxuA8JI+KKSynUt2HT/fRSfdlpQRhg1xtCqVStcLhe7d++mtLS05pz3foRURESCS8FQREQkxJWVlVFSUkJeXh55eXk4nc6aKRlqy+LF0UyYkEZWVgIul+W23rO5zfscLb6dzfYT/8i2gTdR1L8/Rf37B63O5ORkWrduXbMeHh5Ofn4+bdu2ZceOHYSHhwetLhER2ZeCoYiISIgrLi4GYPv27fj9fiD4I4sC+HyQlZXAW2+lsWRJDHFxXkYPeJ3Ltz9P4pc/4I2JYcdll7HzgguCXrfH46FVq1b7bIuIiOCEE04A0GijIiK1TMFQREQkhFhrKSoqwhiD0+kkJyenJhjuCYUQCFLBUlLi4KOPUnjnnUZs3RpO66a7ueeeTQwfnkeXhycQVbSTTXffTd6wYfijo4NW7x4xMTHk5OTUzD8oIiJHn4KhiIhICCkoKGD9+vUAtTrdBATeH3z33UZMm5ZCaamTgcevYVLHZ+m1cCKreo+jMrI1G++/H19sLNTS+30ul4v4+HhycnJq5fgiIvL7KBiKiIiEkL17AmsjFFoL338fw9tvp/HFF/E4HHBz70xu43nafDsTs9JHYUZGTXlfQkJQ6o2IiKBJkyYUFhaSlpZGUVEROTk5dOvWDYfDwcqVK4NSj4iIHB4FQxERkRCxaNGifR4XDabKSsOnnybx9tuNWL06ioSEKq6+ejuXDF7HwD8MxTqd7Bw5ktxLLsHTvHnQ62/cuDFJSUkkJSUBgaknIiIicARhJFMRETlyCoYiIiIhwFpbK6Fw1y4XkyenMmVKKgUFYfRqs5lZAx+lh/97Nt74JBDBmmefpbRLF/xBHOCldevW5OfnExERQW5uLhEREfv87na7a0KiiIjUPQVDERGREBDsULhiRRQTJzZi9uxEfD7D5T3n85fI0Ry34EOc6yspPvlkHKWl+KOj2d27d1DrBoiOjiY5OZmSkhJ27dqF2+0Oeh0iIhI8CoYiIiIhoKKi4oiP4fUGppt4551G/PBDLFFRPkaO3Mntbd7j1H/ehj88nLzBg8m95BIq2rcPQqsPbs9k9DExMfTo0aNW6xIRkSOnYCgiIhICjmT6icJCJ1OnpjB5cirbt4dzXJOdvJ/xLMf9TzjlIwZhKnqypeRWdo0YEbTBZA6kUaNG5ObmAj8HQxERqR8UDEVERELA4YxA+uOPUbz3XiqzZiXh8Tg4t9sPPNhhNCd+PxlnVil5kYPZMGIQNiKCHVdeGbS2RkZGUl5eTlRUFGVlZbjdbuLi4mjRogUej4fCwkINKiMiUs8oGIqIiISAPZPY/xaPxzBnTiLvvZfKsmUxREX5GDFiF49WPUCnD8eBw0H+oEHkXnIJZV271kpbExMTadeuHVVVVZSUlNC4ceOa39q0aVNrI6uKiEjtUTAUEREJAfn5+b/6+/btYUyZksrUqSkUFIRxQvNtfNT/nyTdPoDwlom4vuxKTuq17Dr/fKpSU4PevrS0NBo3boy1FpfLhTGG8PBwYvYbydThcKi3UESkHlIwFBERCQEHepR0z2T0777biC++SMBauPrEL/hzp/9w3MLpOLZUsqG/i7yWwynu25fivn2D1p6IiAgqKipqPiMjI3G59NcGEZFjle7wIiIiIaaszMGMGUm8914j1q2LJD7ey5WjNvPU9+eStHARvogI8oYMYedFF1HeoUNQ646Ojqa0tJR27dqxfPlyIiMjadGiBbGxsUGtR0REQouCoYiISIjYsCGcyZNT+eijFEpLnQxov4qnz51O2j0DCA+3eJ/qyOYhA8kbOhRfLQW1jh07Ul5eTkREBG63m/DwcOLi4mqlLhERCR0KhiIiIiFg9+4q/vjH7lRWGO7qMY0b/C/SctFcWG9YctNMvOFJbL7nnlptg8vlwuFwEB0dDUCnTp30vqCISAOhYCgiIlLH5s2bR1FRPuNu/JAh7/4vMQvW4Y2PZ8fll7PzggvwJiXVSr09evSgqqqK3NxcysvLadSo0T6/h4WF1Uq9IiISehQMRURE6kBZWRmzZs1i0auv0urTT+n+73/jPL09fBrD+qsfoWDQIGx4eK3VHxMTg8PhIDw8nBYtWtRaPSIiUj8oGIqIiBwlhYWFTJ8+nRmTJpHwySdc7fFwLlDhdDKrsJCmycmsGj++1urf81joCSecoEdERURkHwqGIiIiQZadnU1WVhYZGRm0bt2aDz/8kA8++IDczEyu8fl42RjirKWkdWt8d95JxJVXkvrjj7XSlo4dO1JUVMSOHTto27Yt8fHxtVKPiIjUbwqGIiIiQZSdnc3pp59ORUUFxhjcfj8jgX9ERHCSz4ff7cZcdBHceCMx6elgDAButzso9aelpbFjxw4g8I5gbGwsMTExREdHa3RRERE5KAVDERGRIPH5fPzrX/+ivLycjsD11nKty0Wc14tt2RJuuAHHFVdAcvIv9m3SpAk5OTlH3IaoqCgg8A6htRYAYwyJiYlHfGwRETl2KRiKiIgEQVZWFn+57TY6LF1KFtAf8ADF/fvDgw9i+vev6R08kNjY2MMOhpGRkZSXl2OMISEhgdatW+uRUREROSQKhiIiIkdg/fr1vHDddbSeM4dZxhAPlDdvzmdduhB76630Ouec33WcyMjIw6r/pJNOAmDz5s2Eh4fjcDhIPkCPpIiIyK9RMBQRETkMJTk5zL7mGpp/8gnPWIvX5YILLoAbbiCyf38G/Erv4IGYQywP7DPNhKacEBGRI6FgKCIi8ntZi3/+fNbcey/N5s3jPGvZHB9P4Z13knDrrXAEE9H/VjA0xmCtJTU1laSkJCIiInC59J9xEREJDv0XRURE5Lfk58OECZQ+/zzRa9fSDMhMSaHtP/9J12uu+dV3B3+vX5tXMCkpicrKSkpLS0lKSiImJuaI6xMREdmbgqGIiDQoe88xmJ6evs9vfr+fysrKwFJWhvnsMyLffpuoTz/F6fXyI/BefDw9n3ySi//0p6BPEu9wOAgPD6eyshK/30/Tpk1rppqoqqpix44dREdHB7VOERERUDAUEZEG5OuvvyYjI4OqqqqaETyBmjDo8/loC/wRuBJoCeQBrwHjgeUuF59MmcLpp59eK+0zxhATE1MTUBMSEmoGpQkPD6dly5a1Uq+IiIiCoYiINAi7d+/mxhtvpKqqqmZb+/btOfXUU4kxhh7r13Py0qW03rABvzFsPv54Mvv1Y1xuLu9OnYrf78dpLd9++22tBcPu3btjraWyspLt27cTERFRK/WIiIjsT8FQRESOeUuXLmXkyJGsXr2asLAw/H4/7rAwxl17Ld0WLIB334Xdu6FdO3j8cRxXXEGr5s1pBURlZzN15kw8Hg9ut5uMjIxaa6cxBmMMkZGRtGnTptbqERER2Z+CoYiIHLOstbz22mvcfPPNJCYm8tlnnxFdXEzxiy9y6ooVRF1/PURFwUUXwdVXQ9++vxhIJj09nczMzIO+lygiInIsUDAUEZFjUmlpKTfffDNvvPEGgzMyePvii4l/4gmYNQv8fujTBx56KBAKY2N/9Vjp6ekKhCIickxTMBQRkWPOihUrGHnBBSSvWMH3PXvSY+FCTFYWtGgB994LV1wBxx1X180UEREJGQqGIiJyTJn67LMsu/dePvb5aA2wahWMHBkIgxkZEOQpJkRERI4FCoYiIlL/FRXhmTiRDX/7G+du385woKpPH7juOjj/fNCE8CIiIr9KwVBEROqnqiqYPZtdo0cTP3cubp8PC3zavz+njx9PeNu2dd1CERGRekPBUERE6g9r4dtvYcKEwBQTO3digFeBN4CLn36aP991Vx03UkREpP5RMBQRkdC3Zg289VYgEK5Zg9flIjslhWcdDj72+6kCnE4n53o8dd1SERGReumw38A3xrQwxnxmjFlhjFlujLm9enuSMWa2MWZ19Wdi9XZjjHneGLPGGLPEGNNzr2NdWV1+tTHmyiM/LRERqfd27oQXXsCeeip06ID/kUf4Zvt2rgKSvV7+GBWFe+RIHOHhOJ3OWp98XkRE5Fh2JD2GXuAua+1CY0ws8L0xZjbwRyDTWvuEMeZe4F7gL8BgoEP10ht4CehtjEkCHgZ6Abb6ONOstQVH0DYREamHvpk7l+1jxtBnwwaSFizA4fez3OXideA9Y2h14okMHTqU7GHDOP744zHGkJ2drcnnRUREjtBhB0NrbQ6QU/19tzFmBdAMGAFkVBd7A8giEAxHAG9aay0w3xiTYIxpUl12trU2H6A6XJ4NvH24bRMRkfrDVlaS//bb5D7/PF0WLaI3sInAe4MfRkfT6pxzGDp0KPcOHkxKSsov9tfk8yIiIkcuKO8YGmNaAz2Ab4C06tCItTbHGNOoulgzYPNeu22p3naw7SIiUs/t35tXWVnJihUrWLJoEaXTp9MqO5vTtm8n2VosMIHA/xX8Evjj1VfzxUsv4Xa76/YkREREGoAjDobGmBhgCnCHtbbYGHPQogfYZn9l+4Hqug64DiA1NZWsrKxDbq9IbSspKdG1KSHpaF+bS5cu5a677sLr9WKMoXFaGm137OBCv5+LgMZAqcPBN02b8tNJJ7G6VSteeOUVvF4vYWFh9OzZk6+//vqotVfqlu6dEqp0bUpDcUTB0BgTRiAUvmWtfb968w5jTJPq3sImQG719i1Ai712bw5sq96esd/2rAPVZ60dA4wBOO6446wGGZBQtKd3RCTUHK1rc9u2bbz22ms8/fTTVFVV0R24xFou27mT5n4/vrAwSgcMwHfNNUQPG8bAyEgGVu87ctQovS/YQOneKaFK16Y0FIcdDE2ga3AcsMJa++xeP00DrgSeqP78cK/ttxhj3iEw+ExRdXj8FPjHntFLgTOB+w63XSIicvT5fD4+/fRTxowZw/Tp0znO5+PJJk3oV1hIJwKjle0++WS48UacI0YQFxd3wOPofUEREZG6cSQ9hn2Ay4GlxpgfqrfdTyAQvmeMuYbA+AEXVv82AxgCrAHKgKsArLX5xpjHgAXV5f62ZyAaEREJbVu2bGH8+PGMHTuW2M2b+WNUFK8kJJCWlwfbt1PUsyeftGhB8vXXc/LgwXXdXBERETmIIxmV9EsO/H4gwOkHKG+Bmw9yrPHA+MNti4iI1L49A8n07duXwsJCxowZw/qPP+YCa5kXHU0rwJaXY04+GS66CM4/n/jGjTm7rhsuIiIivykoo5KKiEj982vz/3m9Xnbv3l2zZGdnc8stt+DxeOhkLSOBp1wuOlmLNQbTqxdceCHmggugceO6OSERERE5bAqGIiIN0KRJkxg1ahQ+nw+Hw0GbNm3w+/3s3r2bkpISKioq9infhcCEtCOBboAfID0dLrooEAabNDnq5yAiIiLBo2AoItKAWGsZO3YsN998Mz6fDwC/34/b7aZnz57ExsYGlpgY2hUW0nnlStr+8ANxOTn4ga+BP4eFcenkyfQaPrxOz0VERESCR8FQRKSW/dojm0dTfn4+I0aM4KOPPqJXr14sW7aMqqoq3G4348aNI/2UUyA7G6ZMgXHjYNMmcDphwAD4619Z2KwZ85Yt48KMDHpp5FAREZFjioKhiEgt+vrrr+nfvz8+n4/w8HDmzp1bJ+Fw2rRpXHPNNZSXlzN69GhuvfVWvvnmG77IzGRYXByd33wTzj8ftm8HtxvOPBMefRSGDYPkZAB6Ab2GDTvqbRcREZHap2AoIlKLnn76abxeLwAVFRVceuml/OMf/+D8888nPDy81usvKSnhzjvvZOzYsbRv356pU6fSpU0bmD6d9A8+IH3aNMjPh6goGDIELrgg8HmQeQZFRETk2KRgKCJSS1atWsWMGTNwOBwYY3A4HFRVVXHppZeSkpLC1VdfzXXXXUe7du1qpf7s7Gwuu+wy1q9fz99uv53zwsLo8uCD8OmnUF4O8fGBHsELLoCzzoLIyFpph4iIiIQ+R103QETkWOTxeLj00kuJiYnhgw8+4LHHHuPzzz9n06ZNzJo1i379+vHMM8/Qvn17zjrrLKZOnVrTs3ikqqqqeOihh7isTx+uLCggv3t3HvrPf+j69NPw3Xdw9dUwezbs3An//S+ce65CoYiISAOnHkMRkVrw4IMPsnDhQqZOncrw4cMZvtcInoMGDWLQoEFs3bqVcePG8eqrr3LeeefRtGlT/vSnP9GjRw9+/PHHQx+sxlrWv/8+s2+9lZE5OTwGUFAAzZvD/ffzXfPm9LruOjAm6OcrIiIi9ZuCoYhIkM2ZM4ennnqK66+/nhEjRhy0XLNmzfjrX//K/fffz4wZM3j55Zd59NFHa353Op2MGjWKzp07Ex8fX7MkJCT8vB4RQezChewYOxbnjBm0KS/nGiC/c2e45hoYMQKqH1UtycpSKBQREZEDUjAUEQmiXbt2ccUVV9CpUyeeffbZ37WPy+Wq6VW85557eOaZZ7DW4vP5mDhxIn6/f5/yqcAQYBhwJoF3AuKB2cA0Y7j4v//lzD/8IajnJSIiIsc2vWMoIhIk1lquvfZa8vLyePvtt4mKijrkY5x//vlERETgdDqJjIxk3rx5lJWWsjMzk9w77mB3167sMIbXgbMTE1mXns79J5xACnAu8IbDwfebNgX5zERERORYpx5DEZEgGTNmDB9++CHPPPMMJ5544mEdIz09nczMTObNns2w2FiOnzgRRo0ick/YO/nkwPyCQ4cSfeKJdDeGsuxsRp9+Ok6PB7fbTUZGRvBOSkRERBoEBUMRkSBYsWIFd955J2eeeSZ33HHH4R1k0yaYOZP0mTNJnzMHSksD8wsOGgR//WtgfsEmTX6x254wmZWVdegD1oiIiIigYCgicsQqKysZNWoU0dHRvP766zgcv/Mp/aoq+OormDEjsCxfHtjeqhVccUVgjsEBAyAi4jcPlZ6erkAoIiIih03BUETkCN1///0sXryYjz76iCYH6NHbx9at8MkngSA4ezbs3g1hYdCvH1x1VaBXsFMnjR4qIiIiR5WCoYjIEZg1axbPPvssN910E0OHDv1lgaoqmD8fZs4MhMHFiwPbmzeHUaNg8GA4/XSIjT26DRcRERHZi4KhiMhh2rlzJ1deeSWdO3fm6aef/vmHtWth1iz49FOYOzfQK+h0Qt++8OSTgTDYtat6BUVERCRkKBiKSMjKzs4O2QFVrLVcffXVFBQUMPv994mcNevnMLh2baBQ69Zw6aVw5pmBXsH4+Dpts4iIiMjBKBiKSEj6+OOPGTFiBNZawsPDyczMDJ1w6Pcz+b77OGH6dF5q147m/fqB1wvR0YHBYu64IxAGO3RQr6CIiIjUCwqGIhJyVq1axeWXX47P5wMCo35mZWXVSjCcNm0aWVlZDBo0iDOdURQIAAAZ2ElEQVTPPBOn03ngguvXw5w5kJmJd9YsLiwo4ELAxsfD3XfDWWfBaaeB2x30NoqIiIjUNgVDEQkpX3/9NcOGDavpKaysrMTv97NmzRqstZgg9cBZa7nrrrt47rnnAGo+U1JSSElJoV18PP39fnrv3k3X7dtJKiwEoCQujhkeD3Oiovj7N9/QqGvXoLRHREREpC4pGIpIyPjggw+49NJLadGiBTNnziQ3N5fMzEzmz5/P+PHjsdbyyiuvEBYWdkT1FBcXc+211zJp0qSabTHGcFO3bgwyhnYbN9Jm5UoAiowhC5gNzAFWFRcD4Ha7uWr3bhodUUtEREREQoOCoYiEhBdeeIHbbruN3r1789FHHwV67dq1Iz09HWstjz76KI8++ig5OTlMmjSJmJiYw6pnyZIljBw5ki1r1/LMuedSOn06/b1eTrUW95IlgUdB+/QJDBZzxhnEn3QSw51OMoqLiXjsMVY/9xx+vx+fz1drj7eKiIiIHG0KhiJSp/x+P/feey9PPfUUI0aMYOLEiURFRe1TxhjDI488QvPmzbnhhhvo378/H3/8MY0bN/79FVVVMf3hh1nw1FOMczg4LSwM59SpWGPY3rQpOwcMoNmVVwZC4f71A/Hx8VxwwQW8+OKLeDwe3G43GRkZR/4HICIiIhICFAxFpM5UVlZy1VVX8fbbb3PTTTfx/PPPH3zwF+Daa6+ladOmXHjhhaSnpzNz5kw6dep04MJeL3z/PXz2Gb45c6j6/HOGer0MBao6d8Y5aBBkZGD69aNJUtLvam96ejqZmZkhO4WGiIiIyOFSMBSROlFYWMh5551HVlYWTzzxBP/7v//7uwaWGTJkCJ9//jnnnHMOp512GtOmTaNv375QVQWLFsHnn8Nnn8G8eVBSAsD68HA+8XpJHjmSi154gbC0tMNud3p6ugKhiIiIHHMUDEXkqNu8eTNDhgxh1apVTJgwgT/84Q+HtH+vXr2Y/9lnPHTWWWRmZNChSxfS1q6F0tJAgeOPh8sv5+vwcP4wZgwlUVFM+PBDzjrrrFo4GxEREZH6z1HXDRCRhmXJkiWkp6ezadMmPvnkk98fCouKYOZMuO8+6NuXNieeyIQtW3jI5yNnyRJ+6NkT3nsPcnKoXLSIW51O+oweTdMTT2TRokUKhSIiIiK/Qj2GInJUZGdnM27cON5++20SExP58ssv6dat28F32L4dvvoq8EjoF1/A4sXg94PLBb16wZ13wv/8D5U9e/K3W27hgw8+YFTz5jT79ls+/vhjVqxYwZ///GeeeOKJI57eQkRERORYp2AoIrXK4/EwduxYbr/9drxeL8YYJk6cuG8o9Pngxx8DQfCrr+Drr2HdusBvkZFw6qnw0EPQrx/07g3R0TW7RgKTJk3i4osv5u23367Z/s9//pN77733KJ2liIiISP2mYCgiQWWtZenSpcyZM4c5c+bw+eefU1ZWVvO7w+Fg9aJFEBPzcwjMzobqieNp1CgwZcRNN8Fpp8FJJwXmFvwVTqeTnj178v7772Otxel0Yq2tzdMUEREROaYoGIrIYcnOzq6ZtqFZs2ZkZmYye/ZsMjMzyc3NBaBTp05cfdVVdImKInv0aE72eunj93Pi3/8e6CU0Brp0gVGjAmHwtNOgbdvA9kM0YMAAIiIiNMegiIiIyGFQMBSRfULe3lMx+P1+ysvLKS0tpaSkhNLSUkpLS1mwYAF33303VVVVADW9c40bN2Z4RgYjW7fmVGuJX7kyMCDMzp3cAHjcbsq7dcMMHhwIgqeeCgkJQTkHzTEoIiIicvgUDEUaMGstL730Erfddhs+nw9jDI0bN8bn81FSUrLPI6AH4gA6A9d168ZlHTqQsGoVZtIk2PMYZ6dOMGRIIACeeirurl1xu2rvtqM5BkVEREQOj4KhSAO0YcMG3nzzTd544w3W7RnkhUBQTEtLo3fv3sTExBAdHf3zEhVFalkZTbdtwy5YQPGcOZxkLbEAS5fC1q2BgWEuuigQBE8+GRIT6+wcRUREROT3UzAUaSBKS0uZMmUKr7/+Op999hnGGAYOHMhll13GU089VfNu3osvvhjodcvJge++gwUL4LPPAt937QocLCyMkuOOY1VCAnGDBtHh8suhffvDejdQREREROqegqHIMcxay7x583j99deZNGkSJSUltGvXjscee4zLL7+cVq1aATA0PZ21773HaWFhtHzyyUAI3Lo1cBCHIzBAzLBhgV7Ak0+Gbt2ICQ/npDo8NxEREREJHgVDkWNMdnY2H3zwAQUFBcydO5d169YRExPDxRdfzB+vvJI+bdtiFi2CN96ARYtg4UJO3rSJk/ccoGNH6N8/EAB79YIePfaZN1BEREREjj0KhiLHCL/fz+jRo7nnnnvw+/0ADO/alfG33EJ6RATuZctg5EionkoCY6BDh8AUEbfcAj17BuYMDNIooSIiIiJSfygYitRzRUVFvDlmDJnPP0/ili08A3QHegAJy5bBsmXgdAYeBx0yJNAD2LMndO8OsbF13HoRERERCQUKhiL1ibWwbRssWcKOWbPYNH06MWvXcqO13FpdpAxYBkx2Ohlw1120u/BC6NoVIiLqsOEiIiIiEsoUDEVCVWkp/PhjoMdvyRJYvBi7ZAkmLw+ANKASyGvRgp2DBtHk7LOhe3cW5+aSNW8eGRkZtNOcfiIiIiLyOygYitS1igpYtSoQAJctg+XLA5/r19cU8bhcbEtKYn5FBV8BOamp9LnxRi675RZ6pKbuc7j0jh1J79v3KJ+EiIiIiNRnCoYiR0tlJT9Mnsz6jz+mV2QkLYqKsMuWwerVmOrBYvxOJ0WNG7MjOZlNffuyoLyc/y5cyGqvF39uLj179uTBBx9k2LBhuFz611dEREREgkN/sxQJMmdJCcyfDytXwooVsGIFdsUK7Lp1nOj3cyLgA9Y4HCz1+1kGNctqn4+qrVth61bcbjfh4eHstjZwXKeTkSNHct5559XdyYmIiIjIMUnBUOqF7OxssrKyyMjIIP033pvz+XwUFRUxd+5cli1bxllnnfWb+xwyvx82b4bVq/cJgKxcyf/k5PzcFpeLrZGRfF9ezjK/nxXACmAlcEKvXpx55pkkJyczLCWFq1JSSNlriY6OZv78+Zx++ul4PB7cbjcZGRnBPQ8RERERERQMpR7Izs5m4MCBeDweXC4XN9xwA3FxceTn5/9iKSgooLCwEFvdywbw6KOP0r17d3r27EmHDh1o3759zRK713QNvwif1gbm/Pvpp0AA/Omnn7+vWRN4N3CP+Hj8nTqx44QTmJuaSnZhIZ9u2sR6r5emCQkMHjWKNm3aMPVvf6sJeaNHj/7NwJqenk5mZubvDsUiIiIiIodDwVBCWn5+Pvfddx8V1SHM4/Hw/PPP43A4SExMJCkpiaSkJFJSUujYsWPN+vfff8+MGTOw1mKMoaioiE8++YTXXnttn+OnNWrEKa1a0cHpJH/BAtr4fGx0OGjduDEpBQWElZfXlPW7XFQ0a0ZFy5Z4Ro7E06YNK6qqePeHH1hbWsrCRYsoKSnB5XLRr18/rr/1VgYPHkznzp0xxgDQv3//Qw556enpCoQiIiIiUqsUDCUkFRcXM3r0aJ555hmKi4txOp1Ya3G73Xz44YecccYZOByOg+6fnZ3N3Llz8Xg8RIeF8f4TT9AjPp7KH3+keNEiqlatImzTJuLy8gjPza3Zzw9s8PtZvG0bq4GfqpfVwCavF9/GjbBxI8ybt099xhiGDx/OVVddRVhYGEOGDDlguxTyRERERCQUKRhKSCkrK+PFF1/kiSeeIC8vj/POO49HH32UkpKSg/e0eb2wZUsgsG3YABs2kL5xIznHHw/r1hFXXIy55BIAwoHUyEho2xZ69w58tmvHispKLnngAVZXVWHdbsaNG0f/rl05raoKr9dbs1Tttf7uu+8yceJE/H4/DoeD3r17M2LECLKyso72H5uIiIiIyBFRMJSQUFlZydixY3n88cfJycnh7LPP5rHHHqNXr15QXg5bt5Leu3dggJeZM/cJgWzdCj7fvgds2pT4Vq3gnHOgXbuaAEi7dtC4MVQ/2rnH8cDLp512SI95JicnM2XKFA0MIyIiIiL1noKh1ClvaSkf/Oc/TH7uOcJyc3mkdWuGnnsuTX0++NOfAj2Bu3btu5PDAc2aQevW0K9f4LN1a2jVKvDZogWEhx9yWw71MU8NDCMiIiIixwoFQzliB5xKwucLjOi5bVtgycnZ53vJTz/h27iR+MpKLgQu3HOwDRuguBiaNw8EvN69A597ltatA7+FhdXJue5P7wyKiIiIyLGg3gbD/Px8srOz9Zfyo8laKCqCHTsCoW/HDn768kvmPP88jawlD1gRFUUjr5cEjwfnfrv7gcKwMHY4nayvqGAzsBXoPXIkQ66/HtOiRSD0RUcf9VMTEREREWnI6m0w3LVrFwMGDGD69OmcccYZdd2c+mlP0MvLCzyuuedz1659wt8+nx7PPofoCDwA7CIQ8nL8flYnJJDndpMfHk5eeHjge0QEBW431uVi48aNrF27FgCn08ljPXtyjv4ZioiIiIjUmXobDCEwYMmgQYNISUmhXbt2tG3blnbt2tUsbdu2pUmTJnzzzTfH9ntgfn/g8cuCAigsDHzu/T0//+fAt3cIzMsLjOh5IGFhkJYWWBo1gm7dAp9paZTFxfHq1KmM+/hjSE1lbVERlT4fbrebzMxMBv7Gn3F2djann366Bm0REREREQkRIRMMjTFnA/8GnMBYa+0Tv7WP2+3m2muvxev1snbtWubPn8+7776L3+/fp0xVVRXWWsLCwnjttde49NJLayYcr1PWQmUl7N4dCHbFxQf+vve2AwXAoqLAsQ7G5YKUFEhODnx26rTv+t7f93zGx/9i5E6A6dOnc+ONN7Jt2zZuv/NOHnvsMZYsWXJIwVuDtoiIiIiIhJaQCIbGGCfwf8AgYAuwwBgzzVr748H2SUlJYdq0ab8IFVVVVTWPKq5bt46JEyfy5Zdf1vx22WWXceedd5KRkUFGRgYDBgygU8eOGK830HtWVfXzZ2UlVFQc/HP/beXlUFYGpaWBZc/3X9t2sB67ff+AIC4OYmMDn4mJ0KQJdO4MCQmB9cTEA37/9qefyPz2WzIGDDiiAJabm8vtt9/OO++8Q9euXZkyZQqnnHIKcHgDsGjQFhERERGR0BESwRA4BVhjrV0HYIx5BxgBHDQYtti9m/Sbbgo8RmltzWeY3097a2nv94Pfz1Xl5Ww3Boe1hAFxkZGBHrbJk3FNmoQLCGrfYVhYYPCUqKjA557vMTGBRzH33hYV9XPY2zv47b8tOvqAvXcHUlxczJYtW9iyZQtZ777LU089ha/6Mc85c+bQt2/fQzoday0TJkzgjjvuoKSkhL/97W/85S9/we12H86fjoiIiIiIhCBjf+0RxKPVCGNGAmdba6+tXr8c6G2tvWW/ctcB1wF0drtP+rxXLzAG63AA1HzicGCNqfmtsKiI/IIC4lNTSUhJwbpc+J1OdldUkLNzJ9tyc9m8fTtFZWVUAc7wcMo8HiqspcrpZNRVV9G8XTv8bndgCQur+W7322ZdtZe1lyxZwhdffEFaWhqRkZHs3LnzF0tZWdlB9w8LC+PEE0+kW7dudO3alU6dOhEZGXnQ8tu3b+fZZ59lwYIFdOnShbvvvpvWrVvXwpkdW0pKSoiJianrZoj8gq5NCWW6PiVU6dqUUDZgwIDvrbW9gnGsUAmGFwJn7RcMT7HW3nqwfY477ji7atWqoLXBWsvatWvJysrihRdeYPHixTW/xcfHc9555zFw4EAGDBhA8+bNg1bvryksLCQ7O5uvvvqKGTNmsGjRon1+N8bQpEkTmjdvXrO0aNGi5ntubi6XXXYZHo8Hp9PJkCFDWLNmDcuXL8dai8vlokePHvTt25c+ffrQp08f1q9fz9y5c9m1axevvvoqAP/85z+56aabcDr3n4BCDmTPu5MioUbXpoQyXZ8SqnRtSigzxgQtGIbKo6RbgBZ7rTcHth3NBhhjaN++Pe3bt6dLly41o2Y6HA569uzJRx99xOuvvw5Ax44dGThwIAMHDiQjI4PU1NQDT/J+CKy1bNiwga+++qpmWbZsGdZanE4naWlpGGOw1uJwOLj77rv5+9//TthvTPR+oEFeCgoKyM7O5ssvv+Srr77ipZde4rnnnqv5c9jzPwt69+7Nu+++S6tWrQ75fEREREREpP4IlWC4AOhgjGlDYDq8S4BL66oxBxo10+/3s3TpUubOnUtmZiYTJkzg5ZdfBqBdu3Zs3LgRv9+Py+Xi4Ycf5vjjj8fhcGCMweFwHHBZvnw5n376KeXl5axcuZJt2wJZODY2lvT0dC688EL69OnDKaecwtKlS/eZ4uHcc8/9zVC451z2D6qJiYkMGTKEIUOGAIFpPxYuXMhjjz3GzJkzAXA4HAwfPlyhUERERESkAQiJYGit9RpjbgE+JTBdxXhr7fK6bNP+gcrhcNC9e3e6d+/OnXfeSVVVFd9//z1z585l/PjxeKtHF/V4PDzwwAOHVJcxhjPOOIMHHniAPn360LVr1188tlmbUzyEh4eTnp7OQw89RFZWVk34HDBgQNDqEBERERGR0BUSwRDAWjsDmFHX7fi9wsLCOPXUUzn11FMZMGBATW9eWFgYr776Kt26dcPv99cs1tp91t98803GjRuH3+/H4XAwYMAAbrrppl+ts7aneND8giIiIiIiDVPIBMP67HACldPpZMKECTW9c6HyUrPmFxQRERERaXgUDIPkUAOVeudERERERCRUKBjWIfXOiYiIiIhIKHDUdQNERERERESkbikYioiIiIiINHAKhiIiIiIiIg2cgqGIiIiIiEgDp2AoIiIiIiLSwCkYioiIiIiINHAKhiIiIiIiIg2cgqGIiIiIiEgDp2AoIiIiIiLSwCkYioiIiIiINHAKhiIiIiIiIg2cgqGIiIiIiEgDp2AoIiIiIiLSwBlrbV234bAYY3YDq+q6HSIHkALsqutGiByArk0JZbo+JVTp2pRQdpy1NjYYB3IF4yB1ZJW1tlddN0Jkf8aY73RtSijStSmhTNenhCpdmxLKjDHfBetYepRURERERESkgVMwFBERERERaeDqczAcU9cNEDkIXZsSqnRtSijT9SmhStemhLKgXZ/1dvAZERERERERCY763GMoIiIiIiIiQaBgKCIiIiIi0sDVu2BojDnbGLPKGLPGGHNvXbdHGhZjTAtjzGfGmBXGmOXGmNurtycZY2YbY1ZXfyZWbzfGmOerr9clxpiedXsGcqwzxjiNMYuMMdOr19sYY76pvjbfNca4q7eHV6+vqf69dV22W459xpgEY8xkY8zK6ntouu6dEgqMMXdW/zd9mTHmbWNMhO6dUleMMeONMbnGmGV7bTvke6Ux5srq8quNMVf+nrrrVTA0xjiB/wMGA52BUcaYznXbKmlgvMBd1trjgVOBm6uvwXuBTGttByCzeh0C12qH6uU64KWj32RpYG4HVuy1/iTwXPW1WQBcU739GqDAWtseeK66nEht+jfwibW2E9CdwHWqe6fUKWNMM+A2oJe1tivgBC5B906pO68DZ++37ZDulcaYJOBhoDdwCvDwnjD5a+pVMCRwYmusteustR7gHWBEHbdJGhBrbY61dmH1990E/mLTjMB1+EZ1sTeAc6u/jwDetAHzgQRjTJOj3GxpIIwxzYFzgLHV6wYYCEyuLrL/tbnnmp0MnF5dXiTojDFxQD9gHIC11mOtLUT3TgkNLiDSGOMCooAcdO+UOmKt/QLI32/zod4rzwJmW2vzrbUFwGx+GTZ/ob4Fw2bA5r3Wt1RvEznqqh8f6QF8A6RZa3MgEB6BRtXFdM3K0TQa+F/AX72eDBRaa73V63tffzXXZvXvRdXlRWpDW2An8Fr1o85jjTHR6N4pdcxauxV4GthEIBAWAd+je6eElkO9Vx7WPbS+BcMD/R8ZzbchR50xJgaYAtxhrS3+taIH2KZrVoLOGDMUyLXWfr/35gMUtb/jN5FgcwE9gZestT2AUn5+FOpAdH3KUVH9eN0IoA3QFIgm8Hje/nTvlFB0sOvxsK7T+hYMtwAt9lpvDmyro7ZIA2WMCSMQCt+y1v5/O/fzYmMUBnD8+6Qm2fiRpSQlW8sJiwnNQrIaUUQjf4INO4vZWtliYaGkidnZsLASNQvFbohZGGrKRsnisTjnZVB0dd13rvP9bG73vG/ds3h6zn3e9zxnvg6vdNuc6uf7Om7MalQOAMcj4jVlm/0hyhvELXV7FPwYf99is17fzK9bV6RhWQaWM/NJ/X6XUiiaO9W3I8CrzPyQmV+AeWA/5k6tL4Pmyr/KoeNWGD4F9tSToiYozcELPc9JDal9BNeBl5l5dc2lBaA78ekccH/N+Nl6atQk8LHbCiANU2ZeyswdmbmLkhsfZuZp4BEwU2/7OTa7mJ2p9/vUW/9EZr4D3kbE3jp0GHiBuVP9ewNMRsSmusZ3sWnu1HoyaK58AExHxNb6Vny6jv1WjFssR8RRylPwDcCNzJzreUpqSEQcBB4Dz/nex3WZ0md4B9hJWWROZOZqXWSuURp+PwGzmfls5BNXUyJiCriYmcciYjflDeI2YBE4k5mfI2IjcIvSJ7sKnMrMpb7mrP9fROyjHIw0ASwBs5QH1OZO9SoirgAnKSePLwIXKP1Y5k6NXETcBqaA7cAK5XTRewyYKyPiPOU/KsBcZt7842+PW2EoSZIkSRqucdtKKkmSJEkaMgtDSZIkSWqchaEkSZIkNc7CUJIkSZIaZ2EoSZIkSY2zMJQkSZKkxlkYSpIkSVLjvgKqgb9s4TZgOAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "## Projections visualization \n",
    "\n",
    "# Plot the simulated path of Y(t) and the calculated path of X(t) from the last inspection date to\n",
    "# 750 that have been used to generate the learning data set\n",
    "fig = plt.figure(figsize=(15, 7))\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "ax.plot(vTime, vObs, '-k.')\n",
    "ax.plot(vTime, vState, '-r')\n",
    "# plot \n",
    "# i) Y(t) simulated with the estimated parameters of X(t)\n",
    "ax.plot(vTimeProj, mEst, color=[0.8, 0.8, 0.8])\n",
    "# ii) X(t) calculated with estimated parameters\n",
    "ax.plot(vTimeProj, vStateEst, 'b')\n",
    "# iii) X(t) calculated with the true parameters from the last inspection date to 750\n",
    "ax.plot(vTimeProj, vStateProj, '--r');\n",
    "ax.plot([0, tEnd], [L, L], 'r:', LineWidth=2)\n",
    "ax.set_xlim([0, tEnd])\n",
    "ax.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We calculate the hitting time of the failure level "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7f8540430850>]"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3AAAAGbCAYAAACF9nK/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3df9BldX0f8PcnECk/GjCw2lZolgy4M6QkNCpEpyaLoFlEJFNXwUQCjc5alcyQYGcWR02KjKKVSpqiZScaCCaFDabN6q61ovOknWgMkKCAilkICQuJKWDIrKIE/faPcxkfH59l73P3LmfPntdr5jv33nPPc57P8+F6xvd+v+fcaq0FAACAfd8P9F0AAAAA0xHgAAAABkKAAwAAGAgBDgAAYCAEOAAAgIE4sO8CljriiCPacccd13cZo/T1r389hx56aN9ljJb+90fv+6P3/dH7/uh9f/S+P3q/MrfeeuuDrbVVy723zwW4Zz7zmbnlllv6LmOUFhYWsnbt2r7LGC3974/e90fv+6P3/dH7/uh9f/R+Zarqr3b1niWUAAAAAyHAAQAADIQABwAAMBACHAAAwEAIcAAAAAMhwAEAAAyEAAcAADAQAhwAAMBACHAAAAADMVWAq6p1VXVXVW2vqo3LvP/TVfVnVfV4Va1f8t75VfUXk3H+vAoHAAAYm90GuKo6IMlVSc5IckKSV1fVCUt2++skFyT5vSU/+8NJfi3JKUlOTvJrVfX0PS8bAABgfKaZgTs5yfbW2j2ttceSXJ/k7MU7tNbuba19Icl3lvzszyb5ZGvt4dba15J8Msm6OdQNAAAwOgdOsc+zkty36PWOdDNq01juZ5+1dKeq2pBkQ5KsWrUqCwsLUx6eedq5c6fe90j/+6P3/dH7/uh9f/S+P3rfH72fn2kCXC2zrU15/Kl+trW2KcmmJFmzZk1bu3btlIdnnhYWFqL3/dH//uh9f/S+P3rfH73vj973R+/nZ5ollDuSHLPo9dFJHpjy+HvyswAAACwyTYC7OcnxVXVsVT0tyblJtkx5/E8keUlVPX1y85KXTLYBAACwQrsNcK21x5NcmC54fSnJ5tbanVV1aVW9PEmq6nlVtSPJK5NcXVV3Tn724STvSBcCb05y6WQbAABAf9au7cbATHMNXFpr25JsW7Lt7Yue35xueeRyP/uhJB/agxoBAADIlF/kDQAAQP8EOAAAgIEQ4AAAAAZiqmvgAAAA9itXXtl3BTMR4AAAgPE56aS+K5iJJZQAAMD43HRTNwbGDBwAADA+l13WPZ5+er91rJAZOAAAgIEQ4AAAAAZCgAMAABgIAQ4AAGAg3MQEAAAYn6uv7ruCmQhwAADA+KxZ03cFM7GEEgAAGJ+PfrQbA2MGDgAAGJ8rrugezzqr3zpWyAwcAADAQJiBAwAYqNUbt87tWPdefubcjgXsPWbgAAAABkKAAwAAGAhLKAEAgPG57rq+K5iJAAcAAIzPMcf0XcFMLKEEAADG54YbujEwZuAAAIDx+cAHusdzzum3jhUyAwcAADAQAhwAAMBACHAAAAADIcABAAAMhJuYAAAA43PjjX1XMBMBDgAAGJ+jjuq7gplYQgkAAIzPNdd0Y2AEOAAAYHwEOAAAAPYmAQ4AAGAgBDgAAICBEOAAAAAGwtcIAAAA47NtW98VzESAAwAAxueQQ/quYCaWUAIAAOPz/vd3Y2AEOAAAYHw2b+7GwAhwAAAAAyHAAQAADIQABwAAMBACHAAAwED4GgEAAGB8Fhb6rmAmZuAAAAAGQoADAADG573v7cbACHAAAMD4fOxj3RgYAQ4AAGAgBDgAAICBEOAAAAAGwtcIAAAA43PwwX1XMBMBDgAAGJ+Pf7zvCmYiwAEA+6XVG7fO7Vj3Xn7m3I4FsCdcAwcAAIzPO97RjYER4AAAgPH51Ke6MTACHAAAwEC4Bg4AANcMwkCYgQMAABgIM3AAAMD4HHlk3xXMRIADAADG5yMf6buCmVhCCQAAMBACHAAAMD6XXNKNgbGEEgAAGJ/PfrbvCmYiwAEAPIVuv/+RXDDHW/YD4zLVEsqqWldVd1XV9qrauMz7B1XVDZP3P1dVqyfbf7Cqrq2q26vqS1U1vDlKAACAfcRuA1xVHZDkqiRnJDkhyaur6oQlu702yddaa8cleV+Sd0+2vzLJQa21E5M8J8nrnwh3AAAArMw0M3AnJ9neWruntfZYkuuTnL1kn7OTXDt5fmOS06qqkrQkh1bVgUkOTvJYkn+YS+UAAACzOvrobgxMtdaefIeq9UnWtdZeN3l9XpJTWmsXLtrnjsk+Oyav705ySpJHklyX5LQkhyT5ldbapmV+x4YkG5Jk1apVz9m8efMc/jRWaufOnTnssMP6LmO09L8/et8fve/PGHp/+/2PzO1YJz7r8Lkd6+8efiRffXRuh9snzbNf8zSGz/2+Su9X5tRTT721tfbc5d6b5iYmtcy2palvV/ucnOTbSf5Fkqcn+b9VdVNr7Z7v2bELdZuSZM2aNW3t2rVTlMW8LSwsRO/7o//90fv+6H1/xtD7ed4o5N5fWDu3Y/3m7/5hrrh9/76P3Dz7NU9j+Nzvq/R+fqZZQrkjyTGLXh+d5IFd7TNZLnl4koeT/HyS/9Va+8fW2t8l+eMkyyZJAACAp8xFF3VjYKYJcDcnOb6qjq2qpyU5N8mWJftsSXL+5Pn6JJ9u3drMv07youocmuSnknx5PqUDAADM6LbbujEwuw1wrbXHk1yY5BNJvpRkc2vtzqq6tKpePtntg0mOrKrtSX41yRNfNXBVksOS3JEuCP52a+0Lc/4bAAAARmGqBdittW1Jti3Z9vZFz7+Z7isDlv7czuW2AwAAsHJTfZE3AAAA/du/b4EEAACwnGc/u+8KZiLAAQAA47Pp+76eehAsoQQAABgIAQ4AABifDRu6MTCWUAIAAOPzla/0XcFMzMABAAAMhAAHAAAwEAIcAADAQLgGDgBgN1Zv3Dq3Y1184twOBeyJk07qu4KZCHAAAMD4XHll3xXMxBJKAACAgRDgAACA8XnNa7oxMJZQAgAA47NjR98VzMQMHAAAwEAIcAAAAAMhwAEAAAyEa+AAAIDxef7z+65gJgIcAAAwPu96V98VzMQSSgAAgIEQ4AAAgPF5xSu6MTCWUAIAAOPz0EN9VzATM3AAAAADYQYOANgjqzdundux7r38zLkdC2B/ZAYOAABgIMzAAQAA43PaaX1XMBMBDgAAGJ+3va3vCmZiCSUAAMBACHAAAMD4nHFGNwbGEkoAAGB8Hn207wpmYgYOAABgIAQ4AACAgRDgAAAABsI1cAAAwPi87GV9VzATAQ4AABifN7+57wpmYgklAADAQAhwAADA+Kxd242BEeAAAAAGwjVwAADss1Zv3Dq3Y12z7tC5HQv6YgYOAABgIAQ4AACAgbCEEgAAGJ9XvarvCmYiwAEAAOPzxjf2XcFMLKEEAADG5xvf6MbAmIEDAADG56Uv7R4XFnotY6XMwAEAAAyEAAcAADAQAhwAAMBACHAAAAAD4SYmAMA+Y/XGrX2XAIzFBRf0XcFMBDgAAGB8BhrgLKEEAADG58EHuzEwZuAAAIDxWb++e/Q9cAAAAOwNAhwAAMBACHAAAAADIcABAAAMhJuYAAAA4/OGN/RdwUwEOAAAYHzOOafvCmZiCSUAADA+993XjYExAwcAAIzPeed1j74HDgAAgL1BgAMAABgIAQ4AAGAgXAMHAMBcrd64te8SYL81VYCrqnVJfiPJAUl+q7V2+ZL3D0ryO0mek+ShJOe01u6dvPfjSa5O8kNJvpPkea21b87rDwAAAFixiy/uu4KZ7DbAVdUBSa5K8uIkO5LcXFVbWmtfXLTba5N8rbV2XFWdm+TdSc6pqgOTfDjJea21z1fVkUn+ce5/BQAAwEqcdVbfFcxkmmvgTk6yvbV2T2vtsSTXJzl7yT5nJ7l28vzGJKdVVSV5SZIvtNY+nySttYdaa9+eT+kAAAAzuuuubgxMtdaefIeq9UnWtdZeN3l9XpJTWmsXLtrnjsk+Oyav705ySpLXpFtW+Ywkq5Jc31p7zzK/Y0OSDUmyatWq52zevHkOfxortXPnzhx22GF9lzFa+t8fve+P3vdnnr2//f5H5nKcsXjmwclXH+27inE69vADnHN6si+e70+66KIkyW1XXtlzJd/v1FNPvbW19tzl3pvmGrhaZtvS1LerfQ5M8m+SPC/JN5J8qqpuba196nt2bG1Tkk1JsmbNmrZ27dopymLeFhYWovf90f/+6H1/9L4/8+z9BW5YsSIXn/h4rrjdfeT6cM26Q51zerJPnu+POCJJ9r26dmOaJZQ7khyz6PXRSR7Y1T6T694OT/LwZPsftdYebK19I8m2JD+5p0UDAACM0TQB7uYkx1fVsVX1tCTnJtmyZJ8tSc6fPF+f5NOtW5v5iSQ/XlWHTILdzyT5YgAAAFix3c7ft9Yer6oL04WxA5J8qLV2Z1VdmuSW1tqWJB9Mcl1VbU8383bu5Ge/VlX/OV0IbEm2tdasswAAAJjBVAuwW2vb0i1/XLzt7YuefzPJK3fxsx9O91UCAAAA+4a3vrXvCmbiCloAAGB8Tj+97wpmMs01cAAAAPuX227rxsCYgQMAAMZn8j1wWVjotYyVMgMHAAAwEAIcAADAQAhwAAAAAyHAAQAADISbmAAAAOPzznf2XcFMBDgAAGB8XvCCviuYiSWUAADA+HzmM90YGDNwAADA+LzlLd2j74EDAABgbxDgAAAABkKAAwAAGAgBDgAAYCDcxAQAABifK6/su4KZCHAAAMD4nHRS3xXMxBJKAABgfG66qRsDYwYOAAAYn8su6x5PP73fOlbIDBwAAMBACHAAAAADIcABAAAMhAAHAAAwEG5iAgAAjM/VV/ddwUwEOAAAYHzWrOm7gplYQgkAAIzPRz/ajYExAwcAAIzPFVd0j2ed1W8dKyTAAQAwCrff/0gu2Lh1Lse69/Iz53IcWClLKAEAAAZCgAMAABgIAQ4AAGAgXAMHAACMz3XX9V3BTAQ4AABgfI45pu8KZmIJJQAAMD433NCNgTEDBwAAjM8HPtA9nnNOv3WskBk4AACAgRDgAAAABkKAAwAAGAgBDgAAYCDcxAQAABifG2/su4KZCHAAAMD4HHVU3xXMxBJKAABgfK65phsDI8ABAADjI8ABAACwNwlwAAAAAyHAAQAADIQABwAAMBC+RgAAABifbdv6rmAmAhwAADA+hxzSdwUzsYQSAAAYn/e/vxsDI8ABAADjs3lzNwZGgAMAABgIAQ4AAGAgBDgAAICBEOAAAAAGwtcIAAAA47Ow0HcFMzEDBwAAMBACHAAAMD7vfW83BkaAAwAAxudjH+vGwAhwAAAAAyHAAQAADIQABwAAMBC+RgAAABifgw/uu4KZCHAAAMD4fPzjfVcwk6mWUFbVuqq6q6q2V9XGZd4/qKpumLz/uapaveT9f1lVO6vqzfMpGwAAYHx2G+Cq6oAkVyU5I8kJSV5dVScs2e21Sb7WWjsuyfuSvHvJ++9LMsyICwAA7H/e8Y5uDMw0M3AnJ9neWruntfZYkuuTnL1kn7OTXDt5fmOS06qqkqSqfi7JPUnunE/JAAAAe+hTn+rGwFRr7cl3qFqfZF1r7XWT1+clOaW1duGife6Y7LNj8vruJKckeTTJTUlenOTNSXa21r7v686rakOSDUmyatWq52zevHkOfxortXPnzhx22GF9lzFa+t8fve+P3vdnnr2//f5H5nKcsXjmwclXH+27inGaZ+9PfNbh8znQSOyL5/uTLrooSXLblVf2XMn3O/XUU29trT13ufemuYlJLbNtaerb1T7/Mcn7Wms7JxNyy2qtbUqyKUnWrFnT1q5dO0VZzNvCwkL0vj/63x+974/e92eevb9g49a5HGcsLj7x8Vxxu/vI9WGevb/3F9bO5ThjsU+e7484Ikn2vbp2Y5pP8I4kxyx6fXSSB3axz46qOjDJ4UkeTjcLt76q3pPkiCTfqapvttb+6x5XDgAAMDLTBLibkxxfVccmuT/JuUl+fsk+W5Kcn+SzSdYn+XTr1ma+8IkdqurX0y2hFN4AAIB+HXlk3xXMZLcBrrX2eFVdmOQTSQ5I8qHW2p1VdWmSW1prW5J8MMl1VbU93czbuXuzaAAAgD3ykY/0XcFMploE3FrblmTbkm1vX/T8m0leuZtj/PoM9QEAADAx1Rd5AwAA7FcuuaQbA+MWSAAwEKvneLfHa9YdOrdjAQzSZz/bdwUzMQMHAAAwEAIcAADAQAhwAAAAA+EaOAAAYHyOPrrvCmYiwAEAAOPz4Q/3XcFMLKEEAAAYCAEOAAAYn4su6sbAWEIJAACMz2239V3BTMzAAQAADIQABwAAMBCWUALACN1+/yO5YOPWvsuAwVo9x//93Hv5mXM7Fvs/AQ4AABifZz+77wpmIsABAADjs2lT3xXMxDVwAAAAAyHAAQAA47NhQzcGxhJKAABgfL7ylb4rmIkZOAAAgIEQ4AAAAAZCgAMAABgI18ABAADjc9JJfVcwEwEOAAAYnyuv7LuCmVhCCQAAMBACHAAAMD6veU03BsYSSgAAYHx27Oi7gpmYgQMAABgIAQ4AAGAgBDgAAICBcA0cAAAwPs9/ft8VzESAAwAAxudd7+q7gplYQgkAADAQAhwAADA+r3hFNwbGEkoAAGB8Hnqo7wpmIsABwBKrN26d27HuvfzMuR0LACyhBAAAGAgBDgAAYCAsoQQAAMbntNP6rmAmAhwAADA+b3tb3xXMxBJKAACAgRDgAACA8TnjjG4MjCWUAADA+Dz6aN8VzMQMHAAAwEAIcAAAAAMhwAEAAAyEa+AAAIDxednL+q5gJgIcAAAwPm9+c98VzMQSSgAAgIEQ4AAAgPFZu7YbA2MJJQAA9Gj1xq1zO9a9l585t2OxbzIDBwAAMBACHAAAwEAIcAAAAAPhGjgAAGB8XvWqviuYiQAHAACMzxvf2HcFM7GEEgAAGJ9vfKMbA2MGDgAAGJ+XvrR7XFjotYyVMgMHAAAwEAIcAADAQAhwAAAAA+EaOADYi1Zv3Np3CQDsRwQ4AABgfC64oO8KZiLAAQAA4zPQAOcaOAAAYHwefLAbAzNVgKuqdVV1V1Vtr6qNy7x/UFXdMHn/c1W1erL9xVV1a1XdPnl80XzLBwAAmMH69d0YmN0GuKo6IMlVSc5IckKSV1fVCUt2e22Sr7XWjkvyviTvnmx/MMlZrbUTk5yf5Lp5FQ4AADA208zAnZxke2vtntbaY0muT3L2kn3OTnLt5PmNSU6rqmqt/Xlr7YHJ9juT/JOqOmgehQMAAIzNNDcxeVaS+xa93pHklF3t01p7vKoeSXJkuhm4J7wiyZ+31r619BdU1YYkG5Jk1apVWVhYmLZ+5mjnzp163yP974/e92df7f3FJz7edwl73TMPHsffuS/S+/6Moff74jk12TfP9yf9/d8nSW7bx+ranWkCXC2zra1kn6r6sXTLKl+y3C9orW1KsilJ1qxZ09auXTtFWczbwsJC9L4/+t8fve/Pvtr7C0bw3W0Xn/h4rrjdzaj7oPf9GUPv7/2FtX2XsKx98nx/xBFJsu/VtRvTfIJ3JDlm0eujkzywi312VNWBSQ5P8nCSVNXRSf5Hkl9srd29xxUDAADsqTe8oe8KZjJNgLs5yfFVdWyS+5Ocm+Tnl+yzJd1NSj6bZH2ST7fWWlUdkWRrkktaa388v7IBAAD2wDnn9F3BTHZ7E5PW2uNJLkzyiSRfSrK5tXZnVV1aVS+f7PbBJEdW1fYkv5rkia8auDDJcUneVlW3TcYz5v5XAAAArMR993VjYKZaBNxa25Zk25Jtb1/0/JtJXrnMz12W5LI9rBEAAGC+zjuvexzYTUym+iJvAAAA+ifAAQAADMT+fR9VAEZj9Qhu/Q8AZuAAAAAGwgwcAAAwPhdf3HcFMxHgAACA8TnrrL4rmIkllAAAwPjcdVc3BsYMHAAAMD6vf3336HvgAAAA2BsEOAAAgIEQ4AAAAAZCgAMAABgINzEBAADG561v7buCmQhwAACwn1i9cevcjnXv5WfO7Vj7pNNP77uCmVhCCQAAjM9tt3VjYMzAAQAA43PRRd3jwL4HToADoDe33/9ILpjjch8A2N9ZQgkAADAQAhwAAMBACHAAAAAD4Ro4AABgfN75zr4rmIkABwAAjM8LXtB3BTOxhBIAABifz3ymGwNjBg4AABift7ylexzY98CZgQMAABgIAQ4AAGAgBDgAAICBEOAAAAAGwk1MAACA8bnyyr4rmIkABwAAjM9JJ/VdwUwsoQQAAMbnppu6MTBm4AAAgPG57LLu8fTT+61jhczAAQAADIQABwAAMBACHAAAwEC4Bg6AFVm9cevcjnXxiXM7FABzNs/z/TXrDp3bscZOgAMAAMbn6qv7rmAmAhwAADA+a9b0XcFMXAMHAACMz0c/2o2BMQMHAACMzxVXdI9nndVvHStkBg4AAGAgBDgAAICBEOAAAAAGQoADAAAYCDcxAQAAxue66/quYCYCHAAAMD7HHNN3BTOxhBIAABifG27oxsCYgQMAAMbnAx/oHs85p986VsgMHAAAwEAIcAAAAAMhwAEAAAyEAAcAADAQbmICAACMz4039l3BTAQ4AABgfI46qu8KZmIJJQAAMD7XXNONgRHgAACA8RHgAAAA2JtcAwcM0uqNW+d2rGvWHTq3YwEA7E1m4AAAAAZCgAMAABgISyiB0bv9/kdywRyXZAIAA7BtW98VzESAAwAAxueQQ/quYCaWUAIAAOPz/vd3Y2AEOAAAYHw2b+7GwFhCOaV53rL83svPnNuxWLkx/LfcV//GedYFADBGU83AVdW6qrqrqrZX1cZl3j+oqm6YvP+5qlq96L1LJtvvqqqfnV/pAAAA47LbAFdVByS5KskZSU5I8uqqOmHJbq9N8rXW2nFJ3pfk3ZOfPSHJuUl+LMm6JO+fHA8AAIAVmmYG7uQk21tr97TWHktyfZKzl+xzdpJrJ89vTHJaVdVk+/WttW+11v4yyfbJ8QAAAFihaq09+Q5V65Osa629bvL6vCSntNYuXLTPHZN9dkxe353klCS/nuRPWmsfnmz/YJKPt9ZuXPI7NiTZMHn5r5Lcsed/GjM4KsmDfRcxYvrfH73vj973R+/7o/f90fv+6P3K/EhrbdVyb0xzE5NaZtvS1Lerfab52bTWNiXZlCRVdUtr7blT1MWc6X2/9L8/et8fve+P3vdH7/uj9/3R+/mZZgnljiTHLHp9dJIHdrVPVR2Y5PAkD0/5swAAAExhmgB3c5Ljq+rYqnpaupuSbFmyz5Yk50+er0/y6datzdyS5NzJXSqPTXJ8kj+dT+kAAADjstsllK21x6vqwiSfSHJAkg+11u6sqkuT3NJa25Lkg0muq6rt6Wbezp387J1VtTnJF5M8nuRNrbVv7+ZXbpr9z2EP6X2/9L8/et8fve+P3vdH7/uj9/3R+znZ7U1MAAAA2DdM9UXeAAAA9E+AAwAAGIheAlxVHVFVN1bVl6vqS1X1/Mn2X66qu6rqzqp6z2Tb6qp6tKpum4z/1kfN+4vlel9VNyzq771Vddui/S+pqu2T/y4/22ftQ7eS3vvcz9cuen9SVf3JpL+3VNXJk32rqv7L5HP/har6yb7rH7IV9n5tVT2y6HP/9r7rH7Jd9P4nquqzVXV7VX20qn5o0f7O93Oykt47389PVa1Z1MfbquofquqiqvrhqvpkVf3F5PHpk/2d7+dohv4758+qtfaUjyTXJnnd5PnTkhyR5NQkNyU5aLL9GZPH1Unu6KPO/XEs1/sl71+R5O2T5yck+XySg5Icm+TuJAf0/TcMdayw9z73e7n3Sf53kjMm216aZGHR84+n+x7Ln0ryub7rH/JYYe/XJvlY3zXvL2MXvb85yc9Mtv1SkndMnjvf99d75/u989/ggCR/m+RHkrwnycbJ9o1J3j157nzfb/+d82ccT/kM3ORfnH463Z0r01p7rLX290nekOTy1tq3Jtv/7qmubX/3JL1/4v1K8qok/32y6ewk17fWvtVa+8sk25Oc/NRWvX+YoffMyZP0viV5Yvbh8Hz3OyrPTvI7rfMnSY6oqn/+FJe9X5ih98zJk/R+TZL/M9ntk0leMXnufD8nM/SeveO0JHe31v4q3ef72sn2a5P83OS58/3eM03/mVEfSyh/NMn/S/LbVfXnVfVbVXVokmcneWFVfa6q/qiqnrfoZ46d7PtHVfXCHmreX+yq9094YZKvttb+YvL6WUnuW/T+jsk2Vm6lvU987udlV72/KMl/qqr7krw3ySWT/X3u52elvU+S51fV56vq41X1Yz3UvL/YVe/vSPLyyT6vTHLM5LnP/fystPeJ8/3ecG6++4+iz2yt/U2STB6fMdnuc7/3TNP/xDl/Jn0EuAOT/GSSD7TW/nWSr6ebTj0wydPTTWH/hySbJ7MSf5PkX072/dUkv7d4zT4rsqveP+HV+d4ZoFrmGL53YjYr7b3P/fzsqvdvSPIrrbVjkvxKJv9aHp/7eVpp7/8syY+01n4iyW8m+Z9Pfcn7jV31/peSvKmqbk3yT5M8Ntnf535+Vtp75/s5q6qnpQvLv7+7XZfZ5nO/h1bQf+f8GfUR4HYk2dFa+9zk9Y3pTnQ7kvzBZBr7T5N8J8lRk+UcDyVJa+3WdOvyn91D3fuDXfU+VXVgkn+b5IYl+y/+F8KjY6nTrFbUe5/7udpV789P8geTbb+f7y4X87mfnxX1vrX2D621nZPn25L8YFUd9dSWvN9YtvettS+31l7SWntOun80unvR/j7387Gi3jvf7xVnJPmz1tpXJ6+/+sTSyMnjE5fp+NzvHVP13zl/dk95gGut/W2S+6pqzWTTaUm+mC51vyhJqurZ6S76fbCqVlXVAZPtP5rk+CT3PNV17w+epPdJcnqSL7fWdiz6kS1Jzq2qg6rq2HS9/9OnrOD9yEp773M/P0/S+weS/Mxk24uSPLF8dUuSX5zcneynkjzyxNIPVmalva+qfzZZeZHq7kz5A0keekqL3k/sqvdV9YwkqaofSPLWJE/c8dD5fk5W2nvn+71i6aqWLe0k3HoAAAEQSURBVOn+4SiTxz9ctN35fv6m6r9z/uwO7On3/nKS351Msd6T5N+lW2Lwoaq6I92ygvNba62qfjrJpVX1eJJvJ/n3rbWHe6p7f7Bc75PvXaucJGmt3VlVm9P9H67Hk7yptfbtp7LY/czUvU93AbzP/fws1/s/TPIbkxnQbybZMNl3W7o7k21P8o18978Ts1lJ79cnecPkc/9oknNba5YzzW653v9iVb1p8v4fJPntxPl+L5i693G+n6uqOiTJi5O8ftHmy9NdmvPaJH+d7hrExPl+7lbYf+f8GZU+AQAADEMvX+QNAADAyglwAAAAAyHAAQAADIQABwAAMBACHAAAwEAIcAAAAAMhwAEAAAzE/wdaYEyvCljR1wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "## Hitting time calculation\n",
    "\n",
    "# Calculate the time upon which the failure level L is reached.\n",
    "# This can be done with: \n",
    "# i) Y(t) simulated with the estimated parameters of the model \n",
    "mTimeProj = np.tile(np.asmatrix(vTimeProj).T, (1, nRnd))\n",
    "mTimeProj[mEst<L] = np.Inf # if Y(t) is below the failure level L, it is put to infinite\n",
    "vTimeL = np.min(mTimeProj, axis=0).T\n",
    "# ii) the true model X(t) used to simulate the data.\n",
    "vTimePorjTmp = np.asmatrix(vTimeProj).T\n",
    "vTimePorjTmp[vStateProj<L] = np.inf\n",
    "timeL = np.min(vTimeProj);\n",
    "\n",
    "fig = plt.figure(figsize=(15, 7))\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "ax.hist(vTimeL, bins=range(int(np.floor(np.min(vTimeL))), int(np.ceil(np.max(vTimeL)))), density=True)\n",
    "#ax.hist(vTimeL, bins=range(int(np.floor(np.min(vTimeL))), int(np.ceil(np.max(vTimeL)))), density=True)\n",
    "ax.grid(True)\n",
    "ylim = ax.get_ylim()\n",
    "ax.plot([timeL, timeL], ylim, 'r--')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}