diff --git a/Training_Data/Particle_Tracking_Training_Data.py b/Training_Data/Particle_Tracking_Training_Data.py new file mode 100644 index 0000000..7cf767c --- /dev/null +++ b/Training_Data/Particle_Tracking_Training_Data.py @@ -0,0 +1,141 @@ +import numpy as np +from numpy import pi +import tensorflow as tf + +class Particle_Tracking_Training_Data(tf.Module): + def __init__(self, Nt, rings=True): + self.Nt = int(Nt) + self.Ny = self.Nx = 256 + self.d = 3 + ximg = [[[i, j] for i in np.arange(self.Ny)] + for j in np.arange(self.Nx)] + self.ximg = np.float32(ximg) + + x = np.arange(self.Nx) - self.Nx//2 + y = np.arange(self.Ny) - self.Ny//2 + X0, Y0 = np.meshgrid(x, y) + self.X = np.float32(X0) + self.Y = np.float32(Y0) + + if rings: + self.ring_indicator = 1. + else: + self.ring_indicator = 0. + + self._gen_video = tf.function( + input_signature=( + tf.TensorSpec( + shape=[self.Ny, self.Nx, self.Nt, None], dtype=tf.float32), + tf.TensorSpec(shape=[self.Nt, None], dtype=tf.float32), + tf.TensorSpec(shape=[], dtype=tf.float32), + tf.TensorSpec(shape=[], dtype=tf.float32), + tf.TensorSpec(shape=[], dtype=tf.float32),) + )(self._gen_video) + + self._gen_labels = tf.function( + input_signature=( + tf.TensorSpec( + shape=[self.Ny, self.Nx, self.Nt, None], dtype=tf.float32),) + )(self._gen_labels) + + def __call__(self, kappa, a, IbackLevel, Nparticles, sigma_motion): + """a: spot radius scale factor (1.5-4 is a reasonable range) + kappa: noise level (around 0.1 or so) + IbackLevel: intensity level of the random background relative to maximum (must be between zero and one) + Nparticles: the number of particles in the video (larger numbers means slower run time) + sigma_motion: the standard deviation of the random Brownian motion per video frame""" + ## random brownian motion paths + ## Nt, Nparticles, 3 + xi = self._sample_motion(Nparticles, sigma_motion) + + #### translate track positions to img coords + ## Ny, Nx, Nt, Np, 2 + XALL = (self.ximg[:, :, None, None, :] + - xi[None, None, :, :, :2]) + ## Ny, Nx, Nt, Np + r = tf.math.sqrt(XALL[..., 0]**2 + XALL[..., 1]**2) + z = xi[..., 2] + + ### generate video + I = self._gen_video(r, z, kappa, a, IbackLevel) + + ### generate labels + labels = self._gen_labels(r) + + return I, labels, xi + + @staticmethod + def rand(n): + return tf.random.uniform([n], dtype=tf.float32) + + @tf.function( + input_signature=( + tf.TensorSpec(shape=[], dtype=tf.int32), + tf.TensorSpec(shape=[], dtype=tf.float32),)) + def _sample_motion(self, Nparticles, sigma_motion): + #### boundaries + b_lower = tf.constant( + [-10, -10, -30.], tf.float32) + b_upper = tf.constant( + [self.Nx+10, self.Ny+10, 30.], tf.float32) + #### uniform random initial possitions + U = tf.random.uniform( + [1, Nparticles, self.d], + dtype=tf.float32) + X0 = b_lower + (b_upper - b_lower)*U + #### normal increments + dX = tf.random.normal( + [self.Nt, Nparticles, self.d], + stddev=sigma_motion, + dtype=tf.float32) + #### unbounded Brownian motion + X = X0 + tf.math.cumsum(dX, axis=0) + #### reflected brownian motion + ## note that this is imperfect, + ## if increments are very large it wont work + X = tf.math.abs(X - b_lower) + b_lower + X = -tf.math.abs(b_upper - X) + b_upper + return X + + def _gen_video(self, r, z, kappa, a, IbackLevel): + uw = (0.5 + self.rand(1))/2. + un = tf.floor(3*self.rand(1)) + uampRing = 0.2 + 0.8*self.rand(1) + ufade = 15 + 10*self.rand(1) + rmax = ufade*(un/uw)**(2./3.) + ufadeMax = 0.85 + fade = (1. - ufadeMax*tf.abs(tf.tanh(z/ufade))) + core = tf.exp(-(r**2/(8.*a))**2) + ring = fade*(tf.exp(-(r - z)**4/(a)**4) + + 0.5*uampRing*tf.cast(r<z, tf.float32)) + I = tf.transpose( + tf.reduce_sum( + fade*(core + self.ring_indicator*ring), + axis=3), + [2, 0, 1]) # Nt, Ny, Nx + I += IbackLevel*tf.sin( + self.rand(1)*6*pi/512*tf.sqrt( + self.rand(1)*(self.X - self.rand(1)*512)**2 + + self.rand(1)*(self.Y - self.rand(1)*512)**2)) + I += tf.random.normal( + [self.Nt, self.Ny, self.Nx], + stddev=kappa, + dtype=tf.float32) + Imin = tf.reduce_min(I) + Imax = tf.reduce_max(I) + I = (I - Imin)/(Imax - Imin) + I = tf.round(I*tf.maximum(256., tf.round(2**16*self.rand(1)))) + return I + + def _gen_labels(self, r): + R_detect = 3. + ## (Ny, Nx, Nt) + detectors = tf.reduce_sum( + tf.cast(r[::2, ::2, :, :] < R_detect, tf.int32), + axis=3) + ## (Nt, Ny, Nx) + P = tf.transpose( + tf.cast(detectors > 0, tf.int32), [2, 0, 1]) + ## (Nt, Ny, Nx, 2) + labels = tf.stack([1-P, P], 3) + return labels diff --git a/Training_Data/__init__.py b/Training_Data/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/Training_Data/__pycache__/Particle_Tracking_Training_Data.cpython-37.pyc b/Training_Data/__pycache__/Particle_Tracking_Training_Data.cpython-37.pyc new file mode 100644 index 0000000..a739296 Binary files /dev/null and b/Training_Data/__pycache__/Particle_Tracking_Training_Data.cpython-37.pyc differ diff --git a/Training_Data/__pycache__/__init__.cpython-37.pyc b/Training_Data/__pycache__/__init__.cpython-37.pyc new file mode 100644 index 0000000..b53a61b Binary files /dev/null and b/Training_Data/__pycache__/__init__.cpython-37.pyc differ diff --git a/particle_detection.ipynb b/particle_detection.ipynb new file mode 100644 index 0000000..6e25e81 --- /dev/null +++ b/particle_detection.ipynb @@ -0,0 +1,564 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Populating the interactive namespace from numpy and matplotlib\n" + ] + } + ], + "source": [ + "%pylab inline\n", + "%config InlineBackend.figure_format = 'retina'\n", + "from ipywidgets import interact\n", + "import tensorflow as tf\n", + "\n", + "from Training_Data.Particle_Tracking_Training_Data import Particle_Tracking_Training_Data\n", + "\n", + "from tensorflow.data import Dataset\n", + "import keras\n", + "import keras_metrics\n", + "\n", + "from keras.models import Model\n", + "from keras.layers import Activation, Input, LeakyReLU, ConvLSTM2D, Conv2D\n", + "from tensorflow.keras.optimizers import RMSprop\n", + "#from keras.callbacks import EarlyStopping\n", + "from keras.models import Sequential" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "#!pip install keras==2.6.*" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Procedurally generated training data\n", + "The code below demonstrates how to generate training videos and labels. The function also returns the ground truth particle tracks, which might also be useful depending on your goals.\n", + "\n", + "Note that the training generator is a Tensorflow Module and can be easily incorperated into a Tensorflow neural network. Alternatively, you could simply save a large set of data and use another machine learning framework.\n", + "\n", + "Note that the image dimension is fixed at 256x256. The labels are downsampled to 128x128 in the image dimensions. There are two classes (a particle is detected or not detected) per label so the label shape is 128x128x2. Hence, the neural network output should be 128x128x2.\n", + "\n", + "### Our paper\n", + "https://www.pnas.org/content/115/36/9026.short" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "Nt = 20 ## number of frames for each video\n", + "kappa = 0.2 ## standard deviation of background noise added to image 0.1\n", + "a = 1. ## scale factor for the size of particle spots (not true size of particles) 3.0 make random(2.,8.)\n", + "IbackLevel = 0.5 ## relative intensity of randomly generated background pattern; in (0, 1) 0.15\n", + "Nparticles = 10 ## the number of particles (more => slower)\n", + "sigma_motion = 3. ## the standard deviation for particle brownian motion; should be in (0, 10)\n", + "\n", + "## you might consider randomizing some of these parameters when training a neural net\n", + "\n", + "pt = Particle_Tracking_Training_Data(Nt) ## create object instance\n", + "## you can 'call' the object as many times as you want\n", + "## in this example, we only generate one training example\n", + "vid, labels, tracks = pt(kappa, a, IbackLevel, Nparticles, sigma_motion) " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Visualizing training videos and labels" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "80cfe8962449407a9a29342f870c7d87", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "interactive(children=(IntSlider(value=0, description='t', max=19), Checkbox(value=True, description='show_trac…" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "@interact(t=(0, Nt-1, 1))\n", + "def plotfn(t=0, show_tracks=True):\n", + " fig = figure(1, [14, 7])\n", + " fig.add_subplot(121)\n", + " imshow(vid[t], origin='lower')\n", + " if show_tracks:\n", + " plot(tracks[t, :, 0], tracks[t, :, 1], 'rx')\n", + " xlim(-10, 265)\n", + " ylim(-10, 265)\n", + " \n", + " fig.add_subplot(122)\n", + " imshow(vid[t], origin='lower')\n", + " imshow(labels[t, ..., 1], origin='lower')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Particle tracks" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA5UAAANtCAYAAAAAaV+nAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAWJQAAFiUBSVIk8AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3xV9f3H8fc3m4QkjCQQIEjYIwwJSxygoCC4RbG2SrWtWle1Q22rv9rhbK114l5Vq4jVFkQUURDBssLeK0BYSUgIITv3fn9/3MtNQvYBcjNez8cjj3vP93y/53wujwdw3znf8z3GWisAAAAAAJwI8HcBAAAAAICmi1AJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcCzI3wXg9DPG7JIUJSnVz6UAAAAAaJy6STpqrU2s70BCZcsQ1apVq3b9+vVr5+9CAOBkuFwupaenV2o3xqhNmzYKCwur9zGttcrPz1d+fr5KS0tr7Nu+fXuFhITU+xwAADR2mzZtUkFBgaOxhMqWIbVfv37tVq5c6e86AOCkvf7669q7d2+V+0aMGKGLLrpIQUE1//dmrVVaWppWrlyp9evXVxkmY2NjlZCQoJSUFElSWFiY7rvvPgUEcOcIAKD5SU5OVkpKSqqTsYRKAECTkpSU5AuVERERCg4O1pEjRyRJy5Yt0969ezVlyhS1b9++0tj8/HytXbtWK1euVEZGRqX9QUFBGjBggJKTk5WQkKBt27Zp8+bNys/PV/fu3QmUAABUgVAJAGhS+vfvr7lz58paq7y8PP385z/XN998o82bN0uSDhw4oJdfflnXX3+9unXrJmut9uzZo5UrV2rDhg1yuVyVjtmhQwclJydr4MCBatWqla+9d+/e+vWvf61Dhw412OcDAKCpIVQCAJqUyMhIdevWTbt27ZIk7dy5U1OnTtWyZcv05ZdfyuVyqbi4WF9//bVGjRqlr7/+WpmZmZWOExwcrKSkJCUnJ6tz584yxlR5voCAAMXHx5/WzwQAQFNGqAQANDlJSUm+ULl+/XqdddZZGjx4sL777jvl5uZKkqKiovTRRx/JWlthbOvWrTV27FglJSU5WtgHAABUxM0hAIAmp1+/fr77G/ft26esrCx99dVXvkAZHh6uiy66SK1bt6409tixY1q2bJlSUlJ07NixBq0bAIDmiFAJAGhywsPD1bNnT9/2Z599phUrVvi2J02apKioKN1yyy0aM2aMoqOjK4xPT0/Xl19+qaeeekrvv/++Nm3aVOvjRAAAQNWY/goAaJKSkpK0detWSdKOHTt87X369NGAAQMkee6/PP/88zVmzBjt3r1bq1at0saNG30B0lqrrVu3auvWrQoPD9fAgQM1ZMgQ7qEEAKAeCJUAgCappKSkUltoaKgmT55cadGdgIAAJSYmKjExUZMmTdLGjRu1evVq7dmzx9cnPz9fS5cu1dKlS9WhQwedeeaZGjx4cIXVYAEAQGVMfwUANDn79u3TnDlzKrVPmDBBUVFRNY4NCwtT//79FR8fr7Zt21bZ59ChQ5o7d65eeeUVpsUCAFALrlQCAJqU3Nxcffjhh5WeN5mYmKgzzzyzTsf48ssvlZKSUmu/7OxsFRQUKDIy0lGtAAC0BIRKAECTUFpaqr179+rtt9+utC84OFiXXXZZtc+aLC83N1dr1qypcp8xRuHh4WrVqpUiIiI0ePBgAiUAALUgVAIAGr0VK1Zo7ty51U5FHTduXLVTWU+0fPly31XO+Ph4TZ482Rckw8LC6hRMAQBAGUIlAKDRW7JkSZWBMjg4WMOHD9eIESPqdJzi4mItX77ct3322WerS5cup6xOAABaIkIlAKDRy8rKqtR2yy23qEOHDgoMDKzTMfLy8rRw4UIVFBRIkqKjo9WvX79TWicAAC0RoRIA0KgdOnSoUtuDDz6ooKC6/Rd28OBBLV26VOvWratwtXPUqFF1DqQAAKB6hEoAQKOVl5enf/3rXxXajDG1hkGXy6UtW7Zo6dKl2r17d6X9nTp10tChQ09prQAAtFSESgBAo1RaWqoZM2boyJEjFdqttcrMzFRsbGylMfn5+UpJSdGyZct09OjRSvvj4+M1cuRIDRgwQMHBwaetdgAAWhJCJQCg0bHWas6cOVVeZZSkXbt2VQiVBw8e1LJly7R27dpKC/oYY9S/f3+NHDlSCQkJrO4KAMApRqgE0DKs/1jK3CYN/5kU0d7f1aAWy5YtU0pKim973LhxCg4O1ty5cyVJqampSk5O1tatW7V06VKlpqZWOkZ4eLiSk5M1bNgwRUdHN1TpAAC0OIRKAM3fzoXSzJs979d9JE2bJUV18m9NqNb27dt94VGSkpKS1KFDBy1atMjXtnHjRu3bt085OTmVxnfs2FGjRo1iiisAAA2EUAmg+fvf9LL3h7dLb07yBMs2Cf6rCVXKysrSzJkzZa31te3bt0/r16+v1Ld8oGSKKwAA/kOoBNC8ZadKW+ee0LbLGyz/K7VL9EtZqNqqVatUWFhYoS07O7va/kxxBQDA/wiVAJq35a9LspXbc/ZIb032XLFs36PBy0Jlbre7ynsjq9KxY0eNHDlSSUlJTHEFAMDPCJUAmq/ifCnlner3H90nvXmxJ1jG9qm8vyBb+uqPnqudEx+X4vqetlIh7d27V3v37q12vzFG/fr108iRI9W1a1emuAIA0EgQKgE0X+s/lgqP1Nzn2KGyqbAdBpS1p2+S/vUDz1RZSfr+OenyF05frVCbNm0UGBgol8slSYqLi9OwYcO0fft2dejQgSmuAAA0UgH+LgAATgtrpWWv1K1vfqb01iXSgTWe7c2fSa+NLwuUkpSfdeprRAXR0dGaMGGCbzs9PV2hoaG6/vrrNW7cOAIlAACNFKESQPOUtlw6uLbqfa07Vm4ryJLeulR671rpg+ul4mOntz5Uafjw4UpKSvJtz549W+np6X6sCAAA1IZQCaDpKsyR9i6T3K7K+6q7StlrgjThkar3FeVI2744dfWh3owxuvTSS9W+fXtJUklJiWbMmKGioiI/VwYAAKpDqATQdL05SXr9QundqyVXaVl77iFpw6dVjxl9pzTgKil+SFlbTBWL9JxoyxxpxRvSnv9JBbXcp4mTEhoaqqlTp/pWdc3MzNSsWbMqPLsSAAA0HoRKAE2Tq1Q6tMHzfuc30vw/lu1LeVtyl1Qe06ar1GmoZIx0Ybn+mVvqds7Z90pvTJCeOEP6+wDp3SnSlw9Jq//luR+zpMD550EFcXFxuuSSS3zb69ev14oVK/xYEQAAqA6rvwJomgKDpPjB0oHVnu0lz0oJI6R+l0prZ1Q95sge6bHO0o3/ldZ8eHLnP5rm+dk+r6zNBEjtuktx/aS4/p6fjgN5DqZDgwcP1p49e7Ry5UpJ0ty5c9WpUyd17tzZz5UBAIDyuFIJoOnqPrbi9qe3S4d3SGG1rBL6zmXSmvcrtw+4SgqNqn7cwGulDklSQHDV+61bOrxd2jRLWviE9NE06bmh0qsXeEJsKfcF1tfEiRPVsaNnYSWXy6UZM2YoPz/fz1UBAIDyCJUAmq4e51fcLjoqzbhRmvKGNODK+h/v6D5p2E3VnGuc1GWYNP5h6fcHpDuWSVPelM67T+p7idSuhyRT9dh9K6VPbpGeHiB9/RcpZ1/9a2uhgoODde211yo0NFSSlJOTo08//VRut9vPlQEAgOMMCx9UzxjTXtKVkiZLGiips6RiSeskvSnpTWutu1z/bpJ2VTpQmQ+ttddVc65pku6Q1F+SS9IqSX+z1s4+BZ9j5dChQ4cen0IGNBslhZ77G0sLK7YPvl7av0rK2FT1OBPgmSZ71p1S6zjpuWFl92BGxku5B2o+7+UvSGf+qHJ7cb7n/sz0TZ77PdM3SqnfSa7iE84f6Dn/iFukM0Z77vFEjTZt2qQPPyybsjxu3Dide+65fqwIAIDmJTk5WSkpKSnW2uT6juWeyppdI2m6pAOSvpG0R1IHSVdJek3SxcaYa2zlZL5GUlVLT66v6iTGmL9J+pWkNEmvSgqRdJ2kWcaYu6y1z5+CzwI0P8FhnlC24+uK7VVNbZWk4Ahp6A3SyNukdoll7cN/Ki2d7nlfW6CUpKLcqttDwqVOZ3p+jjuW4Vk4aMUbniuhkmRd0sZPPT8dkqQRP/NMrQ0Jr/3cLVR8fHyF7fnz5ysiIkJDhw71U0UAAOA4rlTWwBhzgaQISZ+dcEWyo6RlkhIkTbHWfuxt7ybPlcq3rbU/ruM5RktaLGmHpOHW2uxyx1rpPX9fa23qSXwOrlSi+Vr8rDTvoZr7RHaSRt4qJf9YatWm8v68w9KzZ3qeU1mb/pdLV7/hWSioPlyl0pbPpKWvSLu/q7w/rI006ufS2Afqd9wWYtWqVfrPf/5Tqb13794aP3684uLi/FAVAADNB1cqTxNr7dfVtB80xrwk6RFJYyV9fBKnuc37+sjxQOk9R6ox5gVJD0m6SdIfTuIcQPPVYUDN+8f/URp1uxQUUn2fiPbSjZ94Vo0NjZSWvy4VZFXs03O8NOFRKbYOz7SsSmCQJ5D2v9wzNXbZK57zlXgXnSk8IuVlOjt2C5Camlpl+9atW7Vt2zYNGTJEY8eOVXR0LYs0AQCAU45Q6dzxh+CVVrGvkzHmVkntJR2W9L21dm01x7nA+zq3in2fyxMqLxChEqja/pTq943/o3TOPXU7Tudkz4/keZblBz8o23ftP6X+lzmv8UQdBkiXPuNZ9GfVe9LyV6XsVM89lqjSmjVrKmy3b99ehw8fliRZa7Vq1SqtW7dOI0eO1DnnnKNWrVr5o0wAAFokQqUDxpggSTd6N6sKgxd6f8qPWSBpmrV2T7m2CHkW/zlmra3qRq5t3tfedayruvmtfesyHmiSEsdI+kvZ9tjfeaa65h92/nzI7mOl6AQpZ6809MZTGyjLa9VWGn2n50rqvhVSbJ3+qrc4OTkVpyUPGzZMl1xyiQ4cOKCvvvpKO3bskCSVlpZq8eLFWrlypc4991yNGDFCwcHVPP4FAACcMjxSxJnHJSVJmmOt/aJce76kP0tKltTW+zNGnkV+xkqa7w2Sxx2fp1XdjVzH26u4CQyAJClhhPTDj6W2iVLnYZ5Fb1q1cR4oJc+CObd+K02bJU166tTVWp2AAM/nQJUKCyuu7nvhhZ7f2cXHx+uGG27QjTfeWGEhn8LCQs2bN0/PPfecVq9ezeNHAAA4zViop56MMXdLekbSZklnW2uzahly/Mrmd5JGSrrHWvuMt72TpH2S9llru1QxLlieR5gUWWvDTqJmFuoB0GS53W7NmzdPBw4c0IQJEyqtBHu8z8aNGzV//nxlZ2dX2BcXF6fx48erV69eMjy+BQCAKrFQTwMxxtwhT6DcKGlcXQKlJFlrS40xr8kTKs/zHkMquxJZ3coStV3JBIBmLyAgQBMmTKi1T1JSkvr27auVK1dq4cKFys/3LIKUnp6u999/X926ddP48ePVpUul3+EBAICTwPTXOjLG3CPpeXmeNXm+tfZgPQ+R4X31TX+11ubJc6WytTGm8q/epV7e1631PBcAtEhBQUEaOXKkfvGLX2jMmDEV7qlMTU3Va6+9phkzZvgW+QEAACePUFkHxpj7JT0tabU8gTLdwWFGeV93ntB+/LElE6sYc/EJfQAAdRAaGqrzzz9fd999t4YPH66AgLL/7jZu3KgXXnhBs2fPVm5urh+rBACgeSBU1sIY85A8C/OslGfKa7UPkjPGjDTGVHoYnjHmAkn3ejffPWH3S97X3xtj2pYb003SHZKKJL3ptH4AaMkiIyM1efJk3XHHHerfv7+v3e12a8WKFXr22Wf1zTffqKioyI9VAgDQtHFPZQ2MMdMk/UmSS9IiSXdXschDqrX2Le/7JyQN8D4+JM3bNkhlz6J8yFq7pPxga+0SY8zfJf1S0lpjzExJIZKmSmon6S5rbeop/FgA0OK0b99e1157rdLS0vTVV18pNTVVklRSUqKFCxdq+fLlGjNmjJKTkxUUxH+NAADUB6u/1sAY87CkP9TSbaG1dqy3/08kXSnP40ZiJAVLOiTpe0nPW2sX1XCuaZLulNRfkltSiqS/Wmtnn9ynYPVXACjPWqvt27dr3rx5Sk+veDdDTEyMbrvtNoIlAKDFYfXX08Ra+7Ckh+vR/3VJrzs819uS3nYyFgBQd8YY9erVSz169NDatWv1zTffKCfHs8h2t27dCJQAANQT/3MCAFqkgIAADRkyRAMGDNDy5cu1dOlSjRkzxt9lAQDQ5BAqAQAtWnBwsEaPHq2RI0cqMDDQ3+UAANDksPorAAASgRIAAIcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAAMcIlQAAAAAAxwiVAAAAAADHCJUAAAAAWozsjz7SgT88rJL0dH+X0mwE+bsAAAAAAGgIuV9/rYMP/Z8kqWjrVp3x/nsyxvi5qqaPK5UAAAAAmj1bWqr0vz3l2y5YtUr5//ufHytqPgiVAAAAAJq9Ix//W8U7d1Zoy3xxup+qaV6Y/goAAACgWXPn5yvj+ecqtecvX6785csVPny4JKlgzRod/fJLBbWPUcgZXRWckKCQrl0VEBbW0CU3KYRKAAAAAM1a1ttvy5WRKUkKiotT+IgROjp7tiQpc/p0dR0+XMcWLVLa7XfIlpRUGGtCQ9Xuhh8p5vbbFRAe3uC1NwVMfwUAAADQbJVmZenwa6/7tmPuulOxd98lBQZKkvKWfK9ji77Tvl/9ulKglCRbVKTDr72unZdcqtwFCxqq7CaFUAkAAOBHbrfV3qx8f5cBNFuZL06XOy9PkhTSs4faXHmlQrp2VfQll/j6ZP3zHXX5x9My3mmukRddpIizz1ZQp3hfn5L9+5V228+VdtfdKjl4sGE/RCNHqAQAAPCjv325RZOeWaSFWzP8XQrQ7FhrdWTmTN923L33ygR57gBsf+utkvdxInnfLlJAZJS6vvqKenz1lbo8+4y6vv6aes6fr/hHHlFgmza+Y+TOm6edkyYr6513ZEtLG/YDNVKESgAAAD/5eGWaXlywQ7lFpbrpzWV6e0mqv0sCmhVjjILi4sq2g4N970O7Jypq0iTfdub06QofPlwhXTpXGN/m6qvU/fM5ir7qKl+7Oz9fhx59TKnXTlXBunWn+VM0foRKAAAAP+kbH6n4aM90O7eV/vDfDXro0/Uqdbn9XBnQfESOG+d7n/vV/Ar7Ym67VZIUEBmpsH79ZK2t8hhBbduq06OP6Ix/vqOQHj187YUbNyr12qk6+Oe/yJWbexqqbxoIlQAAAH4yoFO0/nPH2RrcJdrX9s//7dZNby1XTkHlBUMA1MwWFytn9mcqWLvW1xY5vlyo/OZrWXfZL21Ce/VS5388rZ5fz1fsXXfKeKfDVid8+HB1/+Tfir3nHpnQUO9JrbLfe087J03W0c8/rzaYNmeESgAAAD+KiwrTh7eepcmDyhYEWbQtU1e9uFi7D+f5sTKgaXEXFGjPz27R/l//WqnXTlXGs8/KWqtWQ4YosF07SZIrI1OF5QKnJEVNnKjAyMg6n8eEhCjmtlvVffYsRZx7rq+9NCND++79pfbecquK9+49NR+qiSBUAgAA+FlYcKCeu+5M3T2ul69tR0aernhhsZbuPOzHyoCmwV1UpLQ771L+0qW+tswXp2v/b+6TLS1V6wvO97Xnzp9f1SHqLSQhQQmvvKzO/3haQbGxvva8RYu048KLlP73p2WLi0/JuRo7QiUAAEAjEBBg9MsLe+uZ64YoJMjzFS07v0Q/en2pPlrRsq56APVhi4u17xf3KG/x4kr7js6erT033azwocm+thPvqzwZxhhFTZyo7nM+U+SECRX2HX7lFW0eNFhH5849ZedrrAiVAAAAjcjlQzrrXz8bpZjWIZKkEpfVb2au1eOfb5bb3fLu1QJqYq3Vvvvv17EFC3xtMbffrjbXTfVtF6SkKOOZZ3zbxbt2qWjnzlNaR2BkpBRQ9f2Y++65V9snTGjW91oSKgEAABqZ5DPa6tM7zlbfjmX3eb20cIdue3el8oqqfi5eYYlLy1OztGxXVrP+8gqUV7R1m3I/L7sSGBgTo6hJF6vjH/6guPvv9z2HsvTQoQrjTuXVSknKXbCgQh0nKtm9R3lLlpzSczYmhEoAAIBGqEvbcM38+Whd0LfsGXtfbjyka176XgdyCpRfXKrvtmXq719u0dSXv1ffh+bqmpe+17Uvf6+3eN4lWojg+I4K6lS2yJUrM1M7L79Ch/7yiKKvuFxdnntWJiys0rjc+V+dshrceXk6+Kc/+bajL79cfTesV4ff/c7XFhAVpZCEhFN2zsbG8Jus5s8Ys3Lo0KFDV65c6e9SAABAPbncVo/N2aTXvtvla2sdGqTCEpdKq5kO2zEqTP/73bgq9wHNTWlWljKeeVZHPvpIKve4kIDoaMXecYfCBiYp7e675crIrDCu58KFCu4Qd+Lh6u3QY48r6+23JUmBbdqo++dzFNS2rSTP9NyStDQFtGqloJiYkz7X6ZScnKyUlJQUa21y7b0r4kolAABAIxYYYPTgJf312FUDFeS9Z+tYUWm1gVKS0nMLVVjiaqgSAb8KatdO8X98WImf/Fvho0b52t05OTr06KM68OBDir3zLoX26lVhXN533530uQvWb1DWP//p2+7w2wd8gVLyLOQTkpDQ6APlySJUAgAANAE/GNFV7/xkhKJbBfva+naM1I1nnaHnrz9T391/vjpEeR7G7rbSp6v2+atUwC/C+vRR1zffUJfnn1Nw166+9uIdO3TwD39QQGSkgjt39rWbsNCTOp8tLdWB/3vId3U0YvRZirrsMmW++qp2/+gGpT/1lAo3bmwR9zgz/bUFYPorAADNR1ZesTYfPKp+HaPUNiKkwr5Xvt2hR+dsliR1j4nQvF+OUWA1K1ICzZm7uFjZ77yjzOkvyZ2XV2FfQGSkYu64Xe2mTZMxzv9+HH7jTaU/+aQkyYSGqvus/6pg9Wrtv+/+Cv1CzjhDkZMuVtTFFyusd2/H5zvdmP4KAADQQrSLCNHoHjGVAqXkuZoZGRYkSdqZmad5Gw82dHlAoxAQEqL2P/2penwxV22umeJbBVaS3Lm5Ojz9JWW//75sadWrKdemOC1NGc8959uOufMOmcBAHfzTnyv33b1bh6e/pF2XXa4dl1yijOdfOOWPNPE3QmUNjDHtjTE/NcZ8YozZbowpMMbkGGO+M8b8xBhT5Z+fMWa0MWaOMSbLGJNvjFlrjLnHGBNYw7kuMcYs8B7/mDFmqTFm2un7dAAAoLmJDAvWj0ad4duevnBni5h6B1QnKCZG8X/+sxI/nqnwYcN87a6cHB3681+068ordWzx4nod01qrg3/8k2xBgSQptE8ftZs2Tfsf+K3cx45JkoITEhR16aUKCA+vMLZ4+w5lPv+8dk6arP333y/rah73PhMqa3aNpFcljZS0VNI/JH0sKUnSa5JmmBOumRtjLpf0raTzJH0i6QVJIZKelvRBVScxxtwpaZb3uO96z9lJ0lvGmL+d8k8FAACarZvO7qaQIM9XvDV7j2jpriw/VwT4X1j//ur6z3fU+ZlnKtxXWbRtu/b+5Kfa+/PbVZyaWqdjHZ0zR3mLFnk2jFH8n/6o7H++q/zlyz1tAQHq9OQT6vzXJ9VryWJ1fvYZhQ0cWOk4Of/5r1xHjpzsR2sUCJU12yrpMkldrLU/tNb+1lp7s6S+kvZKulrSVcc7G2Oi5AmELkljrbU/sdb+RtIQSd9LmmKMua78CYwx3ST9TVKWpGHW2justfdKGiRph6RfGWPOOr0fEwAANBdxkWG6emgX3/ZLC3f4sRqg8TDGKGrCReo+5zPF3nuvTLmriMe++UY7Lr1Mhx5/Qq6jR2s8Tvb7//K9D+3ZU5KU8Y9/+Nra33qLws88U5JUtGOHjs7+TIXr11c6TtSllyqw3EqxTRmhsgbW2q+ttbOste4T2g9Kesm7ObbcrimSYiV9YK1dUa5/oaQHvZs/P+E0N0sKlfS8tTa13JhsSY96N287uU8CAABakp+dm+i7hZSv/1EAACAASURBVGzBlgxtOlDzl2SgJQkIDVXMrbeox9zPFX3llWU7SkqU9dZb2jFhorI/+FDW7a5yfGj37r73Rdu2KXXqdbIlJZKksAEDFHv77cpPSdGeW25R6tVTlPvll1K5aejho0ap61tvqdOTT8gENI84FuTvApqwEu9r+bt7L/C+zq2i/7eS8iWNNsaEWmuL6jDm8xP61MgYU93yrn3rMh4AADQP3WNba+KAjvp8vWehnle+3amnpw7xc1VA4xIcF6dOjz2qttdfr0OPPaaClBRJkis7WwcffliuI0cUc9utlcZ1+P3vZEtLlfPJJ5X2RV9xhfbc/JOyqbDltB47VjG33apWQ5rf38XmEY0bmDEmSNKN3s3yYbCP93XriWOstaWSdskT5LvXccwBSXmSuhhjwk/cDwAAUJ3bxvTwvf/vmv1Ky873YzVA49VqYJLOeO9ddf77UwqIivK1F23bVmX/gLAwxT/6iCLGnFdp36FHHqkYKI1R5MUTlfjpJ0p4aXqzDJQSVyqdelyeRXXmWGu/KNce7X3NqWbc8fY29RwT4e1X4/8G1T1TxnsFc2hNYwEAQPMyOKGNzureXt/vPCyX2+q1Rbv08GUD/F0W0CgZYxQ2aJDc3hVdJSnq0kuq7OsuKFD6355S3sJvqz9gUJCiL71U7X/2M4V2TzzV5TY6XKmsJ2PM3ZJ+JWmzpBvqO9z7Wp+1vZ2MAQAA0K1jyiZHfbh8r7Lziv1YDdC4ZT73nOS9N7LV0KFqPWZMpT4Fq1dr15VXKfu998oaAwMVkugJjiYkRG2v/4F6fjFXnR57tEUESokrlfVijLlD0jOSNkoaZ609cY3u41cbo1W1qBP6HX8f4x1zuIYx3GEPAADqZUzvWPWLj9KmA0dVUOLSy9/u1AMXs9QCcKLCrVuV899Zvu24X96r8k8OtMXFynjxRR1+5VWp3AI+rceOVfyf/6TA9u1VuGGjgrt0VlAzWdG1PrhSWUfGmHskPS9pvaTzvSvAnmiL97V3FeODJCXKs7DPzjqOiZdn6muatZYbIQAAQL0YY3Tn+T19228u3qX9RwpqGAG0TBnPPOtboTXivHMVPmyYb1/hlq3aNfU6HX7pZV+gDIiIUPwjf1GX6S8qKDZWJiBArQYmtchAKREq68QYc7+kpyWtlidQplfT9Wvv68Qq9p0nKVzSknIrv9Y25uIT+gAAANTLxUkdNaiLZxJVUalbT8+rtDYg0KIVrF6tY/Pn+7bj7rlHkmRdLmW++qpSp0xR0aZNvv3hI0Yo8T//UZurr65wNbMlI1TWwhjzkDwL86yUZ8prZg3dZ0rKlHSdMcb36w1jTJikv3g3p58w5k1JRZLuNMZ0KzemraTfeTdfEgAAgAMBAabClNePU9K05WCuHysCGg9rrdKf/odvO2rSxQrr31/Fu3dr9w03KuOpv/ueQWlCQ9Xhtw+o61tvKqRLZ3+V3ChxT2UNjDHTJP1JkkvSIkl3V/HbiFRr7VuSZK09aoz5mTzhcoEx5gNJWZIuk+fRITMlfVh+sLV2lzHmN5KelbTCGPOhpGJJUyR1kfSUtfb70/MJAQBASzC6R4zG9onVgi0Zclvpybmb9fqPh/u7LMDv8pYsUf7SpZ6NwEDF3HWXsj/4QIeeeFK23EqwYUlJ6vTE4wrt0aOaI7VshMqaHV+uKVDSPdX0WSjpreMb1tpPjTFjJP1e0tWSwiRtl/RLSc9aayut4mqtfc4Ykyrp1/I8/zJAnsWAHrTWvn1KPgkAAGjR7p/YVwu3Zshaaf7mdC3deVgju7f3d1mA31hrlVHuKmXEOWfr0COPKu+778o6BQUp5ue3KeaWW2SCg/1QZdNAqKyBtfZhSQ87GLdY0qR6jpklaVatHQEAABzoFx+lK8/srH+n7JMkPfb5Zn1y+2juCUOLlfvlPBWuX+/bPvG5kyE9e6jT40+oVRLPd60N91QCAAC0ED8Y0dX3fvXeI1qbllNDb6D5sqWlynjmmap3GqN2N92kxI8/JlDWEaESAACgBdiwP0d3vJfi2w4MMIoMY9IaWqac//xHxTt3VrkvJDFR7abdqIDQ0AauqukiVAIAADRz327N0LUvfa/0XM9TzUICA/T01CHqHtvaz5UBDc9dVKQDv3+w2v3FO3cqczoPX6gPfj0FAADQjM1cmaYHPl6rUrdnrcDIsCC9csMwndWDRXrQMh354INa+4T17dMAlTQfhEoAAIBmyFqr57/erqfmbfW1dYoO01s3j1DvDpF+rAzwH3dhoTJfernGPp3++qSiL720gSpqHgiVAAAAzdBn6w5UCJR9O0bqrZtGqGN0mB+rAvyreM8eubKzK7SFJSUp+qorFT15sgKjo/1UWdNGqAQAAGiG9mWXPbg9KMDonZtHKC6KQImWLbRXL0VfdZUK161VxDnnKvrKKxTWu7e/y2ryCJUAAADN0NXJXTR94Q4dyS9RqdvqlW936sFL+vu7LMCvjDHq9Ogj/i6j2WH1VwAAgGYopnWoHpxcFiLfWLxLa/Ye8WNFAJorQiUAAEAzdfXQzjqnZ4wkyW2l+z9eqxKX289VAWhuCJUAAADNlDFGj145UGHBnq98mw/m6pVvq37gOwA4RagEAABoxrq2D9cvLyxbiOSZ+du0KzNPkrQj45imTF+ibg98pp++vVwHcgqqOwwAVItQCQAA0MzdfHaikjpHSZKKS9164OO1crut3l6SqhW7PY9X+GpTus567GuVMj0WQD0RKgEAAJq5oMAAPX7VIAUGGEnS0l1ZmrFir4Z2bVup7+UvLNb6fTkNXSKAJoxQCQAA0AIkdY7Wz87t7tt+ZM4mvbF4V6V+G/Yf1eUvLNZjczapoNjVkCUCaKIIlQAAAC3EPeN76Yz24ZKk3MJSrU0ruyI5IrGdQoM8Xw1dbquXv92pic98qyXbM/1SK4Cmg1AJAADQQoQFB+qxqwZWav/jZQM049az9MU95+ms7u197bsP5+v615bqvplrlJNf0pClAmhCCJUAAAAtyOgeMZo6LMG3ffPZiZo2upskqVtMhN7/2Ug9efUgRYUF+frMWJGmcX9fqM/WHpC1tqFLBtDIBdXeBQAAAM3Jn64YoMTYCEWGBem64V0r7DPG6NrhCRrbJ1YPz9qgOesOSpIyjxXpjvdTdGH/Dvrz5UnqGB3mj9IBNEJcqQQAAGhhQoMCdduYHvrhyDN8K8KeKC4qTC/+MFkv35CsDlGhvvZ5Gw/pwr8v1Lv/2y23m6uWAAiVAAAAqMGEAR0175djdP3IsiuauUWlevDT9brulf9pR8YxP1YHoDEgVAIAAKBGUWHBevTKgfrwllHqHhPha1+WmqVrX/peh48V+bE6AP5GqAQAAECdjOzeXnN+ca5uPOsMX9vhvGJl5RX7sSoA/kaoBAAAQJ2FBQeqxFV2L2VS5yj1iG3tx4oA+BuhEgAAAHW2M+OYZqzY69u+b0JfBQQYfbh8j37777Xam5Xvx+oA+AOPFAEAAECdPTVvq1zeVV9H92ivc3vFaOuhXN3/8TpJUkZukV6bNtyfJQJoYFypBAAAQJ2sS8vRZ2sP+Lbvm9hXxhit3J3ta1u8/bBKXG5/lAfATwiVAAAAqJMnv9jsez9xQEcNSWgjSdp04KivvaDEpQ37j1YaC6D5IlQCAACgVku2Z2rRtkxJUoCRfj2ht2/f5gO5Ffou35XVoLUB8C9CJQAAAGpkrdUTc8uuUl6TnKCecZG+fZsOVrwyuZRQCbQoLNQDAACAGn2x4aDWpOVIkkKCAjRpULxeW7RTPeNaq2dca+UWllbov2J3ltxuq4AA449yATQwQiUAAACqVepy669fbPFtl7jcmvbGMkmSMdK943tXGnMkv0Tb0o+pT8fIBqsTgP8w/RUAAADV+jglTTsy8nzb1qrC+1cX7axy3LJUpsACLQWhEgAAANWauTKtUlv5Wa3lp772LXdlksV6gJaDUAkAAIBqDe7ieWxIUIDRBX3j9PTUwfru/gtU1e2SN57Vzfd+2a4s2fKXNQE0W9xTCQAAgGr9fnI/TRnWRR2jwtQmPMTXPiShjVL2HPFtBxjp8iGd9OicTTpWVKqDRwuVll2ghHbh/igbQAPiSiUAAACqZYxR345RFQKlJJ3XO7bCdmJMhCJCgzT0jLa+tmVMgQVaBEIlAAAA6u3EUNkvPkqSNKIboRJoaQiVAAAAqLfBXdooulWwb9sXKhPb+9qWswIs0CIQKgEAAFBvgQFG4/rG+baTvdNeB3WJVkig5yvmzsw8ZeQW+aU+AA2HhXoAAADgyAMX91VIUIC6x0ZoZGI7SVJYcKCGJLTxPadyeWqWJg2M92eZAE4zrlQCAADAkbioMD1+9SDdcl4PGVP2jJHhidxXCbQkhEoAAACcMtZalbrLnk+Zsifbj9UAaAhMfwUAAMApkVNQovtnrtXcDQd9ba1D+boJNHf8LQcAAMBJW5t2RHe8n6K9WQW+tqTOUXpyyiA/VgWgIRAqAQAA4Ji1Vm8tSdWjczapxFU27XXaWWfod5P7KTQo0I/VAWgIhEoAAAA4klNQovtmrtEXGw752iJDg/TElEGs+Aq0IIRKAAAA1NvqvUd05/spSssum+46sHO0nr/+TJ3RPsKPlQFoaIRKAAAA1Jm1Vm8uTtVjn1ec7vrj0d3020l9me4KtECESgAAANRJTn6JfjNzjb7cWG66a1iQ/jplkCYmMd0VaKl4TiUAvzi875i+/2S7Mvbk+rsUAGix9hzO1y3vrFBOQUmtfVfvPaJJzy6qECgHd4nWZ3edS6AEWjiuVAJocK4St2Y9u1p5OcVa83Warrj3THXsHu3vsgCgRdmwP0fT3liuzGNFyn57uf75k5EKC6566mpqZp6ue+V7FZa4fW03n52oBy7uq5AgrlEALR3/CgBocDtWpSsvp1iSJ2B+9uJaHUnP93NVZay12rL0oLYsPShrbe0DAKAJ2p5+TJnHiiRJy1Ozdef7q1TqclfZd/PBo75A2To0SC/fkKz/u7Q/gRKAJEIlAD/YsGh/he3CYyWa9dwaFeQWy+VyK/9osZ8q81g2a5e+enOjvnpzo3atzvRrLQBwulw+pLMenNzPt/3VpkP63SfrqvxlWruIUN/7nnGtNWFAxwapEUDTQKgE0KCy9udp/7YjkiQTYBQU7Pln6GhGgd74zXd69Z5vteC9zf4sUSvmpPre71yT4b9CAOA0++m53XXbmB6+7Rkr0vTkF1sq9YuLLAuVGblFDVIbgKaDUAmgQe1Yle5737l3G134kwGSKdvvKnFr97rDKi4o9UN1qnSVtNvAGO1am6msA3l+qQcATrf7J/bRNcldfNvTF+zQ69/tqtAnLqpiqOTWAADlESoBNKiomFa+9/u3H1G7+Aide23vCn3cbqtda/0z7XTn6opXJjctOaA5L67VzCdWKDer0C81AcDpZIzRY1cN1Ph+cb62P8/eqE9X7fNth4cEqXWoZ33HYpdbR/JrXy0WQMtBqATQoHoP76C4blGSJHep1cJ/bdHAsZ015MKuFfp99eZGf5SnlC92V9jes+GwJKmk0KVUb9DNzSrUztUZysthChiA5iEoMEDP/WCohp3R1tf264/WaMGWstkl5afApjMFFkA5hEoADcoEGI29vo+Md8pr2uZsbV+RrtFX9lBs10hfv7DWwQ1eW3FhqXIPV3818tsPtmrmEyv0zu+W6POX1unjJ1bKupkCBqB5aBUSqNenDVefDp5/i0vdVj9/N0Wr9mRLkmIrhEpmbgAoQ6gE0OBiu0Zq4Niy+3e++2ibiotcuvo3yerS1/Nb8oR+7Rq8rj0bsmrtc2jXUd/73KxCFeYxBQxA8xEdHqy3bx6hzm08tyoUlLh081vLtT09V3FRYb5+6Ue5UgmgDKESgF+MvKy7wqNDJHkWx1n6350KDA7QZb8Yoh8/frYuvLl/g9d04v2UtWkVFaLQ8KDTVA0A+EfH6DC985MRahfh+Tc6O79EN76+TO5yMzOY/gqgPL4NAfCLkFZBOueaXvrytQ2SpPUL0tTvrHjFdo1URJuyKValJS598eoG7d2UJVlJxrtYrJFkjO+9bwFZYxQSFqiE/u3UMzlOXfq0VUBg7b8/yz9arF31eHyIMdKY63rX6dgA0NT0iG2tN388XD949X/KL3Zpf06h9q874NvP9FcA5fFtCIDf9EyOU0J/zzRXa6UF722u8JtwScrYnavUtZlylbjlKnXLVeJWaYlbpcVulRa5VFLkUkmhS8XHfwpKdSy7SJsWH9CsZ9fozfsW65v3Nmvv5iy5Xe5qa1k9b49Kiyvvj4oJU3Rsq4qNRhp/U3/1GBpXqT8ANBeDE9ropR8lKzjQVNrHlUoA5REqAfiNMUbnXddbgUGef4rSd+dq46J9FfoEhQae1DkK80q0cdF+/fcfq/XWA4u14P0tStuSXSG8FuQWa93CtCrHT7xloM6/oW+Ftgtu6KfeIzqeVF0A0BSc1ztWf7tmsG9xteMOHCnwT0EAGiWmvwLwqzZx4Ro68Qwtn+150Pb3n+5U4pBYRUR7psBGti1bGCIoNFA3//UcyarswdvWMyvW897KWin7YL52rEzX9pR05R0p+216QW6JNny7Txu+3adWUSHqeWaseg6LU+q6w1VepRx9VU/Fdo2UtVZ9R3VU2tZsjbqsu/qMij8tfxYA0BhdPqSzsvOK9fCsskc9pew5ov1HCtSpTasaRgJoKQiVAPxu6ISu2rrsoHLSC1RcUKrFM7frop8MkCSFRgQpKDjAM+W1yCV3qVuh4TU/biS+R7Tie0Tr7Ck9dXBnjrZ7A2Z+TrGvT8HRYq1buE/rFu6r8hhd+rbVkPEJkjxXVMf9uOEXDgKAxuLHZycqPbdILy7Y4Wsb/fjX2vHoJAUGVJ4eC6BlYforAL8LCg7UmB/08W1vW35Iezd6Hu9hjFHrdmVXK49l1/0+HhNgFN+zjc6d2lvTHjtbV/7qTA0c01mtokJqHBfWOljjb+ovwxclAPD5zYQ+6hEbUaHtule+91M1ABoTQiWARiGhXzv1Gt7Bt73wX1tUWuKSJLVuW7YabG6WsxUHAwKMOvVqq/N+0Ec/fvxsXXHvmRpwbqcq+15wYz/f9FsAgIcxRl/cc16FtsN5xdX0BtCSECoBNBrnXNPL99zHnIwCrfx8t6SKobI+Vyqr4ypxK/tgnvZtPVJpX9KYzkocFHPS5wCA5igoMEAb/zRBMa1DFRhgdOmgqn85B6Bl4Z5KAI1GeFSIRl3RQwvf3yJJSvlit3qP6KDW5RbrKb/wTn3lZhVq/cI0bVi0X0X5pZX2t42P0NlX93R8fABoCcJDgrTiwfEqLHEpLPjkVugG0DwQKgE0KgPO6aTN3x/QoV1H5XZZLXhvS4VpscccTH89tOuo1szfo+0pGbInPAfzuC5922rM9X0UFMIXJACoCwIlgOMIlQAaFRNgNPaHfTXj0eWybqv9244oIrpsYZ3cekx/tdZq9nNrtMe76E95UTFhGnRBgvqNjldIGP8UAgAAOMU3KQCNTkyX1ho8LkGr5+2RJG1bme7bdyy77lcqjTFq2zGiQqjs3LuNBl2QoG6DYhTA6q4AAAAnjVAJoFEacUmitq88pGNZRVK5GavHsotkrZUxdQuEA8/vovWL9qlXcpwGXZCg2K6Rp6liAACAlonVXwE0SsGhgTrvuj6V2l0lbhXmldT5ONGxrXTTk+do3I/7EygBAABOA0IlgEYrcVCMEgdXfrxHfR8rEtqKSRkAAACnC6ESQKN27tTeCg6tuMLgqXhWJQAAAE4NQiWARi2yXZhGXJpYoc3JY0UAAABwehAqATR6g87vopiE1r7tovxSP1YDAACA8giVABq9gMAAXXBjP7WKDFZoeJDOSGrv75IAAADgxeoVAJqE2IRITXv0bFlZBQUH1j4AAAAADYJQCaDJCAxmcgUAAEBjwzc0AAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOEyloYY6YYY54zxiwyxhw1xlhjzLvV9O3m3V/dzwc1nGeaMWaZMeaYMSbHGLPAGHPJ6ftkAAAAAHDygvxdQBPwoKTBko5JSpPUtw5j1kj6tIr29VV1Nsb8TdKvvMd/VVKIpOskzTLG3GWtfd5B3QAAAABw2hEqa3evPGFvu6Qxkr6pw5jV1tqH63JwY8xoeQLlDknDrbXZ3va/Slop6W/GmNnW2tT6lw4AAAAApxfTX2thrf3GWrvNWmtP0ylu874+cjxQes+bKukFSaGSbjpN5wYAAACAk0KoPD06GWNuNcb8zvs6qIa+F3hf51ax7/MT+gAAAABAo8L019PjQu+PjzFmgaRp1to95doiJHWWdMxae6CK42zzvvauy0mNMSur2VWX+0ABAAAAoN64Unlq5Uv6s6RkSW29P8fvwxwrab43SB4X7X3NqeZ4x9vbnPJKAQAAAOAU4ErlKWStTZf0fyc0f2uMuUjSd5JGSvqppGfqe+g6nj+5qnbvFcyh9TwnAAAAANSKK5UNwFpbKuk17+Z55XYdvxIZrarVdiUTAAAAAPyKUNlwMryvvumv1to8SfsktTbGxFcxppf3detprg0AAAAAHCFUNpxR3tedJ7R/7X2dWMWYi0/oAwAAAACNCqHyFDLGjDTGhFTRfoGke72b756w+yXv6++NMW3Ljekm6Q5JRZLePOXFAgAAAMApwEI9tTDGXCHpCu9mR+/rWcaYt7zvM621v/a+f0LSAO/jQ9K8bYNU9pzJh6y1S8of31q7xBjzd0m/lLTWGDNTUoikqZLaSbrLWpt6Sj8UAAAAAJwihMraDZE07YS27t4fSdot6Xio/KekKyUNl2fqarCkQ5JmSHreWruoqhNYa39ljFkr6U5Jt0hyS0qR9Fdr7exT91EAAAAA4NQiVNbC/j979x1fd1n3f/x9nZU9m3SPdKV70DLaCkgZZasIqCDKfYuI/gRu4eZG9FZERVRuB4J4q7gYIqjcoAhl2pZVinTQvWeSpkmaZo+zrt8fSU9zkrRJvs3JOUlez8eDxznX9R3ncwqPkneu8bX2Hkn3dPPc30r6rcPPeVTSo06uBQAAAIB4YU0lAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMCxmIVKY8wTxphhsbo/AAAAACD+YjlSea2k7caY24wx7hh+DgAAAAAgTmI9/TVd0o8krTPGnBPjzwIAAAAA9LG+WFNpJM2Q9Lox5kljzMg++EwAAAAAQB+IZajcqZZAaVv/MZI+KWmrMeYOY4wnhp8NAAAAAOgDsQyVMyXdLalJx8KlUcuU2B9K+sAYc24MPx8AAAAAEGMxC5XWWr+19l5J0yX9Qx1HLadJetUY85QxZlSs6gAAAAAAxE7M11Raa/dZaz8i6SOS9qpjuLxaLVNi72RKLAAAAAD0L32xUY8kyVr7D7Vs2HOvJL+ip8SmSfq+pA3GmPP7qiYAAAAAwMnps1ApSdbaJmvt3ZJmSXpZHUctp0h62RjzF2PMmL6sDQAAAADQc30aKo+y1u601l6slqmvReoYLj8uabMx5i5jjDceNQIAAAAAuhaXUHmUtfYZtWzYc7+koDpOif2eWqbELolbkQAAAACA44prqJQka22DtfYuSXMkLVPHUctCSUuNMc8YY8bGr1IAAAAAQHtxD5VHWWu3WmvPk/RpSaXqGC4/ppYpsf9tjPHFr1IAAAAAwFEJEyqPstb+SS0b9jwgKaToKbGpkr4jaaMx5uK4FQkAAAAAkJSAoVKSrLV11trbJc2X9JY6jlpOkvQPY8yzxpiCeNUJAAAAAINdQobKo6y1GyR9VNKf1BImpehw+RFJm4wx32RKLAAAAAD0vYQKlabFbGPMF4wxvzPGbJZULulTOhYk24ZLIylF0j1qmRJ7XhzKBgAAAIBByxPPDzfGDJW0oM0/p6rlUSKRU9q8PzpCqXZ9R8+bJOkVY8zjkm611tbEpGgAAAAAQESfhUpjjFfSKYoOkePan9aubU9wbI+kDEl5ig6Xn5F0tjHm09balb1QOgAAAADgOGIWKo0xYyQt1LEAeYqktuseOwuQ7UPk0XOCktZJelstG/e8ba0tNcakS/qqpK+oZYTz6JTYAknLjDE3WGv/2ItfCwAAAADQRixHKvfpxCON7aeytj1eI+ldtQZISauseZA5gwAAIABJREFUtQ3tP8BaWyfpm8aYX0j6H0nXtrmvT9Kjxphma+1fHX8LAAAAAMBx9cX017aPA2nb1zZE7ldLeDw6ErnBWts+dB6XtfagpOuMMX+U9KikIa2f55L0G2PMMmvt4ZP6FgAAAACADvoiVB6dknpUSNIGRU9lLeqVD7J2qTFmgaQ3JI1o7c6Q9HlJP+yNzwAAAAAAHBPrUGkk1UlapWMjkStbp63GhLV2tzHmK5L+rGOjoxeJUAkAAAAAvS6WofIrahmJXGetDcfwczqw1v7VGFOiltFKI2l6X34+AAAAAAwWMQuV1toHY3XvbtoqaWTr++x4FgIAAAAAA5Ur3gXE0JE27/vseZwAAAAAMJgM5FAJAAAAAIixgTyC95ykZEnzJQ2Lcy0AAAAAMCAN2FBprf2jpD9KkjFmRBenAwAAAAAcGBTTX621B+NdAwAAAAAMRIMiVAIAAAAAYoNQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAAAAAHCMUAkAAAAAcIxQCQAAAABwjFAJAEgYoVCTNmy8Re+/f5UaGvbFuxwAANANhEoAQMI4cmSlyspeVHXNWq1899x4lwMAALqBUAkASBgNDbuj2uGwP06VAACA7iJUAgASRnb26VHtqqp/xakSAADQXYRKAEDC8PnyotqHD6+IUyUAAKC7CJUAgITh8w2JapdXvB6nSgAAQHcRKgEAcRcIVMvakFwun1yu5Eh/Y+NedoEFACDBESoBAHFlbVgbNn5Zq9dco4aGPUpOHhV1/PDhZXGqDAAAdAehEgAQV8XFf9KRIytVXb1a7/3rozIm+n9NlUdWxqkyAADQHYRKAEDcNDYe0M5dP4i0R4+6Tmlpk6POCQQq+7osAADQA4RKAEBcWBvWli13KRRqkCSlpk7S+PH/0WGznlCoKR7lAQCAbiJUAgDiorj4SR2pere15dL06ffL7U6Szxv9WJFwuLHviwMAAN1GqAQA9LmWaa8/jLTHjb1RWZlzJHV8VmUoRKgEACCRESoBAH1u795fRKa9SlJ9wy6Vlb+sUKi5k1DJ9FcAABKZJ94FAAAGH5crKapdUfGaKipek8eToaSk4VHHgsEqrV//RdXUbtCUKd9Rft55fVkqAADoAiOVAIA+N2nS1zRx4p1KSyuM6g8Ga1Vfv6PD+eUVr6q5uVSbNt2mhoZ9fVUmAADoBkIlAKDPud1JKhh3kxacsVRnnP6iCsZ9ScnJY7q8LhSq1+bN/6lwONgHVQIAgO4gVAIA4io9fYomTrxDixYu06nzn9GY0f92wvOra9Zq+YrpamjY0+FYTe1GrXrvUq374Ab5/RUxqhgAALRFqAQAJARjjLKy5qqw8JtdjlpaG9LKd8/Xug9uUCBwRFLLhj4bN96qurqtOnx4udasvU7NzeV9UToAAIMaoRIAkHDS0iZGtSdNvFPDhl7W4bzDh5frX+9/XPX1O7V378NqbDy23rK+fofWrL1Wzc2HYl4vAACDGaESAJBwJk64Paqdk7NIM2f+TAvOeLnDuY2N+/Xuqgu1d98vOhxraNit1WuuUVNTScxqBQBgsCNUAgASTijUKMlE2h5PuiQpLW2SZkz/aZfXe725kfeNjfu0dt2/ydpQr9cJAAAIlQCABHSg6FFJNtJuuzZy+PCPaNiwj5zw+kCgMqrd0LBLjY37e7VGAADQglAJAEg4mRkzo9pr1l6jbdvvUTBYL0maUvhtJSeN7Na9PJ4sjRnz70pJKejtMgEAgCRPvAsAAKC9sWO/IK93iHbsvFfBYK0kqajocVVULNO0qd9Xbu4iTZ/+Y61Ze63ajmiOG/clpaYUKDl5hJKSRig5ebjc7tQ4fQsAAAYHQiUAIOEYYzRy5FXKHXKmtm79hg4fXiZJamoq0tp1n9Gokddo0qSvanzBzdqz9yEZ49Xcub9Xbs7COFcOAMDgQ6gEACSs5KThmjP7EZWWPqftO76rYLBaklRc8icdPrxCU6fep3mn/Eleb7bS0wvjXC0AAIMTayoBAAnNGKMRI67QgjNeUn7eBZH+puYSrfvg31Ra+qySk0fEsUIAAAY3QiUAoF9IShqqWbP+VzNmPCCvNyfSX3Lwz3p31UWqqFgWx+oAABi8CJVdMMZcZYx5yBjzpjGmxhhjjTFPdHHNImPMi8aYSmNMgzFmvTHmK8YY9wmuucwYs9wYU22MqTPGrDLGXN/73wgA+i9jjIYPu1wLznhJQ4deEulvbi7VB+s/r02b71AgUB3HCgEAGHwIlV37hqSbJc2VVNzVycaYj0p6Q9LZkp6V9LAkn6SfSnrqONfcLOl5STMlPSHpEUkjJf3BGPOjk/8KADCw+Hx5mjXzIc2c+XN5vbmR/tLSZ/XuqgtVXv5qHKsDAGBwIVR27TZJhZIyJX3pRCcaYzLVEghDks6x1t5grf0vtQTSlZKuMsZ8qt01BZJ+JKlS0qnW2i9ba2+TNFvSLkn/aYxhO0MA6MSwoRdrwRkva9iwyyN9fn+51m/4ojZu+or8/so4VgcAwOBAqOyCtXaZtXaHtdZ2fbaukpQv6Slr7ftt7tGklhFPqWMw/ZykJEk/t9bubXPNEUn3tTa/6LB8ABjwfL5czZzxgGbP+qV8vvxI/6FDz+vdVRfp4MFnFQ4H41ghAAADG48U6V3ntr6+1MmxNyQ1SFpkjEmy1jZ345ql7c45IWPM6uMcmtqd6wGgP8vPv0DZ2adp+457VVr6rCQpEDiszVvu0O49P9O4sTdq5Mir5XL54lwpAAADCyOVvWtK6+v29gestUFJe9QS5Cd085qDkuoljTbGpPZuqQAw8Hi92Zox/UeaM/s3SvINi/Q3NR3Qtu13640356u84nVVV69VMFgXOV5Tu1ElJX9WMFgfj7IBAOjXGKnsXVmtr8fbevBof3YPr0lrPa/hRB9urZ3fWX/rCOa8E10LAANJXt5iLVjwsvYf+IOKih5VIHBEkhQKNWj9+i9Ikk6Z+7hycxepufmQ1qy5RqFQg2pqN2jqlO/Gs3QAAPodRir7lml97c76zJO5BgAGPY8nQxPG36IPLXpDI0d8osNxr69l19j9+3+rUKjld3bFxU/2aY0AAAwEhMredXS0Mes4xzPbndeTa2pOoi4AGLTc7lRNnvwNJSWNiOr3eXMkSXX1HVYfAACAHiBU9q5tra+F7Q8YYzySxksKStrdzWtGqGXqa5G19oRTXwEAx+fxpGlK4bei+urrd0qSGhp2d3YJAADoJkJl7/pn6+tFnRw7W1KqpHfa7Pza1TUXtzsHAOBQdvbpUe3tO76rQKBGTU3FUf3WhvqyLAAA+j1CZe/6q6QKSZ8yxpx6tNMYkyzp3tbm/7a75veSmiXdbIwpaHNNjqSvtzZ/GaN6AWDQCAQqo9r19Tu0afPtHc4LBlltAABAT7D7axeMMR+T9LHW5vDW14XGmD+0vq+w1t4hSdbaGmPMjWoJl8uNMU9JqpT0EbU8OuSvkp5ue39r7R5jzH9JelDS+8aYpyX5JV0labSkH1trV8bq+wHAYHF0B9i2Dh9e1sl51fK2rrcEAABdI1R2ba6k69v1TdCxZ03uk3TH0QPW2ueMMR+W9N+SrpSULGmnpNslPWit7bCLq7X2IWPM3tb7fFYtI8ibJX3DWvtor34bABik/J2Eys4wUgkAQM8QKrtgrb1H0j09vOZtSZf08JrnJT3fk2sAAN0X8B8Llampk9TQsEudPa0pQKgEAKBHWFMJABgU2q6pzM1dpNGjr+v0vGCgutN+AADQOUIlAGBQaDuttajoCTU1Hez0vECQUAkAQE8QKgEAg0J2zgId+99eWBUVr0mS0tImKz//2FOdggGmvwIA0BOESgDAoDAk90ydOv8vys09K6q/vn6HystfirSDjFQCANAjhEoAwKCRlTVXp8z9Q6fh8iimvwIA0DOESgDAoJOVNe+44dLl8sWpKgAA+iceKQIAGLSOhsvq6jXav/93CgSrNWZ0+0cTAwCAEyFUAgAGvayseZo1a168ywAAoF9i+isAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJQAAAADAMUIlAAAAAMAxQiUAAAAAwDFCJRKaDVk1bjqs5t3VsmEb73IAAAAAtOOJdwHA8dhgWIcf36ymbUckSe4sn1LmDlXaKUPlHZ4W5+oAAAAASIRKJCgbDOvwE1sigVKSQtV+1a0oUt2KInmHpyn1lKFKnT9U7nRfHCsFAAAABjemvyLh2FBYh5/cqqatlZE+k+SOOidQWq/qpXt06GdrFKrz93WJAAAAiJGQtXru0BHdu6tExU38nNcfMFKJhGJDYVX+aauaNh+O9GUsHqPM88eqaUeVGtaWqXHTYSkYliSFawMKljU6Hq0M1frVvKda3qGp8gxNlXGZXvkeAAAA6BlrrZZV1ureXSXaXN8kSfp7WZXeOH2qkt2MhSUyQiUShg1ZVT69TY0bjwXK9A+PVuaScTLGKGVqrlKm5ircEFDJve9KLblSnvwUZ58Xtir/1XoFKxojfSlz8pV1UYE8Ockn9V0AAADQfetrG/TdXSV680hdVP/+Jr8eKSrXLeOGxakydAeRHwnBhq0q/7xNjesrIn3pZ45S1kUFMiZ69DBY7Y8ESneWT+4MZ6OU4fpAVKCUpMYPylXz8l5H9wMAAEDPlDYH9P8279OS97dHBcq2IeVn+w6p3B/o++LQbYRKJISm7UfU+EF5pJ2+aKSyLh3fIVBKUqCoNvLeOzqjx59lw1aNWw6r7JcfdHrck+ds5BMAAAA9c9363fq/Q8c2ZnQb6TMjh+hfC6erMLVl5lhdKKwf7i6NV4noBqa/IiGkTM1V1iXjVf3iHqUtGKGsyyd0Giglyd8mVPpGp3f7M8L+kBpWH1Ld2yUdRijbqlt5UOGGoFLnDZV3VPpx6wAAAMDJqQ6GIu8npCTp0VnjNTmtJUzeM2mkrl2/W5L0x4OH9e+j8zQjnV/+JyJCJRJGxtmj5R2ZpqQJ2ScMcv7iY1MjfN0YqQzVB1T7RpHqV5XKNgW7PD9cH1DdOyWqe6dEnqGpSp03VKlzh8qTndS9LwIAAIBuub1gmG7bekCSdKDdTq/nDsnU4twMLauslZX0rR3F+svcifzCPwEx/RUJJXlSzgl3YLXBsAIH6yNt36iuRyorfr9RdSuKuhUo2wuWNajmpb0q/cF7qnxqq2wg3ON7AAAAoHOfHJ6r07PSJEkBa/W17UWy1kaO3zNplNytPxq+VVWnVw7XxKNMdIFQiX4lUFovhVr+onHnJsuV6u3ymmDZsamurgxvh2dedlfDunI17+cvMgAAgN7iMkY/KBwdFRyfK6uKHJ+SlqzPjsyLtO/ZWSx/mF/yJxpCJfoVf1Hbqa/dW0/pK8iMvA/XBmSbQyc4+/jcQ5K7NTIKAACA7puenqLPj8qPtL+1s1g1bdZa3lEwXJmeltiyp9Gv3xdXdLgH4os1lehXojfp6Xo9pQ1bBUrqOj3mG5OhpEnZsv6QjNctk+SSy+eWSXLLtL66jr73ueXJTT7h1FwAAAA4c8f44fpbWZVK/QGV+YO6f89B3Tt5tCRpiM+j28cN1z27SiRJP9l7SGOTfVqQna4cL3EmEfBvAf1KoIcjlVV/26lwXfRzjVzpXmVdNF6p84YSEgEAABJAhsetb08eqZs27ZMk/a6oQp8cnqtZGamSpM+NztOjJRXa0+hXdTCkf9+4V0bStLRkLcxO18LsdC3ITleej3gTD/ypo98I+0MKlLVu0mMk78gTh0p/Ua3qV0U/0yj9zFHKPH+sXMn8pw8AAJBIPpKfrSdzKrXiSK3Cku7aXqTn502Wyxj5XC59Z9IofWbDnsj5VtLm+iZtrm/Sb1unxBamJmthdpoWZqdrUXa6hiZ1vf8GTh4/WaPfCBysl1rXZXvyUk4YDMPNQZX9fF1UX97nZyp5Uk4sSwQAAIBDxhjdVzhKi9/bJr+1Wl3ToD8drNSnRw6RJF2Ql6W/nTJJL1fU6J2qOm2oazi6f2PE9oYmbW9o0qMlhyVJE1OSWkcyW4LmyGRfX3+tQYFQiX6j7XpKSWrafkRJE7JkPB33m6r6266otmdYKoESAAAgwU1MTdaXxw7VT/cdkiTdu6tEF+VlaUjrtNYzstN1RnbLbLW6YEjvVddrZVWdVlbVaV1tg4LtQuauxmbtamzWEwdbQmZBii8yXXZhdrrGEDJ7BaES/Ua4/tjayGB5oyp+t1Em2aOUablKmTFESYU5cvncalhbpoY1ZVHXDrlmal+XCwAAAAduHTdMzxw6ov1Nfh0JhvS93SX6ydSxHc5L97h17pBMnTukZaf/+lBIq6sbIiFzTU2D/DY6Ze5t9GtvY6X+dLBSkjQ62asbRuXrS2OHxv6LDWCESvQbqXOHqmFduUKVTZE+2xRsCZFry2S8LnnyUlqmybaRMidf3uFpfV0uAAAAHEhxu3Tv5FH6bOv6yScPVuqaEUN0WtaJf55Lc7t1dm6Gzs5teUJAYyisNTX1WlnVMpq5uqZeTeHokFnUFNC3d5Xoo0OzmRp7EnhOJfoN79BUDb/jVOXfNFvpZ46SOzsp6rgNhDsESknKvGBcX5UIAACAXrAkL0sX52VF2l/ddkDBdoGwKylulz6Uk6E7xg/XA9PG6pPDc9XZvv8TU5Ii02vhDH966FeMyyhpfJaSxmcp69LxCpTUq+GDctW9UdTp+amnDpM3L6WPqwQAAMDJ+s7kUVpeWavGcFib65v0++IK3Tgmv0f3KGny67u7SvT38qoOm/oM93n1+dF5+uyoPCW5GGs7Gfzpod8yxsgGQmraWtnpcXdusrIYpQQAAOiXxiT7dHvBsEj7h3sOqrQ5cIIrOrpp0z49WxYdKKemJetnU8fqvYXTdPO4Ycr0uHur5EGLkUr0S+GGgKpf2qv696KfQ+kdla6cKybJk58i43HJuPm9CQAAQH9105h8/bm0UjsamlUXCuuencX65YyCbl+/t6k58n5RdrpuHjtUi3MzZExnE2HhFKES/U7DhnJV/W2XwnVtflPlcSnrwgKlLxop4+76LwkbtgqWN8i/v1aBikYlF+YoeWJ2DKsGAABAT/lcLn2/cLSuWtfyuLjnyqp07YjayGY8XUltM631x1PGaHxq0gnOhlOESvQr/qJaVT65VWq/TjsYVvVLe1T3VpFc6T65071yZfjkTvfJle5VuCEgV5JH4eag/Ptr5T9QK9scilxe91axRn5roVw+pj8AAAAkkjNzMnTlsBw9c+iIJOlr24v0z9OndGsdZEqbWWsN4XDMahzsCJXoV2zIdgyUR4WsQtV+har96tls+5ZrbXNIIlQCAAAknG9NHKlXKqpVGwprV2Oz/nd/mb5SMLzL61LbhMrGEKEyVlhwhn4laVymhlw3TWkLRyhlVp58BZny5KXIJPU8DLpS2/xOxWPkSvP2YqUAAADoLUOTvLprwohI+4F9h7SvsfkEV7RoO/21gVAZM4xUot9JmZmnlJl5HfptIKRQbUChOr/CtQH5i2tV+88D0dfOzlPKjDz5xmbINod06IE1kiRPdrKMiwXbAAAAierfRuXpD8UV2tHQrKaw1e+LK3TPpFEnvCZq+iuhMmYIlRgwjNctT65bntxkSVLKjCHKOGu0/MW18o3N7LBesrHNo0jc2SzaBgAASGRraxq0p83o5HBf17PMoqa/sqYyZgiVGNBcKR4lT8rp9Fio6thfSu4sQiUAAECiqgwEddOmvQq27q0xLzNVnxvdceZae6mMVPYJ1lRi0IoKlYxUAgAAJKSwtbpl834VN7dsxZjtcetXMwrk687ur6yp7BOESgxaoaqmyHtPDqESAAAgET28v0yvV9ZE2g9OG6sxyb5uXcv0175BqMSgZMNW/gO1kTYjlQAAAInn3ao6/WDPwUj7/40ZqiV5Wd2+numvfYNQiUGpafsRBQ+3jFSaZLd8YzLjXBEAAADaKvcH9MVN+xRqXUd5WmaavtbmsSLdwfTXvkGoxKBU93Zx5H3aacPlcvCcSwAAAMTG0XWUpf6WdZS5Xrd+NWOcvD18BBzTX/sGoRKDTuBQvZp3VLU0jJS+cGR8CwIAAECUh/eXafmRY0uVHpo2TiO7uY6yLaa/9g1CJQadurdLIu9Tpg+JPNcSAAAAieHp0mPPEy9I8WleZqqj+7Sd/vpedb32tXnOJXoPoRKDSqg+oPo1ZZF2+pmj4lgNAAAAOrNkyLHNePY2+nXh+9u1ua6xx/cpTEuOBJ4DTX5dvHq73quq66UqcRShEoNK/XulUrBl6oN3VLp8BWzQAwAAkGi+MXGEbi8YFmnvb/Lr0tU79Peyqh7dpzAtWQ9PHyefaVmLWRkI6ap1u/RMm5FQnDxCJQYNGwqrfuWxqa/pi0bKmJ4t9gYAAEDsuYzRneNH6PczC5TWui6yMRzWFzbt1fd2lShkbbfvdcWwHD1zyiQN8XokSX5r9eUt+3X/noOyPbgPjo9QiUHDX1ynUH3LDmKudK9S5+THuSIAAACcyMX52XpxfqEmpBx7pvhD+8t03frdqgoEu32f07LS9OL8ySpMPbaXxk/2HtKXNu9TExv4nDRCJQaNpLGZGnHX6co8f6wyFo+R8fCfPwAAQKKbkpaspfMn67zcY8uWllXWatGqLdpa3/11luNSkvSP+ZN1Tk5GpO+5sipduW6nylsfXQJn+Kkag4o7w6fM88cp40Ns0AMAANBfZHk9emz2eH1l3LF1lpWBkM55b5se2neo2/fJ9Lj1xOwJun7kkEjf6poGXbx6u7Y42AgILQiVAAAAABKe2xjdNWGEfjJ1TFT/93Yf1PBl61QTDHXrPh6X0Q8KR+u7k0ZFwlBRU0BXr9ul+m7eA9EIlQAAAAD6jWtHDOm0v/DNDXqmtLJbm+8YY3TjmHw9Omt8JBBVBIKq6ME6TRxDqAQAAADQr7SdBtvWl7fs18fW7tTG2oZu3WdBdrqObtPjNUajkny9VOHgQqgEAAAA0K80ho+/Y+uq6noteX+77tpepMYudnbd0dAUeT8+JUkeF4+bc4JQCQAAAKBfKfd3nKaa63XL05oJw5L+UFyhz23co+YTBNCGUFgTUpLkkjQ5Lem45+HEPPEuAAAAAAB6oqw5+hEgczJStHR+oXY1NOubO4q1/EitpJZHj3xx0z79ekaBvJ2MQp6Zk6F3FkxTczis2iDPq3SKkUoAAAAA/cpbVXVR7e9NHi2XMZqclqw/zZmg29qsuVxaUa1btuxT6AQb+CS5XMrzMd7mFKESAAAAQL8RDEeHw3NyMnRqVlqkbYzRneOH60tj8iN9z5VV6bat+xXuxs6w6DlCJQAAAIB+47fF5VHtH7d7bqXUEizvnjhSnxuVF+n7c+kR3bW9qFuPHEHPECoBAAAA9AuVgaC+tbMkqm9UcuePATHG6N7Jo/TpEbmRvsdKDutbO0sIlr2MUAkAAACgX7h/T2lUe3Jqxx1b2z5GxGWM7p8yRlcOy4n0/bqoXL8uKu9wHZxjNSoAAACAhLelrlGPFVdE9e1oaNaVa3fqSCCoI8GQjgSCagpb7T17tpLdLgXDVn8vr9Kmusao6547VKWbxgzty/IHNEIlAAAAgIRmrdU3dhSrs4d+vN1uJ1hJ2lTXqM31jXp4f5n2NvqjjiW5jL7QZhMfnDxCJQAAAICE9kJ5dafh8XguXbOjQ1+a26XrR+bpi2PyNTTJ25vlDXqESgAAAAAJqzEU1j27ih1fn+1x6/Oj83XD6DzleIk/scCfKgAAAICE9csDZSpqCvTomgy3SzMzUrRkSJY+M3KI0j3uGFUHiVAZE8aYvZLGHefwIWvt8E6uWSTpG5IWSEqWtFPS7yQ9ZK0NxahUAAAAIKH936EjJzye5/VoVkaKZmekamZ6imZlpGhssk8uY/qoQhAqY6da0gOd9HeYDG6M+aikZyQ1SXpaUqWkyyX9VNKHJF0duzIBAACAxDU/M007GpolSaOTvZqdnqqZGSmalZ6iWRmpGubzyBAg44pQGTtV1tp7ujrJGJMp6RFJIUnnWGvfb+3/pqR/SrrKGPMpa+1TsSwWAAAASEQ/mTpGtxUMU6bHzZrIBOWKdwHQVZLyJT11NFBKkrW2SS3TYSXpS/EoDAAAAIg3lzEal5JEoExg/JuJnSRjzHWSxkqql7Re0hudrI88t/X1pU7u8YakBkmLjDFJ1trmmFULAAAAAA4QKmNnuKTH2/XtMcb8u7V2RZu+Ka2v29vfwFobNMbskTRD0gRJW070gcaY1cc3uwNiAAAgAElEQVQ5NLV7JQMAAABAzzD9NTZ+L+k8tQTLNEmzJP1KUoGkpcaYOW3OzWp9rT7OvY72Z/d+mQAAAABwchipjAFr7bfbdW2U9EVjTJ2k/5R0j6Qrunm7o1tZ2W587vxOb9Aygjmvm58HAAAAAN3GSGXf+mXr69lt+o6ORGapc5ntzgMAAACAhEGo7Ftlra9pbfq2tb4Wtj/ZGOORNF5SUNLu2JYGAAAAAD1HqOxbC1tf2wbEf7a+XtTJ+WdLSpX0Dju/AgAAAEhEhMpeZoyZYYzJ7aR/nKSftzafaHPor5IqJH3KGHNqm/OTJd3b2vzfGJULAAAAACeFjXp639WS7jLGLJO0R1KtpImSLpWULOlFST86erK1tsYYc6NawuVyY8xTkiolfUQtjxv5q6Sn+/QbAAAAAEA3ESp73zK1hMFT1DLdNU1SlaS31PLcysettVE7uVprnzPGfFjSf0u6Ui3hc6ek2yU92P58AAAAAEgUhMpeZq1dIWmFg+velnRJ71cEAAAAALHDmkoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAAAAAIBjhEoAAAAAgGOESgAAAACAY4RKAACAAejQ7p0KBQPxLgPAIOCJdwEAAADoXU11dfrj12+XJylJ42bN0eW3f00ulzveZQEYoBipBAAAGGD2b/pA1oYVaGpUXeVhAiWAmCJUAgAADDB7P1gTeT9u9rw4VgJgMCBUAgAADCDWWu1bvzbSLphzShyrATAYECoBAAAGkCMHS1RTXiZJ8qWkaMTkqXGuCMBAR6gEAAAYQNpOfR0zY47cHvZlBBBbhEoAAIABZN/6Y6GyYDZTXwHEHqESAABggAgFAzqwaUOknZSWJhsOx7EiAIMBoRIAAGCAqK86okBzU6T94kM/0m9uvVFvPfWYDhcdiGNlAAYyJtkDAAAMEOk5Q1Qwd772rlsd6aspP6RVz/5Zq579s4aOn6jpZy3WlEVnKz0nN46VAhhICJUAAAADhMvt1se/+i0Vb9usLW8t1/aVb6mpvi5yvGzPLpXt2aUVj/9OY2fN0bQzz9Hk0xfKl5Iax6oB9HfGWhvvGhBjxpjV8+bNm7d69equTwYAAANGMBDQnnXva+uby7VrzXsKBQIdzvH4kjTx1DM0/azFGj93voyL1VHAYDR//nytWbNmjbV2fk+vZaQSAABggPJ4vZp82kJNPm2hmurrtGPVO9ry5jId2LJRah1YCPqbte2dN7TtnTc0snCaLv7y7coePiLOlQPoTwiVAAAAg0ByWrpmnbtEs85dopqKcm19e4W2vLlMFQf2Rc4p2b5Fj915iz78mc9p9vkXyxgTx4oB9BeESgAAgEEmMy9fp3/0Kp3+0atUvm+PNi5/Tete/ofCoZACzU167Te/0M73V2nJTbcoIzcv3uUCSHBMmgcAABjE8seN1+Lrb9S19/5YQ0aPjfTvXbdaj91xs7a+vSKO1QHoDwiVAAAA0LAJk3Td9x/Q/MuukFqnvTbV1+mFB/9Hzz/wQzXW1sS5QgCJilAJAAAASZLH59M5n7lBn7j7PmXmD4v0b1/5ph6948vavfZfcawOQKIiVAIAACDKmOmzdP3/PKRZ5y6J9NVXHdGzP/i2Xv31z+VvaoxjdQASDaESAAAAHfhSUrXkplv1sTvvVmpWdqR//esv6bE7b1HR1k1xrA5AIiFUAgAA4Lgmzj9d1//oYRWe8aFIX/WhUj19z11a8cTvFPT741gdgERAqAQAAMAJpWZm6bLb7tIlt9yhpLS0lk5r9f7z/6c/fv02le3dHd8CAcQVoRIAAABdMsZo2pnn6Pr/eVjjZp8S6a84sE9//Pptevf/nlY4FIpjhQDihVAJAACAbssYkqcrv/4dnfe5L8njS5IkhUMhvf3043rq7jtVWVIc5woB9DVCJQAAQA9Ulx1SKBiMdxlxZYzR3Asv1Wfvf1AjJk+J9B/cuU2Pf/VWrX3pedlwOI4VAuhLhEoAAIBuCAWDWvrzH+s3t9yg39xyA+sIJeWMGKVPfft+nfmpz8rl9kiSgv5m/fP3v9Jf77tbNRXlca4QQF8gVAIAAHQhGAjo+Z9+X5vfXCZJqqs8rMe/eqv++N+3KxQMxLm6+HK53Trjik/o0/f9RHljxkX6929Yp8f+62ZtW/lWHKsD0BcIlQAAAO2EQyGV79+rTSte16uP/Fw/u+4K7Xp/VYfzSndu18u/fDAOFSaeoQUTdO33fqyRhdMifc0N9XrxoR+ppqIsjpUBiDVPvAsAAABIFLWHK/S3H31PFQf2KhTo3gjkljeXadKpZ6hwwZkxri5xNdRUa+OyV/XBqy+qpjw6QLrcbrk93jhVBqAvECoBAABapWRmqXzf7h4/GuMfD9yvC25q0KzFS2JUWWIq3bld6155QVvfeaNDCDfGpQnzT9eCj39Sadk5caoQQF8gVAIAALTyeL0aMnqsyvftUUZevoYWTFD+uPGqLC5S3thxmnbmYmUNHSZjjGoqyvXX731TR0qKZG1Yr/zyQTXX1+vUy66I99eIqaDfr20r39S6V15Q6c7tHY4nZ2Rq9rlLNOeCS5SZPzQOFQLoa4RKAACANi699U6lZmUpJSPzhOdl5uXrU9/+oZ65726V7dklSVrx+G/lb2zUoquv7YtS+1RNeZk+ePVFbfjnK2qsrelwfPjEyZp74WWasvAseXy+OFQIIF4IlQAAAG0MGT2m2+emZmbpE3ffp2d/+B0Vb90kSVr51yc16bQFGlowIVYldmr9ay/p0O6dWnDlp5QxJK9X7mmt1f4NH2jty//Q7tXvydroZ0+6vV5NWXiW5l54qUZMmnKcuwAY6AiVAAAAJyEpNU1Xfv3beua+u1W8dbMkad+GdX0aKisO7NOrj/xcklS6e4c+fd9P5HK5u329DYcVDAYU9PsV9Dcr2NysPevWaN0rL+hISVGH8zPy8jV3yaWaufgCpWZm9dr3ANA/ESoBAABOkjcpWdPPPjcSKou3btJpl3+8zz7/6PTbo+83vP6K5lxwcdQ5TfV12vr2G9r2zhuqraw4FiD9/m7vdDtu9imau+RSTZh/Wo9CK4CBjVAJAADQC0ZNnRF5X7xti2w4LOPqm0eCVx0qjWq/9dRjKlzwISWlpWn/hg+0cflr2vmvld0Oj235UlI145zzNHfJpcodObq3SgYwgBAqAQAAekHuyNFKzshUU22NmmprVFlS3KP1mSej+tDBqHZTXa1+8flrlTEkX7WHy7t1D7fXK4/PJ48vSR6fT2lZOZp+9mJNO2uxfMkpsSgbwABBqAQAAOgFxhiNmjJdu95/V1LLFNi+CpVVZYc67W8fKIcWTNSMc85XwZx58iYnRQKkx+Pts1FVAAMPoRIAAKCXjJ4aHSpnn39Rn3xu+5HKtpLTMzT9rMWacc75fb4jLYDBgVAJAADQS6LXVW7uk8+01qqxtva4x8+69nrNPq9vwi2AwYl5DgAAAL1k6PgJ8viSJEnVZYdUW1kR8880xmjamedIalnXeda1/6YZHz4vcvztp59Qc0N9zOsAMHgxUgkAANBL3B6vRkyeogOb1kuSirdu1tRFZ8f8cy/80n/ovBu+KI8vScYY+ZsatW/DOtVVHlZDdZXe+cuTWnz9jTGvA8DgxEglAABALxo1dXrk/dHnVsaaMUbepGQZYyRJweZmjZw8NXJ87UvPq/o4m/kAwMlipBIAAKAXjZrSJlT20bpKSbLhsPZtWKcNr7+sne+vUjgUjDpWXXZIWUOH9Vk9AAYPQiUAAEAvGlk4Vca4ZG1Y5fv2qLmhXkmpab12/8baGm345yva9s6bSsvO1tnXfU473ntHG5e9qprysg7ne5NTdMqFl2r09Bmd3A0ATh6hEgAAoBf5UlKVXzBeZXt2SdaqZPtWjZ87/6TvW75vj9a+9Ly2vLlcwYA/0r9n3epOzx9ROFWzzl2iKQvPki855aQ/HwCOh1AJAADQy0ZNnd4SKiXtW79WY2fOltvjjRy31iocCiroDygU8CsY8CsUCCgYCCjk9ysYDCjkDygY8KuxtkabVrymos0bu/zc5PQMTT/7XM1afIHyxhbE6usBQBRCJQAAQC8bPXWG1i59XpK0+oXntPqF55Scli4rGwmLvWnszNmade6FmnTaQnl8vl69NwB0hVAJAADQy0ZNnSGX261wKBTpa6qvO6l7GpdLNhzu0P+5n/1aOcNHntS9AeBkECoBAAB6WVp2ji648Wate+WF1mdFVsva6EDocnvk8Xnl9vrk9nrl8XrlaX3v9vpa2j6f3B6v8sYWaNZ5S/TSwz/R/o3rI/e44cHfKHvY8L7+egAQhVAJAAAQAzMXX6CZiy+QJIXDITXV1cnlcsvt88rj8cq4ev648LaB0u3xECgBJARCJQAAQIy5XG6lZmad1D0aaqqj2mdec/1J3Q8AekvPf0UGAACAPrdj1TtR7WlnnhOfQgCgHUIlAABAP7DlreVR7bTsnPgUAvQTh+oPadXBVbLWxruUAY9QCQAA0A8ULjgz8v5Dn7gujpUAia+isUJX/P0Kff6Vz+uhtQ/Fu5wBjzWVAAAA/cCcCy6WDYeUnJ6hGR8+L97lAAnt/vfuV62/VpL0yIZHdOu8W+Nc0cBGqAQAAOgH3B6P5l/6sXiXAfQLS/cujXcJgwrTXwEAAAAMWEvGLYl3CQMeoRIAAADAgBEIBaLaV0+5Ok6VDB6ESgAAAAADxnul70W1Tx12apwqGTwIlQAAAAAGjL9s/0tU2+NiG5lYI1QCAAAAGDBe3/96vEsYdAiVAAAAAAakyyZcFu8SBgVCJQAAAIABobq5Oqr9sUk8hqcvECoBAAAADAgbKzZGtWflzYpTJYMLoRIAAADAgNB+pBJ9g1AJAAAAYEBoCDZEtd8sfjNOlQwuhEoAAAAAA8LWyq1R7df2vRanSgYXQiUAAACAAaF9qHyj6A01h5rjVM3gQagEAAAA0O+FwiFtP7I9qq8h2KCVJSvjVNHg4Yl3AQAAAABwsvbX7ldjsLFD/zM7nlGaN01el1det7fl1eWVx+WJvG/b7zZuGWPi8A36L0IlAAAAgH6vxl/Taf/yA8u1/MDyk7r3TbNv0s2n3HxS9xjImP4KAAAAoN+bljtNGd6MmNz7V+t/pcONh2Ny74GAkUoAAAAA/Z7P7dPisYv1911/lySdNvw05Sbnqqq5SoFQQMFwUIFwQFsqt/T43otGLlJucm5vlzxgECoBAAAADAgXj784EioP1B7Qb5b8Ri7jUml9qX625mf6x+5/9Pied51+l66dei3rLE+AUAkAAABgQDhjxBnKTspWVXOVSutL9dmln9UH5R+c1D0vGX8JgbILhEoAAAAAA4KRUVVzVaR9okDpMi5dNuEyTc6erKykLGUnZSsrKUsl9SX62ptfkySNyRijnOScmNfd3xEqAQAAAPRrL+x+QXe9eVe3zr33Q/dqScESpXhSOj2+ZcuxNZez82f3Sn0DHaEyQRhjRkv6jqSLJA2RdFDSc5K+ba09Es/aAAAAgET14JoH9ciGR054zldP+6quLLzyuEGyrbajm7PzCJXdQahMAMaYiZLekTRU0t8kbZV0uqT/kHSRMeZD1lr2MAYAAADaeWH3Cyc8flXhVbpu+nXduldDoEFL9yyNtBmp7B6eU5kYfqGWQHmrtfZj1tq7rLXnSvqppCmSvhfX6gAAAIAEdfMpN2tY6jBJ0scnf1wrr1mpxy5+LHL81X2vKhAOHPf6w42H9bedf9Pty2/XGU+eEenPTsrWlJwpsSt8AGGkMs6MMRMkLZG0V9LD7Q5/S9IXJH3GGPOf1tr6Pi4PAAAASGiXT7xcl0+8PKpvTv4cDU8brtL6UlU3V2vVwVU6c9SZkiRrrbYd2aYVB1bojaI3tKFig6xsh/veOOtGed3ePvkO/R2hMv7ObX19xVobbnvAWltrjHlbLaFzgaTX+7o4AAAAoL9xGZcuHHehHt38qCTp7zv/rmA4qBVFLUGyrKGsy3tcM/WaWJc5YBAq4+/omPr24xzfoZZQWaguQqUxZvVxDk11VhoAAADQP108/uJIqFy6d6mW7l3axRXRGKXsPkJl/GW1vlYf5/jR/uw+qAUAAADol4LhoF7f/7p2Vu3Urqpd2nFkh+N7ZfgyerGygY9QmfhM62vHid7tWGvnd3qDlhHMeb1ZFAAAAJBIXMalb779TTUGG3t87ez82ZqWO01Pb3takpTkTurt8gY0QmX8HR2JzDrO8cx25wEAAABox2VcmpQ9SRsqNkT6jIxcxqXCnELNzp8dCY2SNDFroi6dcKkuGn+RJOmmV2+KHGOksmcIlfG3rfW18DjHJ7e+Hm/NJQAAAABJl4y/RPOHzdek7EmalDNJE7ImKMWTEjl+deHVev/Q+zp12KkqzCmUMUbbj2zXTa/epP/f3t0Hy1XXdxx/fxNJAhfyRIuI1F4JgikC7chIeVAIMQgFeWihlamWoWpFlKLIIAOlCkVpqVTKQw3TWhiKFTqU6iiIWASMMmCrBWEaWp5ipEaBRB4SAkng2z/OubhZ7t7k/nL37r1n36+ZM2f2d85v97fky2/vZ8+ec55a+xQAU2Mqp/zmKb16C5OSobL3bq/Xh0bElNYrwEbEdsABwFrg7l4MTpIkSZos3vsb7x1x++5zd2f3ub+89+S9T9zLKbedwnPrngNg2pRpfO6gz7HgDQu6Os6mmdLrAfS7zHwEuBUYBD7Stvk8YAC4xntUSpIkSWNnyeNL+OCtH3wlUA5sNcDiRYsNlAU8UjkxnALcBVwaEQuBpcC+wAKqn72e08OxSZIkSZPCSy+/xJnfOZNtp23LtltVy8BWA8yeMZuDdj6IWdOry5jc/OjNnPPdc9iQGwCYO2Mui9+5mPnbz+/l8CctQ+UEkJmPRMQ+wPnAYcDvACuAS4HzMnNVL8cnSZIkTQZrNqzh1h/fOuy2wZmDfOmIL3HTozdx4T0XkvXNFXYa2IkrF13J4KzBcRxpsxgqJ4jM/AlwUq/HIUmSJE1Wa9Z1PmNs2bPLOOHrJ7D8ueWvtM2bNY8rF13JawdeOx7DayxDpSRJkqRGmDl9Jhe94yJWr1/NmnVrWL1+NU+ufZIbH7oRYKNAudev7MUVC69g9ozZvRpuYxgqJUmSJDXCwFYDHP7Gw1/VPmv6LK564KpXHu/3uv24ZMElbLPVNuM5vMby6q+SJEmSGu203zqNwwersHnUvKO4fOHlBsox5JFKSZIkSY02dcpULjroIi448AKmTZ3W6+E0jkcqJUmSJPUFA2V3GColSZIkScUMlZIkSZKkYoZKSZIkSVIxQ6UkSZIkqZihUpIkSZJUzFApSZIkSSpmqJQkSZIkFTNUSpIkSZKKGSolSZIkScUMlZIkSZKkYoZKSZIkSVIxQ6UkSZIkqZihUpIkSZJULDKz12NQl0XEyq233nru/Pnzez0USZIkSRPQ0qVLWbt27arM3H60fQ2VfSAiHgNmAsvG6SXfXK8fHKfXUzNYNyph3aiEdaMS1o1KTKa6GQSezcw3jrajoVJjLiJ+AJCZb+31WDR5WDcqYd2ohHWjEtaNSvRL3XhOpSRJkiSpmKFSkiRJklTMUClJkiRJKmaolCRJkiQVM1RKkiRJkop59VdJkiRJUjGPVEqSJEmSihkqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqjUmImInSPiHyPipxHxYkQsi4hLImJOr8em3qprITssP+vQZ/+IuDkiVkXE8xHxo4j4WERMHe/xq3si4riIuCwilkTEs3VNXLuJPqOujYg4MiLuiIhnImJ1RNwTESeO/TvSeBhN3UTE4AjzT0bEdSO8zokR8f26Zp6pa+jI7r0zdVNEbB8RH4iIf4uIhyNibf3v+t2IeH9EDPt3sXNOfxtt3fTrnPOaXg9AzRAR84C7gB2ArwIPAm8DTgMOi4gDMnNlD4eo3nsGuGSY9tXtDRFxNPCvwAvA9cAq4N3A54EDgOO7N0yNsz8D9qaqg8eBN4+0c0ltRMRHgcuAlcC1wDrgOODqiNgzM88YqzejcTOquqndB3xlmPYHhts5Ij4HfKJ+/r8HpgHvAb4WEadm5uUF41ZvHQ98AVgB3A4sB14L/C7wD8DhEXF8ZuZQB+ccUVA3tf6aczLTxWWLF+CbQAKntrX/Td2+uNdjdOlpfSwDlm3mvjOBJ4AXgX1a2mdQfXGRwHt6/Z5cxqw2FgBvAgI4uP73vXasagMYpPpjcCUw2NI+B3i47rNfr/87uHS1bgbr7VeP4vn3r/s8DMxpe66VdU0Nbsl7cOlJ3RxCFQintLXvSBUUEvi9lnbnHJeSuunLOcefv2qLRcQuwKFUweGKts2fAtYA74uIgXEemian44BfBa7LzP8caszMF6iOTgB8uBcD09jLzNsz86GsPz03oaQ2/hiYDlyemcta+vwC+Gz98OTC4atHRlk3JYZq4jN1rQy97jKqz7npwEldem11SWZ+OzO/lpkvt7X/DFhcPzy4ZZNzjkrqpsSkn3MMlRoLh9TrW4f5H+454HvANsBvj/fANKFMj4j3RsTZEXFaRCzocD7KUD3dMsy27wDPA/tHxPSujVQTVUltjNTnG237qNl2iogP1XPQhyJirxH2tW76z/p6vaGlzTlHmzJc3QzpqznHcyo1Fnav1//bYftDVEcydwNuG5cRaSLaEfintrbHIuKkzLyzpa1jPWXmhoh4DNgD2AVY2pWRaqIqqY2R+qyIiDXAzhGxTWY+34Uxa+JYVC+viIg7gBMzc3lL2wDwemB1Zq4Y5nkeqte7dWmcGmcR8Rrgj+qHrX/UO+eooxHqZkhfzTkeqdRYmFWvn+mwfah99jiMRRPTVcBCqmA5AOwJXEl1rsA3ImLvln2tJ3VSUhub22dWh+2a/J4H/gJ4K9V5bXOAg6guuHEwcFvb6RnOQf3nL4G3ADdn5jdb2p1zNJJOddOXc46hUuMh6nW3zn3RBJeZ59XnJPw8M5/PzAcy82SqCzltDXx6FE9nPamTktqwnhouM5/IzD/PzB9m5tP18h2qX9DcA+wKfKDkqcd0oOqJiPhTqituPgi8b7Td67VzTp8ZqW76dc4xVGosbOpbt5lt+0lDhk5wf0dLm/WkTkpqY3P7PLsF49IklJkbqG4HAKObgzZ1VEGTRER8BPhb4L+BBZm5qm0X5xy9ymbUzbCaPucYKjUW/qded/qt95vqdadzLtW/nqjXrT8D6VhP9fkLb6Q6If7R7g5NE1BJbYzU53VUtfe45zb1rSfr9StzUGauAf4P2LaukXZ+pjVARHwMuJzqnoEL6it5tnPO0UY2s25G0tg5x1CpsXB7vT40IjaqqYjYjurmwGuBu8d7YJrw9qvXrR/I367Xhw2z/zuoriR8V2a+2M2BaUIqqY2R+hzeto/6z9BVydu/pLJuGiwiPgl8HriXKhg80WFX5xy9YhR1M5LGzjmGSm2xzHwEuJXqoisfadt8HtW3MdfU38Soz0TEHhExd5j2X6f6tg/g2pZNNwBPAe+JiH1a9p8BXFA//EKXhquJraQ2rqK6cflHI2Kwpc8c4Oz64WLUWBGxb0RMG6b9EODj9cNr2zYP1cQ5da0M9Rmk+px7kaq2NMlExLlUF1j5AbAwM58aYXfnHAGjq5t+nXOie/cNVj+JiHnAXcAOwFepLq29L7CA6nD9/pm5sncjVK9ExKeBs6iOaD8GPAfMA44AZgA3A8dm5rqWPsdQfZi/AFwHrAKOorpU+w3A73fxpucaR/W/9TH1wx2Bd1F9g7ukbnsqM89o239UtRERpwKXAiuB64F1VDc13xm4uPX5NTmMpm7qS/jvAdwBPF5v34tf3vPt3MwcCgitr3ExcHrd5wZgGvAHwPbAqZl5eXsfTWwRcSJwNfAScBnDn6O2LDOvbunjnNPnRls3fTvnZKaLy5gswK9RfYuygmoC/THVicxzez02l57WxUHAl6mukPY01Y2CnwS+RXV/p+jQ7wCqwPkLqp9P30/1Dd/UXr8nlzGtj09TXdGu07JsLGoDeDdwJ9WXGmuA/6C6V1jP/xu4dLdugPcDXweWAaupvvFfTvXH/ts38Ton1rWypq6dO4Eje/3+XbpWNwncMUw/55w+XkZbN/0653ikUpIkSZJUzHMqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJEmSJBUzVEqSJEmSihkqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJEmSJBUzVEqSJEmSihkqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJKlBImJGRCyLiGxZno6I7Uf5PHMi4kdtz/NSRJzQrbFLkiYnQ6UkSQ2SmS8AZ7c1zxqmraOIGABuAvZs2/ThzPzylo1QktQ0kZm9HoMkSRpDERHA3cDbWppfBHbLzOWb6DsN+DqwqFjQRVgAAAPJSURBVG3TGZl58ZgOVJLUCB6plCSpYbL6xvgTbc3TgfNH6hcRU4F/5tWB8nwDpSSpE49USpLUUBFxI3BsS9PLwN6Z+cAw+wbwReCktk2XZObHuzdKSdJk55FKSZKa60xgfcvjKcCFHfa9mFcHyi8Cp3dhXJKkBjFUSpLUUJn5MPB3bc1HRsSBrQ0RcS7QfjTyeuBP0p80SZI2wZ+/SpLUYBExF3gYmNPSfFdmHlBv/yhwWVu3m4BjM3M9kiRtgqFSkqSGi4jTqX7e2uoYYDvgGiBa2u8ADq9vTSJJ0iYZKiVJarj6NiFLgV1amn8K7AC8pqXt+8DCzFw9jsOTJE1ynlMpSVLDZeY64Ky25p3YOFDeDxxmoJQkjZZHKiVJ6hMR8T1g/2E2PQS8PTN/Ps5DkiQ1gEcqJUnqH/cN07YKeKeBUpJUylApSVIfiIgzgQ8Ps2k2MHechyNJahBDpSRJDRcRJwN/1WHzFF59ZVhJkjaboVKSpAaLiD8ErmhrXtX2+JCIOGKchiRJahhDpSRJDRURRwNXs/Hn/X3AW4DH23b/64iYOk5DkyQ1iKFSkqQGioiFwPVsfNuQh4B3ZeYK4FNtXeYDHxyn4UmSGsRbikiS1DARsR/wLWCgpfknwIGZubzeZyrVUcs9WvZ5Atg1M58br7FKkiY/j1RKktQgEbE3cDMbB8ongUVDgRIgM18Czm7rvgNwVtcHKUlqFI9USpLUEBGxG7CEKhwOeQZYkJn/1aHPEuDAlqa1wO6Z+ZOuDVSS1CgeqZQkqQEi4g3Av7NxoFwLHNkpUNY+2fZ4a+AzYzw8SVKDeaRSkqRJLiJ2pDpCuWtL83rgqMy8ZTP6fwU4uqUpgX0y84djOlBJUiMZKiVJmsQiYg5wJ7BnS/PLwAmZ+S+b+RzzgfuB1luK3JGZC8ZsoJKkxvLnr5IkTVIRsS1wCxsHSoAPbW6gBMjMpVT3s2x1cEQctWUjlCT1A49USpI0CUXEDKqrvLYfTTwjMy8ueL7XU93HcuuW5geBPTNzQ/FAJUmNZ6iUJEmSJBXz56+SJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJEmSJBUzVEqSJEmSihkqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJEmSJBUzVEqSJEmSihkqJUmSJEnFDJWSJEmSpGKGSkmSJElSMUOlJEmSJKmYoVKSJEmSVMxQKUmSJEkqZqiUJEmSJBUzVEqSJEmSiv0/Cb/AnULl2RYAAAAASUVORK5CYII=\n", + "text/plain": [ + "<Figure size 504x504 with 1 Axes>" + ] + }, + "metadata": { + "image/png": { + "height": 438, + "width": 458 + }, + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "figure(1, [7, 7])\n", + "plot(tracks[..., 0], tracks[..., 1])\n", + "xlim(-10, 265)\n", + "ylim(-10, 265)\n", + "xlabel(r'$x$', fontsize=24)\n", + "ylabel(r'$y$', fontsize=24);" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Basic CNN" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "train_size = 200\n", + "\n", + "train_vid, train_labels, train_tracks = [],[],[]\n", + "for i in range(train_size): \n", + " vid, labels, tracks = pt_single(kappa, a, IbackLevel, Nparticles, sigma_motion) \n", + " train_vid.append(vid[:,::2,::2]) # downsample video to Ntx128x128\n", + " train_labels.append(labels)\n", + " train_tracks.append(tracks)\n", + "train_vid = tf.convert_to_tensor(train_vid) \n", + "train_labels = tf.convert_to_tensor(train_labels) \n", + "train_tracks = tf.convert_to_tensor(train_tracks)\n", + "\n", + "train_vid = tf.transpose(train_vid, perm=[0,2,3,1])\n", + "train_labels = tf.squeeze(train_labels)#[:,:,:,1]" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "val_size = 50\n", + "\n", + "val_vid, val_labels, val_tracks = [],[],[]\n", + "for i in range(val_size): \n", + " vid, labels, tracks = pt_single(kappa, a, IbackLevel, Nparticles, sigma_motion) \n", + " val_vid.append(vid[:,::2,::2]) # downsample video to Ntx128x128\n", + " val_labels.append(labels)\n", + " val_tracks.append(tracks)\n", + "val_vid = tf.convert_to_tensor(val_vid) \n", + "val_labels = tf.convert_to_tensor(val_labels) \n", + "val_tracks = tf.convert_to_tensor(val_tracks)\n", + "\n", + "val_vid = tf.transpose(val_vid, perm=[0,2,3,1])\n", + "val_labels = tf.squeeze(val_labels)#[:,:,:,1]" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [], + "source": [ + "import tensorflow.python.keras.backend as K\n", + "\n", + "def f1_metric(y_true, y_pred):\n", + " true_positives = K.sum(K.round(K.clip(y_true * y_pred, 0, 1)))\n", + " possible_positives = K.sum(K.round(K.clip(y_true, 0, 1)))\n", + " predicted_positives = K.sum(K.round(K.clip(y_pred, 0, 1)))\n", + " precision = true_positives / (predicted_positives + K.epsilon())\n", + " recall = true_positives / (possible_positives + K.epsilon())\n", + " f1_val = 2*(precision*recall)/(precision+recall+K.epsilon())\n", + " return f1_val" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Epoch 1/5\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "20/20 [==============================] - ETA: 0s - loss: 6.6574 - accuracy: 0.9962 - f1_metric: 0.9962 - precision: 0.9964 - recall: 0.9998WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "20/20 [==============================] - 6s 279ms/step - loss: 6.6574 - accuracy: 0.9962 - f1_metric: 0.9962 - precision: 0.9964 - recall: 0.9998 - val_loss: 3.6102 - val_accuracy: 0.9959 - val_f1_metric: 0.9959 - val_precision: 0.9964 - val_recall: 0.9998\n", + "Epoch 2/5\n", + "20/20 [==============================] - 5s 270ms/step - loss: 2.4552 - accuracy: 0.9958 - f1_metric: 0.9958 - precision: 0.9964 - recall: 0.9997 - val_loss: 2.0370 - val_accuracy: 0.9956 - val_f1_metric: 0.9956 - val_precision: 0.9964 - val_recall: 0.9996\n", + "Epoch 3/5\n", + "20/20 [==============================] - 5s 270ms/step - loss: 1.4854 - accuracy: 0.9958 - f1_metric: 0.9958 - precision: 0.9964 - recall: 0.9995 - val_loss: 1.3767 - val_accuracy: 0.9958 - val_f1_metric: 0.9959 - val_precision: 0.9964 - val_recall: 0.9995\n", + "Epoch 4/5\n", + "20/20 [==============================] - 6s 281ms/step - loss: 1.0705 - accuracy: 0.9959 - f1_metric: 0.9959 - precision: 0.9964 - recall: 0.9995 - val_loss: 1.1987 - val_accuracy: 0.9964 - val_f1_metric: 0.9964 - val_precision: 0.9964 - val_recall: 0.9995\n", + "Epoch 5/5\n", + "20/20 [==============================] - 6s 286ms/step - loss: 1.0566 - accuracy: 0.9957 - f1_metric: 0.9957 - precision: 0.9965 - recall: 0.9994 - val_loss: 1.8673 - val_accuracy: 0.9964 - val_f1_metric: 0.9964 - val_precision: 0.9965 - val_recall: 0.9994\n" + ] + }, + { + "data": { + "text/plain": [ + "<keras.callbacks.History at 0x26580485348>" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "total_items = 200\n", + "batch_size = 10\n", + "epochs = 5\n", + "\n", + "model = Sequential()\n", + "model.add(Input((128,128,1)))\n", + "model.add(Conv2D(12, 3, padding='same', activation=LeakyReLU(alpha=0.01), strides=(1,1)))\n", + "model.add(Conv2D(32, 3, padding='same', activation=LeakyReLU(alpha=0.01), strides=(1,1)))\n", + "model.add(Conv2D(2, 3, padding='same', strides=(1,1)))\n", + "model.add(Activation('softmax'))\n", + "\n", + "model.compile(optimizer='adam',\n", + " loss='binary_crossentropy',\n", + " metrics=['accuracy', f1_metric, keras_metrics.precision(), keras_metrics.recall()])\n", + "\n", + "num_batches = int(total_items/batch_size)\n", + "train_data_generator = single_input_fn(total_items, epochs)\n", + "model.fit(train_vid, train_labels, steps_per_epoch=num_batches, epochs=epochs, verbose=1, \n", + " validation_data=(val_vid, val_labels))#, class_weight={0: 0.1, 1: 0.9})#, callbacks=[tensorboard_callback])" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Model: \"sequential_1\"\n", + "_________________________________________________________________\n", + "Layer (type) Output Shape Param # \n", + "=================================================================\n", + "conv2d_3 (Conv2D) (None, 128, 128, 12) 120 \n", + "_________________________________________________________________\n", + "conv2d_4 (Conv2D) (None, 128, 128, 32) 3488 \n", + "_________________________________________________________________\n", + "conv2d_5 (Conv2D) (None, 128, 128, 2) 578 \n", + "_________________________________________________________________\n", + "activation_1 (Activation) (None, 128, 128, 2) 0 \n", + "=================================================================\n", + "Total params: 4,186\n", + "Trainable params: 4,186\n", + "Non-trainable params: 0\n", + "_________________________________________________________________\n" + ] + } + ], + "source": [ + "model.summary()" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [], + "source": [ + "from datetime import datetime\n", + "logdir = \"logs/scalars/\" + datetime.now().strftime(\"%Y%m%d-%H%M%S\")\n", + "tensorboard_callback = keras.callbacks.TensorBoard(log_dir=logdir)" + ] + }, + { + "cell_type": "code", + "execution_count": 162, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "UsageError: Line magic function `%tensorboard` not found.\n" + ] + } + ], + "source": [ + "%tensorboard --logdir logs/scalars" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# ConvLSTM2D" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [], + "source": [ + "train_size = 10\n", + "Nt = 20\n", + "pt = Particle_Tracking_Training_Data(Nt)\n", + "\n", + "train_vid, train_labels, train_tracks = [],[],[]\n", + "for i in range(train_size): \n", + " vid, labels, tracks = pt(kappa, a, IbackLevel, Nparticles, sigma_motion) \n", + " train_vid.append(vid[:,::2,::2]) # downsample video to Ntx128x128\n", + " train_labels.append(labels)\n", + " train_tracks.append(tracks)\n", + "train_vid = tf.convert_to_tensor(train_vid) \n", + "train_labels = tf.convert_to_tensor(train_labels) \n", + "train_tracks = tf.convert_to_tensor(train_tracks)\n", + "\n", + "train_vid = tf.expand_dims(train_vid, 4)\n", + "train_labels = tf.squeeze(train_labels)" + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "metadata": {}, + "outputs": [], + "source": [ + "val_size = 5\n", + "\n", + "val_vid, val_labels, val_tracks = [],[],[]\n", + "for i in range(val_size): \n", + " vid, labels, tracks = pt(kappa, a, IbackLevel, Nparticles, sigma_motion) \n", + " val_vid.append(vid[:,::2,::2]) # downsample video to Ntx128x128\n", + " val_labels.append(labels)\n", + " val_tracks.append(tracks)\n", + "val_vid = tf.convert_to_tensor(val_vid) \n", + "val_labels = tf.convert_to_tensor(val_labels) \n", + "val_tracks = tf.convert_to_tensor(val_tracks)\n", + "\n", + "val_vid = tf.expand_dims(val_vid, 4)\n", + "val_labels = tf.squeeze(val_labels)" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Epoch 1/5\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "2/2 [==============================] - ETA: 0s - loss: 0.6820 - accuracy: 0.6586 - f1_metric: 0.6586 - precision: 0.9965 - recall: 0.4946 WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "WARNING:tensorflow:`add_update` `inputs` kwarg has been deprecated. You no longer need to pass a value to `inputs` as it is being automatically inferred.\n", + "2/2 [==============================] - 69s 34s/step - loss: 0.6820 - accuracy: 0.6586 - f1_metric: 0.6586 - precision: 0.9965 - recall: 0.4946 - val_loss: 0.6239 - val_accuracy: 0.9963 - val_f1_metric: 0.9963 - val_precision: 0.9965 - val_recall: 0.7731\n", + "Epoch 2/5\n", + "2/2 [==============================] - 33s 18s/step - loss: 0.6066 - accuracy: 0.9961 - f1_metric: 0.9961 - precision: 0.9965 - recall: 0.8467 - val_loss: 0.5589 - val_accuracy: 0.9966 - val_f1_metric: 0.9966 - val_precision: 0.9964 - val_recall: 0.8865\n", + "Epoch 3/5\n", + "2/2 [==============================] - 33s 18s/step - loss: 0.5477 - accuracy: 0.9963 - f1_metric: 0.9963 - precision: 0.9964 - recall: 0.9088 - val_loss: 0.5191 - val_accuracy: 0.9966 - val_f1_metric: 0.9966 - val_precision: 0.9964 - val_recall: 0.9243\n", + "Epoch 4/5\n", + "2/2 [==============================] - 32s 18s/step - loss: 0.5121 - accuracy: 0.9963 - f1_metric: 0.9963 - precision: 0.9964 - recall: 0.9350 - val_loss: 0.4916 - val_accuracy: 0.9966 - val_f1_metric: 0.9966 - val_precision: 0.9964 - val_recall: 0.9432\n", + "Epoch 5/5\n", + "2/2 [==============================] - 32s 18s/step - loss: 0.4860 - accuracy: 0.9963 - f1_metric: 0.9963 - precision: 0.9964 - recall: 0.9495 - val_loss: 0.4677 - val_accuracy: 0.9966 - val_f1_metric: 0.9966 - val_precision: 0.9964 - val_recall: 0.9546\n" + ] + }, + { + "data": { + "text/plain": [ + "<keras.callbacks.History at 0x26580a06688>" + ] + }, + "execution_count": 16, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "total_items = train_size\n", + "batch_size = 5\n", + "epochs = 5\n", + "\n", + "model = Sequential()\n", + "model.add(Input((Nt,128,128,1)))\n", + "model.add(ConvLSTM2D(12, kernel_size=(3, 3), padding='same', strides=(1,1), \n", + " return_sequences=True, data_format='channels_last'))\n", + "model.add(ConvLSTM2D(16, kernel_size=(3, 3), padding='same', strides=(1,1), \n", + " return_sequences=True, data_format='channels_last'))\n", + "model.add(ConvLSTM2D(2, kernel_size=(3, 3), padding='same', strides=(1,1), \n", + " return_sequences=True, data_format='channels_last'))\n", + "model.add(Activation('softmax'))\n", + "\n", + "model.compile(optimizer='adam',\n", + " loss='binary_crossentropy',\n", + " metrics=['accuracy', f1_metric, keras_metrics.precision(), keras_metrics.recall()])\n", + "\n", + "num_batches = int(total_items/batch_size)\n", + "train_data_generator = single_input_fn(total_items, epochs)\n", + "model.fit(train_vid, train_labels, steps_per_epoch=num_batches, epochs=epochs, verbose=1, \n", + " validation_data=(val_vid, val_labels))#, callbacks=[tensorboard_callback])" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Model: \"sequential_2\"\n", + "_________________________________________________________________\n", + "Layer (type) Output Shape Param # \n", + "=================================================================\n", + "conv_lst_m2d (ConvLSTM2D) (None, 20, 128, 128, 12) 5664 \n", + "_________________________________________________________________\n", + "conv_lst_m2d_1 (ConvLSTM2D) (None, 20, 128, 128, 16) 16192 \n", + "_________________________________________________________________\n", + "conv_lst_m2d_2 (ConvLSTM2D) (None, 20, 128, 128, 2) 1304 \n", + "_________________________________________________________________\n", + "activation_2 (Activation) (None, 20, 128, 128, 2) 0 \n", + "=================================================================\n", + "Total params: 23,160\n", + "Trainable params: 23,160\n", + "Non-trainable params: 0\n", + "_________________________________________________________________\n" + ] + } + ], + "source": [ + "model.summary()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.6" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +}