extract_modes.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pickle\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython import display\n",
    "from scipy import interpolate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "folder = './dataset/20160/2_poles/noise50/'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open (folder + 'slow_2mat_20148_noise.pkl', 'rb') as fp:\n",
    "    data_X2 = pickle.load(fp) # all values 0 or 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # noise\n",
    "# with open ('./dataset/20160/2_poles/noise25_2/slow_2mat_20160.pkl', 'rb') as fp:\n",
    "#     data_X1 = pickle.load(fp) # all values 0 or 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# data_X = np.concatenate([data_X1, data_X2])\n",
    "data_X = data_X2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(20148, 73, 2, 2)\n"
     ]
    }
   ],
   "source": [
    "print(data_X.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "dataset_size = 20148\n",
    "# dataset_size = 4000\n",
    "assert data_X.shape[0] == dataset_size"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "freqs for the 28148 dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(73,)\n"
     ]
    }
   ],
   "source": [
    "# check plot_dataset.ipynb\n",
    "data_freq = np.array([ 0.59941229,  0.79921638,  0.99902048,  1.19882458,  1.39862867,\n",
    "        1.59843277,  1.79823686,  1.99804096,  2.19784506,  2.39764915,\n",
    "        2.59745325,  2.79725735,  2.99706144,  3.19686554,  3.39666963,\n",
    "        3.59647373,  3.79627783,  3.99608192,  4.19588602,  4.39569011,\n",
    "        4.59549421,  4.79529831,  4.9951024 ,  5.1949065 ,  5.39471059,\n",
    "        5.59451469,  5.79431879,  5.99412288,  6.19392698,  6.39373107,\n",
    "        6.59353517,  6.79333927,  6.99314336,  7.19294746,  7.39275155,\n",
    "        7.59255565,  7.79235975,  7.99216384,  8.19196794,  8.39177204,\n",
    "        8.59157613,  8.79138023,  8.99118432,  9.19098842,  9.39079252,\n",
    "        9.59059661,  9.79040071,  9.9902048 , 10.1900089 , 10.389813  ,\n",
    "       10.58961709, 10.78942119, 10.98922528, 11.18902938, 11.38883348,\n",
    "       11.58863757, 11.78844167, 11.98824576, 12.18804986, 12.38785396,\n",
    "       12.58765805, 12.78746215, 12.98726625, 13.18707034, 13.38687444,\n",
    "       13.58667853, 13.78648263, 13.98628673, 14.18609082, 14.38589492,\n",
    "       14.58569901, 14.78550311, 14.98530721])\n",
    "print(data_freq.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # cut monopole freqs (for stonely)\n",
    "# idx1_start = np.where(data_freq - 0.5 > 0)[0][0]\n",
    "# idx1_end = np.where(data_freq - 6 > 0)[0][0]\n",
    "# # cut dipole freqs for 1st dipole\n",
    "# idx2_start = np.where(data_freq - 2 > 0)[0][0]\n",
    "# idx2_end = np.where(data_freq - 7 > 0)[0][0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(idx1_start, idx1_end, idx2_start, idx2_end)\n",
    "# print(data_freq[idx1_start], data_freq[idx1_end], data_freq[idx2_start], data_freq[idx2_end])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fe253181af0>]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsYAAAHSCAYAAADvxw2lAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAwSElEQVR4nO3df5Dc9X3n+edbEgg7u3FkGOkkgUFea7XBiWSLKUpZZ5XY2DHaGMPGsEWObOE7XGSu2Kyxd+3A2fHumeOOyjo2Wzm8F8rsoqvFdiEvKkR2oxgLciqqojgj4mwCRCchBVtIQROMcXZZC/143x/dLVqjnunvzHy7+/vtfj6qqJ75fr/z1UeNZvo1n35/3p/ITCRJkqRRt2jQA5AkSZKqwGAsSZIkYTCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAGwZNADALjooovysssuG/QwJEmSNOT27t3715k51ulcJYLxZZddxuTk5KCHIUmSpCEXES/MdM5SCkmSJAmDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJQIFgHBHrIuI7bf/9MCJuj4i3RsTjEbG/+bis7WvujIgDEbEvIj7Y27+CJEmStHBdg3Fm7svMd2Xmu4ArgNeA7cAdwK7MXAvsan5ORFwO3Ai8E7ga+HJELO7N8CVJkqRyzLWU4irg+cx8AbgW2No8vhW4rvnxtcDXM/N4Zh4CDgBXljDWWtr7wivc9+QB9r7wyqCHIkmSpFksmeP1NwJfa368IjOPAmTm0YhY3jy+GtjT9jWHm8dGzt4XXuGmr+zh9ZOnOX/JIh762CauuHRZ9y+UJElS3xWeMY6I84EPA9u6XdrhWHa4360RMRkRk1NTU0WHUSt7Dr7M6ydPczrhxMnT7Dn48qCHJEmSpBnMpZRiC/B0Zr7U/PyliFgJ0Hw81jx+GLik7esuBo5Mv1lm3p+Z45k5PjY2NveR98FCyyA2vf1Czl+yiMUB5y1ZxKa3X1jyCCVJklSWuZRS/DJvlFEA7ABuBu5pPj7advyrEfFFYBWwFvj2wofaX2WUQVxx6TIe+tgm9hx8mU1vv9AyCkmSpAorFIwj4s3AB4BfbTt8D/BwRNwCfBe4ASAzn4mIh4FngZPAbZl5qtRR90GnMoj5BNsrLl1mIJYkSaqBQsE4M18DLpx27GUaXSo6XX83cPeCRzdArTKIEydPz1oGsfeFV5wRliRJGgJz7UoxMoqUQdh1QpIkaXgYjGfRrQyirHKLbpyVliRJ6j2D8QIULbdYiDrNShvgJUlSnRmMF6AfXSfKmpXudWgtK8AbriVJ0qAYjBeo110nypiV7kdoLSPA1ylcG+AlSRo+BuMKmC1klTErXSS0dgt63UJrGQG+n+F6IcG2TgG+LGWMtS5/37qMU5JUPoPxgBUJWd1mpbu9kHcLrUXG0C20lhHgy2iRV/SXgIUE27rNji806JUx1rrUytdlnJKk3jAYD9hCQ1bRYD1baC0yhiKhtUhZyUJnx8uYuV7oDHqdZsfLCHpljLVfHVwWqi7jlCT1hsF4wBYasoq+kM8WWouG3oXOCJcxO17GzPVCZ9DrNDteRtArMtaFvmtRloXOjvdrnJKkajIYD9hCQ1YZL+RFx7DQhYb9CmndxlnGDHpdZsfLCLXdxlrGuxZlKGN2vKxxWqcsSfVkMK6AhQTOsl7Ie91dA/ob4ovcZyEz6N1UZXa8jFDbbaxlvGvRGku3XxRmO19WGcRCvxeKlrgYnCWpegzGQ6AfobYMdQnxZYyzKrPj3a7p1zi76RYmi4TNqpRrdHtO+9E5RZI0PwZj9VWdQvxCxlml2fE6jLNbmCxa3lKkPKXXdfLdntN+dE6RJM2PwVjqgVGaHW/dp5e/SBQN8LONo18dOro9p2V1TpEklc9gLPXIqMyOlzWG2cJkncpbWuOdbRHrQjunSJJ6IzJz0GNgfHw8JycnBz0MSUOsNWPcCptV37HQGmNJ6o2I2JuZ4x3PGYwljQrDpiRptmBsKYWkkVGFshFJUnUtGvQAJEmSpCowGEuSJEmMeDDe+8Ir3PfkAfa+8MqghyJJkqQBG9kaYxvoS5Ikqd3Izhh36mkqSZKk0TWyM8Y20Jekzq7fcT37Xtl3zvF1y9bxjQ9/YwAjkqT+GNlgXNZWuJI0bDYs38DBVw9y4vSJM8fOW3Qe71r+rsENSpL6YGSDMdjTVJI6mVg/waMHHj3r2KJYxMSGiQGNSJL6Y6SDsSQNo4WWQoy9eYxr33Et2/dv58TpE5y36Dyue8d1XPSmi3oxXEmqjJFdfCdJw2rD8g2ct+i8s47NtRRiYv0Ei6LxEuFssaRR4YyxJA2ZbqUQRWaUW7PG2/Ztc7ZY0sgwGEtShZTREaJbKUTRxXUT6yd4/gfPO1ssaWQYjCWpT4qE3m6htWhwbp81nl4KUXRx3dibx3jw6gfn9peUpBqzxliS+qRI7W97bW9Le2gtWj/cmjUO4pxSiNa51n1cXCdJDc4YS1JJus3mFpmp7VYGMZdWarOVQsw2oyxJo8oZY0kq4Pod1/PTW3/6nP+u33H9mWu6zeYWnamdrSPEXGZ7W6UQM52baUZZkkaVwViSCiijDGL6NTPN1HYLrWW1UptYP8HGFRudLZakJoOxJBVQJPQWmc0tOlM7W2gta7Z3thllSRpFBmNJKqCMMoj2a7rN1HYLrc72SlL5XHwnSRRrg1ZkwVqRjTHKaINmKzVJKp/BWFLtdQu1ZfQPhuK7wbkxhiTVk6UUkmqv28K4shbOta5baBmEJKmaDMaSaq9bqC1r4VzrOkOvJA0ng7Gk2usWastcOCdJGl7WGEuqtCL1wdB9YVxZC+ckScPLGWNJlVakPhi69/Yto3+wJGm4OWMsaaC6zQi3z/S2zNYfeLZuEEW6RdgGTZJGlzPGkgaq24xw0frg1rWzLYxz4Zyk6fa+8Ar3PXmAvS+8MuihqAIMxpIGqkjHCBfFSeqFvS+8wk1f2cNvfXMfN31lj+FYllJI6q1upRKtGeHt+7dz4vSJjjPCLoqT1At7Dr7M6ydPczrhxMnT7Dn4MldcumzQw9IAFZoxjoifiIhvRMRfRMRzEfEzEfHWiHg8IvY3H5e1XX9nRByIiH0R8cHeDV9S1c11c43Z6oddFCdprmYrldj09gs5f8kiFgect2QRm95+4QBGqCopWkrxb4Cdmfn3gA3Ac8AdwK7MXAvsan5ORFwO3Ai8E7ga+HJELC574JLqYS6ba8zWMcL6YElz1a1U4opLl/HQxzbxyV9Yx0Mf2+RssboH44j4cWAz8ABAZr6emT8ArgW2Ni/bClzX/Pha4OuZeTwzDwEHgCvLHbakupjL5hrOCEuaq9lmhDuVSkx3xaXLuO297zAUCyhWY/x2YAr49xGxAdgLfBxYkZlHATLzaEQsb16/GtjT9vWHm8ckDZmyNt8A26RJOtfeF15hz8GX2fT2CzsG19aM8OsnT3P+kkXnzPq2SiVOnDxtqYQKKVJKsQTYCPzbzHw38N9olk3MIDocy3Muirg1IiYjYnJqaqrQYCVVS1mbb0gaTbPN9hbpGNFtRthSCc1VkRnjw8DhzPyj5uffoBGMX4qIlc3Z4pXAsbbrL2n7+ouBI9Nvmpn3A/cDjI+PnxOcJQ1ePzffkDRaus32FukYUWRG+IpLlxmIVVjXGePM/CvgexGxrnnoKuBZYAdwc/PYzUDr1XEHcGNELI2INcBa4NuljlpSX/Rz8w1Jw6XbxhndZnuLdIxwRlhlK9rH+NeAhyLifOAg8D/RCNUPR8QtwHeBGwAy85mIeJhGeD4J3JaZp0ofuaSeKzIjXKR+WNLwma3+t9tsMHSf7W2F3tlqjFvXGYhVlkLBODO/A4x3OHXVDNffDdw9/2FJqgI335BG00IXvRUpgygSfA296jd3vpNGVJkdJawfloZHkdnebsG3aDcIg6+qpugGH5KGTJkdJawfluplob1/u9X/WvurunLGWBpRdpSQRlMZvX8tg9CwMhhLI6pI/XD7tW6+IdVDt/rgbmUQLnrTKDMYSyPMjhLScCmjGwQYejW6rDGWRpg70kn1s9D6YOt/pZk5YywNqbl0nbB+WKqHMuqDwRlhaSYGY2lIbVi+gYOvHuTE6RNnjs3UdcL6YakeyqoPltSZwVgaUnPpOiGpHqwPlnrLYCwNqbl0nZBUDd06SjgjLPWWwVgaYnadkOqjSEcJcEZY6iWDsVRDRRfWtWaNt+3b5myxVHHd6ocl9Z7t2qQaKrqdMzRmjTeu2OhssVQBs7Va67bNsqTei8wc9BgYHx/PycnJQQ9Dqo2p16bY8sgWjp86fubY0sVL2fmRnc4KSxVVpFSiW42xpIWLiL2ZOd7pnDPGUg21SiRas8YurJOqr+jmG7e99x2GYmlADMZSTU2sn2BRNL6FXVgnVYOlElK9ufhOqikX1knV0q1UwlZrUvUZjKUKcjtnqXq61f8W6SphqzX1QtHXDHVnMJYqyO2cpWopsnCuyK50Ui8Ufc1QdwZjqYLczlnqv9lmhIvOBlsqoV7oNiPsa0Z5XHwnVZBdJ6T+as0I/9Y393HTV/acs3iu6MI5u0qoF7r1rvc1ozwGY6mi7Doh9U+3Vmqt2eBP/sK6Gbdqlnql/fWgZfrrgq8Z5bCUQqoou05I/VOkPtiFc5qPbmUQRRbOtV4Ptu/fzonTJzrOCPuaUQ5njKUKcztnqT+cEVavdCuD6Ha+pciMsK8ZC+eW0JKkoedWy+qFIrO9U69NseWRLRw/dfzM+aWLl7LzIzu56E0XdT3f7q49d7Ft3zb+8bp/zGc3fbZHf6vhN9uW0JZSSH1mv0mpv4q0WpPmo0ibtG5lEEXKJFrsXd97llJIfVb0bTNJ5ei2sE6aryKL4qZfN5/zLa3e9dYP947BWOqzoj9IJZWjaKs1abrrd1zPT2/96XP+u37H9UDxNmmt64KY13n1j8FY6jP7TUrl2/vCK9z35IFz+g+DC+s0f0Xe4Ss629ttYZwL56rBxXfSALQvtphpkYWkYqwhVq8UXRjnorh6cfGdVDH2m5SK69ZRosh2zdJ0ZfUPBhfFDRODsTQg/iCVuisyG1xkcw5puiIdJaDxs/rRA48CM5dKtBbFqf6sMZYGxNXFUndFOkpYQ6z5KLoQ2oVxo8UZY0lSZRWdDXa7Zk3XrVTC/sHqxGAslcwNPKTytGaD3bVOc1WkVKJImQRYKjFKLKWQSuYGHlK5rrh0Gbe99x2GYs1JkVIJyyQ0nTPGUsnaZyBa3MBDmlm3rhNSJ2WVSlgmoXbOGEslcwMPqbhW14nf+uY+bvrKno4bdEidlLX5hguh1c5gLPVA0Z2QpFFXpOuE1ImlEuoFSymkHnADD6mhW5mEPYjViZtvaFAMxlKP+MNYo67I5hx2nVAnbr6hQTEYS3Mwl1Zs/jDWqCu6VbM9iDVd0UXMvjunshmMpTkoOoshyTIJzczNN1RVLr6T5qDoFqKS3KpZMyurowTYVULlMhhLc2ArNulse194hfuePDBjmzU351AndpRQVRmMpTmyFZvUYA9izVfRSYaJ9RNsXLHRn7PqG2uMpTlysYfUUHRxnUZL0UXKdpRQFRWaMY6Iv4yIP4uI70TEZPPYWyPi8YjY33xc1nb9nRFxICL2RcQHezV4aVCcxZDeWFy3OHBxnc4oUj8MlkqomiIzu18U8ZfAeGb+ddux3wS+n5n3RMQdwLLM/PWIuBz4GnAlsAr4FvB3M/PUTPcfHx/PycnJhf1NJEl9120DD42eqdem2PLIFo6fOn7m2NLFS9n5kZ3nhN+p16b41O5P8YWf+4LBWH0TEXszc7zTuYXUGF8LbG1+vBW4ru341zPzeGYeAg7QCMmSpBrptrAOXFync81lkbIdJVQ1RWuME/hmRCTwO5l5P7AiM48CZObRiFjevHY1sKftaw83j0mSaqLIrnUaTUVqiIvUD0tVVHTG+D2ZuRHYAtwWEZtnuTY6HDunXiMibo2IyYiYnJqaKjgMSVI/dFpYJ0GxGmLrh1VXhYJxZh5pPh4DttMojXgpIlYCNB+PNS8/DFzS9uUXA0c63PP+zBzPzPGxsbH5/w0kSaVzYZ1mUnSjIxcpq466llJExI8BizLzb5of/wLweWAHcDNwT/Oxtan5DuCrEfFFGovv1gLf7sHYpdIVbTMkDbvWrnUurNN0RbdrttWa6qhIjfEKYHtEtK7/ambujIg/Bh6OiFuA7wI3AGTmMxHxMPAscBK4bbaOFFKVbFi+gYOvHuTE6RNnjnVqMySNgisuXWYgVkfWEGtYdQ3GmXkQ2NDh+MvAVTN8zd3A3QsendRn7T/sW/yhr2FlqzVNV/RdMzc60rByS2ipzVzaDEl15nbO6qTo5hxgDbGGk8FYmqZ9YYmzxRpWdp1QJ0UX1oE9iDWcDMbSNLYZ0iiw64Q68V0zjbqiG3xII2Vi/QTP/+B5Z4s1tOw6MZrcnEOanTPGUge+RahR4HbOo8fNOaTZGYwlaUjtfeEV7nvygAvrdIabc0izs5RCkoZQq+vE6ydPc/6SRTz0sU3ODMvNOaQuDMaSNIQ6dZ0wGA+3oj2IrSGWZmYphSQNIbtOjJ6iPYitIZZm5oyxJA0hu06Mnrns3GnnHakzg7EkDakrLl1mIB4hReuHW9daQyydy1IKSaohO06oE3fulBbGGWNJqhk7ToymIovrWrPG2/Zts35YmgeDsUZK0VXbUpXZcWI0bVi+gYOvHuTE6RNnjnVaXGf9sDR/llJopBRdtS1VmR0nRlPRzTncuVOaP4OxRkrRFxapylodJz75C+ssoxghrTKJ1i/3sy2ukzQ/BmONFF9YNCyuuHQZt733HYbiEePiOqm3rDHWyHHXJ1XaU/fC6o2wZvMbxw7thhefhp+9fVCjUo8VXf/g4jqpt5wx1shx1ydV2uqNsO2jjTAMjcdtH20c19Cay/qHifUTbFyx0V/qpR6IzBz0GBgfH8/JyclBD0MjZOq1KT61+1N84ee+YDBW9bTC8PgtMPkA3PDg2TPIGjpTr02x5ZEtHD91/MyxpYuXsvMjO/0ZJZUsIvZm5ninc84YayS5aluVtmZzIxTv/s3Go6F46Ln+QaoGg7EkVc2h3Y2Z4s2fbjy2yio01FxYJw2ewViSqqRVRnHDg/C+zzQe22uONbRc/yANnl0pJKlKXnz67JriNZsbn7/4tCUVNVa064S71kmDZTCWpCrp1JJtzWZDcc0V3c65tf5B0mBYSiFJ/fLUveeWRBza3Tiuoeaum1I9GIwlqV/sUTyy7Doh1YPBWJL6pVUvvO2j8MTdbyyys0xiJNh1Qqo+g7Ek9ZM9ikeWXSek6nPxnST10/QexWv+geF4CNh1QhoOzhhLUr/Yo3hobVi+4Uz9cMtsXSecLZaqyWAsSf0yW49i1ZpdJ6ThYCmFhkbRtzKlgbFH8dBq1Q9v37+dE6dP2HVCqilnjDU0ir6VKUm9YNcJqf4MxhoavpWpgXMDj5Fm1wmp/gzGGho20NfAuYHHyJtYP8HGFRv9hVyqqcjMQY+B8fHxnJycHPQwNASmXptiyyNbOH7qOEsXL2XnR3YajNVfrTA8fkujHZsbeEhSpUTE3swc73TOxXcaKq1Z4237tjlbrMFo38Bj86cNxUPAhb3S6LCUQkPHtzI1UNM38LBHce25sFcaHQZjDR0b6Gtg3MBjKLmwVxodBmNJKosbeAwlF/ZKo8PFd5IkdeHCXml4zLb4zhljSZK6sEexNBrsSiFJUgET6yd4/gfPW1ssDTGDsSRppBVtx9Za2CtpeFlKIUkaabZjk9RiMJakIp6699y2a4d2N46r1mzHJqnFYCxJRazeeHZP4lbP4tUbBzkqlcB2bJJaCgfjiFgcEX8SEb/b/PytEfF4ROxvPi5ru/bOiDgQEfsi4oO9GLgk9VWrJ/G2j8ITd7+xkYdbPg+F9lljZ4ul0TWXGeOPA8+1fX4HsCsz1wK7mp8TEZcDNwLvBK4GvhwRi8sZriQN0JrNMH4L7P7NxqOheGjYjk0SFAzGEXEx8IvAV9oOXwtsbX68Fbiu7fjXM/N4Zh4CDgBXljJaSRqkQ7th8gHY/OnGo1s9D5WJ9RNsXLHR2WJphBVt13Yv8Gngb7cdW5GZRwEy82hELG8eXw3sabvucPPYWSLiVuBWgLe97W1zG7Uk9VurprhVPrHmH1hOUQNFW7GB7dgkFZgxjogPAccyc2/Be0aHY+fsO52Z92fmeGaOj42NFby1JA3Ii0+fHYJbNccvPj3IUakLW7FJmosiM8bvAT4cEf8QuAD48Yj4D8BLEbGyOVu8EjjWvP4wcEnb118MHClz0JLUdz97+7nH1mx2trjiJtZP8OiBR8865uI6STPpGowz807gToCI+HngX2Tmr0TEvwZuBu5pPrZ+8uwAvhoRXwRWAWuBb5c+co2cubwlKknwxqK67fu3c+L0CVuxSZrVQvoY3wN8ICL2Ax9ofk5mPgM8DDwL7ARuy8xTCx2o5FuikubDVmySippTMM7MP8jMDzU/fjkzr8rMtc3H77ddd3dm/p3MXJeZv1f2oDWa3J1KPeXOdkPLVmySinLnO9WGu1Opp9zZbqjZik1SEZF5TsOIvhsfH8/JyclBD0M1MPXaFFse2cLxU8dZungpOz+y02Cs8rTC8PgtjT7FtmKTpKETEXszc7zTOWeMVSu+Jaqecmc7SRppBmPVjm+Jqmfc2U6SRlrRne+kynB3KvWEO9vVkm0cJZXJGWNJAne2qynbOI4wO8moBwzGkgSNne2mzwyv2dx5xztVhm0cR5idZNQDBmNJUm3ZxnGEtd7V2fZReOJuS59UCoOxJKnW3NluhNlJRiUzGEuSas02jiPMTjIqmV0pJEm1N7F+gud/8LyzxaPETjLqAXe+kzT8nrq3sSCn/cXy0O5GxwkX11WWrdg0K7+vNU/ufCdptLl6vZZsxTbCirRis5OMesBgLGn4uXq9lmzFNsL8ZVYDYjCWNBpcvV47tmIbYf4yqwExGEsaDa5eryVbsY0wf5nVABiMJQ2/9tXr7/vMGzNRhuPKsxXbCPOXWQ2AwVjS8Hvx6bPfhm29Tfvi04MclQqaWD/BxhUbnS0eJf4yqwGxXZskSaoWW7Gph2Zr1+YGH5KkgbBP8QjrFnw7hd81m60zVs8ZjFUJvkBKo2fD8g0cfPUgJ06fOHPMPsUjotWOrVXi1F46IQ2QNcaqBBv5S6PHPsUjzHZsqiiDsSrBF0hp9NineMTZjk0VZDBWJfgCqQUpsn2sKsk+xSPMdmyqIIOxKsMXSM2b28fWln2KR5Tt2FRRBmNVhi+QmjfrFWvNPsVDpsg7OPYWV0XZx1iVMvXaFJ/a/Sm+8HNfMBhr7p64u1GvuPnTjVkoDYydZkZY+2zw9I4T/rJaWaP0PTtbH2NnjFUpY28e48GrHzQUa+6sV6wUO82MMN/BqSW/ZxsMxpLqz3rFyrHTzIiz40Tt+D3bYDCWVH/WK1aOnWZGnO/g1I7fsw3WGEuSemLqtSm2PLKF46eOs3TxUnZ+ZOfIvcgOpW7bOVtjXElFaohH5XvWGmNJUt/ZaWZIdWuP6Ds4A/HqY4+x/31X8dxPXs7+913Fq489dtb5IjXEfs86YyxJ6iE7zQypVhgev6VRKuFscM+9+thjHPvSvZw8epQlK1ey/BO385Zrrjlz7uhvfI780Y/OXB8XXMDKuz5/5pr22eCWTrPCo/A9O9uMscFYkjRno9TaaaR0K5NoZ3vEvukWfPe/7ypOHjlyztctWbWKtU/sOvP5XXvuYvv+7Zw4fYLzFp3HL639JT676bN9+TtUiaUUkqRS2dppSBXdRdLFdX117Ev3nhWKAfJHP+LYl+4F4OTRox2/bvpxd5jtzmAsSZozWzsNqSI9iG2P2Hfdgu+SlSs7np9+3Bri7pYMegAaDb7tKg2X1gts+9uyvtAOifYexJs/fW7t8GyL66wz7oklK1d2LpVoBt/ln7i9Y6nF8k/cfs7XTKyf4PkfPO8vsTNwxlh94duumren7j13JurQ7sZxDZRvyw6pbmUSP3v7uQF4zeZza5BVmuWfuJ244IKzjrUH37dccw0r7/o8S1atggiWrFp11sK7du4wOzuDsfrCt101b0VrHtV3vi07hCyTqKQiwfct11zD2id28ZPPPcvaJ3Z1DMXqzlIK9YVvu2re2msebQ3VN7O1hmrn27JDxjKJynrLNdcYdvvAdm3qm1HZUUc9YmuovinSE1WS6sp2baoE33bVvNkaqq+6tYaSpGFlMFZfTayfYOOKjb7tquKseZyXbtvDzna+aE9USRo2BmP1lathNWez1Tyqo1YpxMkjRyCTk0eOcPQ3Pncm/HY7X7QnqiQNG4OxpGqzNdScdSuF6Ha+W2soSRpWdqWQpCHTrRSi2/nWArsiXSkkaZgYjCX11lP3NnoOT99S9sWnG7O+3c5rzrrtktXtPNgaStJo6lpKEREXRMS3I+JPI+KZiPjfmsffGhGPR8T+5uOytq+5MyIORMS+iPhgL/8Ckiqu2wYdbuAxZ90W1nUrhbBUQpI669rHOCIC+LHM/K8RcR7wFPBx4JeA72fmPRFxB7AsM389Ii4HvgZcCawCvgX83cw8NdOfYR9jaci1wu5MG3R0Oz9Eim6cMdvXF+kx3O3PWeg4JKmuZutjPKcNPiLizTSC8f8C/D/Az2fm0YhYCfxBZq6LiDsBMvP/bH7N7wP/KjP/cKb7GoylEdBtg44R2MCjjI0z9r/vqs5lEKtWsfaJXaWNVZKG1YI3+IiIxRHxHeAY8Hhm/hGwIjOPAjQflzcvXw18r+3LDzePSRpV3TboGJENPIpunGGPYUkajELBODNPZea7gIuBKyPip2a5PDrd4pyLIm6NiMmImJyamio0WEk11G2DjhHawKNIqLXHsCQNzpz6GGfmD4A/AK4GXmqWUNB8PNa87DBwSduXXQyc875fZt6fmeOZOT42Njb3kUuqh24bdIzQBh5FQq09hiVpcIp0pRiLiJ9ofvwm4P3AXwA7gJubl90MPNr8eAdwY0QsjYg1wFrg2yWPW1JddNugY4Q28CgSaov0GF551+dZsmoVRLBk1ao51ShLkmZWpI/xSmBrRCymEaQfzszfjYg/BB6OiFuA7wI3AGTmMxHxMPAscBK4bbaOFJI0KopsnGGPYUkanDl1pegVu1JIUkMZnSskSTObrSuFO99JUoW4HbMkDY7BWJIqxlIJSRqMOXWlkCRJkoaVwViSJEnCYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJQIFgHBGXRMSTEfFcRDwTER9vHn9rRDweEfubj8vavubOiDgQEfsi4oO9/AtIkiRJZSgyY3wS+OeZ+ZPAJuC2iLgcuAPYlZlrgV3Nz2meuxF4J3A18OWIWNyLwUuSJEll6RqMM/NoZj7d/PhvgOeA1cC1wNbmZVuB65ofXwt8PTOPZ+Yh4ABwZcnjliRJkko1pxrjiLgMeDfwR8CKzDwKjfAMLG9ethr4XtuXHW4ekyRJkiqrcDCOiL8F/Efg9sz84WyXdjiWHe53a0RMRsTk1NRU0WFIkiRJPVEoGEfEeTRC8UOZ+Ujz8EsRsbJ5fiVwrHn8MHBJ25dfDByZfs/MvD8zxzNzfGxsbL7jlyRJkkpRpCtFAA8Az2XmF9tO7QBubn58M/Bo2/EbI2JpRKwB1gLfLm/IkiRJUvmWFLjmPcA/Af4sIr7TPPa/AvcAD0fELcB3gRsAMvOZiHgYeJZGR4vbMvNU2QOXJEmSytQ1GGfmU3SuGwa4aoavuRu4ewHjkiRJkvrKne8kSZIkDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJElAgWAcEf8uIo5FxJ+3HXtrRDweEfubj8vazt0ZEQciYl9EfLBXA5ckSZLKVGTG+EHg6mnH7gB2ZeZaYFfzcyLicuBG4J3Nr/lyRCwubbSSJElSj3QNxpm5G/j+tMPXAlubH28Frms7/vXMPJ6Zh4ADwJXlDFWSJEnqnfnWGK/IzKMAzcflzeOrge+1XXe4eUySJEmqtLIX30WHY9nxwohbI2IyIianpqZKHoYkSZI0N/MNxi9FxEqA5uOx5vHDwCVt110MHOl0g8y8PzPHM3N8bGxsnsOQJEmSyjHfYLwDuLn58c3Ao23Hb4yIpRGxBlgLfHthQ5QkSZJ6b0m3CyLia8DPAxdFxGHgXwL3AA9HxC3Ad4EbADLzmYh4GHgWOAnclpmnejR2SZIkqTRdg3Fm/vIMp66a4fq7gbsXMihJkiSp39z5TpIkScJgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSgB4G44i4OiL2RcSBiLijV3+OJEmSVIYlvbhpRCwG7gM+ABwG/jgidmTms7348+bq+q/9HBf9yRT/4x8kF/4QXv5x+OrPB3/97jG+8cv/b9fzVblHXcYJwFP38uq+Exx7aCcnjx5lycqVLL/pat6y7jz42du7n6/KPeoyzirdA3j1scc49qV737jmE7fzlmuuKe18ne6xYE/dC6s3wprNbxw7tBtefPqs/2+zXrPQ880/oxL/voblHlUZJ9X5XlnwPerynBe8RxW+Z/uRWwatVzPGVwIHMvNgZr4OfB24tkd/1px96OCFTPznZOyHjSdg7Icw8Z+TDx28sND5qtyjLuMEeHXfCY5+8QFOHjkCmZw8coSjX3yAV/edKHS+KveoyzgrdY/HHuPob3zu7Gt+43O8+thjpZyv0z1KsXojbPto40UPGo/bPto4XvSahZ6nQv++huQelRlnRb5XSrlHXZ7zAvdofU/u2/OfuO/JA+zb8586fs/OeL6kcfQjtwxar4LxauB7bZ8fbh6rhJ/51vdZevLsY0tPNo4XOV+Ve9RlnADHHtpJnoqzrslTwbGHdhY6X5V71GWclbrHl+4lf/Sjs6/50Y849qV7Szlfp3uUYs1muOHBxoveE3c3Hm948OyZom7XLPQ8Ffr3NST3qMw4K/K9Uso96vKcF7gHazazb/NvM/Z7v8qpXf87Y7/3q+zb/Ntnfc/Oer6kcfQjtwxaT0opgOhwLM+6IOJW4FaAt73tbT0aRmen/+rYrMe7na/KPeoyToCTR492vKZ1vNv5qtyjLuMcpnvUZZxFrynFms0wfgvs/k3Y/OmzQ3HRaxZ4vi7PeV3uUZdxDtM96jLOlm/993WcOvV+/tmS7fz2yX/Eov++jnVzOF/GOPqRWwatVzPGh4FL2j6/GDjSfkFm3p+Z45k5PjY21qNhdLZk5cpZj3c7X5V71GWcw3SPuoxzmO5Rl3EWvaYUh3bD5AONwDr5wBslD3O5ZoHn6/Kc1+UedRnnMN2jLuNsef+b9vEri7/Fb5/8R9y0+Fu8/0375nS+Cn+Xvv2MXIBeBeM/BtZGxJqIOB+4EdjRoz9rzpbfdDWx+KwJbGJxsvymqwudr8o96jLOYbpHXcZZqXt84nbiggvOvuaCC1j+idtLOV+ne5SiVe97w4Pwvs+8UfLQHly7XbPQ81To39eQ3KMy46zI90op96jLc17gHhzazbrdv8bUlt9h0VWfZWrL77Bu96+d9T076/mK/F0K/V0HrCfBODNPAv8U+H3gOeDhzHymF3/WfLxl3Xms/OQtLFq5ggQWrVzByk/e0lh5WeB8Ve5Rl3G2X7Nk1SqIYMmqVR3vMdP5qtyjLuOs1D2uuYaVd33+7Gvu+vyZleMLPV+ne5Tixac71wO/+HTxaxZ6ngr9+xqSe1RmnBX5XinlHnV5zgvco/U9uW7TL3Lbe9/Buk2/2PF7dsbzJf9deplbBi0ys/tVPTY+Pp6Tk5ODHoYkSZKGXETszczxTufc+U6SJEnCYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCIzBz0GIiIKeCFaYcvAv56AMMZZj6n5fM5LZfPZ/l8Tsvnc1oun8/y+ZzO7tLMHOt0ohLBuJOImMzM8UGPY5j4nJbP57RcPp/l8zktn89puXw+y+dzOn+WUkiSJEkYjCVJkiSg2sH4/kEPYAj5nJbP57RcPp/l8zktn89puXw+y+dzOk+VrTGWJEmS+qnKM8aSJElS31QyGEfE1RGxLyIORMQdgx5P3UXEJRHxZEQ8FxHPRMTHBz2mYRARiyPiTyLidwc9lmEQET8REd+IiL9o/lv9mUGPqc4i4hPN7/c/j4ivRcQFgx5T3UTEv4uIYxHx523H3hoRj0fE/ubjskGOsW5meE7/dfP7/r9ExPaI+IkBDrF2Oj2nbef+RURkRFw0iLHVUeWCcUQsBu4DtgCXA78cEZcPdlS1dxL455n5k8Am4Daf01J8HHhu0IMYIv8G2JmZfw/YgM/tvEXEauCfAeOZ+VPAYuDGwY6qlh4Erp527A5gV2auBXY1P1dxD3Luc/o48FOZuR74/4A7+z2omnuQc59TIuIS4APAd/s9oDqrXDAGrgQOZObBzHwd+Dpw7YDHVGuZeTQzn25+/Dc0AsfqwY6q3iLiYuAXga8MeizDICJ+HNgMPACQma9n5g8GOqj6WwK8KSKWAG8Gjgx4PLWTmbuB7087fC2wtfnxVuC6fo6p7jo9p5n5zcw82fx0D3Bx3wdWYzP8OwX4EvBpwMVkc1DFYLwa+F7b54cxxJUmIi4D3g380YCHUnf30viBc3rA4xgWbwemgH/fLE/5SkT82KAHVVeZ+SLwBRozRUeBVzPzm4Md1dBYkZlHoTHpACwf8HiGzf8M/N6gB1F3EfFh4MXM/NNBj6VuqhiMo8Mxf9spQUT8LeA/Ardn5g8HPZ66iogPAccyc++gxzJElgAbgX+bme8G/hu+RT1vzbrXa4E1wCrgxyLiVwY7Kml2EfEZGqV/Dw16LHUWEW8GPgN8btBjqaMqBuPDwCVtn1+MbwEuWEScRyMUP5SZjwx6PDX3HuDDEfGXNEp93hcR/2GwQ6q9w8DhzGy9k/ENGkFZ8/N+4FBmTmXmCeAR4O8PeEzD4qWIWAnQfDw24PEMhYi4GfgQcFPaR3ah/g6NX4r/tPk6dTHwdET8DwMdVU1UMRj/MbA2ItZExPk0FozsGPCYai0igkbt5nOZ+cVBj6fuMvPOzLw4My+j8e/zicx0Nm4BMvOvgO9FxLrmoauAZwc4pLr7LrApIt7c/P6/ChczlmUHcHPz45uBRwc4lqEQEVcDvw58ODNfG/R46i4z/ywzl2fmZc3XqcPAxubPWXVRuWDcLMD/p8Dv0/hB/nBmPjPYUdXee4B/QmNm8zvN//7hoAclTfNrwEMR8V+AdwH/x2CHU1/NmfdvAE8Df0bjZ707Yc1RRHwN+ENgXUQcjohbgHuAD0TEfhor/u8Z5BjrZobn9P8C/jbwePP16f8e6CBrZobnVPPkzneSJEkSFZwxliRJkgbBYCxJkiRhMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSQD8/1LDdtaMcMJ6AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots()\n",
    "fig.set_size_inches(12,8)\n",
    "n = 0\n",
    "ax.plot(data_freq, data_X[n, :, 1, 1],'.')\n",
    "ax.plot(data_freq, data_X[n, :, 0, 1],'x')\n",
    "ax.plot(data_freq, data_X[n, :, 1, 0],'v')\n",
    "ax.plot(data_freq, data_X[n, :, 0, 0],'o')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "detect curve parts for the whole dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# def detect_non_zeros(k,j):\n",
    "#     ones = np.ones(len(data_freq))\n",
    "#     num_of_nonzero = np.zeros(len(data_freq))\n",
    "#     for i in range(data_X.shape[0]):\n",
    "#         idx_non0 = np.where(data_X[i, :, k, j] > 0)\n",
    "#         num_of_nonzero[idx_non0] += ones[idx_non0]\n",
    "#     idx = np.where(num_of_nonzero > 15000)\n",
    "#     fig, ax = plt.subplots(1,2)\n",
    "#     fig.set_size_inches(8,2)\n",
    "#     ax[0].plot(data_freq, num_of_nonzero, '.')\n",
    "#     ax[1].plot(data_freq[idx], num_of_nonzero[idx], 'x')\n",
    "#     return num_of_nonzero, idx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # first monopole\n",
    "# num_of_nonzero_1m, idx_1m = detect_non_zeros(1,1)\n",
    "# # second monopole\n",
    "# num_of_nonzero_2m, idx_2m = detect_non_zeros(0,1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # first dipole\n",
    "# num_of_nonzero_1d, idx_1d = detect_non_zeros(1,0)\n",
    "# # second dipole\n",
    "# num_of_nonzero_2d, idx_2d = detect_non_zeros(0,0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# calc mu\n",
    "mu_curve = np.zeros([len(data_freq), 2, 2])\n",
    "for i in range(dataset_size):\n",
    "    mu_curve[:, 0, 0] += data_X[i, :, 0, 0]\n",
    "    mu_curve[:, 0, 1] += data_X[i, :, 0, 1]\n",
    "    mu_curve[:, 1, 0] += data_X[i, :, 1, 0]\n",
    "    mu_curve[:, 1, 1] += data_X[i, :, 1, 1]\n",
    "mu_curve /= data_X.shape[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# calc std\n",
    "std_curve = np.zeros([len(data_freq), 2, 2])\n",
    "for i in range(dataset_size):\n",
    "    std_curve[:, 0, 0] += (data_X[i, :, 0, 0] - mu_curve[:, 0, 0])**2/dataset_size\n",
    "    std_curve[:, 0, 1] += (data_X[i, :, 0, 1] - mu_curve[:, 0, 1])**2/dataset_size\n",
    "    std_curve[:, 1, 0] += (data_X[i, :, 1, 0] - mu_curve[:, 1, 0])**2/dataset_size\n",
    "    std_curve[:, 1, 1] += (data_X[i, :, 1, 1] - mu_curve[:, 1, 1])**2/dataset_size\n",
    "std_curve = np.sqrt(std_curve)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 1.0, 'std curve')"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAv0AAAHiCAYAAAB7rnPHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAB2sElEQVR4nO3dfXxU9Zk3/s81mfCkbYxAMIBAEJqihdCQbWnVWKFVuWsElXbpsr3tLq5NV7cFd2u1Pty7Ura6uyrdrd39udqFvRerQqWStj4VbFm2ppawtlVpCgZRIGsohrTKLeTh+v0xc5IzZ848n+f5vHnxSubMyeSbyWTOdb7n+l6XqCqIiIiIiCi6Yn4PgIiIiIiI3MWgn4iIiIgo4hj0ExERERFFHIN+IiIiIqKIY9BPRERERBRxDPqJiIiIiCKOQT8RERFRkoj8tYj8h9/jIHIag34iIiIqCwzoqZwx6CdykIjE/R4DEREFG48V5AcG/RQaIvKaiHxZRH4pIu+IyEMiMklEnhSR34vIj0SkOrnvx0TkkM3XfzzDY48VkXtE5KCI9InIruS2rI+TnDXaIiL/ISK/A/BVEfl/InKmaf8PishvRaQyeftPRWSviPSKyNMiMt3hp4qIqKyJyFdE5HDy2NApIotF5DIAXwXwhyLytoj8IrlvnYj8JLnvswAm5HjspSLyooj8TkReTT5u2jHGfFVBRGaIiIrIKhF5HcAOEXlKRG6wPPYvROSq5OfvF5FnReSt5M/waSefIyo/DPopbK4G8AkA7wPQAuBJJN7EJyDxev5ikY/7DwAWAPgogDMB3ARgKM+vXQpgC4AzAPw9gOeT4zT8EYAtqtovIsuS470KwEQA/wngO0WOmYiILESkHsANAP5AVd8D4FIAr6nqUwD+FsCjqnq6qjYkv+RhAB1IHEfWArgmy2N/CMC/A/gyEu/5zQBeK2B4FwGYkxzTwwA+Y3rscwFMB/ADETkNwLPJfWqS+31LRM4r4HsRpWDQT2HzT6r6pqoeRiJg/pmq/reqngSwFcAHC31AEYkB+FMAX1LVw6o6qKo/TT5mPp5X1e+p6pCq/j+Y3shFRACsSG4DgM8D+Lqq7lXVASQOQPM5209E5JhBAKMBnCsilar6mqq+arejiEwD8AcAblfVk6q6E0BblsdeBeDbqvps8j3/sKr+uoCx/bWqvpM8VmxF6vv/SgCPJ489lyNxovJvqjqgqnsAfBfA8gK+F1EKBv0UNm+aPv9/NrdPL+IxJwAYA8D2oJCHNyy3twD4iIhMRmIWSJE4QQESszjfEJHjInIcwFsABMCUIr83ERGZqOp+AKsB/DWAHhF5JPl+bGcygF5Vfce07WCWhz8bxR8rANPxQlV/D+AHSEwMIflxU/Lz6QA+bBwrkseLlQDOKuF7U5lj0E9R9Q6AccYNEalAIp3Gzm8BvAvgnCIfR1NuqB4H8AyATyOR2vMdVTX2eQPA51X1DNP/sar603x/MCIiyk5VH1bVC5AInhXA3cZdll27AVQn02kM07I89BuwP1YAluMF7AN06/f/DoDPiMhHAIwF8Jzp+/zEcqw4XVW/kGVsRFkx6Keo+g2AMSLyyeQC2tuQuNybRlWHAHwbwL0iMllEKkTkIyIyupDHsXgYwP9GIrf/YdP2fwFwi5GXKSJVIvKpIn9GIiKyEJF6EVmUfA9/F4mrwIPJu98EMCOZ1glVPQhgN4C/EZFRInIBEuvFMnkIwJ8kFwbHRGSKiLw/ed+LAFaISKWINCG/VJwfInFicicSaw2MtWTfB/A+Efls8vEqReQPRGROvs8DkRWDfookVe0D8OcAHgRwGIkZmENZvuSvAPwKwM+RSLm5G0CsiMcxbAMwG8CbqvoL07i2Jh/7kWS1n5cALCnohyMiomxGA7gLiau4/4PEQtivJu/bnPx4TET2JD//IwAfRuK9//8gsVDXlqq+AOBPANwHoA/AT5AI2gHgdiSuAvQC+BukTvhkeryTAB4H8HHz/snUn0uQSPk5kvw57kZ+k05EtmQk64CIiIiIiKKIM/1ERERERBHHoJ+IiIiIKOIY9BMRERERRRyDfiIiIiKiiGPQT0REREQUcXG/BwAAEyZM0BkzZvg9DCKiwOro6PitqmZqMFc2eLwgIsou0/EiEEH/jBkzsHv3br+HQUQUWCJy0O8xBAGPF0RE2WU6XjC9h4iIiIgo4hj0ExERERFFHIN+IiIiIqKIY9BPRERERBRxDPqJiIiIiCKOQT8RERERUcQx6CciIiIiijgG/UREREREEcegn4iIiIgo4hj0ExERERFFHIN+IiIiIqKIY9BPRERERBRxDPqJiIiIiCKOQT8RERERUcQx6CciIiIiijgG/UREREREEZcz6BeRehF50fT/dyKyWkTOFJFnRWRf8mO16WtuEZH9ItIpIpe6+yMUpuNgL+5/bj86Dvb6PRQiIiIi8ki5x4DxXDuoaieA+QAgIhUADgPYCuBmANtV9S4RuTl5+ysici6AFQDOAzAZwI9E5H2qOujOj5C/joO9WPlgO04NDGFUPIZN1y7EgunVub+QiIiIiEKLMWDh6T2LAbyqqgcBLAWwMbl9I4Blyc+XAnhEVU+q6gEA+wF8yIGxlqy96xhODQxhSIH+gSG0dx3ze0hERERE5DLGgIUH/SsAfCf5+SRV7QaA5Mea5PYpAN4wfc2h5DbfLZw5HqPiMVQIUBmPYeHM8WmXesr90g8RUbFEZIyIvCAivxCRl0Xkb5LbQ5kOSkTRYRcDlpuc6T0GERkF4AoAt+Ta1Wab2jzedQCuA4Bp06blO4ySLJhejU3XLkR717HhX7b5Us8dl5+HO7//csqlHwDD+5fbZSAiogKdBLBIVd8WkUoAu0TkSQBXIWTpoEQULdYYsBxjuryDfgBLAOxR1TeTt98UkVpV7RaRWgA9ye2HAJxt+rqpAI5YH0xVHwDwAAA0NTWlnRS4ZcH06uFf9P3P7U+51PPkS90pt7+75xAe33Mo50lAx8Hesn4REREBgKoqgLeTNyuT/xWJtM+PJbdvBPBjAF+BKR0UwAERMdJBn/du1ERULswxYDkqJOj/DEZSewBgG4BrANyV/PiEafvDInIvEjM3swG8UPpQnWdc6ukfGEJlPIYlH6jFz197a/i2AHmdBFgXhgCpJwY8KSCicpEs+NABYBaA+1X1ZyKSkg4qIuZ00HbTl9umg/pxZZiIKGryCvpFZByATwD4vGnzXQAeE5FVAF4H8CkAUNWXReQxAK8AGABwvZ+XarMF3HaXeurPek9K+s939xzKeBJgLALJdmKQT8qQdYw8SSCisEq+388XkTMAbBWRD2TZPa90UL+uDBMRRUleQb+qngAw3rLtGBLVfOz2XwdgXcmjK1E+5Zmsl3qst61rAMwnAcY289UC64lBrpQh60mB3UlCMScGPHEgIj+p6nER+TGAy1BiOigRUTEYC6UqJL0ndOzKMxX6S892EmBsz3ZikCtlyHpSYL1tXE0oZMGxdf98U45KPZHgHxdReRORiQD6kwH/WAAfB3A3IpAOSkThwrr86SId9Ftz9p0oz2S3CCTXiUG2lCHrSYH19sKZ49NOXnJdPbi6cWpRKUfZThRKvb/c/9CIykQtgI3JvP4YgMdU9fsi8jxCkA5KRNHhxMRv1EQ66PerPFMhKUPWkwK72wAKWnCsKCzlKJ+1CblOJHLdnylNiYiiQ1V/CeCDNtsDnw5KRNHixsRv2EU66AeCWZ4p10mB3e1Crh5c3TgVVzdOzfvqQj5rE3KdSOS63y5NiZfaiIiIyA2sy58u8kF/VBR69cDYJ9P9dlcTsq1NyHUiket+uzQlXmojIiIitwRx4tdPkuil4q+mpibdvXu338MgC6cX9hqLaowTAc70E+VPRDpUtcnvcfiNxwsiouwyHS84008Z5TpDLvR+XmojIiIi8geDfvKU9USAC3uJiIiI3Megn3zDGrpERERE3oj5PQAqX3YLe4mIiIjIeQz6yTdGDd0KAWvoEhEREbmI6T3kGy7sJSIiIvIGg37yFWvoEhEREbmP6T1ERERERBHHoJ+IiIiIKOIY9BMRERERRRyDfiIiIiKiiGPQT0REREQUcQz6iYiIiIgijkE/EREREVHEMegnIiIiIoo4Bv0UKB0He3H/c/vRcbDX76EQERERRQY78lJgdBzsxcoH23FqYAij4jFsunYhu/USEREROYAz/RQY7V3HcGpgCEMK9A8Mob3rmN9DIiIiIooEBv0UGAtnjseoeAwVAlTGY1g4c7zfQyIiIiKKBKb3UGAsmF6NTdcuRHvXMSycOR4Lplej42Bv1ttERERElBuDfgqUBdOrh4N5a47/HZefhzu//zJz/omIiIgKxKCfAsua4//kS91pOf/Zgv7l25ajs7czbXt9dT22XLEl5/1EREREURGpnH6We4wWa47/kg/UFpTz31DTgMpYZcq2ylgl5tfMz+t+IiIioqgQVfV7DGhqatLdu3eX9Bgs9xhNRg7/D499Ga+/sz/t/myz8kdPHMWSx5fg5ODJ4W2jK0bjqaufwoSxE2zvt8OZfwoCEelQ1Sa/x+E3J44XRERRlul4EZn0Hrtyjwz6wytT6o1Zrln5ieMmYumspdi6byv6h/pRGavEslnLMGHsBNv7JflvCEN5fw8iCgcWASCicheZoN9IBekfGGK5xwhoqGlAV18X+of6M+4TkxhaG1qHb2c6URCI7f4A0DqvFU/sfwIAhlN9Tg2dyvg9iCh8eCWYiChCQb9duUcKL3MwboghhlgshoGhgeFZ+4M9FXi0az8Wzhxve6JQGavEtPdMQ1dfF8bGx+Lixy5O+15Vo6pwavAUrpx9JRSa8coAEYUTrwQTEUUo6AdSyz1SuORK56mMVWJJ3RI8/drTGMAAYhLDR85ckTJ7d/9nV+AJsZwoSAx3N9+Nr7/wdUw5fQqePPBk2knBRWdfhMNvH0ZrQytUdfhkg7P8RNHAK8FERBEL+ik88snZj0scCsWgDiImMaxZsAZj4mOwuXMzls1ahl8fRsrs3a8Pi20Of/2Z9dhw2QYcPXEUT7/2dMr3MB7XPJu/dNbS4e/BWX6i8OOVYCIiBv3kk3xy9itiFbh0xqVoe7VtOABvndeKV4+/io43O7D/+KM47f0j+//L68A5VecgJolKtNaZ+lwLew3G9+AsP1F08EowEZW7SNXpp/Bondc6HJwbYoghHkuchxoB+erG1Wic1DgcgE8cNxEbLtuAxkmNtjX2m85qwtJZSyGQjAF9ppMCg/E9OMtPREREUcGZfvJEMTn7rQ2tmDB2AjZctiFtf9uFvsmvUdWMM/XGbH8+6Tvs2EtERERRwZl+8oRd99u4xFEhFQBGcuuzzdKbGcG78ZjmVJ1cM/Wt81pTrh4UMmbW7SciIqIwYtBPnrBL56mIVeCTMz+ZEuTnG5BbH7OQSjv5pu/YpiCxog8RERGFENN7yHGZ0mKqRlVhSIdSFtF+ft7nh8tlAiMBeT4KSdUpRr4Lf4mIiIiCjjP95LhMaTEXnX1R2sx8qYtmC7kyUOzjF3M1gYiIiChIGPST4zKlxRSSs58vu5OGjoO9uP+5/eg42OvI4zs9ZiIiIiKvMb2HSpZPZR5rnX23Zsw7DvamdOnddO3Ckmtzs24/ERERhR1n+qlk+VTmsebsuzVj3t51LKVLb3vXsZIfk3X7iYiIKOwY9FPJ8q3M44WFM8djVDyGCgEq4zEsnDnek+9LREREFGRM76GSZapyY63M44UF06ux6dqFaO86hoUzx5ec2kNEREQUBQz6qWCZcvgFAgA5u+m6bcH0agb7RERERCZ5pfeIyBkiskVEfi0ie0XkIyJypog8KyL7kh+rTfvfIiL7RaRTRC51b/jkh0wlOWdWzWSVGyIiIqIAyjen/xsAnlLV9wNoALAXwM0AtqvqbADbk7chIucCWAHgPACXAfiWSHJFJ0VCppKcdzff7WrNfCIiIqJMnCzZHUU5g34ReS+AZgAPAYCqnlLV4wCWAtiY3G0jgGXJz5cCeERVT6rqAQD7AXzI2WGTn4wcfmO238jhrz+znlVuKAXfgImIyAtGye57nunEygfbedyxkU9O/0wARwH8m4g0AOgA8CUAk1S1GwBUtVtEapL7TwHQbvr6Q8ltKUTkOgDXAcC0adOK/gHIffnm8BN1HOwdXkQNwLZngnkfrr0gIiIn2JXs5jEmVT5BfxxAI4C/UNWficg3kEzlyUBstmnaBtUHADwAAE1NTWn3U3A01DSgq68L/UP9w9sqY5WY9p5p6OrrYg4/AUhvjHZ141TbnglON08jIiIySnb3DwyxZHcG+eT0HwJwSFV/lry9BYmTgDdFpBYAkh97TPufbfr6qQCOODNc8gNz+CkTc/qOdZZFgbSeCW40TyMiIjJKdt94ST0nlDLIOdOvqv8jIm+ISL2qdgJYDOCV5P9rANyV/PhE8ku2AXhYRO4FMBnAbAAvuDF48kamOvxGDj+VJ+vM/h2Xn5cyy3J141Rc3Tg1LZWHMzFEROQGluzOLt86/X8BYJOIjALQBeBPkLhK8JiIrALwOoBPAYCqviwijyFxUjAA4HpVHXR85OSaqOXwM4fcOebn0jpr33vilG1jNPNzzuZpRERE/sgr6FfVFwE02dy1OMP+6wCsK35Y5Kco5fBbZ6NLveSX6YSovroeW67YUspQAynbwlzrzL4RxOd6fjkTQ0RE5D125KU0rfNa8cT+J1K2GTn8X3/h66Ga5Xd6NX+mE6L5NfMdGG2w5FqYm2lmn4iIiIIn3+ZcVEaiVIffWM1vXkxaikyLmsN0IpRNoQtzF0yvxvUXz2LAT0REFHCc6Sdb5tn+MAe1TueQZ1rUHKYToUyKXZhLZBCRswH8O4CzAAwBeEBVvyEifw3gz5Do+QIAX1XVHya/5hYAqwAMAviiqj7t+cCJiMpAqIN+LtB0jxHcbu7cHPqg1ukc8qicEFkVszCXyGIAwF+q6h4ReQ+ADhF5Nnnffar6D+adReRcACsAnIdEtbcficj7WPyBiMh5oQ36nV6gWc4yLU49p+oc1uG3EZUTIutJs11jEy8W3fLkPTqSXdqNTu2/F5G9sOnIbrIUwCOqehLAARHZD+BDAJ53fbBERGUmtEE/2y07J9Pi1KazmnDbwtt8HFlwtc5rxavHXw3tCVGmk2avF+by5D26RGQGgA8C+BmA8wHcICL/G8BuJK4G9CJxQtBu+rJDsDlJEJHrAFwHANOmTXN34EREERXahbxOL9AsZ1FfnOqGieMmhm5Rc7ZFukZnXK8X5rJDbzSJyOkAvgtgtar+DsA/AzgHwHwkrgTcY+xq8+WatkH1AVVtUtWmiRMnujNoIqKIC+1MP5v8OCfKi1MpIdciXb9Omu1SiijcRKQSiYB/k6o+DgCq+qbp/n8F8P3kzUMAzjZ9+VQARzwaKhFRWQlt0A+wyY+Toro4lRLyXaTrNZ68R4uICICHAOxV1XtN22uT+f4AcCWAl5KfbwPwsIjci8RC3tkAXvBwyEREZSPUQT85JyqLU2mEeYGsX4t08xGUcZAjzgfwWQC/EpEXk9u+CuAzIjIfidSd1wB8HgBU9WUReQzAK0hU/rmelXuIiNzBoL/MZKrUU19dj3/++D+HenEqjbBbIMsZdXKbqu6CfZ7+D7N8zToA61wbFBFRBuVWPY5Bf5nJVKlnfs384cWpFE7mNy+7BbLsnEtERJRQjtXjQlu9h4rDSj3RZLx53fNMJ1Y+2I7qcaNY3YqIiCiDcqwex5n+MsNKPdEU1IW6REREQVSO1eMY9Jehcq/UE8UcviAv1CUiIgqacqwex6C/DJVzpZ4o5fBZT17K7c2LiIioFOU2Ocagv0y1zmsty0o9djl8YfyDz3TyEsafhYiIiNzHoD/ispXo3HLFFh9G5K+o5PBF5eSFiIiIvMGgP+KylegsR1FJg4nKyQsRERF5g0F/xJkX7RrKcfGuWVjTYJjDT0RERMVi0B9xLNEZDczhJyIiolKwOVcZMDfkKvdZ/rAqxyYiRERE5BwG/WXAmO0XCGf5Q8rI4WeHXSIiD+xaDxzYmbrtwM7EdqKQYnpPmSjXEp1u86o6EnP4iYg8NKUR2Pw54FMbgLrmRMBv3CYKKQb9ZWLiuInYcNkGv4cROV5WRyr3HH6WnyUiz9Q1JwL8zZ8DmlYBux8aOQEgCikG/RHDwMhbblVHslbqIfsTLADo7O3E3I1zh2/ztU5EjqhrTgT8O/8OaL6JAT+FHnP6I6ahpgGVscqUbeVclz8fHQd7cf9z+9FxsLfgrzXWSxjPuRPVkYxKPfc804mVD7YXNa4oMi9Iz4SvdSJyzIGdiRn+5psSH605/kQhw5n+iGFd/sJkKoVZCPNz7sRzzW67CZmuWhniEodCMaiDw9v4WiciR5hz+OuagboLU28ThRBn+iPGjZnnKHOiFKbT1ZFYqSfB7qpVXOKokAoAQEWsAp+c+cmU1/rY+Fhc/NjFmLtxbsr/5duWez5+Igqxw3tSA3wjx//wHj9HRVQSzvRHkNMzz1FmBNj9A0MlBdhOVkcq10o9uWb2gUSgf+mMS9H2ahuWzVqGz8/7PJ5+7WkAidd689RmPHngSU8WVhNRhF2wOn1bXTNn+SnUONMfQazLnz8jwL7xkvqiUnsMRnUkp57rBdOrcf3Fs8om4AfsZ/YFgljybcq4arW6cTUaJzWitaE17bW+unF1Wt4/T3yJiIgY9EdW67zW4cCIsivHADuI7BbqVsYqEY8lLkgawbv1BMv8Wmd6GxERkT1RVb/HgKamJt29e7ffwyAiD+VK56mMVeKq2VdBodjcuRmfrv80blt4W87HPXriKJY8vgQnB09CIFCkv8eFsayniHSoapPf4/AbjxdERNllOl5wpp8oAEopGxpWuRbqGjP7hV61Mqf8zKyayRK2REREYNBP5Ltyrctvl85jVOQxr0cpZr2EcaJwd/PdzPEnIiICq/cQ+a5c6vLnk85jVOQ5/PbhkgJz40QBAJbOWoqt+7aif6ifOf5EVJxd64EpjanVew7sTJTwtKv0Q5GT6RgWpnRRBv0hF4UXYdB0HOz1tFymU2VDg66hpgFdfV0p5TTNDbaMGfgJYycMB+xOYAlbIirZlMbU5lzm5l1UFuyOYWFLF2V6T8jZ5UWH7UUYJH6k2jhVNjTo8k3ncRpL2BJRVrvWJ4J4swM7E9sNRnOuzZ8Ddqxjd94yZHcMC9tEEmf6Q848i2kI24swSPxKtVkwvTpSwX6mK1BVo6owpEMpqTZOpPPk4mTzNCKKmHxn8euagaZVwM6/A5pvYsBfZowJpDCni3KmP+RYl9xZRqpNhcDVVJuoV+vJdAXqorMvGp4pyVR3343nxunmaUQUIfnO4h/YCex+KBHw734o/eoARZ55tj+ME6yc6Y8A5iw7x0i1cTOn30ghOjUwhFHxWCRTejJdgVqzYA3GxMdgc+dm25PTcnhuiCiAcs3im2f/65qBuguZ4lOGjInWzZ2bMTY+Fhc/dnHaPkFeU8mZ/ghgzrKzrB16nZ55tkshippsV6Cy1d0vh+eGiAIo1yz+4T2pAb5xdeDwHq9HSj4zjmHNU5tDt6aSM/0RwZxldxQ782ytAGS+HcVqPZly+AUCIPUKlLmcplUUnxsiCrh8ZvHtynLWNXOWvwwZx7CjJ47i6deeTrkv6NkWDPojIlsgRcXLd2GvOagHkHKicMfl5+HO77+ccuLgdgqR1zKVMpv2nmno6uvK+wqUF+lVREQpss3iM6inDMK4sJdBP1EW+cw8W68GXN04NeVE4cmXutNOHMzpQ1GQKYf/7ua78fUXvl7QzEfUKhkRUcBxFp+KFLY1lczpJ8oiUw19c56/9WqAAikVgJZ8oNaTikB+ypTDX39mPavmEBFRJIVtTSVn+olysM48W2f277j8vJSrAVc3TsXVjVNTUlTqz3pPpFJWCsnhJyIiiqowranMK+gXkdcA/B7AIIABVW0SkTMBPApgBoDXAHxaVXuT+98CYFVy/y+q6tM2D0sUStaZ/d4Tp2zz0M3BfdRSVpzK4SciIiqWtWiGH8K0prKQmf6LVfW3pts3A9iuqneJyM3J218RkXMBrABwHoDJAH4kIu9T1UHHRl3GMs2wBrkubNTY5flHLajPxckcfiIiokKxr0vhSsnpXwpgY/LzjQCWmbY/oqonVfUAgP0APlTC9yGTTJ1Og1wXNmoy5fmXE+bwExGRn9jXpXD5Bv0K4BkR6RCR65LbJqlqNwAkP9Ykt08B8Ibpaw8lt6UQketEZLeI7D569Ghxoy9D5hbQBuZPe8/awKschb0dORERhZdx1T3KRTKclm96z/mqekREagA8KyK/zrKv2GzTtA2qDwB4AACamprS7id7YawLS+GXKa2salQVTg2eCs1rkOlxRETRwL4uhctrpl9VjyQ/9gDYikS6zpsiUgsAyY89yd0PATjb9OVTARxxasDEGVbyXqa0sovOvgiNkxpD8xpkehwRUXTwqnthcgb9InKaiLzH+BzAJQBeArANwDXJ3a4BYKzq2wZghYiMFpE6ALMBvOD0wMtZ2OrCUvhlSitbs2BNqHL4mR5HRETlKp+Z/kkAdonIL5AI3n+gqk8BuAvAJ0RkH4BPJG9DVV8G8BiAVwA8BeB6Vu5xXuu81lDNsFK4ZVq4G5Zg3xCVn4OIiKhQoup/On1TU5Pu3r3b72EQURZHTxzFkseX4OTgSYyuGI2nrn4qlMFyWH8OEelQ1Sa/x+E3Hi+IiLLLdLxgR14iB0VloWhUFu7aMWb7N3duDvXPQUREVIhS6vQTkUVUFopGZeFuJkyPo46Dvbj/uf3oONjr91CIiDzBmX4iB2XqVBu24DLTz7FmwZpIzIyHqW06OY+dPImoHHGmn8hBUVkoGpWfg8gOO3kSUTniTD+Rw8yz5GGZ5c+Uwy/JXnth+TmI8mF08uwfGGInTyIqG5zpJ3JYGPsoZMrhn1k1M1Q/B1E+jE6eN15Sz9QeIiobDPqJXBC2haKZmlbd3Xx3qH4Oonyxk2d562trw75Fi7F3zrnYt2gx+tra/B4SkeuY3hNwUSkBWW7CtlDUuDqxdd9W9A/1D+fw159ZH6qfg4gol762NnTffgf03XcBAANHjqD79jsAAFUtLX4OjchVDPoDrqGmAV19Xegf6h/eFsYSkBQcmU4kz6k6Z3i2nzn8RN7iBI93eu5bPxzwG/Tdd3Hk5ltw5KavIF5bi5o1q3kCQJHD9J6Ay5R2wYCMipUpf7/prKbQrUUgioqo9PgIg4Hubvs7BgcB1eGZf6b8UNQw6A84lk4kp2U7kQzbWgSiqCibCZ5d64EDO1O3HdiZ2O6ReG1tzn303XfRc9969wdD5CEG/SFgPhhE8iBAnsp2ImmsRQjSSSU7p4aHiJwtIs+JyF4ReVlEvpTcfqaIPCsi+5Ifq01fc4uI7BeRThG51L/R+6tsJnimNAKbPzcS+B/Ymbg9pdGzIdSsWQ0ZMybnfhmvCBCFFIP+EAhjCUgKtrCcSBqdU+95phMrH2xn4B98AwD+UlXnAFgI4HoRORfAzQC2q+psANuTt5G8bwWA8wBcBuBbIlLhy8gDICx/lyWpawY+tSER6O9Yl/j4qQ2J7R6pamlB7do7EZ88GRABKuxfcvlcESAKEwb9IcG0C3JSWE4k2Tk1XFS1W1X3JD//PYC9AKYAWApgY3K3jQCWJT9fCuARVT2pqgcA7AfwIU8HHSBh+bssWV0z0LQK2Pl3iY8eBvyGqpYWzN6xHXP2voLJd309beZfxoxBzZrVno+LyE2s3hMSYSsBScGRrVpP0E8k2Tk1vERkBoAPAvgZgEmq2g0kTgxEpCa52xQA7aYvO5TcVrZa57Xi1eOvBvrvsmQHdgK7HwKab0p8rLvQl8DfYFTp6blvPQa6u1m9hyKLQT9RxGUq+9p0VhNuW3ibjyPLzeic2t51DAtnjmcjpZAQkdMBfBfAalX9nYhk3NVmm9o83nUArgOAadOmOTXMQIr8BI+Rw2+k9NRd6H2Kz671iTUEpu9X9YEqVN3ZAlyw2psxEPmAQT+Ri4JQe7t1Xiue2P9EyrYw5QsvmF7NYD9ERKQSiYB/k6o+ntz8pojUJmf5awH0JLcfAnC26cunAjhifUxVfQDAAwDQ1NSUdlJAIXJ4T2qAb+T4H97jXdBvLCY2xmE+ESGKMOb0E7koCLW3y6YqCPlOElP6DwHYq6r3mu7aBuCa5OfXAHjCtH2FiIwWkToAswG84NV4yQcXrE4P7uuavZ1hD8BiYiI/cKafyEV+zLJnurogyUyKMM3yU+icD+CzAH4lIi8mt30VwF0AHhORVQBeB/ApAFDVl0XkMQCvIFH553pVHfR81BR5fW1t6Tn7xmLi5psY8FNZ4Ew/kYv8mGXPdHVhZtXM6FcFIV+p6i5VFVWdp6rzk/9/qKrHVHWxqs5OfnzL9DXrVPUcVa1X1Sf9HD9FU19bG7pvvwMDR46MdNy97Tb0PfpvI4uJrQ3DcglAkzGiQjHoJ3KZ17W3M3X2vLv57sBX6yEiclrPfeuh776bsk1PnkLP3knAoltHUn0KCfwD0GSMqFAM+olc5nXt7UxXF+rPrA9ct10iIrdl6qw7cDTZ7M+8mDhfXBdAIcScfiIPuFl7mzn8RESZxWtrE6k9NtuH1TUXHrCbm4wVuS7Adq0B+wOQSzjTT+QBo/a2G7PszOEnIsqsZs1qdzruWpuMFbguwHatwe13oK+trbRxEWXAoJ8o5JjDT0SUWVVLC2rX3on45MmACOKTJ6N27Z2lzaiba/sXuS7Adq3Bu++i5771xY+LKAum9wRIEBo5UfgYOfxb921F/1B/Wg4/EVG5q2ppcTZtxoEmYxnXGmTYTlQqzvQHSBAaOVE4eV0hiIiorDnQZCxlTUEe24lKxaA/QDKlaTCAo1y8rhBERESlcW2tQdixB4JrmN4TIJnSNBjARYcTKVyZHuOcqnOYw09EFFB2lXpq197J6j1WRg8EI33KvH6CSsKZ/oBhmka0OZHClekxms5qYh1+IqIAylSpBwBm79iOOXtfwewd2xnwA+yB4CIG/QHDNI1ocyKFi2lgRBQaTNUAwEo9BTP3QGhaxYDfIQz6A6h1XivTNCIqU7fcQk7unHgMIiJPGKkaRuBvpGpMafRzVJ5jpZ4CldgDgewxpz+AjEZOFE2t81rxxP4nAOQ3Q8+Ou0QUWuZUjaZViQCuDFM18uoKXK52rU+cBBqviQM7gUdWAuddleiBUHdhXik+fnc3DkPZdQb9RB4zZuo3d27G2PhYXPzYxQU/RmWsEtPeMw1dfV2c5SeiYDOnajTfVHYBP5Co1NN9+x0pKT6s1JNkXbj70ncT2+denfh4eA9wwY2pPRAO7ExuXw1gZM2E8fya10x4Ffg31DSgq68L/UP9w9uCVnad6T1EPjBSuJqnNqctyrWKSxwVUpGyjR13iSg0mKrhTlfgqLAu3N3bBqzYNBLgT2kEdt07khKWTBHr6+zHvkWLsXfOuThy8y2+r5kIw3o7zvQT+cBI4Tp64iiefu3prPtWxCpw6YxL8eSBJ9lxl4jCxVxusa4571SNKHK8K3CUZLsaZJMi1jf+C+j+p4dHAv3BQduH9XLNRBjKrnOmn8hHdotyz6k6J22R7urG1SzlSkThc3hPaoBvBHCH9/g5KvKTXUWnn34TeP6bma8GWar59Gx6Km1m347XayaCXnadQT+Rz6xvEnc33532psFSrkQUShesTp/Rr2sezsUOsr62tuH0kX2LFqOvrS2S39Nz1opOP/0m8MxtwMW3JhbuGrP65sDfkiI20J2+KDpNPA49ccLT5zLox2oG/UQ+s75J1J9Zb/umwVKuRETeyNRMy83A0Y/v6QtrDv9z64BLvgZ89IbU+42rQeYUseRJQXyc2j92RQUgAjnjDIgIBo8f9/y5DPKxmkE/UQBY3yTs3jSMdQBBmzkgIooaP5pplVUDL3O6zkduGAn4zfcbV4NsUsRq/nwVZFTqslQZMwaT7/o65ux9BRXjxkH7+1Pu9+q5DPKxmkE/UQBY3ySC/KZBRBR0pabJ+NFMq6waeBVS0emC1eh7qS/l94ma96N23d9mrIaU8bk8ciTaqVM5sHoPERERRYYTNdv9aKZVNg28CqzolOn3Wbv2Tszesd32W2R6LgGkpPsA3tXxDwLO9BMREVFkOJEmU7NmNWTMmJRtbjfT8uN7+qLAik7F/D7tnkuryKZOZcGgn4iIiCLDiTQZP5ppRbaBl7VEp5Grv2v9yLYsFZ2K+X1an8tMIpk6lQXTe3y0fNtydPZ2pm2vr67Hliu2+DAiIiKicHMqTcaPZlqRbOBllOg0ZvfN6T15KPb3aX4u9y1abJ/uE4th75xzEa+tRc2a1dF77i040++jhpqG4SZMhspYJebXzPdnQERERCFXNmkyYWEt0VlgR2Ynfp8Z030GB6NdHtWCQb+PzE2ZDEHs4EZERBQWkU2TCTNLR918A37Amd9nWrpPRUXaPuWQ4593eo+IVADYDeCwql4uImcCeBTADACvAfi0qvYm970FwCoAgwC+qKpPOzzuSDCaMm3dtxX9Q/2ojFUGsoMbERFRmEQyTSbMrCU66y7MGvj3tbWh5771GOjuHk69yVSpJ1/m18TeOefa7hP1HP9CZvq/BGCv6fbNALar6mwA25O3ISLnAlgB4DwAlwH4VvKEgWyYZ/s5y09ERESRYtNRF5s/l7E2vxediTOtB4hceVSLvIJ+EZkK4JMAHjRtXgpgY/LzjQCWmbY/oqonVfUAgP0APuTIaCPImO0XCGf5iYiIKFo8KNFZKNsc/3gceuJEpJt35Zvesx7ATQDeY9o2SVW7AUBVu0WkJrl9CoB2036Hktsog9Z5rXj1+Kuc5SciInLarvWJCjLmdJIDOxNBZ4YykeQgu+e4rjljeo8XnYmNNB8jhUiqqoB33sHg8eOJ7xXR5l05Z/pF5HIAParakedj2hVEVZvHvU5EdovI7qNHj+b50NE0cdxEbLhsA2f5iYiInGaUjDTSSYx0kymNfo6KMvAq9aaqpQWzd2zHnL2voGLcOGh/f8r9UVzYm096z/kArhCR1wA8AmCRiPwHgDdFpBYAkh97kvsfAnC26eunAkgrjqqqD6hqk6o2TZw4sYQfgYiIiCiDEktGkrf8KLnqxdWFIMgZ9KvqLao6VVVnILFAd4eq/jGAbQCuSe52DYAnkp9vA7BCREaLSB2A2QBecHzkRERERPkooWQkecuPkqvlsrC3lI68dwF4TERWAXgdwKcAQFVfFpHHALwCYADA9ao6WPJIiYiIiIpRYMlI8pfXJVdr1qxG9+13pCwgljFjcPpFzYluvqbSoWHO8S8o6FfVHwP4cfLzYwAWZ9hvHYB1JY6NiIiIqDTmkpF1zYmAnyk+ZGJd2BuvrcXpFzWjb+v3hk8EorC4lx15iYiIKLoKLBlJ5cm8sHf2ju14+yc7XS8d6rVS0nuIiEJt+bbl6OztTNteX12PLVds8WFEROS4AktGkvfsOvD6PZsexcW9nOknorLVUNOAylhlyrbKWCXm18z3Z0BERGXGiw68xYji4l4G/UQUKh0He3H/c/vRcbC35MdqndeKmKS+DcYkxkZ5REQe8aIDbzH8KB3qNqb3EFFodBzsxcoH23FqYAij4jFsunYhFkyvLvrxJo6biKWzlmLrvq3oH+pHZawSy2YtY6M8IiKPBDWNxm5xbxDSjkrBoJ+IQqO96xhODQxhSIH+gSG0dx0rKegHErP9T+xPtBnhLD8ReYVrihLitbWJ1B6b7X4a/v1cAwAVAHqAt76K+m0bQ/v7YXqPh5ZvW465G+em/V++bbnfQyMKhYUzx2NUPIYKASrjMSycOb7kxzRm+wXCWX4i8kxk1xTtWp8ok2p2YGdiu42gptFE8ffDoN9DUXwBEXlpwfRqbLp2IW68pL7k1B6z1nmtaJzUyFl+InJMX1sb9i1ajL1zzsW+RYvTFqZGdk3RlMZEHwQj8Df6JExptN3djw68+Yji74fpPR4ypxEYwv4CIvLagunVjgX7honjJmLDZRscfUwiKl9GRZpsjZ0iu6bI6IOw+XNA06pEB+QcjdC87sCbjyj+fjjT7yHjBWTM9kfhBURERESp8q1IY55NjtQkYF1zIuDf+XeJjyHtiRC13w+Dfo9F7QVEREREqfKtSBPZNUUHdiZm+JtvSny05viHRNR+Pwz6PRa1FxARERGlKqSxU+TWFBk5/J/aACy6dSTVJwSBv10fmCj9fhj0+yBKLyAiIiK/5Vo067VCKtIYa4oiMwl4eE9qDr+R4394j5+jysnoA3PPM51Y+WD7cOAfpd8PF/L6gIsGiYiInJHPolmvWRs79VbF8e/NJ/Ffb30V2PjV4f0iWZP/gtXp2+qaA5/X70YfmKBh0E9EREShlW3RrJ8VYcwVada2r8UL+7YCQ/3D95dzye6+trbAdbo1+sD0Dww50gcmiM3XmN5D4VBgsw8iIioP+S6a9VMUa74Xy7gyM3DkCKA6fGXG75Qsp/vABLE3E4N+CocCm30QEVF5KGTRrF9YsntEvuVM/bBgejWuv3iWI2k9QTzRY9BP4WBu9rFj3UhlgIDnCBIRkbsKWTTrJ5bsTgjDlRknBPFEj0E/hUdEmn0QEZFzqlpaULv2TsQnTwZEEJ88GbVr7/Q9R9yqmJLdQatK5ARHrsyEJOU3aCd6DPopPCLS7IOIKLJ8CsaqWlowe8d2zNn7Cmbv2B64gN9QSMnuoOa+l8qRKzMhSfkNWm8mBv0UDiFu9kFULkTk2yLSIyIvmbb9tYgcFpEXk///l+m+W0Rkv4h0isil/oyaHBWSYMwvhdR8D3LueykcuTITopTfIPVmYtBP3ih09se6/+E9wAU3jjT3sN7O9XhE5IUNAC6z2X6fqs5P/v8hAIjIuQBWADgv+TXfEpEKz0ZK7ghRMOYoF65wRDn33ZErMyFJ+Q1Scy8G/eSNQmd/rPtPaQR23Tuyv/U2Z5OIfKeqOwG8lefuSwE8oqonVfUAgP0APuTa4Mg7IQnGHOXCFY4wVCXyFVN+C8agn7xR6OxPrv3LdTaJKJxuEJFfJtN/jFp4UwC8YdrnUHIbhV05BmMuHJPCUpXIFz6l/IZ9YTWDfvKOefZn0gfS77deCs01W1SOs0lE4fPPAM4BMB9AN4B7ktvFZl+1ewARuU5EdovI7qNHj7oySHJIocFYSKqw5MXhY1JYqhL54vAe+4lAc8qvw6KwsJpBP3nHPPtz5L+BR1aOvNm3fSlx23wp9KffBJ7/ZubZIhdmk8J+Fk8UNKr6pqoOquoQgH/FSArPIQBnm3adCuBIhsd4QFWbVLVp4sSJ7g6YSlNoMBalhb8uHJPCUpXIcxestp8IvGC1a98yCgur434PgMqEefanrhmouzAR5D+yEvhwK/DS46n7//SbwDO3AZd8DfjoDYn9zV9v93glXk41zuKNP2rjLB4A32iJiiQitapqrDy8EoBR2WcbgIdF5F4AkwHMBvCCD0MkJ9kFXXXN+aVyNq1KBMthTNV04ZhEwRKFhdWc6Sdv2M3+rNgETP5g4lLoh1sTt418yOfWjQT8xv7m2SIXLu1F4SyeyE8i8h0AzwOoF5FDIrIKwN+JyK9E5JcALgawBgBU9WUAjwF4BcBTAK5X1UGfhk5+ikKqpg/pJuStKCysFlXbFEpPNTU16e7du/0eBnnNmBkxz+4c+M/EG3/zTYl8UA/tnXMuYPf3III5e1/xdCxEViLSoapNfo/DbzxeRJDdsSCMgT9FmjUbAEgsrA7iOotMxwum95A/MqX7ACP5kHUXevrGH6+tTSzQsdlOREQuYFoMhYQR2Pfctx4D3d2I19aiZs3qwAX82TC9h/xhvRRqOO8q3zrusjwaEZHH3EiLiVJFIAqUsC+s5kw/+cO62OvwnkROv90bv0ezPVE4iyciCpVCF/7mw6gIZFf4IcCWb1uOzt7OtO311fXYcsUWH0bknL62Nh5bA4BBPwWDE2/8u9Yn3uzNX3NgZ+LEIc8yXlUtLXwjIiIKM5cqArkduDbUNKCrrwv9Q/3D2ypjlZhfM9+x7+EHVsYLDqb3UHS4UO+ZdfuJiELI4YpAXjRmap3XipikhmUxiaG1odWx7+EHVsYLDgb9FB0Ot0GPQvc9IqKy5HCjLC8C14njJmLprKWojFUCSMzyL5u1DBPGTnDse/ghCvXto4JBP0WLg7M7nJ0gIgohcw6/Q4UhvApczbP9oZjlz2PRdBTq20cFg36KFgdndzg7QUQUQi5UBPIqcDVm+wUSjln+PNJqWRkvOLiQl6LD4XrPrNtPRBRCLlQEqlmz2rYxkxuBa+u8Vrx6/NXgz/IDeS2aZmW84GDQT9GRbXYnw5t9tmoMXr7JExFRcHkZuE4cNxEbLtvg+OO6xpxW23yT7fGWlfGCgUE/RUeBszu5yohxdoKIiAwMXDOwptXWXchuygHFoJ/KVraFusYbO9/kiYiIMnA4rZbcxYW8VLa4UJeIiKgELiyaJvdwpp/KFhfqEhERlcCFRdPkHs70U9liGTEiogjIo1Y8EXGmn8oYF+oSEYVTSuW1idWomfNNVN387cQMsznP3CfLty1HZ29n2vb66npsuWKLDyMiYtBPZY4LdYmIwiWt8lrPW+juOw24609R9Yd/Ylsr3msNNQ3o6utC/1D/8LbKWCXm18z3bUxETO8h5/FSKxERucS28trJU+j51RmJWvFNq3zPKW+d14qYpIZYMYmFo+EWRVbOoF9ExojICyLyCxF5WUT+Jrn9TBF5VkT2JT9Wm77mFhHZLyKdInKpmz8ABVAebbmJiIiKkbHy2ltvj9SKt048eWziuIlYOmspKmOVABKz/MtmLcOEsRN8HReVt3xm+k8CWKSqDQDmA7hMRBYCuBnAdlWdDWB78jZE5FwAKwCcB+AyAN8SkQoXxk5BZW7LvWMda/YSEZFjMlVYi9eMBxbdOnL88TnwN8/2c5afgiBn0K8JbydvVib/K4ClADYmt28EsCz5+VIAj6jqSVU9AGA/gA85OWgKAXNb7gBcaiUiomiwrbw2ehRqvnxz4kZAasUbs/0Cic4sP9N3Qy2vnH4RqRCRFwH0AHhWVX8GYJKqdgNA8mNNcvcpAN4wffmh5DYqJ9a23D7PuBARUTRUtbSgdu2diE+eDIggPnkyar/2tdSiDHXN9jXkXdTX1oZ9ixZj75xzsW/RYvS1taF1XisaJzVGZ5Y/z/Rdu+eC/JdX9R5VHQQwX0TOALBVRD6QZXexe4i0nUSuA3AdAEybNi2fYVBYsC03ERG5KGiV19IqCh05gu7b70At7sSGlg3+Ds5J5vTdplW2lZIyPRcAAvU7K0cFlexU1eMi8mMkcvXfFJFaVe0WkVokrgIAiZn9s01fNhVAWttTVX0AwAMA0NTUlHZSQCGWrS03g34iKnOs4R49thWF3n0XPfetj16ga07fbb4p7bheVs9FyORTvWdicoYfIjIWwMcB/BrANgDXJHe7BsATyc+3AVghIqNFpA7AbAAvODxuCrILVqcH9z5canULL1sSUSkaahqGq7oYWMM93DJWFMqwPdRypO+W1XMRMvnk9NcCeE5Efgng50jk9H8fwF0APiEi+wB8InkbqvoygMcAvALgKQDXJ9ODiELPuGw5cOQIoDp82ZKBPxHlizXcoydjRaEM20PLnL6boVJS2TwXIZRP9Z5fquoHVXWeqn5AVe9Mbj+mqotVdXby41umr1mnqueoar2qPunmD0DkpWyXLYmI8sEa7tFjW1FozBjUrFntz4CcYq3Wc3gPcMGNI5WRbColRfa5iAB25CUqAC9bEpETWMPdeR0He3H/c/vRcbDX8+9tW1Fo7Z3hz2G3VuuZ0gjsuje1Wo8lfTeyz0UEFLSQlyhUdq1PvDGZ1xcc2JmcqVhd1EPGa2sTqT0224mI8mXM9m/u3MxZfgd0HOzFygfbcWpgCKPiMWy6diEWTK/2dAxBqyjkiDyq9diJ5HMRAZzpp+jKs55wIXjZkoicErka7j5q7zqGUwNDGFKgf2AI7V3H/B5SdLDZZmRwpp+iq8gZimyMmYue+9ZjoLsb8dpa1KxZzRkNIirYxHETseGyDX4PIxIWzhyPUfEY+geGUBmPYeHM8X4PKTqs1XrqLmTgH1IM+inactQTLgYvWxIR+SiZutkRm4v2rmNYOHM8Fgz9Cj/68It44rRPJW57nNoTGda02AM7gUdWAuddlajW42OzzY6DvSO/b/5+i8Kgn6KNMxRERNEypRH9j1yDfzxxPf5zYA4ujO/Fg+Pux9QVG3F93Sy/RxduRlqsEdS/9N3E9rlXJz761GwzCGs2ooA5/RRdedQTJiKikKlrxg/q/xb3xtbjSxWbcW9sPX5Q/7ec0HGCOS12xzpgbxuwYlPqc+tDs02u2XAGg36KrsN7Ui9B2tQTJiKi8Dl7wWV4RD+BL8W34hH9BM5ecFnJj+lnyc9ACeDCXWPNRoWAazZKwPQeii67mYi65kC8gRERUfEWDP0K88b+GD+fdC2ue/NxVA79CkDx7+1MHzEJYFrsgunV2HTtQub0l4gz/URERBQeydTNyhUb8Qd/eg8qV2wsKnXTPLMflvSRvrY27Fu0GHvnnIt9ixajr62ttAe0dtw1Fu6+vyVwabELplfj+otnMeAvAYN+KiuOv2ESEZUra8AIJG7vWu/u93UgddOY2b/nmU6sfLAd1eNGBT59pK+tDd2335FoEKmKgSNH0H37HaUdx6z9bLIt3KXQY3oPlQ3jDVPffRcAht8wAbAEJxFRoayVXszFE9zkQOqmdWa/98SpwKeP9Ny3fvj4ZdB330XPfeuLP4ZZ+9lkWrgb8bTY5duWo7O3M217fXU9tlyxxYcRuYNBP5UNV94wiYjKlQsNEL1i18xrwfTqQAX71kD0kSMDtukZA93dpX0jF/rZhE1DTQO6+rrQP9Q/vK0yVon5NfP9G5QLmN5DZSPTG2PJb5hEROUqgJVe8mEsDL3xkvrALtptqGlAZaxy+Pax99rvF6+tLe0bWRfuBiB/32ut81oRk9SQOCYxtDa0+jQidzDop7KR6Y2x5DdMIqJyFaKA0VqSM+gLQ62B6MMfE5ysTN1HxoxBzZrVxX8TB/rZRGGt3MRxE7F01tLhk6zKWCWWzVqGCWMn+DwyZzHop7JRs2Y1ZMyYlG0lv2ESEZWrEDVAtC7cDUMtfmsg+sLcMXjxTz6K+OTJgAjikyejdu2dpaWnlrgo2pXFxT4xn2RFcZYfYNBPZaSqpQW1a+909g2TiKhchagBYrElOf1u2GUNRC//wt2YvWM75ux9BbN3bC/9+HXB6vSUrAI67mZbKxc2xkmWQCI5yw9wIS+VmaqWFgb5EWPU2A5qxQ2iyApRA0S7hbu5BKFhlxGIbu7c7Ewgumt9ouqS+Xd0YGfiRC3PQN8samvlWue14tXjr0Zylh9g0E9EIRaEgzIRBV8xHV3trg748f7iaCDqcJnVeG1tIrXHZnsYTRw3ERsu2+D3MFzD9B4iCq2wdNEkIv8VunDXuDrgd8MuIxB1JN3EXGZ1x7rUE4AicK1cuDDoJ6LQCspBmYjc41d1mFxlPUNbtcbBMquBWSvnV3fokGF6DxGFVjGX7IkoPPzupJ6pYZff4yqJtcxq3YUlB/6+/8x+dYcOGc70E1GoBb3WNhEVL6jVYYI6rpxCVGa1IA6nLUUVZ/qJiJKsbe8N9dX12HLFFh9GRFTeSqkO42Zlr1BUrbGr1PPSd4H3t9iXWQ17gGxOW2q+Kfw/jws4009ElGRtew8kOjPOr5nvz4CIylyxndTdbsYVig7vRsqLMYt/YCewtw2Ye3XqfgXU5Q+0EHWH9guDfiKiJGvbeyC6nRmJwqDY6jBuV/YKRdWackp5iWraksMY9BMRJVnb3lfGKiPbmZEoDIqtDuN2ZS8vq9aU1BXYwUo9gRai7tB+ElX1ewxoamrS3bt3+z0MIiIcPXEUSx5fgpODJzG6YjSeuvqpQAT9ItKhqk1+j8NvPF5QvkrJ6S90fY9b64FKbkBozIA3rUqkvER1pp9SZDpecKafiMjEmO0XCGf5iUKslMpeha7vcWs9UElpSkx5IQsG/UQlCm2DFsqodV4rGic1Mpe/QCLybRHpEZGXTNvOFJFnRWRf8mO16b5bRGS/iHSKyKX+jJooXaHre9xaD1RSmhJTXsiCQT9RCYwGLQNHjgCqww1aGPiHm6Nt78vLBgCXWbbdDGC7qs4GsD15GyJyLoAVAM5Lfs23RKTCu6ESZVbo+h631gPl6gqcwtqV1qjIY+5KG5VKPVQUBv1EJQhtgxYiF6jqTgBvWTYvBbAx+flGAMtM2x9R1ZOqegDAfgAf8mKcRPkwz97nM2tf6P75yjtNya5E5+bPJbYTgUE/UUlC0aCFyF+TVLUbAJIfa5LbpwB4w7TfoeQ2okCwW9+TrZKO7+uByqlEJxWFHXmJShCvrU2k9thsJ6KsxGabbTk5EbkOwHUAMG3aNDfHRJSidV4rXj3+KlobWvOqpGPe3xfsSktZcKafqAShaNBC5K83RaQWAJIfe5LbDwE427TfVADpZ9AAVPUBVW1S1aaJEye6OlgKGGueOpC4bc5Td5F5fU8+lXScWA+Uqy5/1vvZlZayYNBPVAIvG7QQhdQ2ANckP78GwBOm7StEZLSI1AGYDeAFH8ZHQRagPHW3G34BI3X573mmEysfbE8L7LPe71GJTlasCy+m9xCVqKqlhUE+EQAR+Q6AjwGYICKHAPwfAHcBeExEVgF4HcCnAEBVXxaRxwC8AmAAwPWqOujLwCm4zHnqRTaYcqpxllFJp9iGX/mwu5pg/j5Z789WotOhNB+jYp1RwMKoWAeAx8EQYNBPRESOUNXPZLhrcYb91wFY596IKBJKzFNvqGlAV18X+of6h7cV2zhrwfRqV4J9g3E1oX9gyPZqQtb77Upx1jU7mtefrWJdOQb9fW1t6LlvPQa6uxGvrUXNmtWBfh6Y3kPlw+fcUCIiKkKJeepuNc5yQ666/Ob7f/ThF7Fg6FepD+DyMY0V60aEsU8PZ/qpfBi5ocblT3P+IxERBY/5fbquGai7sOBSlEYpza37tqJ/qN+xxlluyXU1Yfj+A+d7fkxjxboRhV71cCrNrBSc6afywRrGREThki1PvQBuNc7ylQ/HNFasG1HoVY+Gmobhjs2GYtPMisWgnyItrcrAS30juaFNqxjwExEF2QWr09+n65rt89ez8LxxllfppOb1Dh4c01ixbkSmqxuZtgchzYxBP0WWbb7dbbeh79F/Yw1jIqIy0zqvFY2TGr0JsrwqNepDXf6qlhbM3rEdc/a+gtk7tkc+4M/UF6HQqx7Giacx2+9Hmhlz+imybPPtTp5Cz95JqFp0a1G5oUREFE5G4yxPOFBqNCcH1juUq3zz67N1YTZOdgqp3tM6rxVP7E+0KvEjzYwz/RRZGfPtjibP1ovMDSUiIsrJ7dQbh9Y7lKN88+tzdWEu9KqH52lmFgz6KbLyyrcrIjeUiIgoJ7dTbxxa71CO8s2vd6MLs6dpZhYM+imyWGWAiIh8YU69WXTrSKoP15EFQr759bn6JhT7vTdctsGXkrE5g34ROVtEnhORvSLysoh8Kbn9TBF5VkT2JT9Wm77mFhHZLyKdInKpmz8AUSasMkBERL5g6k3g5VvGdcH0alx/8SxXOzF7JZ+FvAMA/lJV94jIewB0iMizAD4HYLuq3iUiNwO4GcBXRORcACsAnAdgMoAficj7VHXQnR+BKLOqlhYG+URE5C27FJu65pLy+oPQ3ClKjNn+zZ2bA92szUk5Z/pVtVtV9yQ//z2AvQCmAFgKYGNyt40AliU/XwrgEVU9qaoHAOwH8CGHx01ERERR41V9+xAKQnOnqPEzv94PBeX0i8gMAB8E8DMAk1S1G0icGACoSe42BcAbpi87lNxGRERElJlX9e1DKAjNnaLGz/x6P+Qd9IvI6QC+C2C1qv4u264229Tm8a4Tkd0isvvo0aP5DoOIiIiiylzffsc61p03mbhnE5bWfChl8ekNB6ehd/GSka7zbW0+j5KCLK/mXCJSiUTAv0lVH09uflNEalW1W0RqAfQktx8CcLbpy6cCOGJ9TFV9AMADANDU1JR2UkBEREQBsmt9YsbdHIAf2JlYnOpkmUhzffvmmxwP+EObGz+lEa3f/Sc8MakKAHD+SwP48PdfwcBgYq514MgRdN9+BwBwLRvZyqd6jwB4CMBeVb3XdNc2ANckP78GwBOm7StEZLSI1AGYDeAF54ZMRETkrI6Dvbj/uf3oONjr91CCy6vUG5fr2zuZG+/p66auGROv/jcsffsEBMA1O04Bg6nJFfruu+i5b737Y6FQymem/3wAnwXwKxF5MbntqwDuAvCYiKwC8DqATwGAqr4sIo8BeAWJyj/Xs3IPEREFVcfBXqx8sB2nBoYwKh5zrB535JhTb5pWJQJyp1NvzPXt65qBugsdT/FpndeKJ/Y/kbKtmNx4X143dc1onbUcr776KE77vX0Il6kbPVE+1Xt2qaqo6jxVnZ/8/0NVPaaqi1V1dvLjW6avWaeq56hqvao+6e6PQEREVLz2rmM4NTCEIQX6B4bQ3nXM7yEFlzn1pmmV87n2HtS3z7cxUy6+vG4O7MTE/34YG973OcRPs8+MztSNnogdeYmIqKwtnDkeo+IxVAhQGY9h4czxfg8puFxOvcEFq9H3Uh/2LVo8sjj1pT5n1wwg/8ZM2Tj1usk7RcjS5bfmz1dBKlIDf3adp2zyWshLRFQQrxb8ETlgwfRqbLp2Idq7jmHhzPFM7cnEg9SbvrY2dN9+B/TddwG4tzjVicZMTrxuCkoRslwFqVr1ZQBAz0OPYaD3HcRra1GzZjUX8VJGDPqpdAzwyMpY8GccoMzBAlEALZhezWA/l2ypNw4F/T33rR8O+A3G4lSng9nWea149firec3yu1Xxxy5FKOPr0OZ4WrXqy8PBP1EuTO+h0rGZClmx1jZR9FywOv1vuK7Z0cmdTItQ3VicWkhjJre64TK1jLzEmX4qnRcVHSh8XK61TRRloa0lX6J4bS0GjqS19vF9capTFX+swpBa1tfWhp771mOgu5spRCHHmX5yhtsVHSh83F7wRxRhbs0sB13NmtWQMWNStgVhcapTFX/sLJhejesvnhXYgL/79jsSJ2Kqw2ss2Pk3nBj0kzMY4JGZpcrE8JUgvi6I8mKuLmNwYmY56KpaWlC79k7EJ08GRBCfPBm1a+8MxMyyExV/wibbGgsKHwb9VDoGeGTlQa1toihzc2Y56KpaWjB7x3bM2fsKZu/YHoiAHxj5nQjEm9/FrvXpx9EDOxPbTdzsCuzlGgtyH4N+Kh0DPLLyYMEfUdSV48xy0LXOa0XjpEZvfhd5FMkwSn7e80wnVj7Y7njgn2kthd9rLKg4DPqpdAzwiIgc5/nMMuVUSMWfkuVRBc3trsBBXWNBxWHQT0Sly/MyNBEVxtOZZQqeHEUy3C75GeQ1FlQ4luwkotKxGReRK4yZZSpT1iIZdRemBP5elPysamlhkB8RDPqJqHTs1UBEZcSTPgrmyZO65kTAb5Piw27SlC8G/RQpbCLio1zNuHatT1wRMG8/sDOx4Dug6z/KtUESEWXXUNOArr4u9A/1D29zvI9CtiIZnFChIjCnnyIjKE1E+trasG/RYuydcy72LVpcPk1McvVqyKMSRdCUa4MkIsrOkz4KLJJBDmPQT5ERhCYiQTnx8Fw+vRryqEQRNOXaIImIsivnPgoUXkzvocgIQhORbCcekU4zyvcydK4UoIAxDuxb921F/1A/D+xENKx1Xiue2P8EAH8mA5h+SIXiTD9FRsFNRFwoMxmEEw9f5HsZOlcKUACxQRIR2fG7jwLTD6lQDPopMgpuIuJCjjm7F2aRTwpQAPl9YCei4PKzjwLTD6lQDPopMgpuIuJCjjm7F2aRLQUo4NggiYjsWDv0dhzsxf3P7UfHwV5PvjfXFVAhRFX9HgOampp09+7dfg+DytWOdSM55otuLfnhWDaU3CAiHara5Pc4/MbjBQVVx8FerHywHacGhjAqHsOmaxe6Xj//6ImjWPL4EpwcPInRFaPx1NVPMeinjMcLLuSl8paj22Ex2L0wWHgSRlHARZvB1951DKcGhjCkQP/AENq7jrke9Buz/Zs7N3OWn3Ji0E/lK89uhxReRglVo6KSUUIVAAN/ChVPmkFRSRbOHI9R8Rj6B4ZQGY9h4czxhT1AkQ0MW+e14tXjr6LjzQ7M3Tg37X6eGJKBOf1UvkKcY075CULvBiIncNGme5ZvW465G+em/V++bXlBj7NgejU2XbsQN15SX1xqT5HFJYx1BY2TGlnNh7Ji0E/li90OI69sS6hS5HDRpnucLH25YHo1rr94VsaAP+tC3xKLS/DEkHJheg8RRVa8tjbRHdlmO1HY+N0MKqrMz6vBjec3r4W+JTQwdKKZINeORBtn+l3k1CVDosBxobGZU/ra2rBv0WLsnXMuBk+cgFSmzuCxhCqFFXtGuMOrqyh2C33TlNjAsNRmgpFp+BXgY5SfGPS7KDJ/PERWLjQ2c4KxcHfgyBFAFXr8OFQVFWeckV/vBqKA87RnRBkFTm513jan8xgLfSsE9gt9HWhgWOqJYWRShAJ6jPIb03tc5NUlQyLPmXNPm1YlZqQCUPXIbuEuBgYg48ZhTvvzAEauBLCEJ4WRsWjTE0bgZPxtm4PSiHGj9KVdOs+maxeivesYFs4cn57ak624RAHvraVU83EiRSgQAnqM8huDfhdF5o+HyE4JuaduybVwlyU8iQoQoMCp42Bv5mDZIUaw7NTEnDmd50/0CbzRcRzLrvrMyPit5TjtikjUNRf8fBsnhmvb1+Lg7w4WXOb1c4fOwce++S7O7FO8VTWIuptmZtw30GsAAniM8hvTe1zm1iVDIt+VmHvqhkwLdI3tLOFJVCBz4NS0yreAf+WD7bjnmU6sfLDdvvKNA4xg2amJOXM6zysyC5/s/Krr6SbmNU2fvmU7zn95MOX+XHFIX1sb3vna32NCnyIGYEKf4u21f4++tjYA6WsV7QL+wKQxB/AY5TcG/S7jwiuKJAdyT93IF65ZsxoyZkzKNvPCXZbwJCpQAAKnvBbABpC5bv8Xr12FyhUbiy7HmQ/rmqah7jdx3Q8HcdErAiC/Bcq5Jkbs1ipaBWKC04ljVAQx6PeApwuviLzgRGMzBxZamWe19i1aDACoXXsn4pMn2y7czXUlgIhMPAqcrH/HxqyyIecC2ABLqdvv8lUTu4A9fmoQf/jjAQD2wbh15v6UTYljADh15AjmbpyLxzofS0kXAoAYYojHEtnigUljZvNNW8zp94CnC6+IvOBE7mmJ+cKZ8vNr196J2Tu2235NzZrVKV8DsIQnhZMnudQOLSzNJp91NsaMuds5/a6zXjWpu9DRwD/TFcvxfZox26ChpgFdfV3Dgfyx9wITf5f+GMfeO/K5JP8NYQiVsUosqVuCp197GgMYcLzyUdG/c4fWR0QNg34i8k8JC62yXYbOtCjX2N5z33pW76FQswZrgAu51B4ETvn+HS+YXh3eYB9IvWpS15wI+B1O8cnUjLCidhIaJ82wDcatVQYf/pig9UnFaNNk/snKxHaDkd5zaugUYhLDmgVrMCY+xvXKR6X8/gO94NhDDPoptPra2hi8hV0JM1/F5udXtbTwdUKhF5WS0GWzzsaDqyaZrmSedeNfYsNlife8TMGv4YW5Y7DwrA/i/O+/Nnxs/dnl0/HCxBeBZBXCK2dfCYWmBPluVj4y1nGUEvQ7eZLsRSUptzDop1Bi6cUIKHHmK9OsFvPzqRxEpSR02fwde3DVJJ8rmXbBb1ziUCgGdRAxieHyL9yNCTeOvI7OOHEU33x8CYCRE0tVTQnynU5jNtZx9A8MObKOw6mTZKevQHiNC3kplFh6MQJKXGiVrVKPuQsmUVRFoSR0ropbVJiqlhbM3rEdc/a+gtk7tqdNgtl13K2IVeCTMz+ZMe/frgqh0+VNrcyVj5wIrI2fwUhNKvYkOayVpAyc6adQKptLwlGWx8yXNYXr9Iua8fZPdg7frrpyWcrtmjWrsX/eBaGeiYkqEXkNwO8BDAIYUNUmETkTwKMAZgB4DcCnVZVnanlyo4us18p6nc2u9YlqZebZfmvDriJkSz/JdIXo8/M+j8NvH8544uh0+k4+nF7HYZ7tL/Yk2ekrEF5j0E+hVDaXhMuYXQrX8e88Mnz/wJEj6Nv6vZSSnADQ/tx+R3NByVEXq+pvTbdvBrBdVe8SkZuTt7/iz9DCyY9gzGllu87GKFtsXPE0pzwWKZ/0E7vgd8LYCVnTc6JQhdCJk+SwV5Jieg+FEi8JR59dCpeVXUpXmGt6l6GlADYmP98IYJl/Qwknt9MsyEXmssUONezKJ/2knJuGOtE3KaX3Qshwpp9CqawvCZeJfFO1rPuFfSYmwhTAMyKiAP4/VX0AwCRV7QYAVe0WkRpfR+gSlgukjEooW2wn3/STKFwhKkYUrliUgkE/hVZYLgmztGhxMqVw2e1nFfqa3tF0vqoeSQb2z4rIr/P9QhG5DsB1ADBt2jS3xucaT2rqUzg53LAr30mPcg9+yxWDfiIXsbRoYcwnSFJVBamshPb3Z9yfKV3hoapHkh97RGQrgA8BeFNEapOz/LUAejJ87QMAHgCApqYm9WrMTvGipj6vJoSQSw27OOlBmTCnn8hFLC2aP+MEaeDIEUAVevw4VBUVZ5wBiCA+eTLO+MwKxCdPHr5tXcRLwSQip4nIe4zPAVwC4CUA2wBck9ztGgBP2D9CuDlVLjCbhpqG4cc38GpCwJVYtjiolm9bjrkb56b9X75tud9DK3uc6afCuVRmLIpYWjR/tgt3BwYg48ZhTvvz/gyKnDIJwFYRARLHnYdV9SkR+TmAx0RkFYDXAXzKxzG6yolygfk+viGstfvLhgcNu/zAdLbg4kw/Fc4oM3ZgZ+K2cYlySqOfowqkTCVEWVo0HU+QoktVu1S1Ifn/PFVdl9x+TFUXq+rs5Me3/B6rW9yumOLF1QSifNg1AOMJaDDkDPpF5Nsi0iMiL5m2nSkiz4rIvuTHatN9t4jIfhHpFJFL3Ro4+ciFMmNRxdKi+eMJEkWdE+UCcz1+2Dv0UvjxBDS48pnp3wDgMss2o6HKbADbk7chIucCWAHgvOTXfEtEKhwbLQWHucxY0yoG/BlUtbSgdu2dzEPPoK+tDfsWLcbeOedi8MQJSGVqTjJPkChK3K6pX8711ylYeAIaTDlz+lV1p4jMsGxeCuBjyc83AvgxEl0UlwJ4RFVPAjggIvuRqNDAhNyocbjMWJSFpbSo16yVjfT4cSAeR8UZZ2Cwr4/lTYmKYNRf73izA3M3zk27n9V8yAtOdL8l5xW7kDdTQ5UpANpN+x1KbqMocanMGJUXLtwlcp5xNWFt+1oc/N1BLqYk35RrA7Agc7p6j9hss62pHPZmK2UtW5mxMAf9rErkqbAu3GU9dAoDVvMhv7EBWPAUW73nzWQjFVgaqhwCcLZpv6kAbFtqquoDqtqkqk0TJ04schjkiwtWpwf3dc3hD4xZlchTYV24y3roFAZcTElEVsUG/ZkaqmwDsEJERotIHYDZAF4obYhEHmFVItdFYeEuy9FRWHAxJRGZ5VOy8ztILMStF5FDySYqdwH4hIjsA/CJ5G2o6ssAHgPwCoCnAFyvqoNuDZ7Kizlg3LdoMfra2pz/JqxKZG/X+pErIIYDOxPb85RPx90wVDbiDCqFBav5EJFZPtV7PpPhrsUZ9l8HYF0pgyKyslZ6GThyBN233wEAzgaJrEpkz0h9Mq58mBdz5ylKC3fd7q5K5BQupiTyTl9bG3ruW4+B7u5AVqBzeiEvkSvsAkZ991303LfeuT8oViXKzJz61LQqcUJU4PMS1oW7dliOjsLCvJgy5yJ0FjMoS0EPVMPCs8nJEhSb00/kKU8CxmxViajk1KewLtzNxO3uqkROy7kIncUMyo417dIIVF1Jn424bJOTQcGgn0LBk4AxqlWJnGJNfbLm+OdQs2Y1ZMyYlG1BWLhb7FoRt7urEjkt5yJ0j4oZeLI+i4Z1HOzF/c/tR8fB3rT7vAhUl29bjrkb56b9X75tuWPfIwjCcDWb6T0UCjVrVqdcNgOCETCWjSJSn+wuGdeuvTNQl5HDcDmWyClGWtrWfVvRP9RvvwjdfEWv+SZXAn7+zXmn42AvVj7YjlMDQxgVj2HTtQuxYHr18P1eBKoNNQ3o6uuKfKO4eG1t4oqJzfZMvO77wpl+CoWqlhbUrr0T8cmTQ1XpJTIKTH3KdMkYAGbv2I45e1/B7B3bff/9heFyLPkj2+xomOUs41niFb1c+DfnrfauYzg1MIQhBfoHhtDedSzlfi+uopdLmeNirmZ73feFM/0UGlUtLb4HiWXLLsWprjnjLKAnC68dEIbLseS9XLOjYWZehD42PhYXP3Zx2j71M+uxZdGtrhQz8OJvjl2zRyycOR6j4jH0DwyhMh7DwpnjU+734ip6XleYIsA4thVyNdvrztmc6Scix4UlmI7a4mJyRq7Z0bAzFqE3T21On2WUCsyfen7ihgvFDLz4m2PX7BELpldj07ULceMl9bYnr15dRY9yozjzVcGqlpaCrmZ73feFQT8ROS4swXRQFxeTv4zZ0QqB7exo2BmL0Fc3rk5Pu4jFUwMyh4sZePE3Vy7pJPlaML0a1188K+PVqkID1WJEtVGccVXwnmc6sfLB9qLSAb08IWLQT0SOK+rA7kDX30JxrQjZyTU7GhXWWUYAODl4Ehc/drFrVVa8+Jtj1+wRQaqUFMUyx05cFfTyhEhU1bUHz1dTU5Pu3r3b72EQkYMKbvhirRBkvV3mRKRDVZv8HoffeLxw1tETR7Hk8SU4OXgSMcQQi8UwMDSQ9WvCkBtv/rlGV4zGU1c/VXZBv7VSEpCYfOHEhnOMmX5jzUSxkwRHTxzFl3d+Gf9w0T848jrNdLzgTD8ROcI6owQUWKnHoxrhRDTCPMt4+TmXo0Iqsu4fltz4qKaTFIKVktzn1FVBr/q+sHoPkcei2PLcsdrbLtcIJ6J0rfNa8erxV7FmwRqMiY8ZrrISlzgUikEdHN43TLnxxs8VlvE6LSwFFcJuwfTq0KQAcqafyENRbXnu2IySyzXCiSideZbRvKiwIlaBT878pKc5/04q967ZYSmoQN5h0E+BFaQFSE6J6uVWR2aUzDn8i24dSfVh4E/kGWtajLnCTwwxxGOpCQJhSfcpR6xORlZM76FAimqr9qhebi2m/XiabF1/meZD5BlzWsyEsROGm3ldfs7lePq1pzGAkYW+/UP9eLTzUTza+WjWxwzD4t+oKaZZVGTtWg9MaUw9lhzYmTi+OFiSNugY9FNuPvyxhKWja6EcCY4DyJGujgV2/SUidxhpMYZMOf+S/DeEoeF97bYBQGdvJ+ZunDt8mycB3mAn+6QpjZmrw5URBv2Umw9/LFGdEfei5blXrAuSq65chrd/spMzSkQRYz4JaJ3Xiif2PwEAw7n+p4ZODe9rt82O9SQgF/NJwvJty9HZ25l1H6IU5upwTasSa8bKsDocg37KzYc/lqjOiEflcqtd+lXf1u+x/jNRxBk5/5s7N+PK2VdCocMz/5WxyrRtdhWACmVdN9BQ04Cuvi70D/Wn7MerCeHg20kbq8Mx6Kc8efzHEqUZcasoXG6NavoVEeVmzvlX1eGZf6Ocp3lbRawCl864FE8eeLLokwC7dQMCyfo1XGAcXHYnbZ78vqzV4eouLLvAn9V7KD8el1L0olU7ZbFrffrv+MDOxHZEN/2KiHIzl8K0a4KVrQKQtQxoZawS51Sdk/G2QBCzhCqVsUrMrJo5vE9c4mlNxcLUT8BrflfGM5eFNbj++2J1OAAM+ikfPv2xVLW0FNbR1Wk5At9IM9ZxGD+/8RqY0giA9Z+JaETrvFY0TmpMCdrM27KWAZUY7m6+O+PtylhlWplQ6z52JxLl2oU3lyD0ijFeD57+vrJVhysjDPopt3L9Y8kR+EaaeR3HjnWpC7nB+s9EkeDQxIZdEyzrtmwnAfVn1me8feXsK7Fs9rK0ANH6NdYTCc7y2wtKrxjzbL8nv68LVqen8tQ1l1W5ToBBP+WjXP9YcgS+kWdex9G0KuXnZvoVUQR4OLGR7SQg1+1MAWK2EwnO8tsLSmomf1/+4EJeomzKebW/ZR1H34FR6Nn0VErVodk7tvs9SiIqlo9lDK29AHLdNioGmQPETP0EOMufWZAq47n5+2JZV3uc6afA8HtxkS2PFzAHhmUdR9/4L6D73od8zQMlIhdkuaIXJHbrBqzs0owoVZBSM938fTXUNAynhBlY0Ykz/RQQdnXfu2+/AwD8SxkxB751zYnyXuWS4mNZx9Gz6SnoYGqJvDCV6Ow42Iv2rmNYOHM8Fkyv9ns4FFKRfB05UMYw1/NibeRXTG8S66w+FaeqpQUHfvsO+v7pHzHhRC9+O64aVX/xRV/fx92YlTc3kTNwrQeDfgqIQNZ9z7aAOWpB/671iTxe4+e6YPXIgr4LVgcmD7QYHQd7sfLBdpwaGMKoeAybrl0YnYCNPBPJ15EDExu5npdATuiUuf+a0YR7Lr0VQwpUCHDjjHrM93E8btTtN9YMmBvHebV2IMiTA0zvoUAIZFBZTguYI1yis73rGE4NDGFIgf6BIbR3HfN7SBRCbryOOg724v7n9qPjYK8DIyyCA5XZcj0vxVaL8f25ibCFM8djVDyGCgEq4zEsnDk+bR8vn3+7uv1GQ7a5G+cO/1++bXnRj+vVLL9xEnzPM51Y+WB74F6/nOl3EBeOFC9Ii4vKUo4FfWHukGwc4PoHhjIe4Ihyyed1ZJ3hyzbj59aVg3xmGUf2uSZ9n7pmdMTmov25/XnNVOZ6XoqZ0InkVZUAWTC9GpuuXej5azMT66y8JP8NYShlv87eTszdOHf4dq7Yynhc6wJwN9mdBAfptcug30G+tZaOgDAHlaVyIt+1INZUHiAxs394T0qlor6X+tCzavHwuKquXIa3f7LTu3E6JNcBjigfdq8jc4ANICVQuuPy83Dn91/OGDjlGxw4fSKRax+7+43xZhpDtr+vfCd0zI8Z9MApChZMr874nPrx/Jtz8I0FuKeGTmX9GutJgJ1zqs7JuQDcSUGfZGLQ7yAuHCmeETx6GvwGgCf5rtYgf0oj8MhK4LyrgCu+MZLKc8GNwwv6+h79N3S3nwY9eWp4XH1bvxfaWvzZDnBE+TK/jqzB8dWNU1MCpSdf6rZNezEC20zBQaEnEubHzJRqYw7Ic+1jvf+7ew7h8T2Hcp7M2J1ctHcdw/mf/TOM+cbdaRM673z2z3B/8mqC3c8Z5MApisyvOz8CV/Os/JWzr4RCh2f+4xKHQjGog8P7xxBDLBbDwNBAxsesjFWi6awm3LbwNtfHbwj6JBODfgf5uXDEUdlmgl3MZ69qaQllQFkKVxYw5wryDS8/Dpxekwj0L7gR2HXvcEpPz1//EHryLWfHVerPAXjyOiTKhzU4ViAlUFrygVr8/LW3hm9XjxuVNoNud+WgkBMJu4DcPAa772kN6Kz7WB9DgJwnM9muJvxTvAqPfOkrOO3//uvwhM47n/0zrOiqwqnfdNr+nL0nTgU6cIoau6s7fjz/5rr9qjo8iVoRq8ClMy7FkweeHI6tltQtwdOvPY0BZA76/Zp0DfIkE4N+h5ln+0M7y28s6jRyus0VHshRrixgtv7+DOYgf8Um4MB/jjQdGxpI2X/gqP3iI08XVvN1SAFmDZ6vbpyKqxunpgRK9We9J+ss/PUXz8qa8pPrRMIakFuD5UzfM9s+1scAgO/uOZRxDHazwNbH/K8ZTbh+x4rh++9/bj9O/aYz489pPH9BDZyiJp/XpheyNWT7/LzP4+nXngaQiK3WLFiDMfExKZOs094zDa///vVwT7q6jEG/w/xYOOI4j7o0ep7LHkCOLGC2zojXNSdm7h/+Q+AjN9gH+UBqbW7L7zcQC6t96hZa6OuSC/jLU6bL+OZAyRq45kqZKPREAkgNyO2CZbvvmWsf6/3Wn9M8BrvAMFd6SD4/J3mnmEXqXjDP/E8YOyEttrJOst7dfDdW/nDl8O1QTrq6TFTV7zGgqalJd+/e7fcwHHP0xFF8eeeX8Q8X/UM4g37DjnUjQeKiWx19aGsuO5DI8wxrznixinoerEH+gZ32OfqzPg788tHE78+ov920CvjZvyS+bsWmlBn0vvFfQM+mpzDQ3Q2pqgLeeQfaP7Io3bffj4uvQ6tifh9r29cOzzYZKmOVuGr2VY7mkopIh6o2OfaAIRXm40VhlXXyC65y7e/G98yHE+Mi7/hRaapQdrHV2va12Ny5GZ+u/zRuW3hb2u1ylel4waCf7BmBo0szrPsWLbafSZ48GbN3bHfs+4RBwVc8rA11jKAfAD7cmpqjbxfkt30JeOnxkdsA+h76e3R/49+hp0z5kfE4Kk4/HYN9ff5diXH5dWhVzOvy6ImjWPL4EpwcPDm8bXTFaDx19VOOnvQz6E/g8YLIW/c/tx/3PNM50szrknpcf/Esv4cFIP1EIDKTriXKdLxgeg+lc6BLYy6BbMblk5wLmAtN35n3hymLcvFOTyLIN7R8A/jA1SmdhXs2PZUa8APAwABk3DjMaX/e0Z83bx68Dq2KeV1GZgE/EZENu/SfoFypsa4DsN72W9DSmBn0U7psXRodCrYCkTMeRHYVa2LxRID/R4+OzOzvuheY02Kfo//8N4GLbx15DJsgPxFEj3yPQJ6EefA6tCr2dRmJBfxERDas61cA5OzlADifqhY2xZTkdnuNGIN+SmdXDtESJJYqtM243C4jaVexZte9iSDenOZipO8035RI3/nZv4yk6xgz4rXzMgb5VoE8CfPgdWhV7OsyEgv4KfpYBpeKZF7cff9z+7OWjjVOAgpp8ma3f9gD/2JKcrvd5JVBP/kitM24nC4jaZO601dxOXqWX4uBdyoQP01R8+dfQNVHbwDe7csvfSePGXHrJcfTL2pG39bvhe8kzGGlvC7NlSaIAollcMkB1nQfa+lYo+FbIU3erL0aotCFuZgr6G43eWXQT57N/tjltoVu0a7TZSQtB+G+h/4e3d9+BjpYAQAYeEfQ/U8PAwCqjhWfvmNmd8mxb+v3UHXlMrz9k53hOglzQbFN4oKWS0qUxqcyuBQtduk+1tKxAApq8mbXqyHsirmC7vYaMQb95MnsTzG5bYFV15w4YBr59IUcMHMsyu351v+FDkrKl+i776LnWw+hasuGotN3zDJdcnz7JzvDdxJWpKjljhLlrZT3L6KkXL0crNuA7E3eotiroWbNarxx6+2oODVS2W1w1Gic/Oyf4f7n9mf8Od1cI8agnzyZ/Skmty2wDuxMbWxVd2H+z1WmnP3kotyBdybbftnAO+LYgtZALtr1UFBqTlN0BbphWynvX3kKWsUScp9dB+Vimrx59l5snYDbtT5RNGNoYCTDocSMh6qWFhz47Ts49k//iAknevHbcdU48dnrsKarCqd+05nx+OPmGrGYY49E4Wae/Wla5fhBIDKBpvkqyKJbR06WDuy033/X+tT7zDP7O9YlvvaCG4H9PwKab0L8tEHbh4mfeXrqhrrmgt6I+trasG/RYuydcy4Qs/+zD3zlJOtzCSRPmtYX9DB2LeeJnNRQ04DKWGXKNicX4xWt0PevIhhXdQeOHAFUh6/q9rW1OfY9KJwWTK/G9RfPSgnwzbc9ZUzAGa/9WBx45rbER2Dkb2VKY0nfZv6frMDpW3+AH/9zG07f+gN0zb/A9vjTcbAX9z+3Hx0HewEkZvsbJzU6vkaMQX+5yRQ4bftS6uxPiQcBc5C5b9FiVFRV2e4X+EDTylpG8vCeRNB+eM/IPuZA1PrGYi23OevjI4tyF92Kmj//M0hFasM8GTMGp196ecrzWcgB1HoQxmD6iUUoFu3aPZdFvCkbi9AqBJHJHaVgaZ3XipikHl4DUco1Wxlch2S7qksUGOYMhx3rEsfhS76W+GhMyDmU8WA+ubE7/hhXn+95phMrH2xHx8He4TViTleCY3pPCQJ9CTcTu/QSo5urteRjkS94u/x9xOOQykpo/0gZqlAEmlbW2XXz8wmkdrsF7BtpmcttWhblVq36MgCg56HHMND7jm1lnULXQ9gdhAEAFRXA0FB4Lr87lIZmXYTG1B5yWmAbtnlQBjcyV3Up+qzrW8xV8lxa72J3/LGWQHWzchGD/hK4XU/VEfl0c/3AVYnqLw7ljNsGmQMDwBlnID5uXLTyPK2BqLl0JpA+s28tt2mzKLdq1ZeHg38A2LdoccHrIcw5tVC13QdDQ5iz95Xifm6/OLQI0Zprmm1hb6H5yaGcDCDHlWvDtkD2/CCyY13fMqbK9fUuQPrxx67jsVtcC/pF5DIA3wBQAeBBVb3Lre/lF7frqToix8JRNN+UyO20KnD2J58gU/v6MLv9+eJ+jiBLC0QvzNxIy1puM48TrIwzZ0eOYO+cc4evBhjlNqWqCnjnnZSrKnZCeRA2v0k//83Em/RHb0i9v8CFV9kW9hZTdSoUkwEBEtVjRbk2bAtt40UqL+b1LXXNiWPJM7clUnw+ekPJGQ+F8PLqs2imWcBSHlSkAsBvAHwCwCEAPwfwGVW1nVZsamrS3bt35/34y79zESb891H80Y8V438HHHsvsHsW0LQfw7cf/lii7GG2fUq9/fDHBONQgWU/HsT43ymOvRfoO/dMvP83pxKpGdWn4fSGGXj7F68Vd7u2FjUXnQko0LPzrcRMYxGPcXr9mXj7hV8MN3s6/UNz8fYLv8LAO5Js/rQKqHl/WrMmc832bLfzDjInT45mSUjjzcOccnLgP0dm9vf/yL4cap5vJPsWLbadOSuFjBmD2rV3hutKi/W5++k3U9+kD+xMXMG6+NaREwFrRQab27u6jmPXb/4H/zLQgi/E23D++87CjDNH44nTPoXFt6yAHn8nbSgyNo6K6hrbv8lY9Tg8NeMEPrg/+/vTbz84EVs+85O8f3wR6VDVptKfyOAo9FgBuHO8yOe9Hij8eLLw1RjO6BtCxZmn4b0Nde6+1+dzPFl5GdD7Gnp++Oui3uvzuW0E96zeQ4GWT/Weti8BCuCKb9jvU8Ttvp91oee7z+f/d5/8u62qryxoMivT8cKtoP8jAP5aVS9N3r4FAFT163b7F/omvmHdVZj/nb0YPTCyTQGYq5ufEkBiQOVg5n1Kvd0vQEwUFUNOPqrltihEAHX7e1TEoQP2lWOcEMogMx/WQNS8RuLDrSMz+yXMRltnm0siEt6DsF0TuZ9+E3huXfp6iUwnBja39ZnbcNfgH+PBgSW4Nv4kbq74D9w1+Mf414El+P7Wv0r5yxmR/W9KoZCU26l7n4wDL35mDj53qyUdLIuIBv0FHSsAd44XuW47cTzhez1RyOSaaCrwdt89X0T3t5+GDpoW+efxdy8VitobV6Wk/eaS6XjhVnrPFABvmG4fAvBhpx78Iz96C0MDqdusB+ZRCmAw+z6l3q5UAOr0o1puq9hky7jwPdw6CIQ5yMyHtRqG4byrEmlTJTbSAkbSSHLm6OcQ+istdidJ1oVXH70h8Vybr7wYFRne7bO9LZd8DX/1k3vw8UmjMf/NH+C/pqzGda8+iHEVb6Ni3BCGTtgVOcv+NyVpt1ONHki8j8Ems67MuHqsAPI7XuS67cTxJPLv9ca3CGv/FSIru+IROY4n2W73PPxMasAP5PV3r4OCnk1PFRT0Z+JW0G83OZbyY4nIdQCuA4Bp06YV9OBD/9NT9MDIO6EPMvNhDUQP7xmpggSUvCjaUNXSMnwQLSbdJ7I5tZkaDWWryGBzu/LdPvxB8vbYc/4cj+w/iC/FH8eWuReifs8biPefGv6W6TO4xeH7GIA8jhUAjxdhw0o9FBlFHE8y3c7UfDMfTv1NuRX0HwJwtun2VAApUYqqPgDgASBxubaQB89UHYCCI7JBZi4elMSzWyiXJh5HxemnY7CvL7pXWqyXXo2rKhfcmL0iQ47bC8ZUYd7YH+Pnk67F0tGP4/VzLsHA936MwROxxMx/7HTg7RMlDz+UC6mdl/NYAfB4ETZ8bVNk5KrwU8DteNv/xcA7xU0ZOfU35VbQ/3MAs0WkDsBhACsA/JFTD16z8jJ03/sQdDAA+ZGh/B4W8ThEJOdC3KzKIcgMCGu6T6bFdJF//u0aDV1wYyLH/48eta/IkOftyku+hj/46A3AT+twzjO34fW1X0XbuKvw8bGdmPzQDej+6ejs7z85bkuFJhZVkqvHCiDP40Vk3+v9yekvywkfip5cFX4KvF3zR28WndPv1PHClaBfVQdE5AYATyNRhu3bqvqyU49fVV8J3LgK//MfP8Bg95vDlRF+94sDGHzrHVTUTsJZH5sIKPA/Pzlqv0+pt03fw9VqC05U7/ngbLz9fEdKs6dc1RaKqdgQ+SAzQMzpPmXL7qrK0MBIwG/cvuRrGE7qLvL2tKEBXH/BLACzgKPPAFOOFPw3mfLe8cefTLyPlTm3jxVAHseLAt7rCz2eePpeH6DqPWX/3kTRYJ1YKvF4UvXhmQAuLb56jwNcqd5TqEKrMRARlZsoVu8pBo8XRETZZTpe2JWlICIiIiKiCGHQT0REREQUcQz6iYiIiIgijkE/EREREVHEMegnIiIiIoo4Bv1ERERERBHHoJ+IiIiIKOIY9BMRERERRRyDfiIiIiKiiGPQT0REREQUcQz6iYiIiIgijkE/EREREVHEMegnIiIiIoo4Bv1ERERERBHHoJ+IiIiIKOJEVf0eA0TkKICDyZsTAPzWx+HkIwxjBDhOp4VhnGEYI8BxFmO6qk70exB+Mx0vgvS7yYbjdFYYxhmGMQIcp5OCNkbb40Uggn4zEdmtqk1+jyObMIwR4DidFoZxhmGMAMdJpQvL74bjdFYYxhmGMQIcp5PCMEaA6T1ERERERJHHoJ+IiIiIKOKCGPQ/4PcA8hCGMQIcp9PCMM4wjBHgOKl0YfndcJzOCsM4wzBGgON0UhjGGLycfiIiIiIiclYQZ/qJiIiIiMhBgQn6ReQyEekUkf0icrPf47EjImeLyHMisldEXhaRL/k9pkxEpEJE/ltEvu/3WDIRkTNEZIuI/Dr5nH7E7zHZEZE1yd/3SyLyHREZ4/eYAEBEvi0iPSLykmnbmSLyrIjsS36s9nOMyTHZjfPvk7/3X4rIVhE5w8ch2o7RdN9fiYiKyAQ/xkbpeLxwFo8XzuHxwvExBupYkRxTaI8XgQj6RaQCwP0AlgA4F8BnRORcf0dlawDAX6rqHAALAVwf0HECwJcA7PV7EDl8A8BTqvp+AA0I4HhFZAqALwJoUtUPAKgAsMLfUQ3bAOAyy7abAWxX1dkAtidv+20D0sf5LIAPqOo8AL8BcIvXg7LYgPQxQkTOBvAJAK97PSCyx+OFK3i8cACPFyXbgOAfK4AQHy8CEfQD+BCA/arapaqnADwCYKnPY0qjqt2quif5+e+ReNOZ4u+o0onIVACfBPCg32PJRETeC6AZwEMAoKqnVPW4r4PKLA5grIjEAYwDcMTn8QAAVHUngLcsm5cC2Jj8fCOAZV6OyY7dOFX1GVUdSN5sBzDV84GljsfuuQSA+wDcBICLn4KDxwsH8XjhOB4vihSGYwUQ7uNFUIL+KQDeMN0+hAC+OZqJyAwAHwTwM5+HYmc9Ei+8IZ/Hkc1MAEcB/FvysvKDInKa34OyUtXDAP4BiTP3bgB9qvqMv6PKapKqdgOJoANAjc/jycefAnjS70FYicgVAA6r6i/8Hgul4PHCWevB44UjeLxwXSCPFUB4jhdBCfrFZltgz5RE5HQA3wWwWlV/5/d4zETkcgA9qtrh91hyiANoBPDPqvpBAO/A/0uLaZI5jksB1AGYDOA0Efljf0cVHSJyKxJpEJv8HouZiIwDcCuAO/weC6Xh8cIhPF44i8cL9wT1WAGE63gRlKD/EICzTbenIiCXxKxEpBKJN/BNqvq43+OxcT6AK0TkNSQuey8Skf/wd0i2DgE4pKrGzNcWJN7Ug+bjAA6o6lFV7QfwOICP+jymbN4UkVoASH7s8Xk8GYnINQAuB7BSg1c7+BwkDty/SP4tTQWwR0TO8nVUBPB44SQeL5zF44ULAn6sAEJ0vAhK0P9zALNFpE5ERiGx8GWbz2NKIyKCRE7hXlW91+/x2FHVW1R1qqrOQOJ53KGqgZtpUNX/AfCGiNQnNy0G8IqPQ8rkdQALRWRc8ve/GAFcQGayDcA1yc+vAfCEj2PJSEQuA/AVAFeo6gm/x2Olqr9S1RpVnZH8WzoEoDH5uiV/8XjhEB4vHMfjhcOCfqwAwnW8CETQn1ykcQOAp5H4A3lMVV/2d1S2zgfwWSRmQ15M/v9ffg8qxP4CwCYR+SWA+QD+1t/hpEvOLG0BsAfAr5D4mwlE5z0R+Q6A5wHUi8ghEVkF4C4AnxCRfUhUEbjLzzECGcf5TQDvAfBs8u/oXwI4RgogHi/KFo8XJQjD8SIMxwog3McLduQlIiIiIoq4QMz0ExERERGRexj0ExERERFFHIN+IiIiIqKIY9BPRERERBRxDPqJiIiIiCKOQT8RERERUcQx6CciIiIiijgG/UREREREEff/Aygxderv4q8HAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 936x576 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(1,2)\n",
    "fig.set_size_inches(13,8)\n",
    "ax[0].plot(data_freq, mu_curve[:, 1, 1], '.')\n",
    "ax[0].plot(data_freq, mu_curve[:, 0, 1], 'x')\n",
    "ax[0].plot(data_freq, mu_curve[:, 1, 0], 'v')\n",
    "ax[0].plot(data_freq, mu_curve[:, 0, 0], 'o')\n",
    "ax[0].set_title('mu curve')\n",
    "# --------------- #\n",
    "ax[1].plot(data_freq, std_curve[:, 1, 1], '.')\n",
    "ax[1].plot(data_freq, std_curve[:, 0, 1], 'x')\n",
    "ax[1].plot(data_freq, std_curve[:, 1, 0], 'v')\n",
    "ax[1].plot(data_freq, std_curve[:, 0, 0], 'o')\n",
    "ax[1].set_title('std curve')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fe251db32e0>]"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsYAAAHSCAYAAADvxw2lAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAhtUlEQVR4nO3df4zc91kn8Pez9jbE/FCdejcKidukss+QIq9bVlG4qohr6NU9UNbSuZIRIOuuKLdSjiunSpAcf5zuj5wi3YFAOsopKmBLlEapaRWL43rkUhBFos2tS4EmYS920iYmITElKRxIaX587o+dDeN4bc+sZ3Zmdl8vKfrOfOY7s89qnOSdJ8/3863WWgAAYKubGnUBAAAwDgRjAACIYAwAAEkEYwAASCIYAwBAEsEYAACSJNtHXUCS7Nq1q914442jLgMAgE3u1KlTf91am1nrtbEIxjfeeGOWlpZGXQYAAJtcVX39Yq8ZpQAAgAjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAkh6CcVXtq6qvdP31t1X1M1V1TVU9VFVPdI47u95zd1WdrqrlqvrgcH+FjXHq6y/mV37/dE59/cVRlwIAwBBsv9wJrbXlJAeSpKq2JfnLJJ9NcleSh1tr91bVXZ3nP1dVNyc5kuRdSb47yf+uqn/SWnttOL/C8J36+ov58U98Md969fW8ZftUPvlTt+b737Hz8m8EAGBi9DtKcVuSM621rydZSHK8s348yaHO44Uk97fWXm6tPZXkdJJbBlDryHzxyW/kW6++ntdb8sqrr+eLT35j1CUBADBg/QbjI0k+1Xl8bWvtuSTpHGc769cneabrPWc7axPr1ne+LW/ZPpVtlUxvn8qt73zbJc/vd+zCmAYAwOhddpRiVVW9JcntSe6+3KlrrLU1Pu+OJHckydvf/vZeyxiJ73/Hznzyp27NF5/8Rm5959suOUbR79jFes7vpQ4AAPrTczBO8qEkX26tPd95/nxVXddae66qrkvyQmf9bJLdXe+7Icmzb/6w1tp9Se5Lkvn5+QuC87j5/nfs7CmIrjV2can39XP+MEO0wA0AbHX9BOMfyz+OUSTJySRHk9zbOT7Ytf5bVfWLWbn4bm+SR6681MmwOnbxyquv9zR20c/5wwrR67m4UOgGADabnoJxVe1I8oEk/6Zr+d4kD1TVR5I8neTDSdJae7SqHkjyWJJXk9w5rjtSDCOw9TN20e/5wwrR/Xa5hxm6hxm4BXQA4FJ6CsattX9I8rY3rX0jK7tUrHX+PUnuueLqhmiYW7D1OnbR7/nDCtH9drmHFbqHHbjHJaD3Q5gHgI3TzyjFptJvl3RcDCNE99vlHlboHmaXe1wC+up7hnUh56R12/1HBQDjZMsG4367pJOon851v+cOI3QPs8s9LgG9n7A7LmF+WAF92MHfbi8A9GvLBuN+u6Scbxihe5hd7nEJ6P2E3XEJ88MK6MPckWWYnw3A5rVlg3HS/ywwwzesLnc/54/LGMq4hPlhBfRhfe6wPxuAzataG/0WwvPz821paWnUZcCGGId53WHNGK92X1cD6aDHLnr93GF/NgCTq6pOtdbm13xNMAYGaRyC/zh9NgDjRTAGAIBcOhhPbXQxAAAwjgRjAACIYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIkmwfdQEAbG6HTx7O8ovLF6zv27kvJ24/MYKKANamYwzAUM3NzmV6avq8temp6RyYPTCaggAuQscYYJMal07t4v7FPHj6wfPWpmoqi3OLG1YDQC90jAE2qXHp1M7smMnCnoU3apmems6hPYey6+pdG1oHwOXoGANMkH66wMPs1Pbbje6uRbcYGFeCMcAEmZudy5PffDKvvP7KG2sX6wKvdmo/+8Rn88rrrwy0U9tPHd21fHr505esYVzGP4CtSTAGGLFhdoGH1aldTzd6cf9izrx05pLn9Bu4AQbJjDHAiPUzC9zvvO7q+ZUa6FzveuaGZ3bM5NjBY5c8Z3H/Yqbq/H81Gb0ANoqOMcAQjFMXuJdObb81r6eOXgxz/APgcnSMAYZgnLrAvXRq+615PXX0qrtrrFsMbCTBGGAI+h0J6DcMLu5fzHuufc9AQ+N6xhiGUcewAjfA5RilAOhBv2MG/Y4E9LprQ/f5xw4eW9fvcqnP7HeMYRh1JL2PfwAMko4xQA/Wc7OMcegC92tcxhh6Hf8AGCQdY2DLGvbNMsahC9yvfmsG2Ex0jIEta5gXyK0ahy5wvyaxZoBBqNbaqGvI/Px8W1paGnUZwBZz7h/O5UOf+VBefu3lN9au2nZVPvcvP7dm4O0+/1LnATC+qupUa21+rdeMUgCbSj/jEcO+QA6AyWKUAthU+r1IbhIvkANgOARjYFPpdy/eYd0sA4DJIxgDm8p6LpLTBQYgcfEdMOb6vbFG4iI5AC7uUhff9dQxrqq3VtWJqvqLqnq8qn6gqq6pqoeq6onOcWfX+XdX1emqWq6qDw7qFwG2nvXcWMMthQFYj15HKX45yedaa9+TZC7J40nuSvJwa21vkoc7z1NVNyc5kuRdSQ4m+XhVbRt04cDW0O/McPf7jEcA0I/LBuOq+q4kP5jk15Kktfat1tpLSRaSHO+cdjzJoc7jhST3t9Zebq09leR0klsGWzawVaz3xhoukgOgX73sY/zOJOeS/EZVzSU5leSjSa5trT2XJK2156pqtnP+9Um+2PX+s501gCT9zw133465l24xAGPsmUeSr30hufF9ye7x6p32MkqxPcl7kvxqa+3dSf4+nbGJi6g11i64wq+q7qiqpapaOnfuXE/FAptDv3PDZoYBxtwzjyRf+IWV4+XOO3578vl7Vo6XO3+D9RKMzyY521r7Uuf5iawE5eer6rok6Rxf6Dp/d9f7b0jy7Js/tLV2X2ttvrU2PzMzs976gQm0nrlhM8MAG2wYYfdrX0he+1bSXls5fu0Lg635Cl02GLfW/irJM1W1r7N0W5LHkpxMcrSzdjTJg53HJ5McqaqrquqmJHuTjNd/DgAjtZ65YTPDABtoWGH3xvcl296S1LaV443vG3ztV6CXGeMk+ekkn6yqtyR5Msm/ykqofqCqPpLk6SQfTpLW2qNV9UBWwvOrSe5srb028MqBiWZuGGAEep3vXSvsXuz81bD72rcuH3Z335IcPTm2M8Y9BePW2leSrLUR8m0XOf+eJPesvyxg0vR7Qd1q1/jTy582NwywEVa7wKsB9ujJ0YTd3beMXSBe1WvHGOCS5mbn8uQ3n8wrr7/yxtrlbsSxuH8xZ146o1sMsF797PDQTxd4E4XdfgjGwEB0j0asutyIxOrcMADr0E8HOOmvC5xsmrDbj17vfAdwSeu9EQcAb9LrbhD97vCw2gV+/89fPkRvUTrGwMC4oA7gCg1rDnjVFuwC90MwBi7KBXUAAzKM3SDGfIeHSSQYAxflgjqAARhmF1gHeKDMGAMXtZ471LkRB8Cb9DMLbA54pARj4KJcUAdwCb1eJNfv3d5235K872NC8QgYpQAuyQV1AGvoZzzCLPDE0DEGLmm1a1wp3WJgc+u1A5ysb6s0XeCxp2MMXJYL6oBNb9g3y2AiCMawhfS7/doqd6gDNr1+tklLjEdsUoIxbCHr2X4NYKL1un+wm2UQwRi2lO4L6Va5oA7YtFwgR59cfAdbiO3XgC3FBXL0STCGLab7ph26xcBEGtb+wWx5Rilgi1ntGn96+dO6xcDkMR7BEAnGsAXZfg2YWOvZPUIgpkeCMUy49WzBZvs1YGLZP5ghEoxhwtmCDZh4vW6plhiPYKgEY5hwtmADJlq/d5xLjEcwNHalgAlnCzZgovW7pRoMkWAMm4At2ICJZUs1xohgDJvAate4UrrFwGRZnRl+/8/3NkYBQ2TGGDYJW7ABE8vMMGNCMIZNwhZswFjpZ6cJGBOCMYyh9exNDDA21rPTBIwBM8YwhuZm597YZWKVvYmBiWGnCSaUYAxjqHuXiVV2mwAmhp0mmFCCMYwhexMDE81OE0wowRjGlL2JgYm2+5bkfR8TipkogjGMKXsTA8DGsisFjDF7EwPAxhGMYYzZmxgYK/YmZpMTjAGAy7M3MVuAGWMA4PLsTcwWoGMMG8Td7ICJtro38WrH2N7EbEKCMWyQudm5PPnNJ/PK66+8seZudsDEWN2b2Iwxm5hRCtgg7mYHTDx7E7PJCcawQdzNDhhLzzySfOEXVo6wxQnGsIHczQ4YK6s7TXz+npWjcMwWJxjDBnI3O2Cs2GkCzuPiO9hg7mYHjA07TcB5qrU26hoyPz/flpaWRl0GAGw97mbHFlNVp1pr82u9pmMMAFvZ7lsEYugwYwwAABGMAQAgSY/BuKq+VlV/XlVfqaqlzto1VfVQVT3ROe7sOv/uqjpdVctV9cFhFQ8AAIPST8f4n7XWDnQNK9+V5OHW2t4kD3eep6puTnIkybuSHEzy8araNsCaAQBg4K7k4ruFJD/UeXw8yR8k+bnO+v2ttZeTPFVVp5PckuSPr+Bnwdg5fPJwll9cvmB93859OXH7iRFUBABciV47xi3J71XVqaq6o7N2bWvtuSTpHGc769cneabrvWc7a7CpzM3OvXF751XTU9M5MHtgNAUBAFek12D83tbae5J8KMmdVfWDlzi31li7YLPkqrqjqpaqauncuXM9lgHjo/v2zqvc5hkAJldPwbi19mzn+EKSz2ZlNOL5qrouSTrHFzqnn02yu+vtNyR5do3PvK+1Nt9am5+ZmVn/bwAjsnp759Wu8fTUtNs8A8AEu2wwrqpvr6rvXH2c5J8n+WqSk0mOdk47muTBzuOTSY5U1VVVdVOSvUkeGXThMA66u8a6xQAw2XrpGF+b5I+q6k+zEnD/R2vtc0nuTfKBqnoiyQc6z9NaezTJA0keS/K5JHe21l4bRvEwaqtd40rpFgPAhLvsrhSttSeTzK2x/o0kt13kPfckueeKq4MJsLh/MWdeOqNbDAAT7kq2awOy0jU+dvDYqMsAAK6QW0IDAEAEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJLGPMZzn8MnDWX5x+YL1fTv35cTtJ0ZQEQCwUXSMocvc7Fymp6bPW5uems6B2QOjKQgA2DCCMXRZ3L+YqTr/b4upmnK7ZwDYAgRj6DKzYyYLexbe6BpPT03n0J5D2XX1rhFXBgAMm2AMb9LdNdYtBoCtQzCGN1ntGldKtxgAthC7UsAaFvcv5sxLZ3SLAWALEYxhDTM7ZnLs4LFRlwEAbCCjFAAAEMEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkrglNFvA4ZOHs/zi8gXr+3buy4nbT4ygIgBgHOkYs+nNzc5lemr6vLXpqekcmD0wmoIAgLEkGLPpLe5fzFSd/0d9qqayOLc4oooAgHEkGLPpzeyYycKehTe6xtNT0zm051B2Xb1rxJUBAONEMGZL6O4a6xYDAGsRjNkSVrvGldItBgDWZFcKtozF/Ys589IZ3WIAYE2CMVvGzI6ZHDt4bNRlAABjyigFAABEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACRJto+6AOjX4ZOHs/zi8gXr+3buy4nbT4ygIgBgM9AxZuLMzc5lemr6vLXpqekcmD0wmoIAgE1BMGbiLO5fzFSd/0d3qqayOLc4oooAgM2g52BcVduq6k+q6nc6z6+pqoeq6onOcWfXuXdX1emqWq6qDw6jcLaumR0zWdiz8EbXeHpqOof2HMquq3eNuDIAYJL10zH+aJLHu57fleTh1treJA93nqeqbk5yJMm7khxM8vGq2jaYcmFFd9dYtxgAGISegnFV3ZDkR5J8omt5IcnxzuPjSQ51rd/fWnu5tfZUktNJbhlItdCx2jWulG4xADAQvXaMfynJzyZ5vWvt2tbac0nSOc521q9P8kzXeWc7azBQi/sX855r36NbDAAMxGWDcVX9aJIXWmunevzMWmOtrfG5d1TVUlUtnTt3rsePhn80s2Mmxw4e0y0GAAail47xe5PcXlVfS3J/kvdX1W8meb6qrkuSzvGFzvlnk+zuev8NSZ5984e21u5rrc231uZnZmau4FcAAIArd9lg3Fq7u7V2Q2vtxqxcVPf51tpPJDmZ5GjntKNJHuw8PpnkSFVdVVU3Jdmb5JGBVw4AAAN0JXe+uzfJA1X1kSRPJ/lwkrTWHq2qB5I8luTVJHe21l674koBAGCIqrULxn833Pz8fFtaWhp1GQAAbHJVdaq1Nr/Wa+58BwAAEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSXNk+xjAwh08ezvKLyxes79u5LyduPzGCigCArUbHmLEwNzuX6anp89amp6ZzYPbAaAoCALYcwZixsLh/MVN1/h/HqZrK4tziiCoCALYawZixMLNjJgt7Ft7oGk9PTefQnkPZdfWuEVcGAGwVgjFjo7trrFsMAGw0wZixsdo1rpRuMQCw4exKwVhZ3L+YMy+d0S0GADacYMxYmdkxk2MHj426DABgCzJKAQAAEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJOkhGFfVt1XVI1X1p1X1aFX9p876NVX1UFU90Tnu7HrP3VV1uqqWq+qDw/wFAABgEHrpGL+c5P2ttbkkB5IcrKpbk9yV5OHW2t4kD3eep6puTnIkybuSHEzy8araNoTaAQBgYC4bjNuK/9d5Ot35qyVZSHK8s348yaHO44Uk97fWXm6tPZXkdJJbBlk0AAAMWk8zxlW1raq+kuSFJA+11r6U5NrW2nNJ0jnOdk6/PskzXW8/21kDAICx1VMwbq291lo7kOSGJLdU1fdd4vRa6yMuOKnqjqpaqqqlc+fO9VQsAAAMS1+7UrTWXkryB1mZHX6+qq5Lks7xhc5pZ5Ps7nrbDUmeXeOz7mutzbfW5mdmZvqvHAAABqiXXSlmquqtncdXJ/nhJH+R5GSSo53TjiZ5sPP4ZJIjVXVVVd2UZG+SRwZcNwAADNT2Hs65Lsnxzs4SU0keaK39TlX9cZIHquojSZ5O8uEkaa09WlUPJHksyatJ7mytvTac8gEAYDCqtQvGfzfc/Px8W1paGnUZAABsclV1qrU2v9Zr7nwHAAARjAEAIIlgDAAASQRjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkPQTjqtpdVb9fVY9X1aNV9dHO+jVV9VBVPdE57ux6z91Vdbqqlqvqg8P8BQAAYBB66Ri/muRjrbXvTXJrkjur6uYkdyV5uLW2N8nDnefpvHYkybuSHEzy8araNoziAQBgUC4bjFtrz7XWvtx5/HdJHk9yfZKFJMc7px1PcqjzeCHJ/a21l1trTyU5neSWAdcNAAAD1deMcVXdmOTdSb6U5NrW2nPJSnhOMts57fokz3S97WxnDQAAxlbPwbiqviPJbyf5mdba317q1DXW2hqfd0dVLVXV0rlz53otAwAAhqKnYFxV01kJxZ9srX2ms/x8VV3Xef26JC901s8m2d319huSPPvmz2yt3ddam2+tzc/MzKy3fgAAGIhedqWoJL+W5PHW2i92vXQyydHO46NJHuxaP1JVV1XVTUn2JnlkcCUDAMDgbe/hnPcm+ckkf15VX+ms/Yck9yZ5oKo+kuTpJB9Oktbao1X1QJLHsrKjxZ2ttdcGXTgAAAzSZYNxa+2PsvbccJLcdpH33JPkniuoCwAANpQ73wEAQARjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIkmwfdQFsXodPHs7yi8sXrO/buS8nbj8xgooAAC5Ox5ihmZudy/TU9Hlr01PTOTB7YDQFAQBcgmDM0CzuX8xUnf9HbKqmsji3OKKKAAAuTjBmaGZ2zGRhz8IbXePpqekc2nMou67eNeLKAAAuJBgzVN1dY91iAGCcCcYM1WrXuFK6xQDAWLMrBUO3uH8xZ146o1sMAIw1wZihm9kxk2MHj426DACASzJKAQAAEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQJIegnFV/XpVvVBVX+1au6aqHqqqJzrHnV2v3V1Vp6tquao+OKzCAQBgkHrpGB9LcvBNa3clebi1tjfJw53nqaqbkxxJ8q7Oez5eVdsGVi0AAAzJZYNxa+0Pk/zNm5YXkhzvPD6e5FDX+v2ttZdba08lOZ3klsGUCgAAw7PeGeNrW2vPJUnnONtZvz7JM13nne2sAQDAWBv0xXe1xlpb88SqO6pqqaqWzp07N+AyAACgP+sNxs9X1XVJ0jm+0Fk/m2R313k3JHl2rQ9ord3XWptvrc3PzMysswwAABiM9Qbjk0mOdh4fTfJg1/qRqrqqqm5KsjfJI1dWIgAADN/2y51QVZ9K8kNJdlXV2ST/Mcm9SR6oqo8keTrJh5OktfZoVT2Q5LEkrya5s7X22pBqBwCAgblsMG6t/dhFXrrtIuffk+SeKykKAAA2mjvfAQBABGMAAEgiGAMAQBLBGAAAkgjGAACQRDAGAIAkgjEAACQRjAEAIIlgDAAASQRjAABIIhgDAEASwRgAAJIIxgAAkEQwBgCAJIIxAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQJNk+6gKYLIdPHs7yi8sXrO/buS8nbj8xgooAAAZDx5i+zM3OZXpq+ry16anpHJg9MJqCAAAGRDCmL4v7FzNV5/+xmaqpLM4tjqgiAIDBEIzpy8yOmSzsWXijazw9NZ1Dew5l19W7RlwZAMCVEYzpW3fXWLcYANgsBGP6tto1rpRuMQCwadiVgnVZ3L+YMy+d0S0GADYNwZh1mdkxk2MHj426DACAgTFKAQAAEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJHHnuy3v8MnDWX5x+YL1fTv35cTtJ0ZQEQDAaOgYb3Fzs3OZnpo+b216ajoHZg+MpiAAgBEZWjCuqoNVtVxVp6vqrmH9HK7M4v7FTNX5fwymaiqLc4sjqggAYDSGEoyraluSX0nyoSQ3J/mxqrp5GD+LKzOzYyYLexbe6BpPT03n0J5D2XX1rhFXBgCwsYbVMb4lyenW2pOttW8luT/JwpB+Fleou2usWwwAbFXDCsbXJ3mm6/nZzhpjaLVrXCndYgBgyxrWrhS1xlo774SqO5LckSRvf/vbh1QGvVrcv5gzL53RLQYAtqxhdYzPJtnd9fyGJM92n9Bau6+1Nt9am5+ZmRlSGfRqZsdMjh08plsMAGxZwwrG/yfJ3qq6qarekuRIkpND+lkAAHDFhjJK0Vp7tar+bZL/lWRbkl9vrT06jJ8FAACDMLQ737XWfjfJ7w7r8wEAYJDc+Q4AACIYAwBAEsEYAACSCMYAAJBEMAYAgCSCMQAAJBGMAQAgiWAMAABJBGMAAEgiGAMAQBLBGAAAkgjGAACQJKnW2qhrSFWdS/L1UdexSe1K8tejLoIr4jucfL7Dyec7nHy+w8k2yO/vHa21mbVeGItgzPBU1VJrbX7UdbB+vsPJ5zucfL7Dyec7nGwb9f0ZpQAAgAjGAACQRDDeCu4bdQFcMd/h5PMdTj7f4eTzHU62Dfn+zBgDAEB0jAEAIIlgvGlV1e6q+v2qeryqHq2qj466JvpXVduq6k+q6ndGXQv9q6q3VtWJqvqLzt+LPzDqmuhPVf37zj9Dv1pVn6qqbxt1TVxaVf16Vb1QVV/tWrumqh6qqic6x52jrJFLu8h3+F86/yz9s6r6bFW9dRg/WzDevF5N8rHW2vcmuTXJnVV184hron8fTfL4qItg3X45yedaa9+TZC6+y4lSVdcn+XdJ5ltr35dkW5Ijo62KHhxLcvBNa3clebi1tjfJw53njK9jufA7fCjJ97XW9if5v0nuHsYPFow3qdbac621L3ce/11W/oV8/Wiroh9VdUOSH0nyiVHXQv+q6ruS/GCSX0uS1tq3WmsvjbQo1mN7kquranuSHUmeHXE9XEZr7Q+T/M2blheSHO88Pp7k0EbWRH/W+g5ba7/XWnu18/SLSW4Yxs8WjLeAqroxybuTfGnEpdCfX0rys0leH3EdrM87k5xL8hudcZhPVNW3j7ooetda+8sk/zXJ00meS/LN1trvjbYq1una1tpzyUrjKMnsiOvhyvzrJP9zGB8sGG9yVfUdSX47yc+01v521PXQm6r60SQvtNZOjboW1m17kvck+dXW2ruT/H3879uJ0plDXUhyU5LvTvLtVfUTo60Ktraq+vmsjIt+chifLxhvYlU1nZVQ/MnW2mdGXQ99eW+S26vqa0nuT/L+qvrN0ZZEn84mOdtaW/0/NSeyEpSZHD+c5KnW2rnW2itJPpPkn464Jtbn+aq6Lkk6xxdGXA/rUFVHk/xokh9vQ9pvWDDepKqqsjLb+Hhr7RdHXQ/9aa3d3Vq7obV2Y1Yu9vl8a02naoK01v4qyTNVta+zdFuSx0ZYEv17OsmtVbWj88/U2+ICykl1MsnRzuOjSR4cYS2sQ1UdTPJzSW5vrf3DsH6OYLx5vTfJT2al0/iVzl//YtRFwRbz00k+WVV/luRAkv882nLoR6fbfyLJl5P8eVb+nenuaWOuqj6V5I+T7Kuqs1X1kST3JvlAVT2R5AOd54ypi3yH/y3JdyZ5qJNp/vtQfrY73wEAgI4xAAAkEYwBACCJYAwAAEkEYwAASCIYAwBAEsEYAACSCMYAAJBEMAYAgCTJ/wdtFpa7uKQnvQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# (s == start; e == end)\n",
    "# 1 or 2 == mode number\n",
    "# d == dipole; m == monopole\n",
    "# a, b - part\n",
    "# ----------------------------------- #\n",
    "fig, ax = plt.subplots()\n",
    "fig.set_size_inches(12,8)\n",
    "# first monopole\n",
    "s1m = 2; e1m = 26\n",
    "ax.plot(data_freq[s1m:e1m], mu_curve[s1m:e1m, 1, 1],'.', color='tab:blue')\n",
    "# second monopole\n",
    "s2m_a = 42; e2m_a = 57\n",
    "ax.plot(data_freq[s2m_a:e2m_a], mu_curve[s2m_a:e2m_a, 1, 1],'.', color='tab:orange')\n",
    "# s2m_b = 13; e2m_b = 33\n",
    "# ax.plot(data_freq[s2m_b:e2m_b], mu_curve[s2m_b:e2m_b, 0, 1],'x', color='tab:orange')\n",
    "# first dipole\n",
    "s1d = 8; e1d = 35\n",
    "ax.plot(data_freq[s1d:e1d], mu_curve[s1d:e1d, 1, 0],'v', color='tab:green')\n",
    "# second dipole\n",
    "# s2d_a = 14; e2d_a = 37\n",
    "# ax.plot(data_freq[s2d_a:e2d_a], mu_curve[s2d_a:e2d_a, 0, 0], 'o', color='tab:red')\n",
    "# s2d_b = 46; e2d_b = 67\n",
    "# ax.plot(data_freq[s2d_b:e2d_b], mu_curve[s2d_b:e2d_b, 1, 0],'v', color='tab:red')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "dataset construction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "66\n"
     ]
    }
   ],
   "source": [
    "# 3 parts\n",
    "input_shape = \\\n",
    "e1d - s1d + \\\n",
    "e1m - s1m + \\\n",
    "e2m_a - s2m_a\n",
    "print(input_shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # 6 parts\n",
    "# input_shape = \\\n",
    "# e1d - s1d + \\\n",
    "# e1m - s1m + \\\n",
    "# e2d_a - s2d_a + \\\n",
    "# e2d_b - s2d_b + \\\n",
    "# e2m_a - s2m_a + \\\n",
    "# e2m_b - s2m_b\n",
    "# print(input_shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "# flatten\n",
    "data_X_new = np.zeros([data_X.shape[0], input_shape])\n",
    "for i in range(data_X.shape[0]):\n",
    "    # first dipole\n",
    "    flat_arr_temp = data_X[i,s1d:e1d,1,0]\n",
    "    # second dipole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_X[i,s2d_a:e2d_a,0,0])\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_X[i,s2d_b:e2d_b,1,0])\n",
    "    # first monopole\n",
    "    flat_arr_temp = np.append(flat_arr_temp, data_X[i,s1m:e1m,1,1])\n",
    "    # second monopole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_X[i,s2m_b:e2m_b,0,1])\n",
    "    flat_arr_temp = np.append(flat_arr_temp, data_X[i,s2m_a:e2m_a,1,1])\n",
    "    data_X_new[i, :] = flat_arr_temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsYAAAHSCAYAAADvxw2lAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAoh0lEQVR4nO3dcYyc9X3n8fd3dtbQhFZZiEGc2dhYuHXWycV0l43dnKtiEKGXNERRkUDpYZJYqCfulFQ9WqCRqp4giorau1OUVAK3xKgERNNyOKTXA9tE9Uk49m4xSrwswrd4sYFi19lemkawOzvf+2MfkwWvvc/aMzu7O++XhGbmmd/MfPcntP7sb77P74nMRJIkSWp3lVYXIEmSJC0EBmNJkiQJg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCYBqqwsAeP/735+rVq1qdRmSJEla4gYHB/8pM5fP9NyCCMarVq1iYGCg1WVIkiRpiYuI0dM9ZyuFJEmShMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBJYJxRPxSRByY9t+PI+JLEXFhRDwdES8Vt13TXnNXRByKiBcj4uPN/REkSZKkczdrMM7MFzNzfWauB3qBnwKPA3cCuzJzDbCreExE9AA3AeuA64FvRERHc8qXJEmSGmOurRTXAP83M0eBG4DtxfHtwKeL+zcAj2bmW5n5MnAI6G9ArZLUfEf2wZ4/mbqVJLWV6hzH3wQ8Uty/JDNfB8jM1yPi4uL4CmDvtNccLY5J0sJ2ZB/1b/4GTI5DxzIqt34Huv27XpLaRekV44hYBnwK+KvZhs5wLGd4v9siYiAiBo4fP162DElqmlcPPEW9Nk6FOvXaOK8eeKrVJUmS5tFcWil+HfiHzHyjePxGRFwKUNweK44fBbqnve4y4LV3v1lm3p+ZfZnZt3z58rlXLkkN9uxkDxNUqWWFCao8O9nT6pIkSfNoLsH4Zn7WRgGwA9hS3N8CPDHt+E0RcV5EXA6sAWzWk7TgXX7l1Xyu/mX+2+SNfK7+ZS6/8upWlyRJmkeReUqXw6mDIt4DHAFWZ+b/K45dBDwGfAB4BbgxM39UPPcHwOeBGvClzPxfZ3r/vr6+HBgYOJefQ5IaYnB0jL0jJ9iw+iJ6V3bN/gJJ0qISEYOZ2Tfjc2WCcbMZjCVJkjQfzhSMvfKdJEmShMFYkppmcHSMrz9ziMHRsVaXIkkqYa77GEuSShgcHeOz2/YyXquzrFrh4a0bzrpnuXTf85F9cHgPrNp0zvsvD+/fydjQbrp6NrP2qmvP6b0kabEwGEtSE+wdOcF4rU49YaJWZ+/IibMKxoOjY9y37SF68yD37V7HHVtvmfl9juyD7Z96++IkbNlx1uF4eP9OVj55M1dQY2LkAYZ5xHAsqS3YSiFJTbBh9UX0Vw9xe/UJrqoeYsPqi2YcN1u7xcvPPcODlXv4nY6/4sHKPbz83DMzf+DhPeTkW5CT5OT41MrxWXwewNjQbjqpUY06ndQYG9p95h/2TBbxJbZthZniPKiduGIsSU3QW3mJby37SrGC+wSVykbgnSu4ZdotNnYMvR1SyRobO4aAz5zyecPnf4SV9Sqd1JjIDkbP/whr3zWmbHtHV89mJkYegKwxQZWuns1nNwlzucR2A9tAGqH0Sn1Ji7U1pZEtQdJiYDCWpGY4vIdKfQKoQ31iKvS9K/CVabdYsf466ge+Rn1ygkq1kxXrr5vx43b9ZBW7J+7mo/EC+/KDXP2TVacE47LtHWuvupZhHjnnIPfqgae4pDZONerUauO8fuApVswUehvYBlK2H3u2cSdX6jupMcHjfPe5bnpXnvoHSRmNbk0p9TM26A+NvSMnWDc5zEcrL7Bv8oPsHVljMNaSZjCWpGZYtWkq5J0Me6s2nTJkw+qLWFatMFGr01mtzNxu0d0/tdI6S8jZsPoivtaxlgO1X6SzWuGuGd7rZHtHbx5kMNaxYfWvnLb8tVddC+e4svnsZA+foPr2yvOzkz385kwDD++ZmqecnLqd4Y8ImH3VtezqZiNX6ssYG9rNFdPea2xo94xzWybwllrJLrlSX+bzrrngMJ/v/ErxB0KV0Qs+DFwx5zmQFguDsSTNUakVu+7+qZXPMwTa3pVdPLx1Q7n3mmXVr8x7lWnvmJNZViUvv/JqPjf45beD+B2nu8T2qk3UK50wCVQ6qczwR0SZVdeyq5tlxpVdqS+jTGtK2daNMivZZVbqy37e2jefJys1Iut0xCRr33weODXUL9ZWEendDMaSNAdz6j0tE2grL9Fb3QOVTZxTSGUqHJ/xa+4S7R2llViV7F3ZxR1bb2HvyAnuONMqaH0N943fPRWgJ9dxR30Nve8aU2bVtezqZqlxJVfqy7QslGlNKdu6UWYlu8xKfelWkVWbiI7zYHKcOM03H+5ioqXEYCypbTRiVauRvaeN7K0tpUR7R1ll+4dnDetMreDuq13B3ryCjmDG3ucyq65lVzfLjpv1D5s5nFw4W2tK2daNMivZZVbqS7eKlPjmo2yriLQYGIwltYVGrWo1sve0bG9tw5QIOWWV7h8uoUyvdakTAkusbs5p3CxKn1xYQunWjRIr2WVW6ufUKjLLHwgN28VEWgAMxpLaQqNWtRrZe9rIFdzSSrR3lFG6f7iEsr3Ws54QWDb4N+gPhEb+cVC6daMYW6bn/Iwr9XP5vFk0ahcTaSGIzGx1DfT19eXAwECry5C0hJ1cMX67r/ST59AH2cg9dxfY/r1zUfpS1UvU9H7zwTj3vY4lzY+IGMzMvhmfMxhLaheeOa9Ga/c/DqTFyGAsSZIkceZgXJnvYiRJkqSFyGAsSZIkYTCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEtaAgZHx/j6M4cYHB1rdSmSpEWs2uoCJOlcDI6O8dltexmv1VlWrfDw1g30ruxqdVmSpEXIFWNJLTO8fyfPbr+b4f07TztmttXgvSMnWDc5zG9XnuBDk8PsHTnRrHIlSUucK8aSWmJ4/05WPnkzV1BjYuQBhnmEtVdd+44xg6Nj3LftIXrzIPftXscdW285ZTX4mgsO8/nOr9BJjQmqjF7wYeCKefxJJElLhSvGklpibGg3ndSoRp1OaowN7T5lzMvPPcODlXv4nY6/4sHKPbz83DOnjFn75vOcX5l6n/Mrk6x98/n5KF+StAQZjCU1XJmT4bp6NjNBlVpWmKBKV8/mU8Zs7Bh6R3je2DF06hut2kR0nAfRQXQsg1WbGvmjSJLaiK0Ukhqq7Mlwa6+6lmEeYWxoN109m09powBYsf466ge+Rn1ygkq1kxXrrzv1A7v7YcsOOLxnKhR39zfjx5IktQGDsaSGOnky3EcrL7Bv8oPsHVlz2l0i1l51LcwQiN/W3U/l1u/MHnq7+w3EkqRzZjCWNCfD+3eecZW34SfDGXolSfPEYCyptDI7Sax983myUiOyTkecPBnuDKvCkiQtEJ58J6m0MjtJeDKcJGmxcsVYEjB7iwQUO0mMPABZO+1OEp4MJ0larAzGkkq1SEC5nSQA+4IlSYuSwVgSY0O7uaJokSCLFonThN5Zd5KQJGmRssdYWuTKXExjtjFlLrYhSdJSZzCWFrHB0THu2/YQP931x9y37aEZg2+ZMWuvupbRTz7C/tX/kdFPztxGIUnSUmcrhbSIvfzcMzxYuafYM/hxvvtcN70rPzPnMWCLhCRJrhhLi9jGjqF3bJ+2sWPorMZIkqSSwTgi3hcR346I4Yh4ISI2RsSFEfF0RLxU3HZNG39XRByKiBcj4uPNK19ausr0Dq9Yfx2V6jLqdFCpLmPF+uvOaowkSYLIzNkHRWwH9mTmtohYBrwHuBv4UWZ+NSLuBLoy8/cjogd4BOgH/g2wE/jFzJw83fv39fXlwMBAA34caWk42RfcmwcZjHXcsfUWeld2zTz4yL7Z9wwuM0aSpDYQEYOZ2TfTc7P2GEfELwC/CtwKkJnjwHhE3AD8WjFsO/A94PeBG4BHM/Mt4OWIOMRUSH72nH4KqY2U7QsGyu0Z7L7CkiTNqkwrxWrgOPBgRDwXEdsi4r3AJZn5OkBxe3ExfgVwZNrrjxbHJJVkX7AkSfOvTDCuAr8M/FlmXgn8K3DnGcbHDMdO6deIiNsiYiAiBo4fP16qWKld2BcsSdL8KxOMjwJHM/P7xeNvMxWU34iISwGK22PTxndPe/1lwGvvftPMvD8z+zKzb/ny5Wdbv7QozXpiXXc/lVu/Q+WaP6By63dsg5AkaR7MGowz8x+BIxHxS8Wha4AhYAewpTi2BXiiuL8DuCkizouIy4E1wL6GVi0tYmUuuAFMheFNv2soliQtLUf2wZ4/mbpdYMpe4OM/Aw8XO1KMAJ9jKlQ/FhFfAF4BbgTIzIMR8RhT4bkG3H6mHSmkpWZ4/07GhnbT1bN5xivIzenEOkmSlpIj+2D7p2ByHDqWwZYdC2oBqFQwzswDwEzbWlxzmvH3AveefVnS4jS8fycrn7yZK6gxMfIAw5x6eeXpJ9aRJ0+sMxhLktrA4T3k5FtE1snJceLwngUVjL3yndRAY0O737GbxNjQ7lPGeGKdJGkxKnPhqeH9O3l2+90M79858/Pnf4Q361VqWeHNegfD53+kWeWelbKtFFLbGxwdY+/ICTasvui0F9vo6tnMxMgDkDUmqNLVs/nUQcWJdV5wQ5K0WEy/8NR9u2e+8FSZb013/WQVuyfu5qPxAvvyg1z9k1Wsnc8fZBYGY6mEwdExPrttL+O1OsuqFR7eumHGcLz2qmsZ5pEz9hgDXnBDkrRglFn4KXN+zNjQbq6Y1io4NrQb3vXv4IbVF/G1jrUcqP0indUKd62+qGk/19kwGEsl7B05wbrJYT5aeYF9kx9k78ia0/7yWHvVtaf8IpAkaSEqsxIM5c6PKfOtae/KLh7eumHWIN4qBmOphGsuOMznO79S/KVcZfSCDwNXtLosSZLOaLbV4LI7Ja1Yfx31A1+jPjlBpdo54/kxZb817V3ZteAC8UkGY6mEtW8+T1ZqRNbpiEnWvvk84KqwJGnhKtMGWHqnpJLnxyz2b03dlUIqY9UmouM8iA6iY9nULwVJklpotl0iTrYB/nblCT40OczekROnjJnTTkltcOEpV4wlSpx40N0/tQm5O0lIkhaAMr3BpdoA3SnpHQzGantld5xwJwlJ0kJRpje4dBug/769zVYKtb0yXzVJkjRfylxIY3pvcCcne4PfxTbAOXPFWG3PHSckSQtF2e3TyuwSYRvg3BmM1fbccUKStFCU3T6tdG+wbRJzYjCWTn7VNDnuV02SpKYa3r/zjPv8lt4+DQy9TWAw1pI32y8hv2qSJM2H4f07WfnkzVxBjYmRBxjmkVP+XSrVIqGmMRhrSSvzSwjwr25JUtONDe3mimmrwWNDu0+9GIbbp7WUu1JoSRsb2v2Os3bHhna3uiRJUpvq6tnMBFVqWWGCKl09m2ce2AYX0lioXDHWktbVs5mJkQcga2f+JSRJ0jmarXVv7VXXMswjZ27vU0tFZra6Bvr6+nJgYKDVZWiJmrXHWJKkM5j16qj8rHXv7a0/P3ma1j21XEQMZmbfTM+5Yqwlb+1V157awyVJUgll9xUu1T+sBc9grMXryD5PTpAkNVXZfYVt3VsaDMZanI7so/7N34DJcehYNnUGr+FYkjRHjdpX2P7hpcFgrEXp1QNPcUltnGrUqdXGef3AU6wwGEuS5qDR+wrburf4uV2bFqVnJ3veseXNs5M9rS5JkrTIlNrSs9hXuHLNH/jtZBtwxViL0uVXXs3nBr9Mbx5kMNZxx5VXt7okSdIiU7ov2ItAtQ23a9OiVWb7HEmSzsQtPdvPmbZrMxhLkiSpbbiPsSRJaiuuBOtsGIy1MLlHsSTpLJXZbUKaicFYC497FEuSzoFXodPZcrs2LTivHniKem2cCnXqtXFePfBUq0uSJC0gw/t38uz2uxnev3PG57t6Nr9jS0+vQqeyXDHWgvPsZA+foPr29jnPTvbwm60uSpK0IJRpk/AqdDpbBmMtOO5RLEnta7atOMu2SXgVOp0NWym04PSu7OKOrbfwnmt+jzu23uIexZLUJgZHx7hv20P8dNcfc9+2hxgcHTtljG0SaiZXjLUg9a7sMhBLUpt5+blneLByD53UmOBxvvtcN70rP/OOMbZJqJkMxpIkaUHY2DFE57Q2iY0dQ8BnThlnm4SaxVYKzavZziSWJLWvFeuvo1JdRp0OKtVlrFh/XatLUptxxVjzxg3XJUln1N0/tXe9F3hSixiMNW/ccF2S2tdsu028rbvfQKyWsZVC88YziSWpPZXZbUJaCFwx1rzxTGJJak9ldpuQFgKDseaVZxJLUvspu9uE1Gq2UkiSpKZytwktFq4YS5Kk5nK3CS0SBmM1zPD+nfYPS5Jm5m4TWgQMxmoI9yiWJEmLnT3Gaoixod1vn1jRSbFHsSRJ0iJiMFZDuEexJEla7GylUEO4R7EkSVrsDMZqGPcoliRJi1mpVoqIOBwRP4iIAxExUBy7MCKejoiXituuaePviohDEfFiRHy8WcVLkiRJjTKXHuOrM3N9ZvYVj+8EdmXmGmBX8ZiI6AFuAtYB1wPfiIiOBtYsSZIkNdy5nHx3A7C9uL8d+PS0449m5luZ+TJwCHDjQkmSJC1oZYNxAk9FxGBE3FYcuyQzXwcobi8ujq8Ajkx77dHimCRJkrRglT357mOZ+VpEXAw8HRHDZxgbMxzLUwZNBezbAD7wgQ+ULEOSJElqjlIrxpn5WnF7DHicqdaINyLiUoDi9lgx/CjQPe3llwGvzfCe92dmX2b2LV++/Ox/AkmSJKkBZg3GEfHeiPj5k/eB64AfAjuALcWwLcATxf0dwE0RcV5EXA6sAfY1unBJkiSpkcq0UlwCPB4RJ8d/KzP/LiL2A49FxBeAV4AbATLzYEQ8BgwBNeD2zJxsSvWaN4OjY+wdOcGG1RfRu7Jr9hdIkiQtMpF5SvvvvOvr68uBgYFWl6HTGBwd475tD9GbBxmMddyx9RbDsSRJWpQiYnDa9sPv4JXvNKuXn3uGByv30EmNCR7nu89107vyM60uS5IkqaHOZR9jtYmNHUN0UqMadTqpsbFjqNUlSZIkNZzBWLNasf46KtVl1OmgUl3GivXXtbokSZKkhrOVQrPr7qdy63fg8B5YtQm6vZChJElaegzGKqe730AsSZKWNFspJEmSJAzGkiRJEmAwliRJkgCDsSRJkgQYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEmAwbnuDo2N8/ZlDDI6OtboUSZKkljIYt7HB0THu2/YQP931x9y37SHDsSRJamvVVheg1nn5uWd4sHIPndSY4HG++1w3vSs/0+qyJEmSWsIV4za2sWOITmpUo04nNTZ2DLW6JEmSpJYxGLexFeuvo1JdRp0OKtVlrFh/XatLkiRJahlbKdpZdz+VW78Dh/fAqk3Q3d/qiiRJklrGYNzuuvsNxJIkSdhKIUmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMZL2uDoGF9/5hCDo2OtLkWSJGnBq7a6ADXH4OgYn922l/FanWXVCg9v3UDvyq5WlyVJkrRguWK8RO0dOcG6yWF+u/IEH5ocZu/IiVaXJEmStKC5YrxEXXPBYT7f+RU6qTFBldELPgxc0eqyJEmSFixXjJeotW8+z/mVGtWoc35lkrVvPt/qkiRJkhY0g/FStWoT0XEeRAfRsQxWbWp1RZIkSQuarRRLVXc/bNkBh/dMheLu/lZXJEmStKAZjJey7n4DsSRJUkm2UkiSJEkYjCVJkiTAYCxJkiQBBmNJkiQJMBhLkiRJgMFYkiRJAgzGkiRJEjCHYBwRHRHxXEQ8WTy+MCKejoiXituuaWPviohDEfFiRHy8GYVLkiRJjTSXFeMvAi9Me3wnsCsz1wC7isdERA9wE7AOuB74RkR0NKZcSZIkqTlKBeOIuAz4BLBt2uEbgO3F/e3Ap6cdfzQz38rMl4FDgJdfkyRJ0oJWdsX4vwO/B9SnHbskM18HKG4vLo6vAI5MG3e0OCZJkiQtWLMG44j4JHAsMwdLvmfMcCxneN/bImIgIgaOHz9e8q0lSZKk5iizYvwx4FMRcRh4FNgcEX8JvBERlwIUt8eK8UeB7mmvvwx47d1vmpn3Z2ZfZvYtX778HH4ESZIk6dzNGowz867MvCwzVzF1Ut3uzPwtYAewpRi2BXiiuL8DuCkizouIy4E1wL6GV97mBkfH+PozhxgcHWt1KZIkSUtC9Rxe+1XgsYj4AvAKcCNAZh6MiMeAIaAG3J6Zk+dcqd42ODrGZ7ftZbxWZ1m1wsNbN9C7smv2F0qSJOm05hSMM/N7wPeK+yeAa04z7l7g3nOsTaexd+QE47U69YSJWp29IycMxpIkSefoXFaM1SIbVl9Ef/UQvXmQwVjHhtW/0uqSJEmSFj2D8SLUW3mJby37CkyOQ8cTVCobcatoSZKkczOXK99poTi8h0p9ggp1KvUJOLyn1RVJkiQtegbjxWjVJuhYBtExdbtqU6srkiRJWvRspViMuvthy46pleJVm6YeS5Ik6ZwYjBer7n4DsSRJUgPZSiFJkiRhMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYLziDo2N8/ZlDDI6OtboUSZKktlJtdQH6mcHRMT67bS/jtTrLqhUe3rqB3pVdrS5LkiSpLbhivIDsHTnBeK1OPWGiVmfvyIlWlyRJktQ2XDFeQDasvoj+6iF68yCDsY4Nq3+l1SVJkiS1DYPxAtJbeYlvLfsKTI5DxxNUKhuB/laXJUmS1BZspVhIDu+hUp+gQp1KfQIO72l1RZIkSW3DYLyQrNoEHcsgOqZuV21qdUWSJEltw1aKhaS7H7bsmFopXrVp6rEkSZLmhcF4oenuNxBLkiS1gK0UkiRJEgZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxvPryD7Y8ydTt5IkSVpQqq0uoG0c2Uf9m78Bk+PQsYzKrd+B7v5WVyVJkqSCK8bz5NUDT1GvjVOhTr02zqsHnmp1SZIkSZrGYDxPnp3sYYIqtawwQZVnJ3taXZIkSZKmsZVinlx+5dV8bvDL9OZBBmMdd1x5datLkiRJ0jQG43nSu7KLO7bewt6RE9yx+iJ6V3a1uiRJkiRNYzCeR70ruwzEkiRJC5Q9xpIkSRIlgnFEnB8R+yLi+Yg4GBF/VBy/MCKejoiXituuaa+5KyIORcSLEfHxZv4AkiRJUiOUWTF+C9icmR8B1gPXR8QG4E5gV2auAXYVj4mIHuAmYB1wPfCNiOhoQu2SJElSw8wajHPKT4qHncV/CdwAbC+Obwc+Xdy/AXg0M9/KzJeBQ4BXspAkSdKCVqrHOCI6IuIAcAx4OjO/D1ySma8DFLcXF8NXAEemvfxocUySJElasEoF48yczMz1wGVAf0R86AzDY6a3OGVQxG0RMRARA8ePHy9VrCRJktQsc9qVIjP/GfgeU73Db0TEpQDF7bFi2FGge9rLLgNem+G97s/MvszsW758+dwrlyRJkhqozK4UyyPifcX9nwOuBYaBHcCWYtgW4Ini/g7gpog4LyIuB9YA+xpctyRJktRQZS7wcSmwvdhZogI8lplPRsSzwGMR8QXgFeBGgMw8GBGPAUNADbg9MyebU74kSZLUGJF5SvvvvOvr68uBgYFWlyFJkqQlLiIGM7Nvpue88p0kSZKEwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJElAiWAcEd0R8UxEvBARByPii8XxCyPi6Yh4qbjtmvaauyLiUES8GBEfb+YPIEmSJDVCmRXjGvC7mflBYANwe0T0AHcCuzJzDbCreEzx3E3AOuB64BsR0dGM4iVJkqRGmTUYZ+brmfkPxf1/AV4AVgA3ANuLYduBTxf3bwAezcy3MvNl4BDQ3+C6JUmSpIaaU49xRKwCrgS+D1ySma/DVHgGLi6GrQCOTHvZ0eKYJEmStGCVDsYRcQHw18CXMvPHZxo6w7Gc4f1ui4iBiBg4fvx42TIkSZKkpigVjCOik6lQ/HBm/k1x+I2IuLR4/lLgWHH8KNA97eWXAa+9+z0z8/7M7MvMvuXLl59t/ZIkSVJDlNmVIoA/B17IzD+d9tQOYEtxfwvwxLTjN0XEeRFxObAG2Ne4kiVJkqTGq5YY8zHgPwA/iIgDxbG7ga8Cj0XEF4BXgBsBMvNgRDwGDDG1o8XtmTnZ6MIlSZKkRpo1GGfm/2HmvmGAa07zmnuBe8+hLkmSJGleeeU7SZIkCYOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSAIOxJEmSBBiMJUmSJMBgLEmSJAEGY0mSJAkwGEuSJEmAwViSJEkCDMaSJEkSYDCWJEmSgBLBOCL+IiKORcQPpx27MCKejoiXituuac/dFRGHIuLFiPh4swqXJEmSGqnMivE3gevfdexOYFdmrgF2FY+JiB7gJmBd8ZpvRERHw6qVJEmSmmTWYJyZfw/86F2HbwC2F/e3A5+edvzRzHwrM18GDgH9jSlVkiRJap6z7TG+JDNfByhuLy6OrwCOTBt3tDgmSZIkLWiNPvkuZjiWMw6MuC0iBiJi4Pjx4w0uQ5IkSZqbsw3Gb0TEpQDF7bHi+FGge9q4y4DXZnqDzLw/M/sys2/58uVnWYYkSZLUGGcbjHcAW4r7W4Anph2/KSLOi4jLgTXAvnMrUZIkSWq+6mwDIuIR4NeA90fEUeAPga8Cj0XEF4BXgBsBMvNgRDwGDAE14PbMnGxS7ZIkSVLDzBqMM/Pm0zx1zWnG3wvcey5FSZIkSfPNK99JkiRJGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIEGIwlSZIkwGAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIAg7EkSZIENDEYR8T1EfFiRByKiDub9TmSJElSIzQlGEdEB/B14NeBHuDmiOhpxmdJkiRJjVBt0vv2A4cycwQgIh4FbgCGmvR5Z2V4/07GhnbT1bOZtVdd29QxkiRJWtiaFYxXAEemPT4KfLRJn3VWhvfvZOWTN3MFNSZGHmCYR04JtY0aI0mSpIWvWcE4ZjiW7xgQcRtwW/HwJxHxYpNqmdGKi96z5pJlb/1CkCTBG1+/4cevnvjpS80Yo1O8H/inVhfRhpz31nDe559z3hrOe2s473O38nRPNCsYHwW6pz2+DHht+oDMvB+4v0mfX1pEDGRmX6vraCfOeWs4763hvM8/57w1nPfWcN4bq1m7UuwH1kTE5RGxDLgJ2NGkz5IkSZLOWVNWjDOzFhH/CfjfQAfwF5l5sBmfJUmSJDVCs1opyMy/Bf62We/fQC1v52hDznlrOO+t4bzPP+e8NZz31nDeGygyc/ZRkiRJ0hLnJaElSZIk2jgYe8nq+RERfxERxyLih9OOXRgRT0fES8VtVytrXGoiojsinomIFyLiYER8sTjuvDdRRJwfEfsi4vli3v+oOO68N1lEdETEcxHxZPHYOW+yiDgcET+IiAMRMVAcc96bLCLeFxHfjojh4nf8Rue9sdoyGHvJ6nn1TeD6dx27E9iVmWuAXcVjNU4N+N3M/CCwAbi9+P/beW+ut4DNmfkRYD1wfURswHmfD18EXpj22DmfH1dn5vppW4U57833P4C/y8y1wEeY+v/eeW+gtgzGTLtkdWaOAycvWa0Gy8y/B370rsM3ANuL+9uBT89nTUtdZr6emf9Q3P8Xpn5xrsB5b6qc8pPiYWfxX+K8N1VEXAZ8Atg27bBz3hrOexNFxC8Avwr8OUBmjmfmP+O8N1S7BuOZLlm9okW1tKNLMvN1mApxwMUtrmfJiohVwJXA93Hem674Sv8AcAx4OjOd9+b778DvAfVpx5zz5kvgqYgYLK5kC857s60GjgMPFq1D2yLivTjvDdWuwXjWS1ZLi11EXAD8NfClzPxxq+tpB5k5mZnrmbraZ39EfKjFJS1pEfFJ4FhmDra6ljb0scz8ZaZaEm+PiF9tdUFtoAr8MvBnmXkl8K/YNtFw7RqMZ71ktZrqjYi4FKC4PdbiepaciOhkKhQ/nJl/Uxx23udJ8fXm95jqr3fem+djwKci4jBTLXGbI+Ivcc6bLjNfK26PAY8z1aLovDfXUeBo8U0UwLeZCsrOewO1azD2ktWttQPYUtzfAjzRwlqWnIgIpnrQXsjMP532lPPeRBGxPCLeV9z/OeBaYBjnvWky867MvCwzVzH1e3x3Zv4WznlTRcR7I+LnT94HrgN+iPPeVJn5j8CRiPil4tA1wBDOe0O17QU+IuLfM9WbdvKS1fe2tqKlKSIeAX4NeD/wBvCHwP8EHgM+ALwC3JiZ7z5BT2cpIv4dsAf4AT/ru7ybqT5j571JIuLfMnXiSwdTiw6PZeZ/jYiLcN6bLiJ+DfgvmflJ57y5ImI1U6vEMPX1/rcy817nvfkiYj1TJ5ouA0aAz1H8vsF5b4i2DcaSJEnSdO3aSiFJkiS9g8FYkiRJwmAsSZIkAQZjSZIkCTAYS5IkSYDBWJIkSQIMxpIkSRJgMJYkSZIA+P8WRgiGvcDjCwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots()\n",
    "fig.set_size_inches(12, 8)\n",
    "# for i in range(10000,data_X_new.shape[0],2):\n",
    "for i in range(1,2):\n",
    "    ax.cla()\n",
    "    ax.plot(data_X_new[i, :],'.')\n",
    "    ax.plot(data_X_new[i+1, :],'.')\n",
    "#     ax.set_ylim([-1, 1])\n",
    "    ax.set_ylim([0, 750])\n",
    "    display.display(fig)\n",
    "    display.clear_output(wait = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with open(folder + 'dataX_extracted20148_noise.pkl', 'wb') as fp:\n",
    "#     pickle.dump(data_X_new, fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_X_new.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ----------------- #"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### the same with safe dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dataset_size = 300"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "load safe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "safe_extra300_params = np.zeros([dataset_size, 7])\n",
    "safe_extra300_freq = np.zeros(30)\n",
    "safe_extra300_slow = np.zeros([dataset_size, 30, 6])\n",
    "for i in range(dataset_size):\n",
    "    params_read = np.loadtxt('./dataset/safe_extra300/my_output/'+str(i)+'_my_params.txt', delimiter=',')\n",
    "    data_read = np.loadtxt('./dataset/safe_extra300/my_output/'+str(i)+'_out_data.txt', delimiter=',')\n",
    "    safe_extra300_params[i, :] = params_read\n",
    "    safe_extra300_slow[i, :, :] = data_read[:, 1:]/0.3048\n",
    "    if i == 0:\n",
    "        safe_extra300_freq[:] = data_read[:,0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "safe_extra300_slow.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "_ = plt.plot(safe_extra300_slow[0, :, :], '.')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "interpolate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "safe_extra300_slow_interp = np.zeros([dataset_size, len(data_freq), 4])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(dataset_size):\n",
    "    f0 = interpolate.interp1d(safe_extra300_freq, safe_extra300_slow[i, :, 0])\n",
    "    f1 = interpolate.interp1d(safe_extra300_freq, safe_extra300_slow[i, :, 1])\n",
    "    f2 = interpolate.interp1d(safe_extra300_freq, safe_extra300_slow[i, :, 2])\n",
    "    f3 = interpolate.interp1d(safe_extra300_freq, safe_extra300_slow[i, :, 3])\n",
    "    safe_extra300_slow_interp[i, :, 0] = f0(data_freq)\n",
    "    safe_extra300_slow_interp[i, :, 1] = f1(data_freq)\n",
    "    safe_extra300_slow_interp[i, :, 2] = f2(data_freq)\n",
    "    safe_extra300_slow_interp[i, :, 3] = f3(data_freq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots()\n",
    "fig.set_size_inches(12,8)\n",
    "i = 100\n",
    "ax.plot(data_freq, safe_extra300_slow_interp[i, :, 0], '.')\n",
    "ax.plot(data_freq, safe_extra300_slow_interp[i, :, 1], 'x')\n",
    "ax.plot(data_freq, safe_extra300_slow_interp[i, :, 2], 'v')\n",
    "ax.plot(data_freq, safe_extra300_slow_interp[i, :, 3], 'o')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "safe_extra300_slow_new = np.zeros([dataset_size, input_shape])\n",
    "for i in range(dataset_size):\n",
    "    # first dipole\n",
    "    flat_arr_temp = safe_extra300_slow_interp[i,s1d:e1d,2]\n",
    "    # second dipole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_extra300_slow_interp[i,s2d_a:e2d_a,1])\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_extra300_slow_interp[i,s2d_b:e2d_b,1])\n",
    "    # first monopole\n",
    "    flat_arr_temp = np.append(flat_arr_temp, safe_extra300_slow_interp[i,s1m:e1m,0])\n",
    "    # second monopole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_extra300_slow_interp[i,s2m_b:e2m_b,1])\n",
    "    flat_arr_temp = np.append(flat_arr_temp, safe_extra300_slow_interp[i,s2m_a:e2m_a,1])\n",
    "    safe_extra300_slow_new[i, :] = flat_arr_temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "_ = plt.plot(safe_extra300_slow_new[0, :], '.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "_ = plt.plot(data_X_new[500, :], '.')\n",
    "_ = plt.plot(safe_extra300_slow_new[90, :], '.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with open('./dataset/safe_extra300/safe_X_extracetd_modes.pkl', 'wb') as fp:\n",
    "#     pickle.dump(safe_extra300_slow_new, fp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# fig, ax = plt.subplots()\n",
    "# fig.set_size_inches(12, 8)\n",
    "# for i in range(300):\n",
    "#     ax.cla()\n",
    "#     ax.plot(data_X_new[i, :],'.')\n",
    "#     ax.plot(data_X_new[20148+i, :],'.')\n",
    "#     ax.plot(safe_extra300_slow_new[i, :],'.')\n",
    "#     ax.set_ylim([0, 750])\n",
    "#     display.display(fig)\n",
    "#     display.clear_output(wait = True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ----------------------"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### the same with sens21 dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open ('./dataset/sens_full_vti_21/data_sens21.pkl', 'rb') as fp:\n",
    "    data_sens21 = pickle.load(fp)\n",
    "print(data_sens21.keys())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(data_sens21['data_slow'].shape)\n",
    "print(data_sens21['safe_slow'].shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# slowness sorting\n",
    "for i in range(data_sens21['data_slow'].shape[0]):\n",
    "    for j in range(data_sens21['data_slow'].shape[2]):\n",
    "        for k in range(data_sens21['data_slow'].shape[1]):\n",
    "            idx = np.argsort(data_sens21['data_slow'][i,k,:,j])\n",
    "            data_sens21['data_slow'][i,k,:,j] = data_sens21['data_slow'][i,k,idx,j]\n",
    "            #data_amps[i,k,:,j] = data_amps[i,k,idx,j]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(data_sens21['data_freq'], data_sens21['data_slow'][0, :, :, 0], 'o')\n",
    "plt.plot(data_sens21['data_freq'], data_sens21['data_slow'][0, :, :, 1], 'o')\n",
    "plt.plot(data_sens21['safe_freq'], data_sens21['safe_slow'][0, :, :], 'k.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# cut freqs\n",
    "freq_start_idx = np.where(data_sens21['data_freq'] - 0.5 > 0)[0][0]\n",
    "freq_end_idx = np.where(data_sens21['data_freq'] - 15 > 0)[0][0]\n",
    "data_freq_new = data_sens21['data_freq'][freq_start_idx:freq_end_idx]\n",
    "print(len(data_freq_new))\n",
    "data_slow21 = np.zeros([21, len(data_freq_new), 2, 2])\n",
    "for i in range(21):\n",
    "    for j in range(2):\n",
    "        for k in range(2):\n",
    "            data_slow21[i, :, j, k] = data_sens21['data_slow'][i, freq_start_idx:freq_end_idx, j, k]\n",
    "print(data_slow21.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# interpolate safe data\n",
    "safe_sens21 = np.zeros([21, len(data_freq_new), 4])\n",
    "for i in range(21):\n",
    "    f0 = interpolate.interp1d(data_sens21['safe_freq'], data_sens21['safe_slow'][i, :, 0])\n",
    "    f1 = interpolate.interp1d(data_sens21['safe_freq'], data_sens21['safe_slow'][i, :, 1])\n",
    "    f2 = interpolate.interp1d(data_sens21['safe_freq'], data_sens21['safe_slow'][i, :, 2])\n",
    "    f3 = interpolate.interp1d(data_sens21['safe_freq'], data_sens21['safe_slow'][i, :, 3])\n",
    "    safe_sens21[i, :, 0] = f0(data_freq)\n",
    "    safe_sens21[i, :, 1] = f1(data_freq)\n",
    "    safe_sens21[i, :, 2] = f2(data_freq)\n",
    "    safe_sens21[i, :, 3] = f3(data_freq)\n",
    "print(safe_sens21.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(data_freq_new, data_slow21[0, :, :, 0], 'o')\n",
    "plt.plot(data_freq_new, data_slow21[0, :, :, 1], 'o')\n",
    "plt.plot(data_freq_new, safe_sens21[0, :, :], 'k.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# falt dataX 21\n",
    "data_slow21_new = np.zeros([21, input_shape])\n",
    "for i in range(21):\n",
    "    # first dipole\n",
    "    flat_arr_temp = data_slow21[i,s1d:e1d,1,0]\n",
    "    # second dipole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_slow21[i,s2d_a:e2d_a,0,0])\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_slow21[i,s2d_b:e2d_b,1,0])\n",
    "    # first monopole\n",
    "    flat_arr_temp = np.append(flat_arr_temp, data_slow21[i,s1m:e1m,1,1])\n",
    "    # second monopole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, data_slow21[i,s2m_b:e2m_b,0,1])\n",
    "    flat_arr_temp = np.append(flat_arr_temp, data_slow21[i,s2m_a:e2m_a,1,1])\n",
    "    data_slow21_new[i, :] = flat_arr_temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# flat safe 21\n",
    "safe_sens21_new = np.zeros([21, input_shape])\n",
    "for i in range(21):\n",
    "    # first dipole\n",
    "    flat_arr_temp = safe_sens21[i,s1d:e1d,2]\n",
    "    # second dipole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_sens21[i,s2d_a:e2d_a,1])\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_sens21[i,s2d_b:e2d_b,1])\n",
    "    # first monopole\n",
    "    flat_arr_temp = np.append(flat_arr_temp, safe_sens21[i,s1m:e1m,0])\n",
    "    # second monopole\n",
    "#     flat_arr_temp = np.append(flat_arr_temp, safe_sens21[i,s2m_b:e2m_b,1])\n",
    "    flat_arr_temp = np.append(flat_arr_temp, safe_sens21[i,s2m_a:e2m_a,1])\n",
    "    safe_sens21_new[i, :] = flat_arr_temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(data_slow21_new[0, :], 'o')\n",
    "plt.plot(safe_sens21_new[0, :], '.')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_sens21_extracted_flat = {\n",
    "    'data_slow21_new' : data_slow21_new,\n",
    "    'safe_sens21_new' : safe_sens21_new,\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with open('./dataset/sens_full_vti_21/data_sens21_extracted_flat.pkl', 'wb') as fp:\n",
    "#     pickle.dump(data_sens21_extracted_flat, fp)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}