pattern.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Pattern recognition\n",
    "\n",
    "### Preliminar version of image recognition.\n",
    "\n",
    "Packages needed\n",
    ">* [polarTransform](https://github.com/addisonElliott/polarTransform)\n",
    ">* [astropy](http://docs.astropy.org/en/stable/install.html)\n",
    ">* [timedistance](https://github.com/AngelMartinezC/Local-Helioseismology) (By the moment this module is not used)\n",
    "\n",
    "### Sources\n",
    "* [Object detection (OpenCV)](https://www.learnopencv.com/image-recognition-and-object-detection-part1/ \"OpenCV\")\n",
    "* [easy-tensorflow](https://github.com/easy-tensorflow/easy-tensorflow)\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from astropy.io import fits\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from timedistance import td  # not needed by the moment\n",
    "import time\n",
    "import polarTransform"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The next step is to find the greatest value of the cube array for every frame. Then it is selected the frame with the greatest value in order to find is pixel-values. These pixel values will be used as center for the sunquake."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'Max value:  3056.45  ---  check value:  3056.45 '"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cube = fits.getdata(\"Data/HMIDoppler.difference.fits\",0)\n",
    "\n",
    "value = np.array([cube[i].max() for i in range(len(cube))])\n",
    "frame = np.where(value==value.max())[0][0]\n",
    "max_data = np.where(cube[frame]==cube[frame].max())\n",
    "\n",
    "y = max_data[0][0]\n",
    "x = max_data[1][0]\n",
    "\n",
    "\"Max value:  {0:.2f}  ---  check value:  {1:.2f} \".format(value.max(),cube[frame][y,x])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f7815f34e80>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAD3CAYAAADhaQjCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJztvXm07VlVHvrNc869RRVVtAopwViFwxh5JkhjE8sGBRuITzSisYkplWdFJYqiBtCMoCbvBaLBhy8RLUEslQhi8yDGRHwMEU2GSFEWgpQGRISSCoVKU9Lce8856/2x97p8d+5vzjX3ae65G9Y3xhl7n9/+/daaq59zfnOtn7XWMDExMTExUcHWSQswMTExMbE5mIvGxMTExEQZc9GYmJiYmChjLhoTExMTE2XMRWNiYmJiooy5aExMTExMlDEXjYmJiYmJMuaiMTExMTFRxlw0JiYmJibKmIvGxMTExEQZOyctwGFw1VVXtfve977yt9HxKGYmv1eePQr4PI87X5W2kiH7PZLPXx+l2++P0vPPR22VlWkkwyjPKI+jqgP13EHav+eT1YX6rVr3R4Us3YOMhcPU1cXAYca36u/VflIt41ve8pa/bK19dOlmh41eNO5973vje7/3ey+4trW1MJ62t7fPf+eK9xOWma1MNNvb2+fT29/fP/89m0CyAcrX1TNqIlTyqTRYvo5ebsbe3l6av6+rra2t8/Wg5OT7vHx7e3tyEu33dZnVfb7sXr7+rCrH1tbWSjl8fv3ZLI9Ru3gZqgrK9vZ2OplkZevP93LyJ2Nvb0+WN1tM1+2b+/v72N3dvSCP/f39lb7IfYM/1XjzZeE2UuC68mXitFT7VqEm6modjcrhn2XZlPxqzEdjoP+fyXf99df/+WqJa9joRaNXjKqUc+fOyfv9dx7IaqLhdFUDq/sri4v6fdTQahKIBohPwy8AEVjWPjEo+Xp6XG6eQDx4Eeoynzp1aiXf/f398+n4iYnTVhNSa21lceTJjK9lkwD/X7Fc1rE+fHo7OzsXKDpeFl9nUfr9O/dnnkBU+/s6UPXCi6Sa5HmBUpOZz2t/f7883hSq/T5bSFT/HGnvyupWi71SzrwMUZv7fs9yVi1oHiv+2aPC5DQmJiYmJsrYaEsDWKzaEZTJ28Ert9fkz549u5IGaxu8mlf8lEpT4WuZ1sQaDWslHkqTirSdir+VXUeZa4uR5cvPKI2L5fDulb29vRUNSrnP+NnMfaGgXFbKrcP1x9q20tArLo3W2vlycP497XX4iCiPra2ttC2rbs8ObjM1/kaWi9fQR+XgMeg1dOW+9d99XlzP63Ij2ThjKyorBz/TPSLr1EsmE2M0Hg+KaWlMTExMTJSx0ZbG1tYWTp8+HfooK+TQ1tbWii8x0sTX1QAyDVdxAZHv1KenSFtF0DI/kJFxih/o1/lTac+sTY7IU6/tnjt3TvqjvebEBLfiV7g8ngNQRKSytiILTGHk/+bP0XPsBx9ZsFkbjbiyKK3omUzrZU6jWxps+WXpcpqcl++nqt+r9FX7KgtHEfWMUT0yl+bv43Hk+9+IJ+K6UOPIc0cK0fynOJejwEYvGsCiI21vb8uG6RgRX76xqqQjY12yKTIpVadUrgXvrok6sTfd1eK3v78vI4eq5LmHGhRqAcvcNsCHJonexvx7a22FMGytpRFQHVFUT1QWzkNFfEVRUWrAZ+4ahppY1eDPJjufp8e67aombW5TT5hXoRZsNd6iAAYPHh9qUlaKBytuXlniscBp+LEV9YGs36t5IBsLvpy+bCpC0Kd7WEz31MTExMREGRttabTWcObMmfPfGYrIBVY1s5GpPzLtK6F4fF8Uj93vV2awN8+5bErLGVkLSj5lziu3iZdF7YdRLiZ1bYRMMzKzC0J2WXYvp2rrzJ2jtHvWPjM3gtq3MHJ9KBei1yyV5cz9QFnTI8s505BVOpxHt/Ky/FSfG1l+LGumoft0VZrRvZmXQZUnGlv+mtL4AU2YK2sma68spF89O/I8HAbT0piYmJiYKGPjLQ1FKAF6Zyygd9V6v6bacBNpyRn5q1b7zDc+Ir4yjVWF+62jTXrNljVl/vTXeBOl0rxVaORoxyvL0OHbQe0i5vtVeGOmlUcciSpbZqWo9h0FHPg0+Dv3U983tre3pRbrEW1sG/E+GdRmPD/euK6yOog0dAUVwJBxjKNNjb4/8+/c59Ztt4x3YutXnT4x4j4Ud+Qtl8gDcBTY6EUDWO28qmFGMe++MqMOlFU6N5CKWlAuDTUAOph4Ux3bdwg1SQJ6V7d6lici4ML4+9OnT6+Ug90wag+Fiq7xE1F05EUlSi0yv1Vb9++8KGQ7x1U5FBHJ7Zu5hzgPP4DVLnB+RpGdnIZ30VUJ/ZELhJ9TsmSLShRxFclXjWyKUFmY1ATLsiq3rFI6OY1MYczKocZ+FGWl4PviSInMgmAOgumempiYmJgoY+MtDa/BslbuXVHKXI4IWq+Rjc4rysJHFYnpy+DBpnG26zaTCcg1PKUhcb7KxeS1XbbAWHtW+0i8xqY0HxXKCNRM9xHRp6wKZQ0qCyJzU7LMytLldL38e3t7UnP0FpPS+M1sxULMrBUvv6+rnZ0d6crLQjfVbm2uW/9s5EqqaM/qehTo4O+PrA8VcODbfHd3d+iy7L9l7isVPKDKlVkrLLNy72V9N0t3HUxLY2JiYmKijI22NPqOcNZgu/+dd2eyRlDRfCL+4qAajdIE+r0MFaapEGkWShP18ilE1pYqm69TxVUov7CXP5JPXYu0U//77u5u2r7KumTORVkBmSWkLJwR71WpD6W1+7Jz/pxGpO1WCHNliW9vb5csXV/OKF/WgJWVong2n2eXK5Mps0xVeVVbcl9X9ZdZYMpCVBvv1FhQ6bF8Wb0oK2/EkayLjV40lHuDIwv8xHXu3Ln08L0+OHhxUY2auYRUFJN6Rrk0+LsictX/ox3Bigzz4MVK5ds/2VXG+zQ6eF+HIsyzelPyqx25SgHosvOikUUJcXpM/HsZqgQ358f1mAUIqPQUsqM5ImQLA8tSdXH6a9FeCqXIZPslVF3xOPL3M0b9KrvGyKKO+Jrqp2ph8nmoPVV8H6fr61y1gSLlo/R8OY4K0z01MTExMVHGRlsafTVVWsmZM2dWwk0ZmZuouo+gGj44elZpBEyEKg24orVHZwV5gizS1lTZ+r6MLIZdmf1K61X5Kqti5Grk/TVZXWf9gAMnojBY/z9rzJ7k3N/fP2+VjYhNn0cUglpxL4z6WmZBKFeKCnTILCOgHs6bWe9s5WUulygMXY2pSh7q9APV/xjKZTVyvXqM3Nodo3rs4PPaDnqGXCjrkaQyMTExMfERgY22NPb39/HBD34QwKoWoTQ4XrlHxBc/01EhwtmPO/IVK77Baz6RJq+sI59exFX0+/rGMM6XrbNMq1LnUmV+4VOnTklNSnEBCipAILPQMg2T5WYNU4VEVnzPChF/ocqRbXAc+etVH+/IODMFxRMpq7aq3asd61G+qmyZ1cT9KnuXvU+Xv0d8g9fm+SVg2QY9lTf3g2qgg2rDfi06UTnDURHgHRu9aDCBpiYLJrb7/coMrUQO+e8eTMBn5B//n5my677TO3o3sJJFkec9HXapRO43/lQDb7SDtiOKzFHHWmfIopgiE1+Run3R4HeTZ8oIB1CoiUFFCfnyRgu7mlT8BKP6DddB5maNJv7MVcZktW8jdhNzHso1qNrIg11C2cKk+rhqtxGUC4wXcV+noz1BPI79PMRYVz4FNa8pheyoMN1TExMTExNlbLSlsbW1hcsvv1yu+kz4Kg1uZF0o15Hae1Ah0VVstdI6qqFxozDD6P9IdtbWWENSxLqvl93d3ZWAg8gV4A+E3N9fPSJapZMdwx1d4/Iol4vPy3/vcir3hYcKs+b7qpqeav9RyK1yh6j8stDszA0UuY6yfR9ZQISy1CKZozIweE9Q5mJSFlhkLVTIZ3UfsOo+i/qp9zJEMndwP1DPVt5YeVSYlsbExMTERBkbbWl0H7EKeey/82dVo4nyGm2k8/lWNjYxRhrrKJwyS5vh5RppJYrA65rNZZddtmIR7O7urhDlrBGePXt2JY8oWCD6jWVW5+2MCEtVtu57zjZ5Rkd3Kz94lFcEZfkpSyNLZ0TuZsEUI0Sab0dV864EdvB5T+yvz8LFOc9sDKr7MkTcUEZis/WtzlzLUJ1fqqHh1XyrOLZFw8x+GsCXALiztfbJy2v3AfAiANcAeAuAr2qtvcsWJXw2gMcCeD+Ab2it3TLKo7XVN/f1htze3j4fHVQlBDndrCNGBGSHItsrURqKyPW/+99Grje1cHbwxM7vxQDGURrZblnlvuABryJzWF41ALJBoNyGbMKrelN5ZO4BrtvKAjs6gqQa4JC5mhRJHsmm3GtR8ATnwdE/KnBDPcuyqL1S2QQXucP6Z3VR9jIrhWykCGZ5KKJZBTVwH89cnNWFQt0bLZIq0OAoSPHjdE/9DIAvdteeCuDlrbVPAPDy5f8A8BgAn7D8uwHAc45RromJiYmJA+LYLI3W2ivN7Bp3+XEAHrn8fhOAVwB4yvL6z7bFkvh7ZnYvM7u6tXZHlofZ4iU0ETHNRxsvZbrgWX+/T9uVR96j3BaKkI7CCoHVMFZ1T4QspFHJEu2q9Zotu/vUfcrC4XOpKgcbRvWXtY3SyC677LKV3zKiPCJAs3dSK8tFYaSNV0lndb/S5H16StNU6Y1cFUp75nQza4H7C4cv8yfLOiLMVcgy56esWtWW67pnMu0/IuD9s2arx68zRtaHH9+qT2aBOVF6h8HFJsLv3xeC5ef9ltcfAOBtdN/ty2srMLMbzOxmM7v5rrvuOlZhJyYmJiYuxKVChKvlT6oFrbUbAdwIANdee22LVnL2H2c+UX6+qtUHcp3/VH56zzccJLzW58Xfo01iXnblp+cTXpnI8+mxP5+1Ok/++Wd6Hp7cVZaGCt3033u6/tW0I81MkcuqHOpZtfmL68dbC9HJqMrKU7+NtNf+qXz3qv9VtGxlSSp+yj8D5Dwaf4/OB+vgfJVlmo0frkvfnzi828tUkUXlodKr8jXqWqWN1PwSBUsoPnPd4AeFi71ovKO7nczsagB3Lq/fDuBj6b4HAnj7KLH9/X28733vG5pdijhUEQjKTTRqyAqBplwaTMxlaUSTgDqyoSOLHOKFU5W338+TfESe82eEERHZ5e5BCzx5K/mUTH2xUsEIPMhUfasJseoe435VnUDU4ucXznXIbF+O6P51SdiqG4sVALUrXl1TUEES/l61W1vJrII4RotMlVivBFVEiFzIBwW3kVogjgsX2z31UgDXL79fD+AldP2f2gKfAeA9Iz5jYmJiYuLi4zhDbn8BC9L7o8zsdgBPB/AMAL9oZk8A8FYAX7m8/dexCLd9ExYht99YyYNNZ8oXQHzWUeamUS4IRpVcU6u+0virWp/SWL02PrIqMq2zukucXUIjYk65Obwbiw+083KzXPxsRgL7skTlHZGOygpVddDrp1tJXg5VzxX5lIyRG0v1jZFW3VE56M/L6mXjZysWGtdfR7T/Qp2n5GXg8vJbGH1whBoz6+zU92VT10ZjS7WvqvuRq89D9QPO46jfEX6c0VNfE/z0KHFvA/DE45JlYmJiYuJocKkQ4QeCmeH06dMXaHWZds/o96vXfPLvDE868vdM640slwrZqQjVkT9Vha9mGs2IE2JZFJnt8+VQS9YW1et3vSz8u8LIT59pU5kVp+pA1aPiSHZ3dyVvlllbjIrGO7KOst3kaixUw5lZQ1eE8KgNfNnYqlDcB8vnXwGrPAtM7qoNpyyTmiMyqPpQG3fXTYPLodpNPRudIJGhYmEfBBu9aGxtbeH06dMANJGm7lf7B/j3/lmZ0EeEpTreItuFqwhuHhQZaRvJk3UUHpSVnbEqcojlU7IoV8pBJk4VDJClp+qI3ReqLVX9ZnsA1t3ZPpJvtOPf56GiYdTEymmrwA8/mY7kZKiJbbQ7Xb07XSk3oz0HPl+/L4tlWIe4zhb7kXvU5zEKBoiiBdeRj79XFeDDYB5YODExMTFRxkZbGgBWQjY7dnd30xV2RFRl5xB1RKa7f4ZDMkcmsSdZ1R4KNvFHx4YrYk7tiVBp+HwjbVwR5l6WkYskcxNFhK9KZ11UyXZlbSmCUbnyqm6x0bNKvqwcmTWjtHt+JtPQWWYVLMDBEr6PqWCAvb096VKrvIhMWe5qt3aURrWe13WFMpQVlY2Fg1iw6r7jwrQ0JiYmJibK2HhLo2sUnVRV2hVrDkrjyjaq9fuz36LvVeuCNQsvC5PFmUbIyKweFTYbWRyK7KxoTSMNLPMvq/Sq5L1ZfDoAoI9pV/WYWQjZGWMe3mrc2tpaOTtpa2v1uP4RV8b/K01e3afeid7RZYnOFvNtzlYKWwMVyyXjixjrvAs768cZ55dZF/6+6hH1qq8pa0dZsFVUeJqIJzoKfmOjF429vT38zd/8jYwYYagJbhTZVJkIowbweUTpqsFd7XR+Yhu936DqilIdLJtIovxU/alBq+rAgwdAtKfE56FcaqPd34qQziYuvr/iIotcfh38PhJFevtFiN1YHcp9wuS9eudJVi+q/Nwn/ZH6/tlsAeb0vMxcB6PJUQWPZP0vG+/8DvfMrTiK0OLF2V+LTnlQC7G/r7q48H1q1/5hMN1TExMTExNlbLSlYWbnQ2691sRm/0E04FG+0bPqPrXCq/OKslDJCKyZVZGRcD5dvr/qioq0e685RtaKCkLwllS0H8LnoUKClWXKZ21xub31oTDa9as0TP5UriNVDh/brwIEIkLVa9IcVp4R6yrt0Y5/hi9HpKFn7kKGsoiVBZaNYbbA1NjPjjJXfVdd436lvAdeFmXtqzlMuTOVfKr+2D1/GExLY2JiYmKijI22NDrUCqo2zShfp3otp9JslQa8zvub/cofnWu0LvnH/ytSTWnZGdGs8qgS5gqKH1CamdJiWfPuvvOMPI3I2Ip2qtLk8FClUat69s/7z2wnsbJq1I76jjNnzpQJUO8vj45u71Aa/4hEzeqKy6gsjcrmUp9OhNFubdU3lOVXzY+Rjbdu5ak+pDiXaFxlO8dHpPtHPBFuZituhey9DOo+YDXKRR1smB0PwXlUFw1FEkYN3aHcDZn7JCJyMzeRGjwsQzagVIw9m+mZ6axk4Le+9UWDo+TU4pfJwvDEsFqsFKrHnfDvnG7WD7J9PUq+aEJXC0mfsNjF5euFZeGjPipkdjQZKTdlxb2i+jMreJxf5jocTZKVsbfOwpkpRhVXZ0U+tWh4F5iqn6NaNKZ7amJiYmKijI22NIAPaR9Ki1Uag1/t1QquyKaMHGOss5J7LSIy55kE888q7ZgtCVVOpTFmpHz/bW9vb6Ue1O7lKOTRv6tdkcDKdXTq1KmVHf8jV5k6PJG1sYoVwGXi/1X9KTeH73/KqojO7vLpcigoy6xeqOUtSa4DlY4qK8ukXB4Z+FnOr+eptOLKuOG+wXXAfceXLXMlRuPJP6ve6hjJ579HO++VJVl5Vs1vHJ7Mn5kb6zCYlsbExMTERBkbb2l0TciHI/IqXdUmvQbE30dcxUjL8fJFhGCFzFbyKK2d60WF9jHvoPiVbDNZdjbW7u7uBXxEv6bqJbOsePe0Op04CwbIuJlqyKNP2+fBMvvvEb/CIa++vEqzZfjd5FyWLJzTf+/3ZUEFXJ4sZDlDxLMoq7ayaS2q5/4sW5QVojzKU/UrlW+myWd9SJVXtYGyPkbPZBxNlM66mJbGxMTExEQZG29pABdqdRkvwVod+0GzaKiOkSYw0ry87zLyeSt4n22k6WVHSvA1r5mpIzRYK2ZNOfOhj3zUSgtSrxzt39X5Wyo99T9r+RkfMjoy5KBlU3UaHfvh20Pl31o7Xw+Zb3x0nExmkRwkSkhB9edMCx9Z8XxfZnVnnBBzR1ze0dEy/f+M28rCrTntjDcZodoeWQj5UUROAYNFw8yeXEjjfa21nzwSadYEE23ZJK9MtWy3LN+XuVL4fxXKGJmI/f8qwag6WB/87Cbi7x1qkKmJMBvI7Hbq5eQQ2CwMUr34iKFcBMrFpOrAT5iRS83nyy61Sjgp8KGj95UbUO2QVpMGuz2zZ1mekUtNBSH4OmDFKDuHSNXzaPd3ReHyUP1eLXbepcvPKCWjt5GSOXIHV9xs0YLn6z5S0nwbVvdfVCd5bt/qGwUPg1GNfS+AKwFclfx993EKODExMTFx6WDknvq51toPZTeY2d2PUJ610FpbMekVGabM0Q7WuJS7ht0cKr0O1rZ7Oso05l2hXA5OwyM7E8nXB3+a2QqRr97zHQUNePcQWxodp0+flq4ZpYV5bahCfvryKhddZsFEJGaH0p5VMIB6dpSv2hGeacCZa4vdrYoI7+CNkMq9osJNOS0lX6bFqjrgfCthn1VXk3INqvZVHoPonKlMfgXuz4rMVhZY1hf5WV+/I9kyT4ZK5yhIcGBgabTW/sUogco9ExMTExMfHigR4Wb2JADPB3AXgOcCeCiAp7bWXnaMspXQV3dlQahXPmbhgyNStGPdsNlog4/HSBNQGguHXGaafBY2y2kyf5FtvOunC58+fXqFV+L2UMcaZGVmjSuzwEaaaPYq1ogjUTJmFs6ID8nITsXDqPOwRmSr7xPcD9gK99blyP/O+fdne5tHZGvG0UVWeSaLx8hCZAu/GhJcsVKj/nKYkOVs3hilp8Z5xi1l/N5BUI2e+qbW2rPN7IsAfDSAb8RiETnxRaM3ih8MrbWVmHj+nZ/34AZUBwvyZ6WxKmXo9/s8IrO156fKpiYGllO95VANFHZpAQuzv9cpuwD8YYLVPSaRq0JFufgFRF3j7152QL+tLhtEkWtLKQX+PkX8j5SRzL3H7imlBHG5qzvvVT9VfbyD2zfr2yN3q3Lh+HEZRSZm7VbtBwqZq1QpMjw+svHLsitlM3NjRe64zMXu72Xs7u4eyaJR3afRc3osgOe31l5L19aGmX2Xmf2Rmb3ezH7BzO5mZtea2avM7I1m9iIzO33Q9CcmJiYmjgdVS+M1ZvYyANcCeJqZXQWgxmI6mNkDAHwHgAe31j5gZr8I4KuxWJB+tLX2QjP7CQBPAPCcSprKnGPNYhTq5u9jRHsEOjKiT2lc61ofEfGqXAHZrmkfKuvzUnsystBhtlYyszozyUfEf6atqTptbfWIbaVRKy0xOircl41dAd4FF8maBRwobZJ/71DvJlfaKbdl5uZQryHN3JU97Ug+9Uxkqal29c+oHf/87Lquxr29PekSUmVQ1/0ZWqocIzelspIVlCuU7/eyKKuWz1dTgRiHQdXSeAKApwL41Nba+wGcxsJFdVDsALjczHYAXAHgDgCfD+CXlr/fBODLDpH+xMTExMQxYLS572Hu0oMO6xNrrf2Fmf0IgLcC+AAWvMhrALy7tdaZu9sBPKCSnlqJAb0ZKuIMlNbutQNFOrLGqvyVmX+bobgF5ccdaScZ6cyyXHbZZRf8psJw2XfPG/n4/RO+DpT/O7NCVHtwORRG78mo+P397x7KuuP7lXy+vyhLVyHamZ2R/0pT5/yyDZ2K5+D0K2TsSFNmotbXs5k+R0xZASPewsuZ7XJX/ZT7RmYRr2NV+Px53lDvbuF0M8uef1PyeR5TWWqjdqti5J7698vPuwF4OIA/xILL+PsAXgXgs9bN0MzuDeBxWLi63g3gxQAeI26VpTOzGwDcAAD3uc99cObMmQvcEkx+Z/sbVKOrAaXcEuse6czEYdY5uWNXTG6fR0dGHKrDCRuR42fPng3zjUzkbEd4ddHga9nizHXvF9jIlVfZdb6/v79SV6O2VSRrli9HM6l+oMrLMvhTAFQbjSY9Tle5NKr9zpct2ltSSYMPs6xOxgpR+/fn1GKvIuGiSLWRTOtMyr6u1D4wBu+3yvZAZXkdVuHvGO3T+LzW2ucB+HMAD2+tPaK19nAsQm7fdMA8Hw3gz1pr72ytnQPwKwA+E8C9lu4qAHgggLcHMt24lOMRV1555QFFmJiYmJg4CKpE+N9trb2u/9Nae72ZfcoB83wrgM8wsyuwcE89CsDNAH4LwOMBvBDA9QBeMkpoa2sLl19++QXX1I5cT/ICOvyy/767u7uibWxvb6+kF+WrkJGEXg7+XiXXFLa3t9NziHjHt3JtZUejs5aTlX1EmCqteN3zgDIo8i+yMrOACbWHgrXUChFZfXe1IuXZmlaWbkY0Rxq4snCyfsDPeZdfZE17sLuGUX3JWeZOzLRtJsKVZcAhutnYVM+O+rOy9iuut+pvgLaYFI7C2qgS4beZ2XPN7JFm9rlm9lMAbjtIhq21V2FBeN8C4HVLGW4E8BQATzazNwG4L4DnHST9iYmJiYnjQ9XS+EYA3wrgScv/X4liOKxCa+3pAJ7uLr8ZwKetm5b377I25sNM1YuAWBNQ5C4jW6UzHzVrV0eRB/v9O5S2y+dHKQuMNSpFEvoztHZ2dqTvOdNuFCeUaYZqY5tPx6enZMrqUbUHY1SOTKsbEcSKqPeWQdQ3vHyKBOayKUsok5XvUyclM1TARoc64l3lqUKRfV14mZXGr9ows3rUKQqZdaHaiKFOYKhwp5VyVDwUPu8o36oVP0Jp0WitfRDAjy7/Lhns7+/jgx/8IABNfHnwPWoLPt+XTQyK6MtisM1WY+fVAWqjSYVl8Oa8yjdarFR9eKj0VDSMlxUY14saNJl8ypxXrgWGimar7kPIiM29vb2VSSJK1z+7s7OTllO5ELntfbAHg9tZHd2euTd4sVLuTAV1AoNX3FQ6Kr3R7vmREqEmQu9y5jplOb2rNpo31N6hfq0frxK5tpTCk+3QV6i4o/m+yH12FKiePXUdgB8A8HH8TGvtQUcixcTExMTERqDqnnoegO/CYj/F8b/lYw2MdhUzopU204JGLhJldSitya/2TCqr+PJIhmqZgLqLQVk4Ub0qMlcRpZkLaqTxqHpWZKI/6C+Tm9Nl+TICVrU9uz35/ef+mei96z49/91fy2SOXB9Ke/YyR/0sqz/lruHfqmSsL2O0B6qS70jW0YkOI829Q+1pqZY30/izYAqVZnU8RScXHIW1UV003tNa+6+Hzm1iYmJiYqNRXTR+y8x+GIs+rB1qAAAgAElEQVQ9FWf6xdbaLcci1SGgfJ2sqXi/68hfrs78YSi/a6Yx8P2edFRElSJ8mSjNNFG1wXF/f3+F4FZli7S/iibKqG6ayngOpTVF/veM6xm9REjtivfty6HXfne8l0WhEjQwknNE3vpyqDxUOlE9V3gYxWOpgA3la1eENN+nxkxm2fO1o9CsFX+W1cUII2KdZVd1oMZ0x2iMXUxL49OXn4+gaw3A5x9agkOgE8zK3FOTASPbW8BmuiLSVJx8FvmgiDSe+Pk3lY46lpll9WCTvOLuUvJFE85oIoryVBEy0X6EbOJSA4DrTE1SFYKbZVAKhQpC4MW+qngo8lzViydKs4gkTk8tGpGC4mWq7sweBTpk+1KU4sFtlCky/D/Xo+rPKl//u5q81fiLFDcvC9cf5+UDGJRrC9CuQ5Vepqhme4Gqrr0RqtFTn3ckuU1MTExMbDSq0VP3xGJfxecsL/02gB9qrb3nuASrwMzOh7t1KMKSw9pUGpkVwOlWQwC9xhVZOv35zEWjoMKEvaz9U2kqXlNWLrpuxflyqjy6FpxpOeoFMLzzXmmarHmfOnVqJX8+SLGnl2mEnK5yX3jLIDqryBPwbO2qfhARlV6WEVRb+vSU/KpO1bUIvv+OCPjMElZuyqgOMkK/Y7RvhtPKXFaqrjpUeC2HXo+sQJa156Fcv5kng8vRkVlMfI0/j2KfRjWFn8biVa9ftfx7LxZv7puYmJiY+AhCldP4+NbaV9D/P2hmtx6HQOuiE7YVn3y0M7ZjFM420rj9/SOLpLoze0SWRVDaBmtr/UTbc+fOrZDKW1tb8nW5HWpn++g1s5WNhqzV8WcmA/vrlcavoNpG7Z739csWbNcwz549Ky1TVQcKWZ9UJLD/jcH9qnrWFSOrF4V163mk7SrrSZ37pcLU1UmwKt3IQotkV5a4Ooo+4s88V8GWyWjO8Yi8DL4PqTEThfWui+qi8QEz+6zW2u8uBbsOi8MGTxS8I1x1TtWAanLsYDeLcgWoBWI0MPt92QDl+zg/QMf7R0SkR7TIqT0ZyiXTFxWus+4m4mt9cVEDmQe0mtDV7mVlkqs4eb9YRYuzvxZBtYefgKPd1dnEP0I1Oilzv/ACqhbv7B0gjOp9apLysvjvEUaRUOo3nrSzgIjRfq1M+VKBDPycCi7I3vHCykZlX5SqF+Wii6AUnqNAddH4VgA3LbkNAHgXgG84EgkmJiYmJjYG1eipWwE8xMzusfz/vccq1ZrY3t5eMf1YI/QaATAOdfOINNYKkctWikqDUXVBqHIoZIepsRabkXBs9nfSeaTteE1KkZ2KwGP5/L0eSj5vdXCbs7ZWJW1VnqqulLwjt1SHeuFXdE8kV9bXgHzXfpYuyzXKP3MJZQECmXs2yi+ztlX7srWvvBEsZ+bizNxsyu0VtYW3cKKTCdRbDtW9fpyzG7XaD6soEeFm9n+Z2b1aa+9trb3XzO5tZv/mSCSYmJiYmNgYVN1Tj2mtfV//p7X2LjN7LIB/eTxi1dB93Xt7eythkOoVnBy6qTbaMEHl/X+8wme7YFlTybSlkUbD96nvPb3MZx9pwErj8ukqAq+/WlelGeWhdlf3z1OnTq2ETI+Oq87O0VFtxLKNzvryVs9IG+f0Mp4j8vd3ZBsm1QZHnz/nxf05yzey7CqcjAoQUHXP5ehQGnB0eoOCspK91h5xZVVLqAIeH2zpVrmeanm9rJFnQZVNzS9HYW1UF41tM7ustXZmKdDlAC47dO6HhJmtxO9n5iP/XyWTO1prK8dvjCYDftZ/VwNPdeIoPzWpeCj3T7Vz8jNqR3MWPaWCEHhSyRa8qA18vlFZfOz8aKCo/qJ2UitwnaqJXz3r74vy8C4wnghHLjW1aHCET5QXf1/X3Qp8qM6zxV71e/5dtcfIRat+U4tahtGknAWPqDHIfSmLSON8lWKpds+r71F5/PfqGxIzVBeNnwfwcjN7PoAG4JsA3HTo3CcmJiYmNgpVIvzfmdkfAng0AAPwr1trv3GskhXQNTxlivEqrTTbHqqrzExFOip3ktJOlbXAsmRa4mgfSWYGq2ciQtCDwzlHx24rmbzVE71gKrOOOCZfaXtKW/PWTM+bEbVH5nbyWj7Lx22UhXUq8LOZC2fk3hsdoqmsMtUPMpdM5PrysnCe/gDHiEDOysby+p3WUR/K3IoVdxujej+P/cx6Y3eccpmPXNPZ3iGeh3wd8Dji+6q71zNULQ1g8U7w3dba/2dmV5jZVa21uw4twcTExMTExqB69tQ3A7gBwH0AfDyABwD4CQCPOj7Rxtja2sIVV1wREl9K4/Fa7t7e3oovVvEhkbWgVn2Wr0OFzo12jPf/fdqRX7hybXt7e2Wn99bWlgzZ8/KdOnVq5RprhNmZTUpDUvexlcJ5qfZYl6/hulX9oCNrZ2WFjojkjFhnRJp5/8yCECKLyiO7FnEbntuKxkIFyipU/IXKg9stO5aeteysbyirPzrZWD3ToSwIhrJqM8+ICrCobAvg+/b3V19cVnm+gqql8UQAnwbgVcuM32hm9zt07odEaw1nzpxZyxxVBGNHlgZP8gcZeJxO/y07JE09m/0euT7U4PZvM1MLEz/br+3u7qZuuI7oTXqKHI/KNbo+mnzUgsP5V97zXZGr/5bFxHM9Zm47NYmpBUfdNyJIfb1wYAJHsCni1aenDoZUC+eI5Of/s53ofH+XvwfAKJdVFsk2gnoLY3XBYfmiez0yxZZPW+BTF7JIKZ+Gx1EQ4dUDC8+01s72f8xsBwtCfGJiYmLiIwhVS+O3zez7AFxuZl8A4NsA/OfjE6uG7hpRJqoilZVGaGYr7prIFdCh8lPupAr57eEJaYWI/PQEbaSFZZpoFm7aWluxUlQwgHLLKY21BzL4cig3ViV+njXW/skvNMrIWD5VgNPz30eksQoGUHstlIap+i7nr+pPhUB7jFwmXC/emo72w6iwXl+OkZU8cpUoF53KT9VpNvYUAc9ur+w0A+XG4nb2bqTIdZnVAVtqft9Z1WLisaAstcOgmsJTAbwTwOsA/DMAv44T3tg3MTExMXHxUQ253QfwUwB+yszuA+CB7RCMipndC8BzAXwycH7fx58AeBGAawC8BcBXtdbeNZAr1Nwif2qmgXSL49SpUysnhe7u7p5fuaNNgj5dvqZe2aq010p66ndFrnF6aqd8FjQQneyZ8ReK86mStwqscfmyRbuwff1FFkqFsOT8lOaYtd+pU6ekJpqR7JkskZXs+26V01O8zogwVVaw4tGywIBqv1ccnQokUBaYegEb35PxcmxxZuMoshqyVzNniO7z8u/s7KxYcltbq+9iH3keDoNq9NQrAHzp8v5bAbzTzH67tfbkA+b7bAD/rbX2eDM7DeAKAN8H4OWttWeY2VOxsG6eMkqoHyHiBw0Psg4+RsS7lfjZ/jz/zm4sZaIq9w8Pbr9z3T/j813XdFdlitxYUVr8jCqHGkTcOftR6qO0+X8/MUSuGS9/RCCr8mWuLTUBK/mrkSjqyHiecLJJJ3tDXVQOTwyr+uM9FOpAz8z1xnl0qEWB2zJ7lp/Joo6iAAbl8lPPZkea+Prm+yL3jeqTqu/2oIJRn8zGt5KPr6nd5EpOr8RFrul1UXVP3bMtTrb9RwCe31p7OBYb/daGLU7K/RwAzwOA1trZ1tq7ATwOH9plfhOALztI+hMTExMTx4cqEb5jZldj8arX7z9kng/Cgh95vpk9BMBrADwJwP1ba3cAQGvtDiuE9PZ9GmfPnpUmtiKRvMtAuXV4NVcWhCLMR6t+JD/fxxomP5tpKpEWHpVtpG0olxDnmbkClMtlRP5lZCd/Zppg5lbi/SaZFafqfh3S0afNcqp3sft+GCGzLkYhrfyptHt/v3JzrON2UntffD9QeyhUXwN0aLjq415WZaUocHnZra1cxIrgrgYDeItTWfH8ndtDzU1e5ihMWLmmj+JFTFVL44cA/AaAN7XWXm1mDwLwxgPmuQPgYQCe01p7KID3YeGKKsHMbjCzm83s5rvumhvSJyYmJi4mqkT4iwG8mP5/M4CviJ9IcTuA21trr1r+/0tYLBrvMLOrl1bG1QDuDGS5EcCNAPAxH/Mx7c1vfjOuuOIKXH755YsCLVf106dPn/fzMp/gd3/zSs/hmv6MFiaglA844w9G2rGyDHwajEhjVVqJ0jpH2mvPQ2lcyqc88on78o6040z7roZF9/v4JUycV8WKiCwnZVVWCNfW2nneZ1RnyoL11onaZKfaQ2nPChziOUKlrymeSJ3FxLzYKA+1iVIR5uqlUxn67+xlYDkyYl21UcaPqZNvq1xD1OcyC3HE262Ldc6eOhK01v6Xmb3NzD6xtfYnWBxF8obl3/UAnrH8fMkorXe961148YtfjNtuuw3XXnstAOC6664DAFx99dW4xz3uAQC47LLFKe6qUbmDMXmmdhb7t9aNBvy6HUxNrKOIEbVwqTwUVIet7BOJkB3dzfWcuayq7rNIfkVO+nyjBaeqDCjZ1eTtJx8umyJys37ArhRVV5lSwAR8to+J+z3XrTo9Ievj7I7xbcR1xf0lWzTU/+p+li97p3yWhnKBKURu7cwFy/OMIuqjfEbl4HZTbruOyl6nCi76orHEtwN4gS0ip94M4BuxcJX9opk9AcBbAXzlCck2MTExMRHgRBaNtnjn+CPET2sdgPhJn/RJeOUrX4krr7wSd9658Ga94Q1vAAC89a1vxXve8x4AHzoGnV1Oyn3B2lNfqVVI7iicr0KkMlHl0/fI4rIzspM1GpYzC/XNZB9ZOKocVSI8c4txfkqTzywrDi1V1huXJzPxM6g+otyeUZ1Vd7v7+7k9MotOkcWsnSoLkL/7PUssw8i69Ne4//X8lUto5EZRJLQKWumI3J4Vq3bkWmWZsnrJDvRUaUfeCD93RfVccRsfBGstGmb2WVgcXPj61trLjlSSiYmJiYlLHumiYWa/31r7tOX3b8bitNtfBfB0M3tYa+0ZF0HGEO9+97vxa7/2a7jb3e52Xhvqq/jOzs55TuOe97zn+WtKe1QrstfGI24h09oVRoSX0pCUb9yDtRIOxasSXxWfcnSPIiKzfDOtVJU30vQ8yR1ZNerlQGrneJWM9PkpYnMdi0lpp4pU7hiFkyoNmcnz/ttIa87S89cUBxbVs/K/qzGoLHFldXuM3jOvoKzQLNBBWfaqHNyv2NJQbe4tYlV/nDafUJF5CNRG5sNgZGnwNuYbAHxBa+2dZvYjAH4PC9L6xHDu3Dm8/e1vlxEr6s1V7HZSA487QUYEVyfEzB1zGESx9mrgeVdaNLiz9PhatnCqeHHONyN8eVLzC7GaSKr1qIhNtTBFaVbIWB7cfJS1OrrFp6MUFRUVxflxP624Qrm8/KlI9ExRUNFJozrjRcpDya6i3VTafiLmdLhe1LPZ4gfkbsksuGVnZ0cqQUpmP/FHezei/P01P9cxuF7WVYwURovGlpndGwuS2lpr71wK9j4zO/x7AycmJiYmNgqjReOeWOzYNgDNzP5WW4TMXrm8dqLY3t7GVVdddcG1dc1RFRao7uv3rpO2kiszEUc7fDMNfUSUevfEKD0mJzn/LFRQmeRcz/5ZZc6rNlB5qSOsVRCCOpfMP9M/VSy+b49Io/f1y9p9Fmqs3CH8UiRlubJmm/VJVS+q7rkNlKWhws/9mFFHqEfICN8syGR/f/UNnZyecjGpfNVY4fvVtQrOnTuXejI6uE+q/pf1Z64DZSGuS/YfBOmi0Vq7JvhpH8CXH6kkExMTExOXPNYOuTWzb2ut/TiAPzsGedZGtItV+TrZz8wagHo/tvKnZxoIaxjZyq7INaVVqbT9c/wb38Mhkr5sytJgjZWPh8/CBzOo+tnd3V3xb0f16LU1RfTxOTpVziMjTznkNuOgImvF+9j53CSW3WvtXPeZ7EoTVW0Z9R9f9lOnTq20K3MB1Z33o7OMKvUXhX77dlLWjNKyqzvMR2NNWS7ravL8LFtRyuruGFndnthmDiwqi5flMBhFT/mjzw3A08zsbgDQWnvWoSU4JPxk0iu07wIH8o4YRSL56zx4/FEknD/LkE1cSgZFREbPVsm6ittHuWbU4FFQQQgqP3UEROSqyyaabLJVYHNeQZ0CwFCRJ4pQVbL68ip3HL93XblmoiCFLA/lrlF9UpHOfnJUfShyY6k9QT5dtQiNFuwqic7yKVdPBpWO+u0w4MXDB6aokwsil5XvQ2peU+WIrq2LkZP+BwF8OoArAVy1/Nxefr8qeW5iYmJi4sMQI/fU/wbgWQDuDuAHW2vvN7PrW2s/ePyi1eB3k/ZV+uzZs+mhZSokjt1FXQtTpGimofvvfH9/JrpPkZ6R5lUNnVX3eFdKVQsbaeOZ5s/kbkdkffj2UHXbWkt3BY/M9eo7l5ULMdOK+XDMrBz8rH9BEqPqKusYuUe5LTk8OCrb6Dwl9X/VMlXuwnUDRaJ8ebe5z4Otx0owAKc96i9KPmVdKjdgxXugysFtqSymdd3LI6SWRmvtra21xwP4HwB+08wefyS5TkxMTExsJKpHo7/EzH4TC3fV7ccrUh3dShjtxlZWRUekPatn1W+e5Bqlw/9XT33NLJdMM1Na3enTp0saTbTBzEOR1BFx7fNVVp6y3lS7KatMhfp2eRjK78/5qTwybkGFRiryWZWXrS1OV/UDD/VCo2gsKD99J+uVtlvlETh9xQNGaURQO9F9PhGyuh/lpTgcVRejPpFxpqo92EuSbWZVQSFZmZRVrdroIChHT7XW3g/gew+d4xGjcsTG6PlsdzATr76xokGZEVC8CPkFiTu7Mjc5rWwyydwJioRjE7rXBUdPKRIuGwCq0/OuVeXuynboR2VTO8ezSWUUZaXcST6CTC2I0aJfdQdU9/9USd3KxBUtsGqxz8pbVZoyxc7MVkj5kTtwRJj79h1NmGqxyhRGL7/HqJ+qvRZe5kih6L+P9qpUXbDrIu2tZnbjKIHKPRMTExMTHx4YWRpfZmYfTH43AJ93hPKsDa/hsMaSEdf+hUp8n1qlR5qKMo0VRhZJh9JsFTk5Cqv0VoB/73JFZhUumYX1VklotjTUe9wz60MFJqhQaaXFs2WZkbs+7Yp8I0LVg0OHOV3f/1g77YS5spIjS0i5Qyrun729PenqybT7DkU0r2sh8LUqKR/tJs/Kyi6pzBWqAmO4rTKLi60LRU6rsvn7lMWmrFozWwkCOCqLY7RoVNxRv3MUgkxMTExMXPoYHSNy08US5CDY2tq6YBMfMA6N6ysyv584uy/6P8uvsqKr9JSPsnqmj5mlYXnqZVNVeTOyGFi19ljzUcd9q93SI3kq1gd/z7RoFdbLUBo1W0Te8hv1DZV2dgKusvyi49f9se+RHBVSnrVnfjbTnrmd1w3hzrRstpyzNCJZKs8oyzniNxWf6S0/zkPxdnyP5yyVhajKxPWcnVjM4HnhKKyNk3rd65GhN7zqdArqvspOzCitCjkZTWoqDxVBkaWpyC7uJMpdkw04li8j4xUycpyPEcncK1F6atHIFpIRjsIVxciCBjhPPzlFk2OFtFUE8uhYj3Wjk0ZunZHbqerG8ulxEILqkwwfgMHPqrbMXMQcsOHLE33PFm7lRuXjjHiMKfK+I1MERy7iakBGFbWwjYmJiYmJCaz/ute7t9bed1zCrIvW2gqxyxpzRaOOiEMf967C36I8qhqvv29EUvJvSlvz5eBrnEc1dPPMmTMAxrtWfb7swunaEB8HzRaH19KUi2SEqmuGocpeIYtHL/gZafoZMZzJp6we5bKK0lKumcz9w+Mkc3korZjTUNcr1qXSlFkW5cLpUGWL5PN9J2pfdcaXmksqbqKzZ8/KMZOdqNChwtSBC08iiGSO2mNdlCwNM/tMM3sDgNuW/z/EzH780LlPTExMTGwUqpbGjwL4IgAvBYDW2mvN7HOOTaoimIwaaeb+e19x1Xk/SoupaoQqnUhb8xpFpC0o7VRtDlL5KwIv2yzFPmWvxUb1kpF/3brgMFK2PnxAwjrlUESul4+hdtVmGqGqK2XRRecaZeG+I+1PWb9ePmVJKi1bPRtp6qpsamx5bkG1B3/n8lbCon06Xk6WV/XTirXKfKYCl8e3m+rjI4s2G7/nzp2TvKLnLrk91ObcUbkv9o7wtzmBcjv8IqBXajTwq0dsqArPiNIonZ6vvz9K1+ebHR/Av7fWpCyVQQasRuyMYr+rCycvxMp053d1dDn6d57EKx2bF9PRbmSPaHBX08tcRzypsKz8yTJEfUotFpnMDO/+idIYKRwRook6k68apcb3e7mi4ILM5RzJ72VmN7TK1z+rfleTvLpvf39/RVlS7m+Gj5iKylQ9EuYwqC4abzOzzwTQzOw0gO/A0lU1MTExMfGRg+qi8S0Ang3gAVgcWPgyAE88TMZmtg3gZgB/0Vr7EjO7FsALAdwHwC0Avr61djZLo2seIzJ65LaoaHXqwD0FRYZF+WdyVw5D4/tYk8/SZxmyuPpoR7PSMJVrKwuD5Df49WtsavNeAi+zl8Nfy97MmLVHZAUoDVhpdcrtqSwML99IO83ObGJERG/2XGaZHsR9lllPowAKby1sb2+vtE1mpWfl8vkry6CD81WyZt4G5ZZT+zQ4FJ77kLJwMjJbuVE5D++t4PFxGJSI8NbaX7bWvq61dv/W2v1aa/+ktfZXh8z7SbjQWnkmgB9trX0CgHcBeMIh05+YmJiYOGKULI2lFfDtAK7hZ1prX3qQTM3sgQD+IYD/E8CTbbGMfz6Ar13echOAHwDwnEE6qXXAZ9bwM4DmDxSxpAg85TdU2lCGSIPzMnHamc8zIvQzrZg3Bymir1KWEYmZacDA6sYstSFxVFeRVdn/V1pYlf/J8mWZ1Qm/Xa7+AqoonFOVo0Mda63CkzMtPCpDhaOLkGm7Ssse9aXM2lcb+pTmzX04s+IVx6jy6IisvYp1p4IVlEz+e0dvf8VlsEWeWecZr3QQVN1T/y+A5wH4zwCOIuf/G8C/wIdeGXtfAO9urfWZ73YsXGFDdHMyM0mjoxrWQUQMZ/HqWYcYuY7YheN/U4jcROq+rF7UZJ/Fjys3Fueh3AK8N8M/qwaoIic5jw6OQMnISTWxR+SuT08t5uoASY4MU/WsjlfJXC583Ssq0X2jBTubRPn/yntfmMgdvWpALSSqH/hrGbkc/V5xKSs5fRpVpbBSDpWf2i+mIqWyvPxv2Vg4DKqLxgdbaz926NwAmNmXALiztfYaM3tkvyxulTViZjcAuAEA7nvf+x6FSBMTExMTRVQXjWeb2dOxIMDP9IuttVsOkOd1AL7UzB4L4G4A7oGF5XEvM9tZWhsPBPB29XBr7UYANwLANddc07IQtMzkHZHK/lrkHuD8omcjZKZ7RsapZ0YapLpXHbSmdhtzjLhKIysn15nSfJRrQRG4qizq2Ois3lhbVG67bM9PJkskm0pPWQkji6D/r8hnn67C6PwolQdfG+1y9zIoVwlbP5k8o7b3Y4a1cd7npPbXqN8qJH8UYOFl52Pko/O0PLjdMpL6IOGzVVfjuqguGn8PwNdjwTv0krfl/2uhtfY0AE8DgKWl8T2tta8zsxcDeDwWEVTXA3jJumlPTExMTBwvqovGlwN4UBuEwB4STwHwQjP7NwD+AAsOZQgOcQUuXNUzDbi6CUZtQFKaFOevzqlR+Sp4rUT5OqO8vZyRRq1ITPVsh/LdK+0pI+FUmZV8avetIpAVn3Tq1KmSppxtKuvl9FDlYUvMa5jRWWWjsFqfh3+e8+WyqhdLqX7A91QsnCoiHquyI5z7pL/uZevlzHgizo/7UhacwcjCVzMrT3EfKlRfnRYw4huyQBEva5beUVgf1UXjtQDuBeDOQ+dIaK29AsArlt/fDODT1k2jV44a1B2Z6c4dkScQb46OCLXR0QksL3+q8lRkzgaZWjiVLKMOlC1W1c7HaXg3l7rG4LwqBK2qU7UIKYVitAO5WheKuOQ2UG9OrCoAPj0lw2i/yeiUAgUVhHBQ+bhOM8VILSRqjCk5DxIllPU/PuGgQ7m/1eKys7NT7jtVuTMXZ/RWyqNEddG4P4A/NrNX40JO40AhtxMTExMTm4nqovH0Y5XiEPDmYOY64vuyVZ+tD+ViikxTf+0gZaloByPLirWOTBvmulIktdKkKntG+BmlgSpTmw82VFaFagdlGUS/+fTU/pr+DO/vUbJU9jeM2nFEwHqLJcrXpxdZoSqvzEKraqmVfqjyBi60xjL3WSZLZE0ryyUjmlXdZ9q7kknt/Fdli8qR1SFbPVm6mdt41IeqKC0arbXfPnROExMTExMbj3TRMLPfba19lpndBVywb8IAtNbaPY5VugGY3FLkVUZsZ5YGay/VlVn5YjPfb0YW+vsUWZudSdOhtJ2IsFQno3ptvDV9uq7PIyKu/a5pDkvtGtRBiP91+ZURuH6UpeHz3d/flzvy1TOZ5ZKdEaSsX04vCyGP6kdZ2BkybXxkGWSkvJKJn2FkdaCOKFdpZLyJekadg6VOlFD1sru7m1oEKkBgFCyTtavygmRlOwhGlsbdAaC1dtXgvhNDb5RscI8aTXX2qimbEe8qj/5sRFj5CbPqBlIkHE8+HOnjO77q2MrdxnsyWObKQsKTI5v9/shnnnxHE0jVLaH6hiqvmrx9esoFpiYLJfOoPHwtO0yy6k7K3GgqQIDbKFN41G+qX41cOGri5XT9gsgupmxsjcY+K0NZXY6Uvoo7KYqi8/fz9ZFSUB1v2bXDYKRGHw/9PjExMTGxkRhZGvczsydHP7bWnnXE8qwFMzt/GJzXcqK4dR+mWdX+olVapaOerboCsmO8R+6f7HycTOMCaiYvu6eUWc1pKe3Ua0hqv4SyZiJLzKc3csP4NFgGtjS4HD6NSKP3L5OK7vVuDtaU1X4Yltm7XyJLN4PSXDPylJFZvSNrn+EDDqJ7vAzK8mMLJ8uL66q6w30UTKo8pk8AAB+lSURBVFGx+HisjvqzvzZyjWduXLa6K+VZB6NFYxvAlYA8G2piYmJi4iMMo0XjjtbaD10USQ4I5X8EtL9UaVcKkT9QEcNVC6KymYy1F5XvyOpRZF3lnCzetTriDjJNOdOgOI8sTJnPD1L8AafPYbpe1hEnoPzlWf2N6tGTmJE1qGRRvyttUwVnKC12ZBX157zVqOogs5Z8ObitfdmyvhulV+EMVB9R/VRB1aNCZJlm1j7LUuEgqhj1oYx7OyqMFo1L2sIwM5w6deqCjtEnHbXrd2TGZcQhk7ZRmv3/kStolIbKv4KKS2i04KhrXAcVwlVNXNGuWp5ke1pq53iWp4ruqpr2Pg2WmfuGmgh5QfHHW0Qka7aQZIu4cgOqOuA01cTfr6k8eKJWe3j4WV8eNbFGsnpZoj5ZPbrfyzLKVyF7X4WaN/h+7hueqFeKoJrQR4Q5X/PXVYQlK6VqMT0MRjX6qEPnMDExMTHxYYPU0mit/fXFEuQg6NoXuweUW0KF+WVmJj+v3F9q5VYa6MjkHRFZ/bNKbFfcRK2t7pdQ2q4i7zkdtgZG50b5fEfajtK4fHlHIcv8v3JVZJqZuk8FIfD93o2l9vqMNEeu00xrrlp+PQ1ldXPfVdaZ6hvV3ekMPy5ZG8/elhe5QjOtXbVRNra4zbP0/DNRGdWYjl4A5+tFWQZReqoOKu64yAuyLg6/p3xiYmJi4iMG1bOnLkm01sIdlxEyDSkjOyNfotohqrSWSH6GyiMiBKM0ulxelii8taeR8QgjgtFra1tbW+l7nTMSPQoQ8GlE/JLXmlWdKutNWY2qr0THpnuZ2QrILBeWhbXxCs/FPIcqL2unqj0ynsHL6++vnMjK8o04g6w9qlB1qqwjtrBU/2MZIllG4clVzoX7WtbvOtR5corv4nSyOewgmJbGxMTExEQZG21pAKvHiChNmTUMFWHkoULYFKJNWD5fhcwPH8mQhdGp+yPfeI8w6xsjWQsb+YK7dcJWitfSlD9/e3v7vNan3jXheRZf7szCyHz8akOYaqMoykVFkI3OdPLp+eucL5eZ6yxr3w62IBS/kmmY3JYjrdhbQoeJwOH647ZUR+Co9soih1QebD1WNO8oL2+lcB+O+FP/m+IqRjyRss4r4fv8THU+q2LjFw1f0dyxs4WBn6+4L0bE8Ghir7qYVCy+z0tdU+lzmKvqsCod7sSZy0BN/Ko8Cuosob6QKZcLsHo2FQ8Uduv4xefs2dWXTSo3ljp7itMZuej8s2qfRs9bffo8/DW+N5NTuVx4gcgmn6jPZS5L5WrMxh3XiwrrVROcSkftueK2z4IpOP+qsqTaTfVzH3Lb5eHfWIFiZC4wzjdbwDitigvxIJjuqYmJiYmJMjba0jBb7ApWmpkKYfPf+//qtaEVM05pmiNNXl1Tv2UkHOcXPRPJojQgLkd/HWmUb4VAjuog0+45/f6s0rR482bXLFnz8u6ayJ1UdVN6La7qSuQgiY4R8a/kVG2pNP+sv7amw539xsqq60IRzQzVD7JNl0q+kbWq7stI79EcMTrVOqsbZe3z/f7kAuV+jAInskAHZbFzGVVo/VFYH9PSmJiYmJgo48PC0mBCLdPa+fu6Z1NFvuWMxBxpfx1q85KXzUORXJUjNFi7ynzju7u7K++64HszHzBrQyNrTxG+Pm0mfFUbcXkyyyUjuP0zPT2lhau+4a2E0TEO3G8yC7GXZ2dnp2QJRNZMR3ZUCf+eYcQ7ZNwff+e8+ncmyUeBBv0zI7i5DpTmrfp4ZtH5MnqZMmK9Q3GNEZmt0lPcjL+f72MrrropM8NGLxrc6HytQx2D7s8IUlCTKHdOhdEBb2pyVNe8DCNCWi0+WSevnJ8FLCYpFdvv3RIqioSJSJbTL7BKdn/dp8dydFnUG9v4PuV+VJO3yj9b1Bh+0usKjce6EUg8MWSuNMZo0vb3dUTtppApKJkLJIro8+ly3qP0vJtNEc1cth4UwOfJReNWfe/I6lQpKNmCyAqKv87pqTwjpVntWar2uwzTPTUxMTExUcZGWxrAqkbEK33VvPVp8S7zEUnoXT1VLVCdksnPK5Mzcmn1z4rJyxpw5urh3zv5vLOzI7UXr2221i7Q5vpn5prxZeU8qqfrsvwq9HrkvlNWh0+Xv2dEdGQZ9LqsWpKKxOSyqSAE9YpfVdbs/KOo7NF9o/GRIbIGM7Kd0/djYW9vL3UxqROVVftyetlOeS9PlJ4K7GBUzipTeUTnnPlrR+GaAqalMTExMTGxBi66pWFmHwvgZwH8LQD7AG5srT3bzO4D4EUArgHwFgBf1Vp7VzHNEunMZBg/q+B91FG+UV6RfJnmxVpxxk9k8o7kiTT2av2pZ5WWw++ViGQdkXr8bEZijpCF3HIZK5bpyOJQRLOyJPuz0aml3uph/7vaPc/5e0sj6sOVUGQFbvPqy64yqA2Jasyw9yALLvBp+/RY81blza5lPEfUf5TFtK6lNrJcFPfhre79ff2Ol3VxEu6pXQDf3Vq7xcyuAvAaM/tNAN8A4OWttWeY2VMBPBXAU7KEuhtEuXL8fYA+TjsyybP0Rh0sc5uMrvmJMNq16vMYxV+rzsTuuAyj+skm/kw+tUCcO3duSML6Z/j+zO2YufJ8mSKoCZbzzQZoFBBRiYZR5YgmZd/f1V4VllmdAqAWRJYpC4SI/u/PVsaCqqtqRJr/7lF94VdVKVRQbrYM0ckUflFWC4WSjxWPo8ZFd0+11u5ord2y/H4XgNsAPADA4wDctLztJgBfdrFlm5iYmJjIcaJEuJldA+ChAF4F4P6ttTuAxcJiZvcrPI/t7e2yy2V7e3sY7tlRNT2V9qdcYNl7k73cPj+lcfmzhPjZTBsZaR+cRkas87XMEmJt1pdXEdzRIZBKK/YkZkRYKjdRFiLNyFwfmTYeaXrK6u3I3HYcRtoPmlTpjtw1qg6yPsEWiXI/sqZesRBH9aJcyVl5+He2kpQW7s8v42eic7ey3zJXHvdNf5/yFHA5lAwdajc512nVSj8MTowIN7MrAfwygO9srb13jeduMLObzezmu+666/gEnJiYmJhYwYlYGmZ2CosF4wWttV9ZXn6HmV29tDKuBnCnera1diOAGwHgmmuuacv01H2hBQJoDSkLZezX/WfFb9i5F0akbSjCbaQJeln42Szf7Nyg0Smco82M/juHMSutWGlhmXxbW6svemIfMEP5q7M2zzS+qA6UP9qDr2X1Gz2rfq/6y9X9qv6U5acsEn/uV6R5K229skmW+YaRVag0eQ8uWxZIMOIxRmVTz6i5hE9mVs9E8vD4UJyg8gr4QJXD4iSipwzA8wDc1lp7Fv30UgDXA3jG8vMllfS6e0qZtx18FIOPk+eJRk0MKjqEG1LFeWdQAy97VkV9VNL2UPKvG9mk3E5RXj66ZkTGchmzvTF8zQ8G3sWuJhr+7hemzB3I90duD7XnQcmeTQLc/yruBrXoqp3F0b6KyuTNB4LyZ/aeeZbFBwhEqLh6OD1uv3UXDbUzm+uxEgihFJRIwczmpkwBGAX4qHlIpcs4CnL8JCyN6wB8PYDXmdmty2vfh8Vi8Ytm9gQAbwXwlScg28TExMREgou+aLTWfhdAtJQ/ap20urYUhS2OVmch28q1LKSQoXbuqnxZs1Va5Ej79+nxs+qsI4UKWRa5ebzGzQfpKXOZ71danUcUkqnK4PeyqJDq6C1uilj34PT4zDKvoatAh0iLVTufVdinf1ZZM1EemZtDjY/qjmF2U6mXOnnrY50w8CyYIts9P7JglLbeP3m8ZB6KyK3px2AUvq8soSw//xxj5AXpUAdSRgEC62LuCJ+YmJiYKGOjz57qBPOITFKalCKlMqJx5M+PTnuNoEJQWXuOdkb7MmWIfOhKu1Y++UpYqvIVq2ey18MeBFxXKvy46itWJ4RmFipr99n9vGtaWWIqWEH1tQyRxsoyRPJl71r3z/i2HPF7VeI6swK4bMqCHAWA+DxUgAWPj1F61XJ4ORUnObIgFFgm3+ajwIjMu3IQbPyi0Se3jOjjSlYmpToGw6cXma2+I3AHyo6P8LJG6fl0ehpqIFejJPykp97VEJGdKkqjctw6u0PYRadi8f0EyIR55qrjiZp/y3ao8/1ePpZZTSqZi4nT4TpQx4Osu2fElw/Ij4dXaYwidTJ3qzp+g12XCmoR4rpQi1/m/lFE8yhfpUBlrlVVV6yA+vklUoyU602ND9Xvq3tavCx87AzLNN1TExMTExMXFRttaQD5wX6VULzIXdPRV26lJfLvPs/ou9Jssr0FEVTZPCGstHZFZp89e7ZspivtKyP12GWhdqwrrdOnwRomh0crTb7iAmjtQ++FHx27XSHlWa4OpU0qLZr7lSK4+T5FmCsLLGtLXz7/qTRgtbs/u6byUNe4vJV3hFe15MgiyQJUuF+psZX1XZZVWZyj8G+PqluMP9XL6KoW7LqYlsbExMTERBkbbWmY6Vdq9t8yoiqzFvrzDKWFqecOEgKocJDNUBlvwqSxP2JbaUqKC4hCCpXP1lsQI82LNf+qtaAI6WoIrYKXWYW0Rs9VyM1o42IlD0ZGmI8CE7JgD+aYWGbf15TFpK6p/KIQ08ziWxfRZsvM49Ch+KnI0ojy9sisad+HGZFlpZ5Rm1QVJ3QU2OhFA4hdTzwAqlv2VWPw9+pu6GyAqvxGE+u6k0qW10g+VVcd0SSl8siQEf9R1Fg2KfNAPXfu3AXpqElvHah+pcrpJym1qAGr7guOouu/8S5sTt9POtGkUo0m8r8p1xaXSZWtisx9HLWRrxeFSPnL9p5kLtYoGEWR41k/qCJTbJWL0+fn5cz2lhyk3RSme2piYmJioowPC0tDmdUKfB9rDmq1XzdOfqTNVkxE5f5R2uJIi8ieYe2ln8Ol3En8jNKKOX9FYiptyGuOI9da1XJRhLk6/nqk1XmoY/cjl0//zqGPPs2ob6idzz4P1UZcRr6valFV7/OuGSZZVdhsVs/8e2Z9RPdXXjGg+lpkxavTAnw9q/6syHYlM7tbfVisx0GtgMg9qnaYT0tjYmJiYuKiYqMtjU6EK/JPaazsF2bt2K/Iysd6GOKSn8+IUJVmZJlkvmmVb1WWjujd1V7OkQ832/WreKfoJUyqDMrn7dt3d3d3RUZlTUX8SrXdfHpRm3oNfWRJcnoH4eF6OSpQsvB31pQzvpAx4v8q1zgNn1/EzXRZVShtJpuyIPhZ/lScWke1v4x+q841yuOhxtFRbO7b6EWjtYYzZ87Iw9eA1QmtLzJ8H5t2/JuKHMrinkcuIf+bl9XLPIr4UmlkA1TJkk3oWRoMRVJGewo6eMHx+zl2d3dlxE3meuP/extyNJYichUxnC0amfkf3Ze5GpVbZ0Si+37KMjC8CzbaB5RFz3Geahwxad+vZYvVunsVot8z5ZDvUe+9z1xW/GxlIR65W7m/ZnuQFEZuw8pY9fll19bFdE9NTExMTJSx8ZZGD7Gsko4dI22jayqZZpv9Fv2eaTGs7SoNqmpVMJTMXrON9iMorU65XzLiM9uprORTWrs6DFLF4nP9qTOCshdmKUtDkaJKW6weeqjyHrnoMheJIupVHrzrPHOpKCs0si59nbI7c3SOk3KbqHpRz2ZBF9xWWZ/ktqy6edV5Y74OVLlV3+Cy+TJn//c0KidYHDXZzpiWxsTExMREGRtvafCf+t2jr878jmN1pLMilqocREYSs2ZTOR22qrFG9ypNuUrCZlq2ylc9m70uV3FRig/Z2dmRVoC3BtVpqayZVd973qE0dGWp+XR82TJEFqeCz0OFFfe8PRRX4fNTMo808I7obLHK+UeqT46CLnjMqrrPoCx21t6z8Zv1IbZ+laXbEZ2NlY3pUV/KOL+jIL8ZG71obG1t4W53uxsAPcFlphqTir5hoqMAfOQGT3AqbbWQqAmHJ8SqiZqZ+MoM5ndOeFM7InIrRzuMDnBjWbLByO2m3i+hFhc/kCKiURG5FdP+7Nmzaw9uRrbbXcmZ9UklX+TmU4tBJiej0v9GY2tdjNy8/Jt3zahjhEZk+yhYQaXtXav7+/srAQcjxcgHD/AzymWl+rhy20V1pvZpHAWme2piYmJiooyNtjTM7LxWpsgwtdIq89ub1eqo5kgzzcjxSGbgQhJznb0OXb4KWFPpu7+VpsIYaS8V+dZxaWQvYcpCUHvb+2u+HIxR2yits2Ktjtw6yvrl+7wmur29veJSG7nKVL5ZGSNrQY0j30+5D2V9UVkGozbiZ72lpF5qFoWNq/7i71PWb9W9o9pSpaPmF+WyYvRnOfzcy+TLodo9cw0fBtPSmJiYmJgoY6MtDWBVM1ZayUh7UJqe0kCUFlbZLKMsHNYO1iWqFKehNJ/o/CVFWGZWWfbsQXza3GaVo8xZs83OA4osyUwGRhZmui7BGGmnGeHLXBNrm/1a9dwlBdVPR+VcB1F7qPGptH9f9+pEApYvs9JHnJWSeTSmFU/kw7v5u+IqVPhtFpIetYUqmzoL7Liw0YtGa22FNFVQDc2/VSafCKqTqAVi3cGoJhOfZ5cV0OSpKm9kkmeLgJoks0iykVuHo6KyhYnBx4L0vA5TpypfNZArpHfUVpV4enVsiuovfG0UDFCJ+hkpVT4vj2xCHykU2UKnXExKzux6RIRX3g7I9azcO3xqRDYu1FjIFiOVnuo3kXvPl2WktBwG0z01MTExMVHGJWVpmNkXA3g2gG0Az22tPSO7v7XFjvDoADq/Uo/OCMo0Uf99JFeU7rppAXmMfZY/x41zWpkLj90i/dkzZ86c/61C2rLl0hFprD5fLhdr6v376dOnV57NiOYor2wPTbYrX1kSyvKrWpcjV+jIpeHLNtqHw//7tlTg+7qWHe2bUf0+C3Rg2f2zo2CUqmtt5BZT6Olkb8NTgQ5d7kyG6De2sKtzzrou1q2trQO7Hy9I59ApHBHMbBvAfwTwGAAPBvA1Zvbgk5VqYmJiYoJxKVkanwbgTa21NwOAmb0QwOMAvCF7yHMISivpUGGGI+000wRG2su61ory51fDcVVoJGs9vINWcR/+2s7Oznn+oIfr8u+Zv56R+e7VfRHh7PNjvoZlUn5oD7ZcRtaEklNZAV5mlXbVQuQ+zZaLIjmVhVjRJrmvqd3L6vRaZcGOzpnyckayq3FZ1bKzvqPqVJHPaowehDOrtHnk5VDl4JMr1PP+mmpLFXBwGFxKi8YDALyN/r8dwKdnD5gZLrvsMlnhEaGqXBrqN0Wu+R3hUcdWE2vWcToUsal2vLKsirDs+fJkn+WrsLOzc8Gb/TivSJaMAFX1rOpFuY6qC3vkZsnkUu/iWJcsVIMcWB2gqmxe6fHPZS4/NUlF/SWDUoyyCYyfUZN3Ru5m+TPUnoyRS0rJnEUIRpO8d99G+aqTBpTSmrmMVHrcX7NIx8zV7r/38qzjFo9wybinAKjSrIxyM7vBzG42s5vvuuuuiyDWxMTExESHHQUxchQws38A4Adaa1+0/P9pANBa+7fJM+8E8D4Af3lRhDwcPgqXvpybICOwGXJugozAlPMosQkyAgs5795a++iDPHwpLRo7AP4ngEcB+AsArwbwta21Pxo8d3Nr7REXQcRDYRPk3AQZgc2QcxNkBKacR4lNkBE4vJyXDKfRWts1s38O4DewCLn96dGCMTExMTFxcXHJLBoA0Fr7dQC/ftJyTExMTExoXEpE+EFx40kLUMQmyLkJMgKbIecmyAhMOY8SmyAjcEg5LxlOY2JiYmLi0seHg6UxMTExMXGRsLGLhpl9sZn9iZm9ycyeetLydJjZx5rZb5nZbWb2R2b2pOX1HzCzvzCzW5d/j70EZH2Lmb1uKc/Ny2v3MbPfNLM3Lj/vfYLyfSLV161m9l4z+85LoS7N7KfN7E4zez1dk3VnC/zYsq/+oZk97ITl/GEz++OlLL9qZvdaXr/GzD5A9foTJyhj2MZm9rRlXf6JmX3RxZAxkfNFJONbzOzW5fWTqsto/jm6vtl3RW7SHxbRVX8K4EEATgN4LYAHn7RcS9muBvCw5fersAgjfjCAHwDwPSctn5P1LQA+yl37dwCeuvz+VADPPGk5qc3/F4CPuxTqEsDnAHgYgNeP6g7AYwH8Vyw2sH4GgFedsJxfCGBn+f2ZJOc1fN8JyyjbeDmWXgvgMgDXLueB7ZOS0/3+7wH8qxOuy2j+ObK+uamWxvlzqlprZwH0c6pOHK21O1prtyy/3wXgNiyOSNkUPA7ATcvvNwH4shOUhfEoAH/aWvvzkxYEAFprrwTw1+5yVHePA/CzbYHfA3AvM7v6pORsrb2stdbf0PV7AB54MWSJENRlhMcBeGFr7Uxr7c8AvAmL+eDYkclpi/M5vgrAL1wMWSIk88+R9c1NXTTUOVWX3MRsZtcAeCiAVy0v/fOlCfjTJ+n2ITQALzOz15jZDctr92+t3QEsOiCA+52YdBfiq3HhgLzU6hKI6+5S7q/fhIWm2XGtmf2Bmf22mX32SQm1hGrjS7UuPxvAO1prb6RrJ1qXbv45sr65qYtG6Zyqk4SZXQnglwF8Z2vtvQCeA+DjAXwKgDuwMGVPGte11h6GxXH0TzSzzzlpgRTM7DSALwXw4uWlS7EuM1yS/dXMvh/ALoAXLC/dAeBvt9YeCuDJAP6Tmd3jhMSL2viSrEsAX4MLlZoTrUsx/4S3imtpfW7qonE7gI+l/x8I4O0nJMsKzOwUFg32gtbarwBAa+0drbW91to+gJ/CRTKpM7TW3r78vBPAr2Ih0zu6ebr8vPPkJDyPxwC4pbX2DuDSrMslorq75PqrmV0P4EsAfF1bOreXLp+/Wn5/DRZ8wd85CfmSNr4U63IHwD8C8KJ+7STrUs0/OMK+uamLxqsBfIKZXbvUQr8awEtPWCYA532bzwNwW2vtWXSd/YRfDuD1/tmLCTO7u5ld1b9jQY6+Hot6vH552/UAXnIyEl6AC7S4S60uCVHdvRTAP11GqnwGgPd0V8FJwBZvyHwKgC9trb2frn+0LV6GBjN7EIBPAPDmE5IxauOXAvhqM7vMzK7FQsbfv9jyOTwawB+31m7vF06qLqP5B0fZNy82u3+EUQKPxSIy4E8BfP9Jy0NyfRYW5t0fArh1+fdYAD8H4HXL6y8FcPUJy/kgLKJQXgvgj3odArgvgJcDeOPy8z4nLOcVAP4KwD3p2onXJRaL2B0AzmGhrT0hqjssXAD/cdlXXwfgEScs55uw8GP3/vkTy3u/YtkXXgvgFgD/+wnKGLYxgO9f1uWfAHjMSdbl8vrPAPgWd+9J1WU0/xxZ35w7wicmJiYmythU99TExMTExAlgLhoTExMTE2XMRWNiYmJiooy5aExMTExMlDEXjYmJiYmJMuaiMbFxMLO95cmhf2RmrzWzJ5vZ1vK3R5jZjyXPXmNmX3vxpF3J+wP9JNTBvc81swcfMJ+3mNlHJb+/wMz+2swef5D0Jz6ycUm97nVioogPtNY+BQDM7H4A/hOAewJ4emvtZgA3J89eA+Brl8+cBP60y56htfZ/HJcArbWvM7OfOa70Jz68MS2NiY1GWxyBcgMWh9uZmT3SzH4NAMzsc+l9Bn+w3AH/DACfvbz2XUvt/3fM7Jbl32cun32kmb3CzH7JFu+eeMFyty3M7FPN7H8srZzfN7OrzGzbFu+pePXykL1/NpJ9mfcfm9lNy2d+ycyuWP72iqXV9HG2eAfCR5nZ1lLWL1ze80+W+d9qZj/ZdyBT+nc3s/+ylPP1ZvaPj7LuJz4yMReNiY1Ha+3NWPRlfyLv9wB44lKz/2wAH8DiXQK/01r7lNbaj2JxBs8XtMXBjf8YALu2HgrgO7F4H8GDAFy3PLbmRQCe1Fp7CBZHSHwAi13M72mtfSqATwXwzctjLkb4RAA3ttb+PoD3Avg2V7Y/x+KdFz8B4LsBvKG19jIz+6SlvNcty7cH4Otc2l8M4O2ttYe01j4ZwH8ryDMxkWIuGhMfLlCndf53AM8ys+8AcK/2oXdIME4B+Ckzex0Wp+gyj/D7rbXb2+LQvFuxcG19IoA7WmuvBoDW2nuX6X4hFmf43IrFUdT3xeK8oRHe1lr778vvP4/FMRAXoLX2XCxeqPMtWCyEwOL9Ig8H8Oplno/CYmFjvA7Ao83smWb22a219xTkmZhIMTmNiY3H8kC4PSyshk/q11trzzCz/4LF2Tu/Z2aPFo9/F4B3AHgIFkrUB+m3M/R9D4vxYtBHRxuAb2+t/caa4vu0VtJeuqz6i5KuBHDXMr+bWmtPCxNu7X+a2cOxKP+/NbOXtdZ+aE35JiYuwLQ0JjYaZvbRWLhu/kNzB6mZ2ce31l7XWnsmFuT438Viwr2KbrsnFpbDPoCvx+K1shn+GMDHmNmnLvO4yhZHY/8GgG+1xbHUMLO/szw9eIS/bWb/YPn9awD8rrjnmVi88+JfYXFMOLA4dO7xy0CA/g7oj+OHzOxjALy/tfbzAH4Ei1eVTkwcCtPSmNhEXL50yZzC4iVCPwfgWeK+7zSzz8PCSngDFm+o2wewa2avxeJ00h8H8Mtm9pUAfgvA+7KMW2tnl4Ty/2Nml2PBZzwawHOxcF/dsiTM34naq3JvA3C9mf0kFieQPod/NLPPxYIjua61tmdmX2Fm39hae76Z/Uss3ry4hcXJq08EwK/D/XsAftjM9pe/f2tBnomJFPOU24mJiwRbvH7z15ak9Mr/F1mWn1nm/UsXO++JzcZ0T01MXDzsAbhnZXPfccLMXgDgc3EhfzMxUcK0NCYmJiYmypiWxsTExMREGXPRmJiYmJgoYy4aExMTExNlzEVjYmJiYqKMuWhMTExMTJQxF42JiYmJiTL+f8BpBEydtpq+AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f782bc6f5f8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Testing the real sunquake with timedistance\n",
    "\n",
    "flare = \"Data/HMIDoppler.difference.fits\"\n",
    "\n",
    "# 192-260 the values of x and y found by eye \n",
    "\n",
    "# image = td(flare, x0=x, y0=y, theta0=230, theta1=240, savfile=None,rad0=0,radius=120)\n",
    "# final = image.tdplot(save=None,colorbar=False,plot=False)\n",
    "\n",
    "t0 = 230*np.pi/180\n",
    "t1 = 240*np.pi/180\n",
    "\n",
    "image_td = []\n",
    "for ij in range(len(cube)):\n",
    "            polarImage, ptSettings = polarTransform.convertToPolarImage(cube[ij],\\\n",
    "                                center=[x, y], initialRadius = 0, \\\n",
    "                                finalRadius = 140, initialAngle = t0, finalAngle = t1)\n",
    "            slices = []\n",
    "            for i in range(len(polarImage[:,0])):\n",
    "                slices.append(np.mean(polarImage[i,:]))\n",
    "            image_td.append(slices)\n",
    "\n",
    "plt.ylabel('Time [45 seconds]')\n",
    "plt.xlabel('Distance [pixels]')\n",
    "plt.imshow(image_td,origin='lower',cmap='Greys_r')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Time is 13.6 s\n"
     ]
    }
   ],
   "source": [
    "def ML(x1,y1):\n",
    "\n",
    "    s = time.time()\n",
    "    b=0\n",
    "    ii=2\n",
    "    images = []\n",
    "    while b < 240:\n",
    "        ii = ii+1\n",
    "        b = ii*10\n",
    "        a = b-10\n",
    "\n",
    "        t0 = a*np.pi/180\n",
    "        t1 = b*np.pi/180\n",
    "\n",
    "        # The next step is a conditional if to break the loop when the condition \n",
    "        # of sunquake is reached\n",
    "\n",
    "        image_td = []\n",
    "        for ij in range(len(cube)):\n",
    "                    polarImage, ptSettings = polarTransform.convertToPolarImage(cube[ij],\\\n",
    "                                        center=[x, y], initialRadius = 0, \\\n",
    "                                        finalRadius = 120, initialAngle = t0, finalAngle = t1)\n",
    "                    slices = []\n",
    "                    for i in range(len(polarImage[:,0])):\n",
    "                        slices.append(np.mean(polarImage[i,:]))\n",
    "                    image_td.append(slices)\n",
    "        images.append(image_td)\n",
    "        \n",
    "    e = time.time()\n",
    "    print(\"\\nTime is {0:.3} s\".format(e-s))\n",
    "    \n",
    "    return images\n",
    "    \n",
    "images = ML(x,y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Optional using timedistance module\n",
    "\n",
    "def using_timedistance_module():\n",
    "\n",
    "    def ML(x1,y1):\n",
    "\n",
    "        s = time.time()\n",
    "        b=0\n",
    "        i=2\n",
    "        images = []\n",
    "        while b < 240:\n",
    "            i = i+1\n",
    "            b = i*10\n",
    "            a = b-10\n",
    "\n",
    "            # The next step is a conditional if to break the loop when the condition \n",
    "            # of sunquake is reached\n",
    "\n",
    "            image = td(flare, x0=x1, y0=y1, theta0=a, theta1=b, savfile=None,rad0=0,radius=None)\n",
    "            final = image.tdplot(save=None,colorbar=False,plot=False,noprint=False,notime=False)\n",
    "            images.append(final)\n",
    "\n",
    "        e = time.time()\n",
    "        print(\"\\nTime is {0:.3} s\".format(e-s))\n",
    "        return images\n",
    "\n",
    "    images = ML(x,y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}