{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Subpixel Accuracy and Uncertainty Estimation\n", "\n", "## How is subpixel resolution possible?\n", "\n", "Tracking resolution can exceed the traditional diffraction-limited resolution of a microscope. If a particle's image spans multiple pixels, we can find its position with subpixel accuracy by taking the average position of these pixels, weighted by brightness. This is referred to as the particle's *centroid*.\n", "\n", "This tutorial will use simulated images to demonstrate how well trackpy's implementation of the Crocker--Grier algorithm determines particle locations under various conditions. It also demonstrates trackpy's uncertainty-estimation capability and introduces the most important sources of error in particle tracking experiments." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "\n", "import matplotlib as mpl\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pandas as pd\n", "\n", "import trackpy as tp\n", "import pims" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Demonstration of trackpy's uncertainty estimation on noisy images\n", "\n", "We'll draw a simulated image with a single Gaussian blob. The class below will make it easy to generate copies of the images with simulated camera noise." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "class SimulatedFrame(object):\n", " \n", " def __init__(self, shape, dtype=np.uint8):\n", " self.image = np.zeros(shape, dtype=dtype)\n", " self._saturation = np.iinfo(dtype).max\n", " self.shape = shape\n", " self.dtype =dtype\n", " \n", " def add_spot(self, pos, amplitude, r, ecc=0):\n", " \"Add a Gaussian spot to the frame.\"\n", " x, y = np.meshgrid(*np.array(list(map(np.arange, self.shape))) - np.asarray(pos))\n", " spot = amplitude*np.exp(-((x/(1 - ecc))**2 + (y*(1 - ecc))**2)/(2*r**2)).T\n", " self.image += np.clip(spot, 0, self._saturation).astype(self.dtype)\n", " \n", " def with_noise(self, noise_level, seed=0):\n", " \"Return a copy with noise.\"\n", " rs = np.random.RandomState(seed)\n", " noise = rs.randint(-noise_level, noise_level, self.shape)\n", " noisy_image = np.clip(self.image + noise, 0, self._saturation).astype(self.dtype)\n", " return noisy_image\n", " \n", " def add_noise(self, noise_level, seed=0):\n", " \"Modify in place with noise.\"\n", " self.image = self.with_noise(noise_level, seed=seed)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The uncertainty in particle location can be estimated from the radius of gyration, the mask size (i.e., `diameter` specified to `locate()`) and the signal-to-noise ratio. The relation was worked out in a [detailed study of particle tracking errors by Savin and Doyle](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1305040/pdf/623.pdf). For a hat-like spot, the uncertainty in location of a hat-like spot is\n", "\n", "$$\\epsilon = \\frac{N}{S}\\frac{l_\\text{noise}}{2\\pi^{1/2}}\\frac{w^2}{a^2}$$\n", "\n", "where $N/S$ is noise-to-signal; $l_n$=1 pixel, the length scale of noise correlation; and $w$ and $a$ are feature size and mask size respectively. (The uncertainty in the location of a Gaussian blob is similar and more computationally intensive.)\n", "\n", "Trackpy performs the computation and gives the results in the `'ep'` column of the results set returned by `locate` or `batch`. This can be used as an estimate of uncertainty." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAEICAYAAABGaK+TAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzsnXl4FdX5+D8nCyELCQlbIGETRCDI\noiwiCFRABIq0Yi2t2mq16A/9iguiUkVERcGFVoVWLVQsIFgEKwURNLIJQRBCCBIgEJaEfQtJWAPn\n98e9ud478waGJCT3hvN5nnmSeefMzJmZd9575rznfY/SWmMwGAyGwCOooitgMBgMhpJhDLjBYDAE\nKMaAGwwGQ4BiDLjBYDAEKMaAGwwGQ4BiDLjBYDAEKAFrwJVS9yilFpXDeXoopbKv4PG/Ukr98Uod\n3+D/GF02lBS/NuBKqa5KqZVKqVyl1FGl1PdKqQ4AWuvpWuvbKqhe9YpeBKXUTqXUQaVUpNf2h5RS\nS5wcS2vdV2s99QpV1TFKqTCl1GSl1C6lVJ5SKlUp1ddSpqdSKkMpdVIp9Z1SqqFl/ylKqRNKqf1K\nqaeEczyvlBqrlKqilJrtvndaKdWjHC6xQjG6XH6Usy7fpJRa7H6mh5RS/1FK1fUqp5RS45RSR9zL\nOKWUKqtr9VsDrpSKBv4HvAfEAQnAy8CZiqyXm37AQq/1YGBYBdWlrAgB9gDdgRjgBeAzpVQjAKVU\nTWAO8CKu57EWmOW1/2jgWqAh8AtghFLqdss5+gML3P+vAO4F9pf5lfgZRpfLnfLU5VjgQ6CRu3we\n8C+vckOAXwFtgNbAAODhMrhGF1prv1yA9sDxi2y/H1jhtX4bsAXIBSYBS4GHvMsCbwHHgCygr9e+\nDwCb3Td/B/Cw17YeQLbl3HOAO93/7wSeA44C1d2yh4AlXuVvBta467YGuNlr2xKvejZ11zsXOAzM\n8irXHFjsPs8W4O5yeAZpwCD3/0OAlV7bIoFTQHP3+l7gNq/trwAzvdZjgYNAsOUc2UCPitY3o8tG\nl0ury+5tNwB5XusrgSFe6w8CKWV1XX7bAge2AueVUlOVUn2VUrHFFXT/os4Gngdq4FKKmy3FOrnl\nNYHxwGSvT5mDwC+BaFwvwASl1A3FnCsU6IZLAYtYi0t5hwvl44D5wLvuur0DzFdK1RAO/wqwCJeC\nJOJqseH+pF0MzABqA4OBSUqplsXUcZJS6ngxS5q0j3CMOkAzYJNblARsKNqutS4AtgNJ7mdT13u7\n+/8kr/U+wLda6/NOzl/JMLp89ehyN6/z2M4lHKtU+K0B11qfALoCGvgIOKSU+tL9MKz0AzZpredo\nrQtxKZj103yX1voj902fiush1XGfa77Wert2sRSX4t1STNW6ARu01nkW+Sjg/5RStSzy/sA2rfW/\ntdaFWutPgQxcn1JWzuH6DKuntT6ttV7hlv8S2Km1/pf7GOuBz4HfSBXUWg/VWlcvZmldzHV5cL/Y\n04GpWusMtzgKV2vKm1ygmnsblu1F27zvwwKuQowuXx26rJRqjevePeMltp4rF4gqq35wvzXgAFrr\nzVrr+7XWiUAroB7wV6FoPVx9XkX7aVyf5t7s99p+0v1vFIC7VZTidkQcx/US1SymWv0QHp7WOh1X\nP+dzQt12WWS7cPWDWhkBKOAHpdQmpdSf3PKGQCfv1gdwDxBfTB1LjFIqCPg3cBZ4zGtTPq5WnTfR\nuD7V873WrduKjtkb377Wqwqjy5Vbl5VSTYGvgGFa6+UXOVc0kO9+rqXGrw24N+5fz49xKb+Vfbg+\n0wCX59d7/WIopcJwtQDeAuporavjUurifiFFpXfzEvBnfBV6Ly6l9aYBkGPdWWu9X2v9Z611PVyO\njkluxdgDLLW0PqK01v+vmGv6h1Iqv5hlk7SPez8FTMbVmhuktT7ntXkTLkdMUdlIoAmu1uIxXM+g\njVf5Nvz8KdkBV6vxUHHnvpowuly5dFm5RrB8A7yitf63pSo+57Icq9T4rQFXSjVXSj2tlEp0r9cH\nfgekCMXnA9crpX6llAoBHsX5L3oVIAw4BBQq13AjcUiXUqoxEKa13ixt11pn4vJmP+4lXgA0U0r9\nXikVopT6LdASVwvHevzfFF0vLgeVBi64yzZTSt2nlAp1Lx2UUi2Kqccj7pdCWi7W//Z3oAUwQGt9\nyrJtLtBKKTVIKVUV16dimtdn6SfAC0qpWKVUc1wv/8fubf1wPSPvaw1zHwegilKqall9VvobRpcr\nry4rpRKAZOB9rfU/hHp8AjyllEpQStUDnvY6VunRfuCllxZcv/yf4fp1L3D//QCI1rLn/nZczqIi\nz/0q4D6prFumgabu/x8FDgDHcX1yzQRe1RbPPa7PsPctx9kJ9PJarw+cxtdz3xX40V23H4GuXtuW\n8LPnfrz7OvNxOVW8vdfX4VKcQ8ARXErTtgzvd0P3PTntPn/Rco9XmV64+jxPuevdyGtbGDAFOOG+\nl095bVsLtBfum7YsjcrqevxpMbpceXUZ15eKtpwn32u7ct+Lo+5lPKDK6lqV+ySVCnc/VTauB/Zd\nGR53AS6lvyqdcSXB7ahbDyToyqhsVxijy/6DP+qy33ahXC5KqT5KqerufsCRuH75pE/U0rAEKLOX\n6CohBnjaXxQ+EDC67Lf4nS5Xmha4Umo08H+4+gF/Ah7XWq+u0EoZDCXA6LLBKZXGgBsMBsPVRqXp\nQjEYDIarjZDyPJlSyjT3DWWK1tovhh6WRreDg4NtsvPnS55xoKyPVxqkulgp62u9cOGCTRaIPQ1O\ndLtcDbjBYLBTvXp1m+zIkSOO9pWGzlerVs0my8uzRsvLhjMkxJlJKCwsdFQuJibGZ12q7/Hjx20y\nqZx0zuhoa0Al5Ofn22Tnzp2zySoDpgvFYDAYAhRjwA0GgyFAMQbcYDAYAhTTB24wlAHBwcFERUX5\nyAoKCmzlQkNDbTKn/d0ScXFxNtmZM/aJfpw6CqtWrWqTSQ5AqT/aifM0N9eaxVU+Z5UqVWyyEydO\n2GTHjh2zySSfgtTPLt27o0eP2mROiYyMtMkkHZCoUcM3pbpUXwnTAjcYDIYAxRhwg8FgCFAqpQHv\n06cPGRkZbNu2jWeffbbU5QwGf8Gpzvbq1Yt169aRmprKU0+5JlVPTEwkOTmZTZs2kZ6ezuOP/5wp\n9vHHH2fjxo388MMPDB069Ipfh6GMKOe0mtb0oWW+BAUF6czMTN24cWMdGhqqU1NTdYsWLUpcziz+\nvZSn/l5sCQ4O1jExMT5LSEiIbQkPD7ctpdHtzp076xo1avgs0dHRevv27bpVq1Y6NjZWp6Wl6RYt\nWuj4+Hjdrl07DeioqCi9ZcsW3aJFC52UlKQ3btyow8PDdUxMjE5OTtatW7fWUVFROioqSkdGRtoW\nqX7BwcG2xXpPpP2qVq1qW6Kjo22L0/tUvXp12yKVi4uLsy2l0UWn90larM8wODjYkW77nROzT58+\nvPHGG4DLGdO5c+fLiqLq2LEjmZmZZGVlATBz5kwGDhzI5s2bS1SuInA6r0EgRpddTViDTG699Vb+\n8pe/AJfW7bCwMJusXbt2Np3t3r07q1at8inXq1cvsrKy2LPHNTPbnDlzGDhwIG+88Qb797tmY8vP\nz2fz5s0kJCRQvXp1Vq9ezalTpzh16hTJycn06dOHN99887KuV3KUSk5LK+Hh4TaZ5Jx0inRPJQer\nU4el9D5K53AaBCU5O0vqyPY7A/7ee+/RrVs3j6J5s2zZMjHKbPjw4Xz77bcAJCQkeBQXIDs7m06d\nOtn2cVrOYCgr3njjDbp27epIt4uMxvPPP09ycjLgXGfr1q1LTs7Ps5zt3buX66+/3qdMw4YNadeu\nHatXryYhIYHXXnuNuLg4Tp06Rb9+/Vi7dm3pLtZQLvidAV+wYAFpaWlMnz6dJ5980mdbt27dKqhW\nBkPp+eabbxzrttQCLysiIyP5/PPPeeKJJ8jLyyMjI4Nx48axaNEiCgoKSE1NrbDcKYbLw68MeOfO\nnVFKUbduXVGBnLTAc3JyqF+/vmdbYmKiT2ukCKflDIayoEOHDpel21IL3KnO7tu3j4SEn+cirlev\nnqdcSEgIn3/+OdOnT2fu3LmeMlOmTGHKlCkAvPbaa2RnZ5fmcg3lhT85Md955x396KOPetarVat2\n2Y6E4OBgvX37dt2oUSOPo6dly5YlLlcRi1LK0VLR9fSHpaKdlxdzYno7xyZNmqRHjBjhWLfDwsJs\ni1OdjYuL01lZWbp169a6Zs2aeuPGjZ5yU6dO1RMmTLDtU6tWLQ3o+vXr682bNxfrcLwSS2xsrG0p\nzfGszyEmJkZ0sDo9ntN3TzqvVM6ps9ORTfUnA96hQwedlpamN2zYoFNSUvQNN9xQogfYt29fvWXL\nFp2ZmalHjhzps23+/Pm6bt26lyx3pZagoKArulxthr6iDXfREhQUpCMiInwWSbfT09P1mjVrdPfu\n3T0jJKQRF8WNVrHqrPc+X3/9tW7WrJlYrkaNGrpfv35aa63T09N1WlqaXr9+ve7Xr59WSully5bp\nTZs26dTUVH3rrbeW6FmEhobaFmsZyVhLo0GKRsB4L07r4dSQSov0w+l0X6ejX6T3tqS6Xa4z8ph8\n4BAUdGWH3kvPszyfcXmj/SQfeHBwsLaGhJ88edJWTgrzlvJXS+lPT506ZZNJ6VSlkHNrqDbIozBK\noytSmgDrdcTGxtrKSKM8zp49a5NJaWIlrClswdloGJB9D1JqAgmnIfySDSgmh/kldbtSBvIYDAbD\n1YAx4AaDwRCgGANuMBgMAYrpAy8jpEgvqa/LicxpH5kkk4aoOZ0jMBD7yv2lD9ypbltTzoLct1ua\nVKcRERE2mdSn7HRaNGnorpTu1UmfulQ3KTLx0KFDjso5Tdcq4fRZOMXpe+sU0wduMBgMlRhjwA0G\ngyFAMQbcYDAYApRShdIrpXYCecB5oFBr3b4sKmUwVDRGtw2BQFnkQvmF1vpwGRzHL5GCDCSZ5MSU\nggKkYAerQ0g6luRwkoI9JJnTORJNAiMbl6Xb1nSi0rN2Ok+mU4el5GCUnNGSw04KMpFSokr1y8vL\nc1Sudu3aPuuSM1VCClCSAqPq1q3rqG6Sc1K6T05TwkrvY2kcllantdPAI9OFYjAYDAFKaQ24BhYp\npX5USg0piwoZDH6C0W2D31PaLpSuWuscpVRtYLFSKkNrvcy7gFv5zQtgCDSMbhv8nlK1wLXWOe6/\nB4G5QEehzIda6/bGCWQIJIxuGwKBErfAlVKRQJDWOs/9/23AmDKrmZ8gOSwlZ43ksJQcR5Jzxhqd\nJjlSJKeJ5NSRMtFJ1yA5Np1GZwZixOblUBLdDg4OtkUKSs406d5Jzi+nUYLSOSSkyE6nkY2SPtas\nWdPROay67fR92rZtm02WmJhok50+fdom27dvn00mIb0D0nsmRZ1KlCZS1FrOqUO0NF0odYC57gcS\nAszQWi8sxfEMBn/B6LYhICixAdda7wDalGFdDAa/wOi2IVAwwwgNBoMhQDEG3GAwGAIUk072Ejh1\nWErOScnREx8fb5NZp7uyTs0F8nRaR44csckOHjxokx0+bA8mlJydkkNIiuz0J/w5nazTacwkh6Xk\n7JOQHNnh4eE2meTsrFWrlk0mpXGVnJMtW7a0ySRdadCggU1mRYo6lKKCpeNnZ2c7Kie9K9K9k56Z\n0yhop5GdTh3PJp2swWAwVGKMATcYDIYAxRhwg8FgCFCMATcYDIYApSzSyVYaJMeRNM+dUyem5LC8\n5pprbDKro8dpRNeePXtsMgnJOenUYSk5cCp7JGZJUErZnM+S40xCcmrFxsbaZMeOHXN0PCllq5Si\nWIowbNasmU0mRSxKjtekpCSbrHr16sXWswjJmSg5e6UIywMHDthk0vU71VnJuS+975KTVbpP0oCE\nssS0wA0GgyFAMQbcYDAYAhRjwA0GgyFAMQbcD7hh1SoihExzFyMiP59u6elXqEYGQ9lw4+rVl63b\nkQUFdDe67QjjxPRCclhKMik6U3I8ShFs9evX91m/MSWFPvPmcdP69SwZNYozMTHExMTY9vOOVgvL\nzaXHpEnEZGcTW706P950k2ebNM+h5BCSot+ka3UaEXi1Oza11qITy4rk1JOcbk6j9SSkZyZFZ0oO\nO6mcFJ3oZG7Xm378kb6LFtElNZXlY8ZwJiZGjEz0dmKG5eZy6yuvEJOTQ5UqVfi+bVvPNkk/pWhn\n6fql+yk5cZ1GhTpNMSsNFpCwXpvTdLKmBV7BbG7VitzERGKys+kxZgxhl5jMNCw3lx5jxhCTnc2h\n2rXZ3KpVOdXUYLg8NjZvzonERKKzs7ll1ChHul1kvPfHxbFBGBFj8MUY8ArmZFQUS0aN8jHiVYRW\nNPga79zERKY/+CAnheFcBoM/UBAZyfIxY3yM+MV0u8h45yYkMOk3vyHfMhmEwY4x4H7AmZgYHyPe\nYcQIm6JbjfeSUaOM8Tb4PWdiYnyM+M0vvCDqtrfxTn7xRWO8HXLV9oE77duV+t2koAhJJg3ilwIg\nYmJiICaGdW+9RYcRI4javZtOzz/P5okTKYyLI+ToUa579VWisrPJb9CAdePHU7V6daKEDIXSOaV+\nR6m+0j2RZE77564mgoODbc9W6k+VfBRljdN+V6kfV/KXSLpy3XXX2WSNGjWyyWJjYyE2lrS//Y12\nTz1F9K5d3DJqFOvfeYdzsbGEHjtGu9dfJzInh4KGDUmfMIHY2FjaCHos1U0KcHMa3CP1bUtI+i7Z\nBQnpHkt+BklXnGBa4H7E2erVWTN+PCcbNyYiK4sWjz5K+PbttHj0UaJ27ya/QQPWjB/PWQfRbQaD\nP3EuNpb177xDQcOGRO7aRbunniIyK8v1d9cuCho2JHXCBM4JEaiG4jEG3M84W706mydO9Bjx1vfc\nQ0RWljHehoDHasQ7PvigMd6lxBhwP6QwLo7MV1/1kW0YOdIYb0PAcy42lk2jRvnIfnrpJWO8S4gx\n4H5IyNGjNH3hBR9Zm7Fji/XgGwyBQuixYySNGeMja/nyy4Q6TNRl8OWqdWJKSMEokgNDytAnyaTA\nDsnpYg3Sue7VV4lw93lvGDmSNmPHErV7NzcMH+4J9gHIFyLcnGYZlBw4xjlZcs6fP+/IESUFd0mO\nYqeZDJ0iBd5I+i4FpElIupeVlWWTeWfpDMvN5dZXXyXSPZJq9ZNP0mnCBGJ27aLV44+z9KWXPLot\nTe22f/9+m0yaQlCaflC6fqmclAVSCniSkM4hOU+l7IbW8zo9p2mB+xFFQwW9HZb5jRqxZvz4ywr2\nMRj8Dc9QQbfxXjZ6NCcaNGDZ6NEe3e7+8stGty8TY8D9BO9x3laH5dnq1W3BPkbRDYGCzzhvt/Eu\nammfiYlh6UsvGSNeQowB9wOsQTrSaBNrsE+PMWMuO0mQwVDeWIN0vI13EZIRN7rtDGPAK5iI/Hxb\nhGVxo02sRvyeyZONohv8lsiCAluEpdV4F2E14nf/4x9Gtx1w1ToxJQeO5EySykkOS8k5KTmidu/e\n7bPePT3dk5hq+n33cfLQITHC0vv4O//wB37/z39S6+BB4lesYJlXQqvDhw/b9pWmrJIcm5IT82rP\nMng5WCMWpQhYyXEmlZMi/UrjZJZ0tk6dOjaZlLVPcmyuXbvWJrNOF3jbtm3E5ORwqHZtZvzxj5w8\nepQEIQpx06ZNnv9TBw7k4ZkzqXPgANUXLybVa5o26R2TokQlnZUcjNI9djptnTQIIEII/5cc26XJ\nNGnlqjXg/sLSVq2IjY1lc6tWjnObnIyKYsZDD1Fn+XIf420w+BOLrr2WxMREMq6/3rFuF0RG8sHg\nwTRdv54lwhybBl+MAfcDvPN5O+VkVJQx3ga/Z13nzpe9T0FkpDHeDjF94AaDwRCgGANuMBgMAYrp\nQvFCchJJzgopwlKKTJMixyQHi9XRIaWElSIsJSepFJkmRXVJEWLGiVk6rLoiOf+k++nUqSVNx+Y0\nPa00LZrkYJOmKJOcfZKzXHLGWnVPSkMrOXallLB79uyxyZwOPHAa2VqtWjWbTHo+kvO0pClhS4Np\ngRsMBkOAYgy4wWAwBCjGgBsMBkOAYgy4wWAwBCjGiXkJJCem5ACUnJgS0r5WB4vkDJLqIUWmSQ4X\nqW6lmSPQIGN1qEnObkkmpZiV5n906rCUnKeS01FCcsRJUZyS41FKk2q9DukapHlipfpK+0r1rVu3\nrk3mdGCA5BSVHJalmRNTeves90CKnpYwLXCDwWAIUIwBNxgMhgDFGHCDwWAIUIwBNxgMhgDlkk5M\npdQU4JfAQa11K7dsNPBnoGjiupFa6wVXqpIVidM5MSWnjuQUdDJfn+Q0keohOUOkNLFOy0kRbJU5\nErOsddt6ryQHluTEdOoQlPaNEfJrlyYiUEq7um/fPptMujbJ8WbVW6cpXHft2mWTSc5ZKXJScoBK\n+u4U6X2Url9CevcknA6CsOKkBf4xcLsgn6C1buteKqXxNlR6PsbotiGAuaQB11ovA+xjmgyGAMfo\ntiHQKU0f+GNKqTSl1BSlVGxxhZRSQ5RSa5VS9ik8DAb/xOi2ISAoqQH/O9AEaAvsA94urqDW+kOt\ndXutdfsSnstgKE+MbhsChhJFYmqtPbkelVIfAf8rsxoFAJJjz2naWcmZIjlJnCCd02lK3KvNYemU\nstRtKSpWmjdRcv5J5SRKqjsgRyc6TbsqOVlzcnJsMmtUpJTuWKqH5LB06ngPF+bdlPaVIkCl1MsS\nUqSshOTslKI9JVvhhBK1wJVS3k/l10B6ic5uMPgZRrcNgYSTYYSfAj2AmkqpbOAloIdSqi2ggZ3A\nw1ewjgbDFcHotiHQcTIK5Xda67pa61CtdaLWerLW+j6t9fVa69Za6zu01vaBouVMdHQ0c+bMYe3a\ntaSlpfHggw+W6fEnT57MgQMH2LhxY7Fl+vTpw6ZNm8jIyGDEiBGOjhsdHc2nn37KihUr+OGHH/jj\nH/9YVlUG4LbbbuOnn35iy5YtF61Tnz59yMjIYNu2bTz77LOXlEPx9yQrK4u0tDTWr1/PmjVryvR6\nypJA0G1v/bgSet27d282bNhAeno6w4cPF8s8/vjjbNy4kfT0dIYNG+azbeLEiWzfvp2UlBRH56tW\nrRqffPIJ3377LcuXL+fee+8t9TVY6dmzJ2vXrmX9+vU8+eSTYpk+ffqwbt06UlNTeeqppzzysLAw\nVq9eTWpqKunp6YwePRqAxMREkpOT2bRpE+np6Tz++ONlXu8SobUutwVXq+aKLA888ICeOnWqZ71q\n1aplevxbbrlFt2vXTm/cuFErpWxLaGiozszM1E2bNtVVq1bVqampunXr1jo0NNRnCQsL81mGDBmi\np0+friMiInRERISOi4vz/H+pxXqssLAw27kyMzN1kyZNdFhYmE5NTdVJSUk6KCjIp+7BwcE6MzNT\nN27cWIeGhurU1FTdokULHRQUJMqle+J9r7KysnSNGjWu2LP2XspTf8tCt6XneLFyjzzyiEc/JL2O\ni4uzLZdTl+3bt+vmzZvratWq6Q0bNvg8X0AnJSXpjRs36vDwcB0cHKwXL16smzRpoqtVq6arVaum\n+/Tpo7t27ao3bdqkg4KCbEtCQoLP8thjj+mZM2d66lqvXj1xqVGjhm2Jjo62LVFRUT5LdHS03rFj\nh27durWuUaOGTktL0x06dPC5piK9btWqlY6NjdVpaWn6xhtv1FFRURrQkZGRGtAhISE6JSVFd+rU\nScfHx+t27dppQEdFRektW7bY7pW0VKlSxbZI721JdbvShNKvW7eO7t27s2bNGkaPHu1xCiQnJ9Or\nVy8AXnnlFd59990SHX/58uUex4V0I9u3b09mZibbt2/nzJkzzJo1i379+nHu3LmLLmvWrKFLly4s\nXbqU5557jhMnTtCkSRMWL17M6dOnOX36NC1atODLL7/0rBctlzp2u3btyMzMZMeOHZw9e5ZZs2Zx\nxx132OreoUMHMjMzycrK4ty5c8ycOZOBAwfSsWNHUS7dE4Mzzp49a1uio6Nty8mTJzl58iSrVq3y\n6Iek11prnn/+eV5//XXvHxMbISEhtqV169Zs27aNjIwM8vLy+PTTT32eL0CLFi1YvXo1p06d4vz5\n8yxdupQ777yTwsJCCgsLWbZsGYcOHUJrzYULF2xLTk6Oz7J27Vo6d+7M4sWLGT58OCdOnCA3N5eE\nhAQ+//xzcnNzyc3NpU2bNsyZMwellGeRriE/P99nadmyJVu3biUtLY0jR44wY8YMevbsSfXq1T1L\nz5492blzJ+np6Rw7dowZM2bQu3dv8vPzUUpx8uRJlFJUqVKF0NBQqlatSl5eHlu3biUyMhKtNZs3\nbyYhIcHnXsXExNgW6XlL98n7OiUnZ3FUinzg0dHRjB8/ntatW1NQUEBycjKpqal88cUXvPTSS4wZ\nM4batWvTrl077rjjDs9+y5YtE0Nxhw8fzrfffntZdUhISCA7O9uznp2dTceOHS9Z73HjxtGuXTsK\nCgr45ptvSE1N5b///S/XXHMNQUFBXLhwgbfeesv2ebtkyZJL1j0hIcFnItjs7Gw6deok1l0q53R/\nK1prFi1ahNaaDz74gI8++uiS+xjsONHrGTNm0Lp1a37/+9979pP0WinFiBEjSE5O9sicPN/09HRe\ne+014uLiOHXqFP369WPt2pINe4+OjubVV1/lpptuoqCggAULFpCWlsa8efPYvHkzjRs39uj8mDFj\nGDVqlGffefPmiSNTnnzySZ931ck11a1b12fEjLVMUFAQa9eupWnTpkyaNMl2vQ0aNKBdu3asXr26\nRPehLKkUBvzhhx/m66+/9gwBWrVqFfHx8YCrlaiU4qmnnqJHjx4+w+y6detWIfUtYsiQIT71TklJ\nIT4+Hq01mzZtIikpiWuvvZbdu3ezfv16n3179OghHlMaRljedO3alb1791KrVi0WL15MRkYGy5cv\nr+hqBRxO9Hro0KHccccdl9RraZIQJ2RkZDBu3DgWLVpEQUEBqampJZ744+GHH+abb77xXM/q1aup\nU6cOgKdV26JFC5o2bUp2djZpaWmefQcMGCDq9pX4Arxw4QI33HADMTExzJkzh5YtW/LTTz8BruGN\n06dP54knnhCHiZY3lcKAt2t+GyksAAAgAElEQVTXjqlTp/qsz5s3D4BWrVpRt25djhw5YksYU5Yt\n8JycHBITEz3riYmJ4rhYa70/+eQTz3rbtm353/9cw45Xr15Nly5deOSRR+jXr59tXyct8JycHOrX\nr3/JOhVXzun+Vvbu3QvAoUOHmDt3Lh07djQGvAQ40evc3FxHei21wJ0+3ylTpjBlyhQAXnvtNZ8v\nzcu9nlmzZnnW27Rpw1dffeVZX7NmDZ07d+bPf/4zgwcP9tnXaQvcyTXt27fPp/ujuOvOzc1lyZIl\n9OrVi59++omQkBCmT5/OrFmzmDt37mVc+RXEHx09l7tMnDhRP/fccxrQ/fr106tWrdJKKR0fH683\nbNigmzdvrhctWqT79OlTqvM0bNjQ5rArWkJCQvT27dt148aNdZUqVXRqaqpu2bKlrZy3g2fSpEn6\n+eef18HBwXrAgAE6JSVFh4SE6KCgIH3HHXfow4cP6zFjxojOoeIW73MFBwfb6pSUlKSVUmK5Ro0a\neZyVLVu2LFZ+sXsSERHhcQZFRETo77//vtT3/WJLRTsvL1e3Q0JCbIvknHOq1999952+6667dGxs\nrGdxel4nzxfQtWrV0oCuX7++3rx5s46JidHh4eGe5brrrtPp6emO3tNRo0bpyMhIfeedd+rVq1fr\nqKgoHRkZqSMjI/Vdd92lDx8+rF9//XVds2ZN2+LEYVvcNVWvXt2z1KhRQ2dlZYnXXatWLV29enWt\nlNLh4eF62bJletCgQToyMlJPnz5dv//++x4np3WJiYmxLVI568CG0NBQ26AIp7rtl0p+uUujRo30\n6tWrdWpqqp47d66uW7euDg8P1ytXrtS9evXS4BoxsXLlyhKfY8aMGXrv3r367Nmzes+ePfpPf/qT\nBvT8+fN13bp1NaD79u2rt2zZojMzM/XIkSNLVO+ibU2bNtU5OTme0QfSyBcn9b5YnZzU/WL7S/ek\ncePGOjU1Vaempur09HRH96E0S0Ub7qIlODjYx4jGxsZe8gf2Yi90bGysbtOmjV67dq3euHGj/t//\n/lfmel2lShV9xx136K1bt+rt27frF198UdSNZcuW6U2bNunU1FR96623OnovLlffrTp/qRFWRcvl\n6LwTfb/++uv1unXr9IYNG/TGjRs996RLly5aa603bNig169fr9evX6/79u3rc15phJH0w1mWuq3c\nhrVccBsdgwPee+891qxZ4+likTzT5fns/BWttXOX/RXE3ZL2kUl5uaV+XCkHtxTmLeXNLg1OJ9wt\nL7x1XgrVl+5dafJ8lzVS+gPpfkrzCUg40e1KM4ywsnDNNdewefNmwsPDffrHDYbKitH5kmNa4AGC\naYHLmBZ4yfG3Frg3pgVuWuAGg8FQqakUwwivBkxr2/+xPiNp/kdpKJw0nthpK01CGl4qncNpKlrp\na0C6toKCgkuWc9piLml61dISG2ufv8Ppl4/0zKSvq9I8WyumBW4wGAwBijHgBoPBEKAYA24wGAwB\nijHgBoPBEKCU9zDCQ8CucjpdTeBwOZ3rSlIZruNKXUNDrXWtK3Dcy8bo9mVTGa4BKli3y9WAlydK\nqbW6EswWXhmuozJcgz9RGe5nZbgGqPjrMF0oBoPBEKAYA24wGAwBSmU24B9WdAXKiMpwHZXhGvyJ\nynA/K8M1QAVfR8D0gSul7gH+qLW+7QqfpwcwTWudeKmyJTz+V8BMrfXUSxY2XBUY3TaUFL9qgSul\nuiqlViqlcpVSR5VS3yulOgBoradfaQW/SL3qKaWy3f/vVEodVEpFem1/SCm1xMmxtNZ9/UHBlVJh\nSqnJSqldSqk8pVSqUqqvpUxPpVSGUuqkUuo7pVRDy/5TlFInlFL7lVJPCed4Xik11iIb5c5n3uty\njhXoGN0uP8pbt933KFMpla+UWqiUqudVTimlximljriXcepyZi2+BH5jwJVS0cD/gPeAOCABeBmo\nmKQIvvQDFnqtBwPDKqguZUUIsAfoDsQALwCfKaUaASilagJzgBdxPY+1wCyv/UcD1wINgV8AI5RS\nt1vO0R9YULSilGoC/AbYZynn5FgBi9HtcqfcdNv9VTMWGOg+VhbwqVe5IcCvgDZAa2AA8HBZXCRA\nhc9k4jVbT3vg+EW23w+s8Fq/DdgC5AKTgKXAQ+5th4EC4CBQ6L6pfb32fQDYDOQBO4CHvbb1ALIt\n554D3On+fyfwHHAUqO6WPQQs8Sp/M7DGXbc1wM1e25Z41bOpu965wGn3ku7e1hzYDlxwy3cC/a7w\nM0gDBrn/HwKs9NoWCZwCmrvX9wK3eW1/Bfgv8B3wk9f9Dcb1QuQAJ9zXdADo5bWvdKyZFa2T/qbb\nbh3YA+QD+4FjRrcrXLdXup9Hqnu5F9eMOk3c+64Ehngd60Egpayuy29a4MBW4LxSaqpSqq9Syp4W\nzI37F3Q28DxQA5ey32wpVgX4CxAGjAcme326HAR+CUTjUvgJSqkbijlXKNANWOwlXotLWYcL5eOA\n+cC77rq9A8xXStUQDv8KsAiIdddnqPsYke7zbcT1QnUGonApulTHSUqp48UsadI+wjHqAM2ATW5R\nErChaLvWugDXS5fkfjZ1vbe7/28GPK21bgmMw6XI17m3fwt8p7VugutlKTpvccdKclLvAKEsdXsc\nLp1+EVcQidHtS1AOup2mtW6rtW6Ly8gDtJLORRnrtt8YcK31CaArrhvzEXBIKfWl++Zb6Qds0lrP\n0VoX4lKo/ZYye7TWH2mtzwNTcT2UOu5zzddab9culuJStFuKqVo3YIPW2pqPcxTwf0opa7RUf2Cb\n1vrfWutCrfWnQAauTycr53B9ptXTWn8DFE0Z/ktcCp0KXNBarwc+x9X9YENrPVRrXb2YpXUx1+XB\n/SJPB6ZqrTPc4ihcrSdvcoFq7m1YtucC4Vrrde71nsA2XN0FVYC+yJ/mxR3LnhM1QLkCur3L6Lbf\n6HYm0EYp1VopFY7r3mmgaHYH67lygaiy6gf3GwMOoLXerLW+X7u85K2AesBfhaL1cH1KFu2ngWzv\nQwG1lVI/KqWGaK1PuuVRAO5WUIrbmXQc10tTs5hq9cOrH9frnOm4+jWfE+pmDanehethWxkBKOAH\npdQmflbihkAn97HHK6XOA38CGhRTxxKjlAoC/g2cBR7z2pSPqxXnTTSuT8d8r3XrtqJj3o7LqKzG\n9ekeAnyplJqCr95d9FiVhTLSbY2rZRyvlBri3m50uxjKSbd34OoGWoOr62m/u2zRM7OeKxrIdz/X\nUuNXBtwb96/lx/z8KeLNPsAzFMr9a+Y9NOo1XL/wfYFHlVLdvMqG4frFfwuoo7WujkuJi/tFFJXc\nzUvAn/FV4L24lNSbBrj6gH3QWu/XWv9Za10Pl2PjVVyt1T24+g8bAqHuZXxxdVRK/cPtAZeWTdI+\n7v0UMBmXMg7SWntn29+Ey/FSVDYSaIKrdXgM1zNo41W+DT9/onbD1QJ53N36jAbOA7WBu3E9q8+U\nUs86OFaloxS63RWXP2ETFr12lzW6/fN+5aXbf8el1+G4upSScDVW0qVzUda6XVad6aVdcDk2ngYS\n3ev1ge+Bj7TF0YOrRZGHy7sbguvX9Rw/O1C8y47G1WrRuBwr1XAZk+64lKYvcBJ4VVscPUBjYIel\nnjvxdcB9BBzB7ejB1Td4HPi9u26/da/X1HZHz2+8rjcJ1y/5FncddwH38bOS34Hr87Us7/k/gBQg\nSthWC9fn3iCgKq5+vxSv7W/gehFj3c9uH66WSSiuz8qvvMrWAOLdSwf3s/pN0XmLO1ZF66Q/63aR\nXrtlRrcrTrer4vohVrh8FQXAWK/tj+Byeibg+oLZBDxSZtdZ0crtdaEJwGe4fs0L3H8/AKKtSu5e\nvx2Xc6jIU7/KrRSR7pu2wv3/SndZDTR17/sorpEQx3F9Ys0sRskfA96/hJLXdyvnEi9ZV+BHd91+\nBLp6bfNW8vHu68zH5UR5np899dcB3wCHcL1E24CFZXi/G7rvyWn3+YuWe7zK9MLVx3nKXe9GXtvC\ngCm4RpYcAJ5yK/En7vX2XmXrev3/pPv59rrYsSpaH/1Qtx/CZfzud+v0Stw/cka3K1S3r8M1wqXA\nXf4nINhru3Lfi6PuZTzuAMqyWAImEvNiuPulsoF7cP26z3VvCgFmaK1fK+FxF+BS8uI+M8sMpdSn\nuF6wmriU5CX3eltcyrgT15Aw6xhqv0Ep1RVYjmvoZtFn4kjgdwTQdfgTXro9HHjWLS6VXruPa3T7\nMvBX3Q5YA66U6oPLiXAKeAZXy+MarfWpi+54eecYAbxXlses7CilmgE3atcIBUMJMLrtn/ijbgey\nAR8N/B8ux8hPuJwKqyu0UgZDGWB02+CUgDXgBoPBcLXjt8MIDQaDwXBxQsrzZEopW3NfCkgKDg6W\n9nV0jnPnzl26UDFUqVLFJjt79qxNJtUvKMj+WyjVRboOJ19BYWFhjvaT6isRGhpqk0nX4PT6z58/\n76h+ZY3Wuswyu5WGsLAwHRER4SM7fvy4rZx03yU9kZ7FhQsXSlHDkmO9LoCTJ08KJS9NbKw9i8Cx\nY8cc7RsTE2OT5eZaAypLR2RkpE0mvQPSM3N6bZINqF69us96QUEBZ86cuaRul6sBBwgJ8T2lZAyk\nGyGVk8jJscUUOCY+Pt4m2717t00mKZJkYPftszujpRfYidGtX7++TVZYWGiT7dy585LHAvlapRdV\nOp70fKQX6dSpq8c/FhERwS9+8Qsf2dy5c23lateubZNJOisZkry8iglObd68uU22bt06oeSl6dmz\np002e/ZsR/t27drVJps/f36J6lEcbdq0scmysrJsMund7tWrl032n//8xyarWrWqTXbbbb7ZhBct\nWnTRehZhulAMBoMhQDEG3GAwGAIUY8ANBoMhQCnXYYSSE7Nu3bq2clL/UmmcH9dcc41NtmPHDkf7\nNmnSxCY7fPiwTVbWzhQr1arZs6uWpk/0uuuus8m2bNlik0nONKvDBeDo0aOOzis97/37rdlS5XNI\nz9tfnJiSbt90002O9pV8GdJ9mjdvnqPjNW3a1CbLzMx0tK/Uzyzp3vr1620y6TlaadasmaN6SO/d\nV199ZZMlJdlTa2/aVPJcUS1btrTJJB+F0/e9QQN7kkXJrybhRLdNC9xgMBgCFGPADQaDIUCplAa8\nT58+ZGRksG3bNp599tnLLhcTE8N//vMfNm/ezMKFC2nXrh1169Zl2rRpLFy4kFWrVvHII4+Ux6UY\nDD441e3XX3+dlJQUcZhdUFAQ69atc9wlY/BfKp0BDwoKYuLEifTt25eWLVvyu9/9TuwTlMq1aNEC\ngL/97W8sXLiQFi1aMGDAADIzMyksLOT111/n9ttvp3fv3jz00ENiP7LBcKW4mM5amTNnDn/605/E\nbcOGDWPz5s1XsqqGcqLcA3msHDlyxGf9tttu45VXXgHgzJkzdO7cGa216MCSnCtJSUlkZmZ6Bt/P\nnDmTG2+8kW+++can3J133snevXsJCgqiQYMGLF68mD/84Q988skn/OIXv+CFF14gMTGRjIwM2zkS\nExPZsWMHTZo0Ye/evQCcPn3aVu7MmTM2mfRjcuDAAZ/1evXq2cpIDkaJuLg4m8ypg1FCik6VAnSk\nQKPs7GybTHJQS0FFUsCTU6d1RRATE0OPHj18ZP/973991vv06cMbb7wB+Op2VFQUVqw6AdCxY0eb\nbv/lL3/hiy++8Ck3e/ZsMjMzadiwIWfPniUzM9MT/ZeQkED//v0ZO3YsL7zwgs2pKOnxjz/+aJNJ\nOnDDDfa5k626JzkdpYCnrVu32mQSksNSesesAYQgOxN/+uknm6x9+/Y2mfRerFy50iarVcs6ragc\ntVzSAMQKN+BWJkyYwC233CJ6tJctW+ZjtItGSLzwwgssWbIEcCnonj2eKQXJzs4WRwPUqVPHx5js\n37+fVq1aUb9+fY4ePco777xDixYtSElJYdiwYT6hw4mJiSQlJYmeeIOhON577z26desm6vbXX3/t\nY8iLIo9fe+01VqxYAci63b1798uqw4QJE3j22WfFxo8h8PA7A75w4ULS0tKYPn06Tz75pM+2bt18\npgC8IkoYEhJCq1atePHFF0lNTWX48OE899xzjBo1CnCFOH/wwQeMHj2a/Pz8SxzNYPiZBQsWFKvb\nffr08VmXhs2Wlv79+3Po0CHWrVt32Ybf4J/4lQG/6aabUEpRt25d8TPDSQs8JyfH53M+MTFR/Dw5\ncOCAz1jb+Ph49u/fz759+9i3bx+pqamA63P0uedck3OHhITw+eef88UXX7Bw4cLSX7DhqqFz584X\n1W0nLXBJt61dkBejS5cuDBgwgL59+1K1alViYmJ48803eeaZZ0p6WYYKxq8M+KBBg9i2bZtHwatV\nq+YTrOKkBb5mzRquvfZaGjVqRE5ODoMHD+bee++1lUtLS6Nhw4YkJiZy4MAB+vfvz2OPPcahQ4fY\nt28f11xzDTt27KBnz56efrHJkyezefNmPvroo7K8bMNVwG9+8xu2bt1arG47aYFLuv3xxx87rsPI\nkSMZOXIkAN27d2fUqFHGeAc45RqJGRwcrK2ZuLz7ljt06MDkyZPRWnPq1CmGDh3qyXomZfCSHC5B\nQUH07duXd955h+DgYP71r38xa9Ysz/aPPvqIv/zlL6xcuZK+ffvy17/+leDgYKZMmcKMGTMAaNGi\nBW+88QZVqlQhIyODBx54gKSkJFasWEFaWponrefIkSM90WGSk0TKWii1mKwZ6g4ePGgr4zRiNTEx\n0SY7dOiQTSY5WEuTNlRKkSkZIadZGyWsGSnPnz/v15GY3lxMtwcOHGgrL2Wo/Oabb7j99tt5++23\nCQoKYurUqbz++uue7V9++SUPP/wwe/bsYcaMGfTo0YOaNWty4MAB3n77bT799OeZwG6++WZeeOEF\n3nzzTZ9zzJw503bevn372mRSVKTkeLaO1JIy+xUUFNhkkmNbck5K74rkAG3btq1NJvkiJJnTyFYp\na6Pk7JUc1JItc6LbfmXAL8blGHArUii99AAaNWpkkzlNz2oMuDHgJd33cgy4Fen9lULz69SpY5NZ\n09+CMeCBZsAr3Thwg8FguFowBtxgMBgCFGPADQaDIUCp8HSyEtJ0UlI/WXR0tE124sSJEtTMhRRd\nJfUnSn3b0r41atSwyaRIr7JE6oeUUlpu377dJnPqA5AiRaX7Lt0nqW9X6qOXSEhI8Fk/ePAgZ8+e\n9Ys+8JiYGN2lSxcfmRSuLt3Pu+66yyYrivD1Ror0k/qxv/vuu4tVtUxwOrertS/bScrZ0uI0VbKE\n9CykKd+c+qScYk2fm52dzenTp00fuMFgMFRWjAE3GAyGAMUYcIPBYAhQShWJqZTaCeQB54FCrbU9\nbZfBEIAY3TYEAqVyYrqVvL3W2j5JpFzedjJpkLxUJykA5vjx405OWyqkdJDWgBKQnTOSQ1FKp2kN\nbpDm2+vUqZNNJs3NKTknW7dubZNJzlTJYSshOY8lp45Tx5GEFBgl1e9KBfKUhW5LSNclpYSQUuda\nA75ADmRxiuR4l5zM1jB/cOVuKQlOHeW9evWyyb7//nubTEpjLNkFp/fJmq4D5PdYCryR9F16Rxs2\nbGiTWfPjHDhwwJGD3nShGAwGQ4BSWgOugUVKqR+VUkOkAkqpIUqptUqptaU8l8FQnhjdNvg9pc1G\n2FVrnaOUqg0sVkplaK2XeRfQWn8IfAilyxdhMJQzRrcNfk+pWuBa6xz334PAXKBjWVTKYKhojG4b\nAoESOzGVUpFAkNY6z/3/YmCM1rrYmQ6qVKmirY4YabIFp9m/JAejlLVQimqTkuo7RTqvlHVMcnZ6\n54AuQrpeK1L0o3QsydkrOYlq1qxpk0kOFym7oTTXpVMkJ9auXbtsMqd6eSWcmFdStyWH4KpVq2wy\naRpAyUkm3bvS0Lt3b5tMiux06vC2RmxKz/Wll16yyb799lubTJrDUspaKL0DknPfqZ2pKJzodmm6\nUOoAc90PKASYcTEFNxgCCKPbhoCgxAZca70DaFOGdTEY/AKj24ZAwQwjNBgMhgDFGHCDwWAIUMp1\nUuPCwkJb6lBpGi+njgSnaUglpOgqa7pSkB12UuScdB3h4eE2Wbt27Wwy69RjkhNGcpzu2LHDJrvx\nxhttMikVr3S8ork+vTl37pxNJiEdT3o+Tqeok7A6qKVp4SqK8+fP25zK0jReUiSmdO+k++50ejsJ\nKcJw2bJlNtnixYttsmHDhtlkUsSmNMGy1RkrTSfmPVdnEdI0c5IT84cffrDJpCngpBTVaWlpNtng\nwYNtskWLFtlkR48etcmcIulFampqiY5lWuAGg8EQoBgDbjAYDAGKMeAGg8EQoBgDbjAYDAFKhc+J\nKTkOpQi2isJptFarVq1ssvbt7SmkpUhRa0pMKf3r0qVLbTLJSSg5Q+rUqWOTSVGX69ats8kkZ6/k\niJPSa0r7Ss4vKTpTQnKAXql0spdLbGys7tmzp49MutZ58+Y5Ol5YWJhN5tRpe/vtt9tkCxc6i0Ny\n6uyUkJxzI0eO9FmXIiKld0cqt3atPWfYtGnTHNVDimKVomIlp6gUyX333XfbZJ999plNJiENgpCi\nqp3otmmBGwwGQ4BiDLjBYDAEKMaAGwwGQ4BiDLgfcGNKChH5+Ze1T2RBAT0zMq5QjQwGQyBQrk7M\nkJAQHRMT4yOTIpokp5sUweUU6XiSg+3EiROOZM2bN7fJpNSx0jyRf/rTn3zWr120iPYff8zJxo3J\nmDSJwrg40tPTbft5z6MXeuwYNwwfTtSuXSy7+27Su3f3bJMcrJJTODk52SaTHJsS0rWWJj1vaVLW\n+osTMzQ0VFtT9ErzpDpFSiebkpJS4uNJSBGgkuP9m2++scmcRmf++te/9vxfb+5c9nbpwtnq1X3K\nSHpX5MgPPnKEmEWLOPq73zFo0CBbOSnqUopklt4BKU1uaXT7t7/9rU02a9YsR+fwvk/gSqd79OhR\n48T0d3Z36sTJxo2JyMqi+dChhFwiRNfbeB+Nj2f7DTeUU00NhpJTb+5cmr37Lu2feYYqDicjDz5y\nhMYPPki9sWOJE8LtDcaAVzhnYmLImDTJx4hXFVr94Gu88xs25L9PPMEpYUiSweBvHOrRg4KGDYna\nvduRES8y3lW3b+d0kybk3nZbOdU0sDAG3A8ojIvzMeJ933zTZsStxnvdW28Z420IGM7FxpI6YQL5\nDRpc0oiH5eb6GO+syZM5X6NGOdc4MCjXbITBwcFER0fbZFZK098tUd3S5wb2qZ4A8h06EjME52HH\njvYpE6+55hqbTApQaN++PcTGsn/GDBLuu4/YzEzufP99j+IGHzlCg/vvJ3LXLgoaNWLjX/8KsbGE\n5ebajrV161abTOr/i4uLs8mkvkOrzwLkZ3ZceBmjoqJsslyhztLzlvwHVr9FafqYy5qwsDCaNGni\nIytN/TZv3myTNWvWzCaTnrdEDcEA9ujRwyaTAlmkzIhSoJGU9fK6667zWd/6wQc0HzqUqKwsOj73\nHBmTJhHkNRVd6LFj3DhuHFV37SKvfn1WvvQSZ0+cgBMnRL+A1Z6A8/5uidL4cqSp8aSpEKXAIOu9\nc1oP0wL3I87XqEHOv//N6SZNqLp9O40ffJCwbdto/OCDRO7cSUGjRqT99a+cs6SfNRgCBevXZvOh\nQwk9dgxwG+9nniGqyHi/9prN4WnwxRhwP+N8jRpkTZ7sMeLX3nknVbdvN8bbUGmwGvEbn3mGSPff\noi5CY7ydYQy4H3K+Rg32vPmmj2zzSy8Z422oNHgb8ahdu+g8ZIjHeP/45pvGeDvEGHA/JPjIEeo/\n84yPrMXLL3s+NQ2GykBhXBzbX3vNR7bxL38xDZXLoFwDeYKDg7V1aiPJQSJNfyQ52CQk59ypU6ds\nMslJJjk/pAxwksNSCmKQAn4kp0uLFi08/4fn5fGbv/+duH37OFq3LosffJDekyd71uc9+SSn3ceQ\nAjHetLTcQc7EJjnEpHssZUmTnI6So1hybDpFmo5Oeo7+EsgjZdqU+NWvfmWTffHFF2VeHytWByvI\n2R1vvfVWm2zDhg022YgRI2wyKTOgNdPkpEmTPP8HHT5M3F13EerliD3TtCl7p03jzkcecXTOd955\nxybr1KmTTTZx4kSbzClOA3Qkh7/ktHeagdVkIwwwwvPyuGviRI+x/vKJJziWkMCXTzzB0bp1idu3\njwETJhQ7TtxgCBS8jfeZpk3ZPX8+Z5o2JSwzk3r33kusw7lYr3aMAfcTiox3jQMHPMa7qKV9Ojqa\neU8+aYy4oVLgbbzPNWvG3mnTOOv+W2TEP9q+3RhxBxgD7gd4G+8jder4GO8iJCMujQM3GPwZq/E+\nOnu2J0jnfI0aHiPe5PRpY8QdYAx4BVPt1Ckf4z370UdtxrsIqxG/ZdQoY8QNAUPMmTM2433BkgCs\nyIhvr1rVGHEHlLsT0xqdJ2X7k4gVPNPHynhUhuSEkKYPk5x9ktOtf//+l9x3YHY2w7Zt41CtWnzy\nwAOcjIoSHbve2fmqnTrFiIULSTh+nM+6dWPZ9dd7tknRpPPnz7fJpHNITlHJuSI5E6XoTGt2PnAe\nZdugQQObbPfu3TaZPzsxe/XqZSsnZfazZqIDmDt3rqPzShn6Pv/8c0f7OkWa8kzSHylS8uGHHwYg\n5OhRkv7v/4jIyiKvfn1WjR3rGSooZaMc9/TTvPHDDzTMz2dXVBTPdexI73vusZX76KOPbDKnWRul\nqfwkx67EgAEDbLJdu3bZZGlpaTaZdD+lDKROdLtcQ+kNdv6bmEizZs34KSmJk0LouUReeDjjb7+d\nLnv3+hhvg8FfqfHdd6LxLo7csDCe69jRY8Rv8aO0Cf6EMeB+wFph2NOlyAsPN8bbEDAccH8pZLZt\n6zhIp8iI37J/P/9r2JC7rmQFAxTTB24wGMqFA4MGXXaEZW5YGP/zmszE4Isx4AaDwRCglKsTMzQ0\nVFsj9qTplJo2bWqTSUfAqyIAABQcSURBVFOFSelaJUeC09SMUspNpxGgUtpIadq2zp0722TW+nlH\nZhbx9ddf22SSQ1ByxEpI6UqlaNL69evbZHv27LHJpOhZySkq3U8pilM6h4Q/OzFLgxTtK6V6jY+P\nt8mk92fFihWOzut0OkPJEbdjxw6brGXLlj7rjz32mK3M008/bZPdfPPNNpmkO5Kz++TJkzbZl19+\naZPdfvvtNpkUdbpv3z6bTHovunbtapN96nAmIet7m5+fT2FhoYnENBgMhsqKMeAGg8EQoBgDbjAY\nDAGKMeAGg8EQoJTrOPCgoCDbPImSE1OK9JNSMErOydLMaSc5AJ06MaV57iQHy08//WSTWefJlO6J\nlNZWcnZKzirrvIQgOyylyDQpElNy9l64cMEmO3jwoE0mRfBJDsvyiLwtS8LDw20peiU9luawlPRE\nuseSY09y2ks64HRggNNIWSlyUErlbHVkS7otpdjduHGjTSa9YytXrrxoPYuQHJZSqmTJYSnNn3v6\n9GmbrDTpk6W0s04wLXCDwWAIUIwBNxgMhgDFGHCDwWAIUC5pwJVSU5RSB5VS6V6y0UqpHKVUqnvp\nd2WraTCUPUa3DYHOJSMxlVLdgHzgE611K7dsNJCvtX7rsk4mRKtJEYySs0KKarQ6RAFq165tk0lz\nQkrlJKebNHek5IgKCrL/FkrRalK61/2WTGvSM7nxxhttMslxKEXDFRQU2GRS5KTkAJbS/UoOVQlp\nHkarwxagoZDrQnLOWXWgsLCQCxculDgSsyx1OzQ0VFsdj5JDUEodK0XFRkRE2GTr1q2zySQHteRM\nkxz0kh47xen7Y02pK6XTleawtEZwgux0nD179kXreTEkR/Ftt91mkx09etQmc5pmWdLjwsJCR/Ur\nkzkxtdbLAPsVGAwBjtFtQ6BTmj7wx5RSae7PUPuYL4MhcDG6bQgISmrA/w40AdoC+4C3iyuolBqi\nlFqrlFpbwnMZDOVJiXRb6s4yGK40JTLgWusDWuvzWusLwEeAPXXaz2U/1Fq311q3L2klDYbyoqS6\nLflADIYrTYkiMZVSdbXWRSFLvwbsYVkOkRyWUuSTFMEnOQOk6DIJaQ5H6byS01FyHEmOUsmBITko\nrU5ByZlodXSCPF+nVE5C2rdq1ao2mRQRKTmUpWhK6RwSkpNVwtrKvRKpkEuq2yEhIbaoQ0lnpbku\nJWe8pHcSkt5JSPOTSvouPe8OHTrYZFJ64549e9pk1uu49957bWWmTZtmk61evdomk5yOTo8nIUVZ\nr1q1yiaT5smsUqWKTSZFNzulbt26PutSxKrEJd8wpdSnQA+gplIqG3gJ6KGUagtoYCfw8GXV1mDw\nA4xuGwKdSxpwrfXvBPHkK1AXg6FcMbptCHQCtuMuOjqa6dOns3TpUlatWsUf/vCHMj1+nz592Lx5\nMytWrODRRx8Vy/zxj39k3rx5fPnll8yYMYOwsDAAwsLCWL16NStWrGDVqlU8//zzjq5n9uzZrF69\nmrS0NB588MEyvR5wXVNGRgbbtm1j2LBhYpkePXqwfPlyvv/+e9vsKb/4xS9Yvnw5y5YtY+jQoYDr\nWr/88ktSU1NJT09n9OjRnvLe9+ftt98WPzsNJSc6Opo5c+awdu3aK6Izt956K6tXr2bKlCncfffd\ntu1Fel707J2+g+Hh4QwbNowxY8YwduzYMq+3t54/++yzYhlr3b31dty4caxZs4aFCxeWab2uBAE7\nK/2gQYPIz8+ne/fugNx/W1KCgoJ4//33ue222zh//jwLFixg0aJFbNu2zVOmdu3a3HffffTv358z\nZ84wZswYBg8ezNSpUzlz5gy33norISEhhISEsHDhQhYvXszatcUPxLnzzjvJy8ujU6dOnD9/vkyv\np+iaJk6cSO/evcnOzmbdunUsXLjQpw81KCiIsWPHMnjwYPbt28eCBQv4+uuv2bp1q2fbb3/7Ww4f\nPsy8efNYvHgx27ZtY/DgwWzZsoWQkBBWrFjBV199xZ49e3zuz4QJE+jfv7+YZc5QMgYNGkReXh7t\n27vGB5T1OzB+/HgGDRpEvXr1ePfdd0lJSWH37t2eMkV6XlBQQEhICGlpaaxZs0YMTPKmQ4cOnD59\nmlGjRgGlC8axopTy0fM1a9bw5Zdf2upkrXuR3h44cIDPP/+cTz75hLffLnYAkt9Qrga8SpUqtjn8\nvBWiCCfOqXXr1vHSSy+RnJzM/Pnzefnll9Fa8/333/Puu++yfPlyRowYQbVq1Xx+XYuQnJjVqlUD\nXAq2c+dODh8+TO3atVmwYAG//OUv+eCDDwCX4ygvLw+tNXv27OHEiRNERET4OGSLHHPh4eEEBQWR\nl5d30ZSRa9eu5cUXXyQlJcVzPQBJSUl8+OGHHgdR27ZtGTt2LP369RNTWhZHx44dyczMJCsrC4Dp\n06fTo0cPUlJSPGWuu+46tmzZ4nEgTZs2jX79+rFz505uuOEGdu3axf79+9m9ezf//ve/6dixI8nJ\nyYBrftKqVasSGRlJvXr1PJGpubm5nDhxAqUUGRkZYmSrFDkopba16oWU+rOiOHPmjM8PPEC/fvYo\nfMmJKTksBw8ebJNZnaJxcXH06tWLjRs3smzZMiZNmoTWmokTJ/LGG2/wzTffMGbMGKKjo3n88cdt\nx5Pm3SxKCdu8eXOOHj1KfHw8x48fZ9GiRdxwww2edMh33XWXz35VqlTh/Pnz7Nmzx2cgQevWrX3K\nLViwgN27dzNt2jSefvppFi1axPTp0z163qVLFwCGDh3KQw89xHPPPefZV4pEtTosb7rpJo4fP07z\n5s1p3rw569evZ8SIETzwwAO2fYve0dDQUEJDQ9FaU1hYyMqVK0lMTASKj5os61FHJdXlgGyBR0dH\nM378eFq3bk1BQQHJycmkpqbyxRdf8NZbbzF8+HBq1qxJq1atuP/++z3e9tmzZ3tCyL098GPHjuX7\n77/3rNerV88nVHb//v20adPGpw579+7lrbfeYvfu3Zw6dYpFixaxePFiz/agoCB+/PFHmjZtysSJ\nE8UJaZ1cz08//cQ111xDUFAQFy5cYNy4cT6fhcuWLfP88HgzfPhwvv32W896QkKCz0iR7OxsWwiz\nVKZ58+aAa/JcbyXz3j8oKIh58+bRsGFDpk2b5pkY9r333iM9PZ3Tp0+TnJxMcnKyOAmv4fKJiIjg\nvvvu49e//jWnTp1iypQpZGRk8O233/Lyyy8zevRoatWqRdu2bX3ybS9dutSjL97h+u+99x5r1qzx\nrNeoUcNnJMShQ4fE8PagoCDGjRtHfHw8H3/8MampqRetd1GD6pZbbuHkyZN88cUXDBw4kP/+978+\nej5kyBA+/PBDn33nzZtnG7FTWFjoo+sJCQk+9T58+LA4aqyo7tZ3tMhwBwoBacAffvhhvv76a0+e\njlWrVnkMw+rVq1FKMWTIEO666y4uXLjgaW17txqkFvjlUL16dQYOHEjjxo05fvw4//nPf7jnnnuY\nPn064Bry1q5dO2JiYpg7dy5JSUls2rTpsq9Ha82mTZto2bIlTZs2Zffu3T4vSbdu3Up1HWXBhQsX\nGDBgANWqVeMf//gHzZo1Y//+/fTv35/WrVuTm5vL1KlTufvuu1m2bFlFV7dS0Lt3b1JTUz2t9w0b\nNniGCi5fvhylFE888QQ9e/b0GX5Z1OUIcgv8crlw4QLPPPMMERERPPLIIzRr1uyiOVbuv/9+kpOT\nPXlN1q5dS3x8vEfPk5KSuPbaazl48KBtSPCAAQNsx3M63K64ulvf0ZJOrFBRBKQBb9euHVOnTvVZ\nnzdvHuD69KtTpw7Hjh2zjS922gLfu3evzwxA8fHxtsREvXr1Iisry6NAc+bM4eabb/YY8CJyc3P5\n7rvvuP3224s14Be7HoCUlBQ6d+7MkCFDGDhwoM++TlvgOTk51K9f37OemJhoS8gjlSkaV75//36f\nsarS/nl5eaxatYpu3bqRk5PDrl27PGNt582bR6dOnYwBLyMaN27MkiVLPOvNmzf3rLdq1Yr4+HiO\nHj1q655x2gI/cuSIz9jxWrVqcejQoWLrc/LkSVatWkX37t0vasCvv/56Zs6c6bP+2WefAS4979Kl\nC0OHDuX111+37eukBZ6Tk+NT75o1a15yVi3vd3TWrFkXLetvBKQBP3bsGO3atePrr7+mX79+REdH\ns3LlSuLj43nvvfd44IEHeOWVV+jRo4ePkjttga9bt44mTZrQsGFDCgsL6d+/P0899ZRPmd27d3PT\nTTcRHh7OqVOn6Nmzp8dJWbNmTc6dO0dubi5Vq1ald+/ejBs37rKvp4iUlBQ+/vhjPvjgA1vgk9MW\n+Jo1a7j22mtp1KgROTk5DB48mN///veXLPP0008DrhZe48aNqV+/Pvv27fPsX3St4PLsd+3alQ8+\n+ICjR4/Svn17z/3p3r0769evd1RXw6XJz8+ncePGgEsHoqKiWL9+PTVr1uS9997jzjvv5K9//St9\n+vTxCbpx2gLfunUr9erVo06dOpw8eZJevXp5/DJFREdHU1hYyMmTJ6lSpQpdu3blH//4x0XrnZub\nS+vWrUlOTqZ3795Uq1bNo+tFej5x4kTR6Dppga9Zs8ZT7yNHjtCtWzfGjx9v2+9y31F/pVwN+Llz\n5xx11hcNx/PG+5f3n//8J//85z8ZPHgwO3fu5K677iIiIoI5c+bwyiuvsH37dv72t78xcuRIli5d\n6njORe/UqUOHDmX27NlUrVqVzz//nB07dhAUFMSHH37I8OHD2bdvH4sXLyYtLY3z58+zYcMGZs+e\nTVxcHM2bN2fixImAq5/ts88+Y/78+cVe74wZM3j//fe577772LNnD4888gh16tTxbD9y5Ahnzpzh\ntddeu6Tj0rsFXcSePXs4f/48jz32GF9//TXBwcFMmTLF45CaP38+Dz30EGfOnOH5559n8eLFBAcH\nM2PGDJ8fkmHDhvGvf/2LkJAQpk2bRk5ODklJSfz973/3udbJk11Dqdu0acOSJUsoLCxk/fr1vPvu\nu2K0mtOUvdYvDadpOcsDrbUtotb7C6gI6Vql1LFWhyjAjz/+6PP/rFmz+PTTT8nK+v/t3UFoVFcU\nxvH/N2lRYpuNA0VqSbFLK9ggXdimuCq0YmgRCl10JbaLFhrQjQVRIi6EtoibQmsaCbQ1i1aDISlF\nKNaNCy0Sk7qMQkRqN9qou+ntIjOSzDtlXjNO3rvX89uYPJPJPeOZY+bed8+dY9euXdy9e5eTJ0+y\nb98+ZmdnGRoa4tixY0xOTrJ79+7M41nTBVNTU48+3rNnD8ePH6enp4fz58+zsLBAtVrl4MGDnDhx\ngu7ubgYHB6lUKlQqFYaHhxkZGVn2eH19fcs+P3LkCGNjYwwMDDA3N8fAwABbtmx59PeNu76snZjN\n7z4BxsfHl31eq9XYv38/Q0NDdHV1cfr0acbHxx8tOk5MTLB3716q1SojIyN0dXUte412d3dz6tQp\n+vv7Wb9+PZcuXWJsbCxzS6G1K9TKR2vRunnHLpBZY7N2hFpa9gN/nCqVSmjeXm1tN25VwBusQ26b\nt6SC3bs376G51gvO6vtrbQe3+ghbmg+1bS5yR48e5cKFC4yOjrZ8rP8q4HlYh9JasVr3c7czd7jS\nAv7gwQNqtdqK+4E/Tlav+56enszXWYu4VgG33iEuLeANVrsCq5DkLeBWv26rcFrtHiYmJjLXmvuf\nWwWtcafKgQMHmJmZ4dy5c0xPT+caR3MBB/s5tnqVWw3IrH+Lxp0xS61GAb937177/cBdcXp7e7l4\n8SJr167NVbydi9HGjRs5e/Ysa9asWbb241qLcg78SXHz5k36+/uLHoZzHTU/P7/sVkeXn/8G7pxz\nkVr1OfDmVqTttGB83C0drcezrlnz8dY1a07Mak25UtbWaWsjgtVi1/pea97VmifMex6iNX9u2bRp\nU67vtRbA85wbuBqsOfC8rLNDrVtDrdbLmzdvzlyzzh21djy3Y+fOnZlr1mvFmgNutn379sw16yzW\nmZlsZ19rzcuan7a+19pNac2BWzXSmse2zuy07qaxfsbDhw+tn+tz4M45lyov4M45Fykv4M45Fykv\n4M45F6lVXcSU9BeQPSiyM6rAyjvdlEcKcXQqht4QQnY1qQCe2/9bCjFAwbm9qgV8NUm6HELYVvQ4\n2pVCHCnEUCYpPJ8pxADFx+FTKM45Fykv4M45F6mUC/jXrb8kCinEkUIMZZLC85lCDFBwHMnOgTvn\nXOpS/g3cOeeS5gXcOecilWQBl3RD0jVJVyVdLno8eUj6VtIdSTNLrh2WdKsex1VJbxc5xlYkvSDp\nV0l/SJqV9Gn9elRxlFWMeQ2e2x0dV4pz4JJuANtCCNFsFJD0BnAfGA0hvFy/dhi4H0L4vMix5SVp\nA7AhhPC7pGeBK8A7wHtEFEdZxZjX4LndSX6gQ0mEEH6T9GLR42hHCOE2cLv+8YKk68DzxY7KFc1z\nu3OSnEIBAvCLpCuSPix6MG36RNJ0/W1o9iDPkqq/YF8BGqfTRhlHyaSU1xBpTpQpt1Mt4K+HEPqA\nt4CP62/hYvQV8BKwlcX//b8odjj5SHoG+BEYDCH8TaRxlFAqeQ2R5kTZcjvJAh5CuFX/8w5wBni1\n2BGtTAjhzxBCLYTwD/ANEcQh6WkWE/y7EMJPEGccZZRKXkOcOVHG3E6ugEtaV19kQNI64E0ge55S\nBOoLJw3vUvI4JAkYBq6HEL5ccj2qOMoopbyG+HKirLmd4iLmc8CZxeebp4DvQwg/Fzuk1iT9AOwA\nqpLmgUPADklbWZz7vAF8VNgA83kN+AC4Julq/dpnwPuRxVFGUeY1eG53UpK3ETrn3JMguSkU55x7\nUngBd865SHkBd865SHkBd865SHkBd865SHkBd865SHkBd865SP0LqVEL2bBU7zkAAAAASUVORK5C\nYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig, axes = plt.subplots(2, 2)\n", "\n", "frame = SimulatedFrame((20, 30))\n", "frame.add_spot((10, 15), 200, 2.5)\n", "\n", "for ax, noise_level in zip(axes.ravel(), [1, 20, 40, 90]):\n", " noisy_copy = frame.with_noise(noise_level)\n", " features = tp.locate(noisy_copy, 13, topn=1, engine='python')\n", " tp.annotate(features, noisy_copy, plot_style=dict(marker='x'), imshow_style=dict(vmin=0, vmax=255), ax=ax)\n", " dx, dy, ep = features[['x', 'y', 'ep']].iloc[0].values - [16, 10, 0]\n", " ax.set(xticks=[5, 15, 25], yticks=[5, 15])\n", " ax.set(title=r'Signal/Noise = {signal}/{noise}'.format(\n", " signal=200, noise=noise_level))\n", " ax.text(0.5, 0.1, r'$\\delta x={dx:.2}$ $\\delta y={dy:.2}$'.format(\n", " dx=abs(dx), dy=abs(dy)), ha='center', color='white', transform=ax.transAxes)\n", " ax.text(0.05, 0.85, r'$\\epsilon={ep:.2}$'.format(ep=abs(ep)), ha='left', color='white',\n", " transform=ax.transAxes)\n", "fig.subplots_adjust()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The effect of `diameter` on precision\n", "\n", "The size of the feature, as specified in the call to the `locate` function, has an important effect of the precision. If a user underestimates the particle's diameter, precision suffers. Thus, it is best to err on the high side, although this a larger diameter comes at some cost in performance.\n", "\n", "As demonstrated below, a too-small `diameter` often biases the location toward the pixel edges." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEOCAYAAABiodtuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xl4FOeV6P/vQQjEYhazCSMLgcHY\n7GCJxQazI8CS8A42OBjHYRyP75PrJPP8PEMmYTyjiZPJnUxyL57EYzuOJ8rEscf2tBYEmMXsRgKL\nHWzAIMQqhFiF0HZ+f3RDhNZuSd3Vap3P8/RDV/VbVackoaOq963ziqpijDHGVNbK6QCMMcYEH0sO\nxhhjqrHkYIwxphpLDsYYY6qx5GCMMaYaSw7GGGOqseRgjDGmGksOxhhjqrHkYIwxpprWTgfQUN27\nd9eYmBinwzDGmKBXUVHB5cuXKSoq4vTp0+dVtUd92zTb5BATE0N2drbTYRhjTNBSVfbu3UtmZibF\nxcVMnDiRKVOmHPdm22abHIwxxtTu8uXLpKen89VXX9GnTx+SkpLo2bOn19tbcjDGmBCiquzcuZPV\nq1dTXl7OzJkzGTt2LK1a+dbFbMnBGGNCxIULF0hNTeXYsWPExMSQmJjInXfe2aB9hVRyKC0tJS8v\nj+LiYqdDMR4RERFERUURHh7udCjGhKyKigq2bdvGunXrCAsLIzExkVGjRiEiDd5nSCWHvLw87rjj\nDmJiYhr1RTFNQ1UpKCggLy+Pfv36OR2OMSHp7NmzuFwuTp06xaBBg5gzZw6dOnVq9H5DKjkUFxdb\nYggiIkK3bt3Iz893OhRjQk5ZWRkbN25k06ZNRERE8MQTTzBkyJAm+/0XUskBsMQQZOz7YUzTy8vL\nw+VykZ+fz/Dhw4mPj6d9+/ZNeoyQSw7GGBOqSkpKWLduHdu2baNTp04888wz3HvvvX45liUHP1u2\nbBkdO3bk8uXLPPzww0yfPt1vx/rnf/5n/u7v/q7J91tYWMgLL7zAkSNHiIiI4N1332Xo0KFNfhxj\nTO2++eYbUlNTKSwsJDY2lunTp9O2bVu/Ha/l1laKjASR6q/ISL8c7vXXX/drYgB3cvBVeXm5V/sd\nOXIku3fv5v333+d73/teQ8IzxjRAcXExLpeL999/HxFh0aJFPPLII35NDNCSk8PZs76t90FycjL3\n3nsvEyZM4NChQwA8//zzfPTRR4A7UcTFxTF06FCWLFmCqgIwefJkXn31VWJjY7n//vvJysri8ccf\nZ+DAgfzoRz+6tf8//OEPjBkzhpEjR/JXf/VXlJeX89prr3H9+nVGjhzJggULam0H0LFjR37wgx8w\nYsQItm7dWu/57N+/n6lTpwJw3333cezYMc42wdfJGFO3gwcPsnz5cnJycnjwwQd56aWXCFRNuZab\nHPxkx44d/OlPfyInJ4eMjAyysrKqtXnllVfIyspi7969XL9+nbS0tFuftWnThuzsbF566SXmzp3L\n8uXL2bt3L++99x4FBQUcOHCADz74gM2bN5OTk0NYWBgpKSm88cYbtGvXjpycHFJSUmptB3Dt2jXG\njh3Lrl27mDBhAq+++iojR46s9nrjjTcAGDFiBB9//DEA27dv5/jx4+Tl5QXgq2lMy3Tt2jU++ugj\nPvjgAzp06MCLL77IjBkzAvq8kPU5NLGNGzfy2GOP3Ro5kJSUVK3NunXr+PnPf05RUREXLlxgyJAh\nJCYm3tZ+2LBhDBkyhN69ewPQv39/Tpw4waZNm9ixYwdxcXEAXL9+vcZ6KWvWrKm1XVhYGE888cSt\ntr/85S/rPKfXXnuN733ve4wcOZJhw4YxatQowsLCfPq6GGPqp6rs2bOHzMxMSkpKmDJlCg899JAj\n/9/8nhxE5F0gATinqtV6McU91vFXwBygCHheVXf6Oy6nFBcX8/LLL5Odnc3dd9/NsmXLbnui++Z9\nxFatWt12T7FVq1aUlZWhqixatIif/vSndR6nrnYRERG3/bC9+uqrrFu3rlq7+fPn89prr9GpUyd+\n97vf3dpvv3796N+/v28nboz5i8hIUs6eZSmQC0QDycDjnTvzyauvEhUVRVJSEj161FtZ228CcVvp\nPWBWHZ/PBgZ6XkuAfw9ATH7z8MMP8+mnn3L9+nWuXLlCamrqbZ/fTATdu3fn6tWrt/ohvDVt2jQ+\n+ugjzp07B7hrqRw/7q7AGx4eTmlpab3tqvrlL39JTk5Otddrr70GwMWLFykpKQHg7bff5uGHH26S\nJzCNaalSzp5lCXAcUM+/S4CPL10iPj6exYsXO5oYIABXDqq6QURi6mgyF3hf3b2y20Ski4j0VtXT\nfg2sV6+aO5979WrUbkePHs28efMYMWIEPXv2vHVb56YuXbrwne98h6FDhxIZGVnt8/oMHjyYf/qn\nf2LmzJlUVFQQHh7O8uXL6du3L0uWLGH48OGMHj2alJSUWtv56sCBAyxatAgRYciQIbzzzjs+78MY\n8xdLcd8mqazIs/7YuHGBD6gGcnOkjF8P4k4OabXcVkoD3lDVTZ7lNcD/p6rVZvIRkSW4EyzR0dEP\nVP1L+MCBA9x///1NHr9pHPu+GHO7ViLU9JtXgAo//04WkR2qGltfu2Y1WklV31LVWFWNdfqSyxhj\nGuLMmTNE1/JZbeudEAyjlU4Cd1dajvKsM8aYkFFWVsaGDRvYvHkzybhvgVS+tdQed6d0sAiG5OAC\nXhGRPwFjgUt+728wxpgAOnHiBC6Xi/PnzzNixAjm9uwJ585VG620oJF9nk0pEENZ/wuYDHQXkTzg\nJ0A4gKr+BsjAPYz1MO5EutjfMRljTCCUlJSwdu1avvjiCzp37syCBQsYMGAAPPooC4AFTgdYh0CM\nVnqmns8V+Gt/x2GMMYF05MgR0tLSuHjxInFxcUybNs3v9ZCaUjDcVjLGmJBx/fp1Vq1aRU5ODt26\ndWPx4sVERwdTV7N3mtVopeZo2bJl/OIXv+DHP/4xn332mV+P1ZCqrN44ePAg48ePp23btvziF7+4\ntb64uJgxY8YwYsQIhgwZwk9+8hO/HN8Yv/BDZeYDBw7w5ptv3qpb9tJLLzXLxAAtPDmkpKQQExND\nq1atiImJuVWYzh+ac8nuO++8k1//+tf88Ic/vG1927ZtWbt2Lbt27SInJ4fMzEy2bdvmcwzGOKIJ\nKzNfvXqVDz/8kD//+c907NiR73znO0ybNo3WrZvvzZkWmxxSUlJYsmQJx48fR1U5fvw4S5YsaZIE\nEWolu28+6V21IqSI0LFjRwBKS0spLS21aUFNi6Kq7Nq1i+XLl3Po0CGmTp3Kiy++eKtgZnPWYpPD\n0qVLKSq6/QH2oqIili5d2qj9hmLJ7rqUl5czcuRIevbsyYwZMxg7dmyjvn7GNBcXL14kJSWFTz/9\nlB49evDSSy8xceLEkKlY3HyveRopNzfXp/XeCsWS3XUJCwsjJyeHixcv8thjj7F3716bQtSENFUl\nKyvrVh/i7NmziYuLC7mr5habHKKjo2usUurvzqPmWLLbG126dGHKlClkZmZacjAh6/z587hcLk6c\nOME999xDQkICXbp0cTosv2ixt5WSk5Nv/XV/U/v27UlObtwD7KFYsrs2+fn5XLx4EXBfmaxevZr7\n7rvPp/MxxjG1PY1cw/ry8nI2btzIb37zG/Lz85k7dy4LFiwI2cQALfjK4Wan7dKlS8nNzSU6Oprk\n5ORb6xsqFEt2nzlzhtjYWC5fvkyrVq34t3/7N/bv38/p06dZtGgR5eXlVFRU8PTTT5OQkODz/o1x\nxJkzXjU7ffo0LpeLM2fOMHjwYGbPnn1rIEYoC0jJbn+IjY3V7Ozbq3pbaejgZN8XE6xSUlJq/QOx\nrKyMzz//nM2bN9O+fXseeeSRkPg59rZkd4u9cjDGtGw3h7PfHLV4czg7wMSJE3G5XBQUFDBy5Ehm\nzpxJu3btnAw34Cw5GGOcExlZ+4yMXt72aajahrO/+uqr/PVf/zVdunRh4cKF3HPPPX6NI1iFXHJQ\n1ZAbUtacNdfbliZAmvApZV/VNmw9Pz+fMWPGMG3aNNq0aeP3OIJVSI1WioiIoKCgwH4hBQlVpaCg\ngIiICKdDMaaa2oat9+nTh9mzZ7foxAAhduUQFRVFXl4e+fn5TodiPCIiIoiKinI6DGOqSU5O5sUX\nX7ztOaP27dvzs5/9zMGogkdIJYfw8HD69evndBjGmCB35coVwsPDmTNnDuvXr6ewsLDJhrOHipBK\nDsYYUxdVJScnh1WrVlFWVsb3v/99PvzwQ1q1Cqk77E3CkoMxxjm9etU+WqmJFRYWkpaWxtGjR4mO\njiYpKYlu3bo1+XFChSUHY4xz/DxcFaCiooKsrCzWrFmDiDBnzhxiY2NtVGM9LDkYY0JWfn4+LpeL\nvLw8BgwYQEJCAp07d3Y6rGbBkoMxJuSUl5ezefNmNmzYQJs2bXjssccYNmyYXS34wJKDMSaknDp1\nCpfLxdmzZxkyZAizZ8+mQ4cOTofV7FhyMKalc7CERVMqLS1l/fr1bN26lQ4dOjBv3jwrId8IlhyM\naekcLGHRVI4fP47L5eLChQuMGjWKmTNn2pP5jWTJwRjTbN24cYPPPvuM7OxsunTpwnPPPUf//v2d\nDisk2JMfxhhHpaSkEBMTQ6tWrYiJiSElJcWr7b7++mvefPNNsrOzGTduHN/97nctMTQhu3Iwxjim\nrjkVaitjUVRUxMqVK9m9ezc9evTg29/+ttXv8oOAzAQnIrOAXwFhwNuq+kaVz6OB3wNdPG1eU9WM\nuvZZ00xwxpgGqGt4p59/P8TExNQ4t3nfvn05duxYlVCU/fv3k5GRQXFxMRMmTGDixIm0bm1/4/oi\naGaCE5EwYDkwA8gDskTEpar7KzX7EfBnVf13ERkMZAAx/o7NGENAS1hUVducClXXX7lyhfT0dA4d\nOsRdd91FUlISvQIQX0sWiJQ7BjisqkcBRORPwFygcnJQoJPnfWfgVADiMsaAo8NVo6Oja7xyuDnX\ngqry5ZdfsmrVKsrLy5kxYwbjxo2zQnkBEIjk0Ac4UWk5Dxhbpc0yYJWI/C+gAzA9AHEZYxyWnJx8\nW58DuOdUSE5OprCwkNTUVL755hv69u1LUlISd955p4PRtizBkn6fAd5T1ShgDvCfIlItNhFZIiLZ\nIpJtE/qYkBIZ6b73X/UVGel0ZH61YMEC3nrrLfr27YuI0LdvX37729/Sv39/3nzzTU6ePElCQgKL\nFi2yxBBgfu+QFpHxwDJVjfcs/y2Aqv60Upt9wCxVPeFZPgqMU9Vzte3XOqRNSHGwUziYnDt3DpfL\nxcmTJxk4cCAJCQl06tSp/g2N14KmQxrIAgaKSD/gJDAfeLZKm1xgGvCeiNwPRAB2aWBMC1FeXs6m\nTZvYsGEDERERPP744wwdOtQK5TnI78lBVctE5BVgJe5hqu+q6j4ReR3IVlUX8APgP0TkVdyd089r\nIMbYGmMcd/LkSVwuF+fOnWPYsGHEx8dbobwgEJABwp5nFjKqrPtxpff7gYcCEYsxJjiUlpaybt06\ntm3bRseOHZk/fz6DBg1yOizjYU+PGGMC7tixY7hcLgoLC3nggQeYPn26FcoLMpYcjAkGDj6IFkjF\nxcWsXr2anTt30rVrV771rW/Rr18/p8MyNbDkYEwwaEbzJjTUV199RVpaGlevXmX8+PFMmTKF8PBw\np8MytbDkYIzxq2vXrpGZmcnevXvp2bMn8+bNo0+fPk6HZerRqOQgIotV9XdNFYwxJnSoKnv37mXF\nihXcuHGDyZMnM2HCBMLCwpwOzXihsVcO/wBYcjDG3Oby5cukp6fz1Vdf0adPH5KSkujZs6fTYRkf\n1Fs+Q0R21/LaA4RWb5kxLVRDJ9ypSlXJzs5m+fLlHD16lJkzZ/LCCy9YYmiGvLly6AXEA4VV1guw\npckjMsYEVEMm3KnJhQsXSE1N5dixY/Tr14/ExES6du3ql5iN/9VbW0lE3gF+p6qbavjsj6patRRG\nQFhtJWOahi8T7tSkoqKCbdu2sW7dOsLCwpg5cyajRo2y0hdBqslqK6nqt+v4zJHEYIxpOt5OuFOT\ns2fP4nK5OHXqFIMGDWLOnDlWKC9E2FBWY1q4+ibcqUlZWRkbN25k06ZNRERE8OSTTzJ48GC7Wggh\nlhyMaeHqmnCnJnl5ebhcLvLz8xk+fDjx8fG0b98+UOGaALHkYEwLd7PTeenSpeTm5hIdHU1ycnK1\nzuiSkpJbhfI6derEs88+y8CBA50I2QSANx3S36/rc1X91yaNyEvWIW1M4Bw9epTU1FQuXrxIbGws\n06dPp23btk6HZRqgKSf7ucPz7yAgDnB5lhOB7Q0LzxjTHBQXF7Nq1Sq+/PJL7rzzTp5//nn69u3r\ndFgmALwZrfQPACKyARitqlc8y8uAdL9GZ4xxzMGDB0lPT+fatWs8+OCDTJ482QrltSC+9Dn0Akoq\nLZdgT0gbE3KuXr1KZmYm+/bto1evXjzzzDPcddddTodlAsyX5PA+sF1EPsH9dPRc4D1/BGWMCTxV\nZc+ePWRmZlJSUsKUKVN46KGHrFBeC+V1clDVZBFZAUzEPc/zYlX90m+RGWMC5tKlS6SlpXH48GGi\noqJISkqiR48eTodlHOTrUNZyoAJ3cqho+nCMMYF0s1DeZ599hqoya9Ys4uLiaNWq3pqcJsR5nRxE\n5HvAd4D/xn1b6Q8i8paq/l9/BWeM8Z+CggJcLhe5ubn079+fhIQEK5RnbvHlyuHbwFhVvQYgIj8D\ntgKWHIxpRioqKti6dSvr16+ndevWJCUlMXLkSCt9YW7jy7Wj4L6tdFO5Z50xpgk01ZwKdTlz5gxv\nv/02n332GQMGDODll1+2CqqmRr5cOfwO+MIzWgngUeCdpg/JmJanqeZUqE1ZWRkbNmxg8+bNtGvX\njqeeeorBgwc3er8mdNVbPuO2xiIPAA95Fjc6OVrJymeYUNLYORXqcuLECVwuF+fPn2fEiBHEx8fT\nrl27Ru3TNF9NWT7jFlXdAexocFTGmBo1Zk6F2pSUlLBmzRq2b99O586dWbBgAQMGDGjw/kzL4sto\npVhgKdDXs50AqqrD/RSbMS1GQ+ZUqMuRI0dITU3l0qVLxMXFMW3aNCuUZ3ziS4d0Cu5+hydwF91L\n8PxbLxGZJSKHROSwiLxWS5unRWS/iOwTkT/6EJcxzV5ycnK1ORHqmlOhNtevX+d//ud/+MMf/kDr\n1q1ZvHgxc+bMscRgfObLbaV8VXXV3+x2IhIGLAdmAHlAloi4VHV/pTYDgb8FHlLVQhHp6etxjGnO\nvJ1ToS4HDhwgIyODa9euMWHCBCZNmkTr1jZli2kYrzukRWQa8AywBrhxc72qflzPduOBZaoa71n+\nW892P63U5ufAV6r6treBW4e0MW5Xr14lIyODAwcOEBkZSVJSEr1793Y6LBOk/NEhvRi4DwjnL6Uz\nFKgzOQB9gBOVlvOAsVXa3AsgIpuBMNzJJNOH2IxpcVSVXbt2sXLlSkpLS5k6dSoPPvigFcozTcKX\n5BCnqoP8GMdAYDIQBWwQkWGqerFyIxFZAiyBhnfUGRMKLl68SFpaGkeOHOHuu+8mKSmJ7t27Ox2W\nCSG+JIctIjK4cl+Bl04Cd1dajvKsqywP+EJVS4FvROQr3Mkiq3IjVX0LeAvct5V8jMOYZk9VycrK\n4rPPPgNg9uzZxMXF2RPOpsn5khzGATki8g3uPgdvh7JmAQNFpB/upDAfeLZKm09x92f8TkS6477N\ndNSH2IwJeefPn8flcnHixAnuueceEhIS6NKli9NhmRDlS3KY1ZADqGqZiLwCrMTdn/Cuqu4TkdeB\nbM8IqJXATBHZj7tm09+oakFDjmdMqCkvL2fLli18/vnnhIeH8+ijjzJ8+HC7WjB+5VP5jGBio5VM\nS3D69GlcLhdnzpxh8ODBzJ49m44dOzodlmnG/FI+wxgTGGVlZaxfv54tW7bQoUMHnn76ae6//36n\nwzItiCUHY4JMbm4uLpeLgoICRo4cycyZM61Qngk4X2or/Qr439pc70MZE+Ru3LjBmjVryMrKokuX\nLixcuJB77rnH6bBMC+XLlcMVwCUi81X1mojEAz9W1Yfq29AYU7fDhw+TlpbGpUuXGDt2LFOnTqVN\nmzZOh2VaMK+Tg6r+SESeBdaLSAlwFaixiJ4xxjtFRUWsWrWKXbt20b17d1544QXuvvvu+jc0xs98\nua00DfgOcA3oDbygqof8FZgxoUxVbxXKu379OhMnTuThhx+2QnkmaPjyk7gU+HtV3SQiw4APROT7\nqrrWT7EZE5KuXLlCRkYGBw8epHfv3ixcuJDIyEinwzLmNr7cVppa6f0eEZkN/DfwoD8CMybUqCo5\nOTmsWrWKsrIypk+fzvjx42nVypdpVYwJjAZfw6rqac+tJmNMPQoLC0lLS+Po0aNER0eTlJREt27d\nnA7LmFo16ganql5vqkCMCUUVFRVs376dtWvXIiLMmTOH2NhYK31hgp71fhnjJ/n5+bhcLvLy8hgw\nYAAJCQl07tzZ6bCM8YolB2OaWHl5OZs3b2bDhg20adOGxx57jGHDhtnVgmlW6k0OIvL9uj5X1X9t\nunCMad5OnTqFy+Xi7NmzDBkyhNmzZ9OhQwenwzLGZ95cOdzh+XcQEAe4PMuJwHZ/BGVMc1NaWsr6\n9evZunUrHTp0YN68edx3331Oh2VMg9WbHFT1HwBEZAMwWlWveJaXAel+jc6YZuDYsWOkpqZy4cIF\nRo0axcyZM4mIiHA6LGMaxZc+h15ASaXlEs86Y1qkGzdusHr1anbs2EHXrl157rnn6N+/v9NhGdMk\nfEkO7wPbReQTz/KjwO+bPiRjgt/XX39NWloaV65cYdy4cUyZMsUK5ZmQ4ssT0skisgKY6Fm1WFW/\n9E9YxgSnoqIiMjMz2bNnDz169OCpp54iKirK6bCMaXK+FN4TYDDQWVVfF5FoERmjqtYpbUKeqrJv\n3z5WrFhBcXExkyZNYsKECVYoz4QsX36y3wQqgKnA67jnd/hv3COYjAlZly9fJiMjg0OHDnHXXXeR\nlJREr17W3WZCmy/JYayqjhaRLwFUtVBE7CarCVmqys6dO1m9ejXl5eXMmDGDcePGWaE80yL4khxK\nRSQMUAAR6YH7SsKYkHPhwgVSU1M5duwYMTExJCYmcueddzodljEB40ty+DXwCdBTRJKBJ4Ef+SUq\nYxxSUVHBF198wdq1awkLCyMhIYHRo0db6QvT4niVHDyd0RuAHcA0QIBHVfWAH2MzJqDOnTuHy+Xi\n5MmT3HvvvTzyyCN06tTJ6bCMcYRXyUFVVUQyVHUYcNDPMRkTUOXl5WzcuJGNGzcSERHB448/ztCh\nQ+1qwbRovtxW2ikicaqa5bdojAmwkydP4nK5OHfuHMOGDSM+Pt4K5RmDj6OVgIUicgy4hvvWkqrq\ncH8EZow/lZaWsm7dOrZt20bHjh2ZP38+gwYNcjosY4KGL8kh3m9RGBNA33zzDampqRQWFvLAAw8w\nffp0K5RnTBW+JIdFtax/vb4NRWQW8CsgDHhbVd+opd0TwEdAnKpm+xCbMfUqLi5m9erV7Ny5k65d\nu7Jo0SJiYmKcDsuYoORLcrhW6X0EkADUO1rJ82zEcmAGkAdkiYhLVfdXaXcH8D3gCx9iMsYrhw4d\nIj09natXrzJ+/HimTJlCeHi402EZE7R8Kbz3fyovi8gvgJVebDoGOKyqRz3b/QmYC+yv0u4fgZ8B\nf+NtTMbU59q1a2RmZrJ371569uzJvHnz6NOnj9NhGRP0GlM1rD3gTTnKPsCJSst5uDu3bxGR0cDd\nqpouIrUmBxFZAiwBiI6O9jlg03KoKnv37mXFihXcuHGDyZMnM2HCBMLCwpwOzZhmwZeqrHvwlM7A\n3XfQAy/6G7zYbyvgX4Hn62urqm8BbwHExsZqPc1NC3Xp0iXS09P5+uuv6dOnD0lJSfTs2dPpsIxp\nVny5ckio9L4MOKuqZV5sdxK4u9JylGfdTXcAQ4H1noeOIgGXiCRZp7TxhaqyY8cOVq9ejaoSHx/P\nmDFjrFCeMQ3gy/+aMcAFVT0OLAb+7LkdVJ8sYKCI9PNUcZ0PuG5+qKqXVLW7qsaoagywDbDEYHxS\nUFDA73//e9LT0+nTpw/f/e53G1RBNSUlhZiYGFq1akVMTAwpKSl+itiY4ObLlcPfq+qHIjIBmA78\nC/DvVOk/qEpVy0TkFdyd12HAu6q6T0ReB7JV1VXX9sbUpaKigm3btrFu3TrCwsJITExk1KhRDSp9\nkZKSwpIlSygqKgLg+PHjLFmyBIAFCxY0adzGBDtR9e7WvYh8qaqjROSnwB5V/ePNdf4NsWaxsbGa\nnW0XFy3Z2bNncblcnDp1ikGDBvHII49wxx13NHh/MTExHD9+vNr6vn37cuzYsUZEakzwEJEdqhpb\nXztfrhxOishvcT+v8DMRaYtvt6WMaRJlZWVs3LiRTZs2ERERwZNPPsngwYMbXSgvNzfXp/XGhDJf\nksPTwCzgF6p6UUR6Y88kmADLy8vD5XKRn5/P8OHDiY+Pp3379k2y7+jo6BqvHGzYtGmJfHkIrgj4\nuNLyaeC0P4IypqqSkhLWrl3LF198QadOnXj22WcZOHBgkx4jOTn5tj4HgPbt25OcnNykxzGmOfDp\nITgR6QoMxF0+AwBV3dDUQRlT2dGjR0lNTeXixYvExsYyffp02rZt2+THudnpvHTpUnJzc4mOjiY5\nOdk6o02L5EuH9Iu4ax9FATnAOGCrqk71X3i1sw7p0FdcXMyqVav48ssvufPOO0lKSqJv375Oh2VM\ns+aPDunvAXHANlWdIiL3Af/c0ACNqcvBgwdJT0/n2rVrPPTQQ0yaNMkK5RkTQL4kh2JVLRYRRKSt\nqh4UEZsdxTSpq1evsmLFCvbv30+vXr145plnuOuuu5wOy5gWx5fkkCciXYBPgdUiUghUH9phTAOo\nKrt372blypWUlJQwZcoUHnroISuUZ4xDfBmt9Jjn7TIRWQd0BjL9EpVpUS5dukRaWhqHDx8mKiqK\npKQkevTo4XRYxrRovlRlFWAB0F9VXxeRaGAksN1fwZnQpqpkZ2fz2WefoarMmjWLuLg4K5RnTBDw\n5bbSm0AFMBV3qe4rwH/j7qTTqc94AAAS2klEQVQ2xicFBQW4XC5yc3Pp378/iYmJdOnSxemwjDEe\nviSHsao6WkS+BFDVQk+VVWO8VlFRwZYtW1i/fj3h4eHMnTuXESNGNLr0hTGmafmSHEo980ErgIj0\nwH0lYYxXzpw5g8vl4vTp09x3333MmTOnUYXyjDH+48vN3V8DnwC9RCQZ2Az81C9RmZBSVlbGmjVr\neOutt7h8+TJPPfUU8+bNqzUx2JwKxjjPl9FKKSKyA5jmWZWkqgf9E5YJFSdOnMDlcnH+/HlGjBhB\nfHw87dq1q7W9zalgTHDwpXxGLLAUiMGdVBRAVYf7K7i6WPmM4FZSUsKaNWvYvn07nTt3JiEhgQED\nBtS7nc2pYIx/+aN8RgruEt17sL4GU4cjR46QmprKpUuXiIuLY9q0aV4XyrM5FYwJDr4kh3yb0tPU\n5fr166xatYqcnBy6devG4sWLfZ4LweZUMCY4+JIcfiIibwNrgBs3V6rqx7VvYlqK/fv3k5GRQVFR\nERMmTGDSpEm0bu1TRXjA5lQwJlj48r93MXAfEM5fbisplSYAMi3P1atXycjI4MCBA0RGRrJgwQJ6\n9+7d4P3ZnArGBAdfOqQPqWrQVGG1DmlnqSq7du1i5cqVlJaWMnnyZMaPH2+F8owJcv7okN4iIoNV\ndX8j4jIh4OLFi6SlpXHkyBGio6NJTEyke/fuTodljGlCviSHcUCOiHyDu89BAHVqKKsJPFVl+/bt\nrFmzBhFh9uzZxMXFWekLY0KQL8lhlt+iMEHv/PnzuFwuTpw4wT333ENCQoIVyjMmhPnyhLRN7NMC\nlZeXs2XLFj7//HPCw8N59NFHGT58eGheLURGwtmz1df36gVnzgQ+HmMc5PtYQ9NinD59GpfLxZkz\nZxg8eDCzZ8+mY8eOToflPzUlhrrWGxPCLDmYakpLS/n888/ZsmULHTp04Omnn+b+++93OixjTAAF\nJDmIyCzgV0AY8LaqvlHl8+8DLwJlQD7wgt3GckZubi4ul4uCggJGjhzJzJkz6yyUZ4wJTX5PDp45\nIJYDM4A8IEtEXFWGxH4JxKpqkYh8F/g5MM/fsZm/uHHjBmvWrCErK4suXbrw3HPP0b9/f6fDMsY4\nJBBXDmOAw6p6FEBE/gTMBW4lB1VdV6n9NmBhAOIyHl9//TVpaWlcvnyZsWPHMnXqVNq0sUn+jGnJ\nAjGTex/gRKXlPM+62nwbWFHTByKyRESyRSQ7Pz+/CUNsmYqKivjkk0/44x//SJs2bXjhhReYNWsW\nbdq0aZkT7vTq5dt6Y0JYUHVIi8hCIBaYVNPnqvoW8Ba4y2cEMLSQoqrs37+fFStWcP36dR5++GEm\nTpx4q1Bei51wx4arGnNLIJLDSeDuSstRnnW3EZHpuCcTmqSqN6p+bprGlStXyMjI4ODBg/Tu3ZuF\nCxcSGRl5W5ulS5feVhUV3FcZS5cuDe3kYIy5JRDJIQsYKCL9cCeF+cCzlRuIyCjgt8AsVT0XgJha\nHFUlJyeHlStXUl5ezvTp0xk/fjytWlW/s+johDv2IJoxQcHvyUFVy0TkFWAl7qGs76rqPhF5Hcj2\nTCD0L0BH4EPPk7e5qprk79haisLCQtLS0jh69Ch9+/YlMTGRbt261dre0Ql37EE0Y4JCQPocVDUD\nyKiy7seV3k8PRBwtTUVFBdu3b2ft2rWICI888ggPPPBAvaUvbMIdY0xQdUibppOfn4/L5SIvL48B\nAwaQkJBA586dvdrWJtwxxng92U+wscl+alZeXs6mTZvYuHEjbdq0YdasWQwbNqz5FMqrK85m+rNq\nTDDxx2Q/JsidOnUKl8vF2bNnGTp0KLNmzaJDhw5Oh2WMaYYsOYSA0tJS1q9fz9atW+nYsSPz589n\n0KCgmdHVN7161T5ayRgTMJYcmrljx46RmprKhQsXGD16NDNmzCAiIsLpsBrOhqsaExQsOTRTN27c\nYPXq1ezYsYOuXbvyrW99i379+jkdljEmRFhyaIa++uor0tPTuXLlCuPGjWPKlClWKM8Y06QsOTQj\nRUVFZGZmsmfPHnr06MFTTz1FVFSU02EZY0KQJYdmQFXZt28fK1asoLi4mEmTJjFx4kTCwsL8c0Ar\nYWFMi2fJIchdvnyZjIwMDh06xF133UVSUhK9/D1yx0pYGNPiWXIIUqrKzp07Wb16NeXl5cyYMYNx\n48bVWCjPGGOamiWHIHThwgVSU1NxuVx8/vnnXLhwwUpYGGMCypJDEKmoqOCLL75g7dq17N69m4yM\nDIqLi4EWNOGOMSYo2D2KIHHu3DneffddVq1aRf/+/dm6deutxHDTzQl3jDHG3+zKwWHl5eVs3LiR\njRs3EhERwRNPPMGQIUNqvToIyIQ7VsLCmBbPkoODTp48icvl4ty5cwwbNoxZs2bRvn17wOEJd2y4\nqjEtnt1WckBpaSkrV67knXfe4fr16zzzzDM8/vjjtxIDuCfcqbwMNuGOMSZw7MohwL755htSU1Mp\nLCzkgQceYPr06TUWyrMJd4wxTrLJfgKkuLiY1atXs3PnTrp27UpSUhIxMTG1b2BPKRtj/MAm+wki\nhw4dIj09natXr/Lggw8yefJkwsPD697InlI2xjjIkoMfXbt2jczMTPbu3UvPnj2ZP38+d911l9Nh\nGWNMvSw5+IGqsnfvXlasWMGNGzeYPHkyEyZM8F+hPGOMaWKWHJrYpUuXSE9P5+uvv6ZPnz4kJSXR\ns2dPp8MyxhifWHJoIqrKjh07WL16NapKfHw8Y8aMsUJ5xphmyZJDEygoKCA1NZXjx4/Tr18/EhMT\n6dq1a+N2ak8pG2McZMmhESoqKti6dSvr168nLCyMxMRERo0ahYg0fuc2XNUY4yBLDg109uxZXC4X\np06dYtCgQTzyyCPccccdTodljDFNIiDJQURmAb8CwoC3VfWNKp+3Bd4HHgAKgHmqeiwQsfmqrKyM\njRs3smnTJtq1a8eTTz7J4MGDm+ZqwRhjgoTfe0tFJAxYDswGBgPPiMjgKs2+DRSq6gDgl8DP/B2X\nVyIjSREhRoRWnn8/CA8nLjGRoUOH8vLLLzNkyBBLDMaYkBOIK4cxwGFVPQogIn8C5gL7K7WZCyzz\nvP8I+H8iIupwbY+Us2dZAhR5lo8DSwCuXWPBY485FpcxxvhbIMZZ9gFOVFrO86yrsY2qlgGXgG4B\niK1OS/lLYripyLPeGGNCWbMahC8iS0QkW0Sy8/Pz/Xqs4uJiaptWJwDT7RhjjKMCkRxOAndXWo7y\nrKuxjYi0Bjrj7pi+jaq+paqxqhrbo0cPP4ULBw8eZPny5dQ2rU4AptsxxhhHBSI5ZAEDRaSfiLQB\n5gOuKm1cwCLP+yeBtU70N1y9epUPP/yQDz74gA4dOpAMtK/Spj1g0+0YY0Kd3zukVbVMRF4BVuIe\nyvququ4TkdeBbFV1Ae8A/ykih4ELuBNIwKgqu3fvJjMzk9LSUqZOncqDDz5I2LJlcPYsS3HfSorG\nnRgW2FPKxpgQ1+In+7l06RJpaWkcPnyYqKgokpKS8OctK2OMcZJN9lMPVSUrK4s1a9agqsyaNYu4\nuDgrlGeMMbTQ5HD+/HlSU1PJzc2lf//+JCYm0qVLF6fDMsaYoNGikkNFRQVbtmxh/fr1hIeHM3fu\nXEaMGGFPOBtjTBUtJjmcOXMGl8vF6dOnuf/++5kzZw4dO3Z0OixjjAlKIZ8cysrK+Pzzz9m8eTPt\n27fnqaeeYvDgqqWdjDHGVBbSySE3NxeXy0VBQQEjRowgPj6edu3aOR2WMcYEvZBMDiUlJaxZs4bt\n27fTuXNnFixYwIABA5wOyxhjmo2QSw5HjhwhNTWVS5cuMWbMGKZOnUrbtm2dDssYY5qVkEkO169f\nZ9WqVeTk5NCtWzcWL15MdLRVQTLGmIYIieSwf/9+MjIyKCoqYsKECUyaNInWrUPi1IwxxhHN+jfo\nlStXWLFiBQcOHCAyMpKFCxcSGRnpdFjGGNPsNdvkUFRUxJtvvklpaSnTpk1j/PjxhIWFOR2WMcaE\nhGabHC5evEjPnj1JTEyke/fuTodjjDEhpdlWZRWRfNzTOgdKd+B8AI8XaHZ+zVconxvY+TW1vqpa\nb+npZpscAk1Esr0pc9tc2fk1X6F8bmDn5xSrT22MMaYaSw7GGGOqseTgvbecDsDP7Pyar1A+N7Dz\nc4T1ORhjjKnGrhyMMcZUY8mhChGZJSKHROSwiLxWw+dtReQDz+dfiEhM4KNsGC/O7fsisl9EdovI\nGhHp60ScDVXf+VVq94SIqIgE3QiRunhzfiLytOd7uE9E/hjoGBvDi5/PaBFZJyJfen5G5zgRZ0OI\nyLsick5E9tbyuYjIrz3nvltERgc6xmpU1V6eFxAGHAH6A22AXcDgKm1eBn7jeT8f+MDpuJvw3KYA\n7T3vv9tczs3b8/O0uwPYAGwDYp2Ou4m/fwOBL4GunuWeTsfdxOf3FvBdz/vBwDGn4/bh/B4GRgN7\na/l8DrACEGAc8IXTMduVw+3GAIdV9aiqlgB/AuZWaTMX+L3n/UfANGkek1DXe26quk5VizyL24Co\nAMfYGN587wD+EfgZUBzI4JqAN+f3HWC5qhYCqOq5AMfYGN6cnwKdPO87A6cCGF+jqOoG4EIdTeYC\n76vbNqCLiPQOTHQ1s+Rwuz7AiUrLeZ51NbZR1TLgEtAtINE1jjfnVtm3cf8l01zUe36eS/W7VTU9\nkIE1EW++f/cC94rIZhHZJiKzAhZd43lzfsuAhSKSB2QA/yswoQWEr/8//a7Z1lYy/iMiC4FYYJLT\nsTQVEWkF/CvwvMOh+FNr3LeWJuO+6tsgIsNU9aKjUTWdZ4D3VPX/iMh44D9FZKiqVjgdWCiyK4fb\nnQTurrQc5VlXYxsRaY378rYgINE1jjfnhohMB5YCSap6I0CxNYX6zu8OYCiwXkSO4b6v62pGndLe\nfP/yAJeqlqrqN8BXuJNFc+DN+X0b+DOAqm4FInDXJQoFXv3/DCRLDrfLAgaKSD8RaYO7w9lVpY0L\nWOR5/ySwVj09SkGu3nMTkVHAb3EnhuZ0vxrqOT9VvaSq3VU1RlVjcPepJKlqtjPh+sybn81PcV81\nICLdcd9mOhrIIBvBm/PLBaYBiMj9uJNDfkCj9B8X8C3PqKVxwCVVPe1kQHZbqRJVLRORV4CVuEdP\nvKuq+0TkdSBbVV3AO7gvZw/j7mCa71zE3vPy3P4F6Ah86Oljz1XVJMeC9oGX59dseXl+K4GZIrIf\nKAf+RlWbw1Wtt+f3A+A/RORV3J3TzzeTP8wQkf/Cnbi7e/pMfgKEA6jqb3D3ocwBDgNFwGJnIv0L\ne0LaGGNMNXZbyRhjTDWWHIwxxlRjycEYY0w1lhyMMcZUY8nBGGNMNZYcjDHGVGPJwRhjTDWWHEyL\nJyLLROSHVdZ1EZGXG7i/27YVkZja6vg7qabzrvRZnfMPmNBnycEEPU9JgUD/rHbBPXdHQ+Kpddtm\n5D2gOVV1NU3MkoNxnIj8vWcGsE0i8l8i8kPPX9uHROR9YC9wt2emur2e1//2bHvbX+WebZd51h8Q\nkf/wzIq2SkTaVWq3VES+EpFNwKAawnoDuEdEckTkX2qIZ2JNx61pW8+6sNpiqbSPdSIyw/P+n0Tk\n/1b5PEZEDorIe57YU0RkuqdE99ciMqZS22pfKy/PG/Bq/gET6pyebcheLfsFxAE5uIuo3QF8DfwQ\niAEqgHGedg8Ae4AOuOs/7QNGedrtrbS/H+Ku+x8DlAEjPev/DCyssq/2uCePOQz8sEpcVfdbNZ4a\nj1vHtjXGUuWYDwPrgQVAOhBWQ0xlwDDcf9jtAN7FPXvYXODTer5W9Z53XV8De7WslxXeM057CPgf\nVS0GikUktdJnx9U9KxbABOATVb0GICIfAxOpXrmzsm9UNcfzfgfuX3Z4tvtEPbPeiYi3Rfkqx+Or\n2mK5RVU3iLvi4feByapaXst+9gCIyD5gjaqqiOyptM/avlataNh5mxbIbiuZYHbNizZl3P5zHFHp\nfeX5KMppfBXiyvHUddya1BuLiAwDegMlqnrFi/1UVFquqGmfxjSUJQfjtM1AoohEiEhHIKGWdhuB\nR0WkvYh0AB7zrDsL9BSRbiLSto7tK9vg2Vc7EbkDSKyhzRXct7lqU9dx69u2GnHPF5yC+/bQVWnc\nFJ+1fa28OW9jAPtLwzhMVbM8tzd24/6Fuwf3vNxV2+0UkfeA7Z5Vb6vqlwCemv/bcc+cddCLY+4U\nkQ+AXcA53BPNVG1T4Ono3Yt7Lu3lVT4vre249W1blYi0Bz4GfqCqB0TkH4GfAZn1nUsd5/ceNX+t\n6jzvSjFVm39AVd9pSDymebL5HIzjRKSjql71/JLcACxR1Z1Ox2VMS2ZXDiYYvCUig3Hft/+9JQZj\nnGdXDsYYY6qxDmljjDHVWHIwxhhTjSUHY4wx1VhyMMYYU40lB2OMMdVYcjDGGFONJQdjjDHVWHIw\nxhhTzf8POY/JrDgZamwAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "good_results = []\n", "bad_results = []\n", "steps = np.linspace(0, 1, num=10, endpoint=True)\n", "\n", "for s in steps:\n", " frame = SimulatedFrame((20, 30))\n", " frame.add_spot((10, 15 + s), 200, 2.5)\n", " feature = tp.locate(frame.image, 13, topn=1, preprocess=False, engine='python')\n", " good_results.append(feature)\n", " feature = tp.locate(frame.image, 9, topn=1, preprocess=False, engine='python')\n", " bad_results.append(feature)\n", " \n", "good_df = pd.DataFrame([pd.concat(good_results)['x'].reset_index(drop=True) - 15, pd.Series(steps, name='true x')]).T\n", "bad_df = pd.DataFrame([pd.concat(bad_results)['x'].reset_index(drop=True) - 15, pd.Series(steps, name='true x')]).T\n", "\n", "fig, ax = plt.subplots()\n", "ax.plot([-0.1, 1.1], [-0.1, 1.1], color='gray')\n", "ax.plot(bad_df['true x'], bad_df['x'], marker='s', lw=0, mfc='r', mec='r', mew=1, label='diameter=9')\n", "ax.plot(good_df['true x'], good_df['x'], marker='o', lw=0, mfc='k', mec='k', mew=1, label='diameter=13')\n", "ax.set(xlabel=r'ground truth $x$ mod 1', ylabel=r'measured $x$ mod 1')\n", "ax.set(xlim=(-0.1, 1.1), ylim=(-0.1, 1.1))\n", "ax.legend(loc='best');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## How tracking errors can qualitatively affect physical interpretation\n", "\n", "Random errors in particle tracking can lead to systematic errors in the derived statistics, such as mean squared displacement, that are qualitatively misleading. The nature of the distortion depends on the particulars of the experiment. For simple Brownian diffusion---a random walk---random errors in locating a particle add a constant offset to the mean squared displacement, $\\langle \\Delta r^2\\rangle$. This and other cases are derived in detail in [the paper by Savin and Doyle](http://dx.doi.org/10.1529/biophysj.104.042457)." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "N = 10000\n", "traj = pd.DataFrame(np.random.randn(N, 2).cumsum(0), columns=['x', 'y'])\n", "noisy_traj = traj + np.random.randn(N, 2)\n", "traj['frame'] = np.arange(N)\n", "noisy_traj['frame'] = np.arange(N)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAERCAYAAACHA/vpAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd0VNX6xvHvTiCEXgLSIaFDCi30\nYgEBpYvoxYqCXL32q15BEURRUbn2dsH2sytNREWQJiAivYXepLdAQighZfbvjxNi6BMyk5lMns9a\nrEVmzpx5E8Y8vmefvbex1iIiIuKOIF8XICIieYdCQ0RE3KbQEBERtyk0RETEbQoNERFxm0JDRETc\nptAQERG3KTRERMRtCg0REXFbAV8X4Glly5a14eHhvi5DRCRPWbp06SFrbblLHRdwoREeHs6SJUt8\nXYaISJ5ijPnLneN0eUpERNym0BAREbcFTGgYY7obY8YkJib6uhQRkYAVMGMa1topwJTY2Nh7zn4u\nNTWVXbt2kZyc7IPKRLwnNDSUKlWqULBgQV+XIvlEwITGxezatYvixYsTHh6OMcbX5Yh4hLWW+Ph4\ndu3aRUREhK/LkXwiYC5PXUxycjJhYWEKDAkoxhjCwsLUQUuuyhehASgwJCDpcy2Zko/C4W1ef5t8\nExoiIgFry2zs+61J/eZ2cLm8+lYKjTwqPDycQ4cO+bqMTMWKFQNgzpw5dOvWzcfViOQTJxPgh4fg\n817sOprO0JQ7sV7uPvPFQLg/sdZirSUoSHktIjmwbgr2x8ewxw8yNq0b3xS9nee6NvX6Jct8Fxoj\npsSxds9Rj56zQaUSDO8eecHnt2/fTufOnWnRogVLly7l559/ZtSoUSxevJiTJ09y4403MmLECMDp\nIO68806mTJlCamoq48aNo169esTHx9OvXz92795Nq1atsNZmnv+1117j448/BmDgwIE88sgjbN++\nnS5dutCyZUsWLFhAs2bNuOuuuxg+fDgHDhzgyy+/pHnz5mfU2bVrV1566SViYmJo3LgxvXv3Ztiw\nYQwbNoyqVavSr18/evbsyZEjR0hNTWXkyJH07Nnzgt/34sWLGTRoEOPHj6dmzZo5+RGLyGlH98Av\ng2HtZDYHRfDvU88R0/wqfriuHsVDvX/rtf53N5ds2rSJf/3rX8TFxVG9enVeeOEFlixZwqpVq/jt\nt99YtWpV5rFly5Zl2bJl3HfffYwePRqAESNG0LZtW+Li4ujduzc7duwAYOnSpXzyySf8+eefLFy4\nkLFjx7J8+XIANm/ezGOPPcb69etZv349X331FfPnz2f06NG8+OKL59TYrl075s2bR2JiIgUKFOD3\n338HYN68ebRv357Q0FAmTZrEsmXLmD17No899tgZ4ZXVggULuPfee5k8ebICQ8QTXC5YNBb7TjNS\n103lldSbGBjyCkMG9OOF3tG5EhiQDzuNi3UE3lS9enVatmyZ+fV3333HmDFjSEtLY+/evaxdu5aY\nmBgAbrjhBgCaNm3KxIkTAZg7d27m37t27Urp0qUBmD9/Pr1796Zo0aKZr503bx49evQgIiKC6Oho\nACIjI+nQoQPGGKKjo9m+ffs5NbZr14633nqLiIgIunbtyq+//sqJEyfYtm0bdevWJTU1laeeeoq5\nc+cSFBTE7t272b9/PxUqVDjjPOvWrWPQoEFMnz6dSpUqefCnKJJPHd4GUx6Gbb+xyDTkyeT+XNO6\nJT93qkPRQrn7azzfhYavnP6lDrBt2zZGjx7N4sWLKV26NP379z/jXvtChQoBEBwcTFpa2mW/5+nz\nAAQFBWV+HRQUdN7zNmvWjCVLllCjRg2uvfZaDh06xNixY2natCkAX375JQcPHmTp0qUULFiQ8PDw\n884RqFixIsnJySxfvlyhIZIT6anw5wfYWS9wymUYkTqAxWV68PpdDWlcrbRPStLlKR84evQoRYsW\npWTJkuzfv5+pU6de8jXt27fnq6++AmDq1KkcOXIEcLqD77//nhMnTnD8+HEmTZpEu3btLquukJAQ\nqlatyrhx42jVqhXt2rVj9OjRtG/fHoDExESuuOIKChYsyOzZs/nrr/OvpFyqVCl++uknhgwZwpw5\ncy6rFpF8b9dSGHM1TB/KAlckV514mcItB/LjQ+18FhigTsMnGjZsSOPGjalXrx5Vq1alTZs2l3zN\n8OHD6devH5GRkbRu3Zpq1aoB0KRJE/r37585qD1w4EAaN2583stP7mjXrh0zZ86kcOHCtGvXjl27\ndmWG0K233kr37t2Jjo4mNjaWevXqXfA85cuX58cff+S6667j448/pkWLFpdVj0i+k5wIs0ZiF43l\naIEw/pPyCBtLX8Wbt8XQokaYr6vDXGggM6+KjY21Z2/CtG7dOurXr++jikS8S5/vAGEtrPoOpg/F\nHj/IuKAujDzZh37to3i0Yx1CCwZ79e2NMUuttbGXOk6dhoiIr8VvcQa6t89ja6H6PHTqIdIrNOKL\nu6OJqVLK19WdQaEhIuIr6WnwxzvYOS+RSkFesvfw9fGrebBTXQa1r0HBYP8bdlZoiIj4wt5VMPl+\n2LeKxYVa8WDirdSpXYfpvaKpFlbE19VdkEJDRCQ3pZ6Eua9if3+Tk8El+U/ao8w3rXimbyQ3NKns\n9ysXKzRERHLLtrnO2MXhrcws1JHHEm+ibXRtZvSMpGyxQpd+vR9QaIiIeNvJIzBtKKz4goTQKjyQ\n+jTrghvz4i1RdI2p6OvqssX/RlnELVoa/fL179+f8ePH+7oMhg0bxowZM3xdhniTtRA3Cd5tiV35\nNd8UupEWCSMpF9OJGf++Ms8FBqjTyHVaGv1Mzz77LOHh4fTv39+t49PT0wkO9u796rnlueee83UJ\n4k2Ht8HPT8DmX9lXpC73nHqQgwXr896dUXSoX97X1V22/Peba+pg+KSrZ/9MHXzRt9y+fTt169bl\njjvuICoqip07d3LfffcRGxtLZGQkw4cPzzw2PDyc4cOH06RJE6Kjo1m/fj0A8fHxdOrUicjISAYO\nHHjO0uhRUVFERUXxxhtvZL5nvXr16N+/P3Xq1OHWW29lxowZtGnThtq1a7No0aJz6uzatWvmaruN\nGzfO/KU2bNgwxo4dy7Fjx+jQoUNmbZMnT77o97148WIaN27Mli1b3PiHubDw8HCefPJJmjRpwrhx\n4xg7dizNmjWjYcOG9OnThxMnTgBOB/HQQw/RunVratSokdlNWGt54IEHqFu3Lh07duTAgQOZ5545\ncyaNGzcmOjqau+++m1OnTmW+55AhQ2jUqBGxsbEsW7aMzp07U7NmTT744INzaty+fTv169fnnnvu\nITIykk6dOnHy5EkAVqxYQcuWLYmJiaF3796ZS8Bk7XgGDx5MgwYNiImJ4fHHHwfg4MGD9OnTh2bN\nmtGsWbPMVYfFz6WlwNxX4b2WpG//nXdDBtD68DNENm3P9H+3z9OBAfkxNHxES6PnTFhYGMuWLeMf\n//gHN9xwA4sXL2blypXUr1+fjz76KPO4vXv3Mn/+fH788UcGD3bCfNKkSWzYsIG1a9fy2WefsWDB\nAgCSk5Pp378/3377LatXryYtLY33338/81zVqlVjxYoVtGvXLvMX/MKFC88I+aw2bdrE/fffT1xc\nHKVKlWLChAkA3HHHHbz88susWrWK6OjozL1TTouPj2fSpEnExcWxatUqhg4dCsDDDz/Mo48+yuLF\ni5kwYQIDBw7M8c9RvGzbPPhfO5g1kjVFW9L62Ct8HdyNzwa0YlSfGErk0vLl3pT/Lk9dN8onb6ul\n0f+2evVqbr/9dgD27dtHSEhIZoc0c+ZMwsLOXV/n5ptvzvz7mjVrGDp0KAkJCRw7dozOnTtnPter\nVy+CgoJo0KAB+/fvz/zZ9evXj+DgYCpVqsQ111wDwIYNG4iIiKBOnToA3Hnnnbz77rs88sgjAPTo\n0QOA6Ohojh07RvHixSlevDiFChUiISGBUqXOnKkbERFBo0aNMv/ttm/fTmJiIgkJCVx55ZWZ79G3\nb98zXleyZElCQ0MZMGAA3bp1yxwTmjFjBmvXrs087ujRoxw7dixz/Ej8yP44mPEsbJrOySKVebrA\nU0w6EMXdbSN4rFMdioQEzq/awPlO/JyWRv9bdHQ0K1asANwf08j68+vfvz/ff/89DRs25NNPPz1j\nJd2s33NO11XL+vM6+2d5vp9f1mOCg4MzL09dSoECBVi0aBEzZ85k/PjxvPPOO8yaNQuXy8XChQsJ\nDQ3N0fchXnTyCMx+ERZ/SHpICSaVGsDT+9oRUSGMSXfG0Kiqfy0B4gm6POUDWho9Z5KSkqhYsSKp\nqal8+eWXlzy+ffv2fPvtt6Snp7N3715mz54NQN26ddm+fTubN28G4PPPP8/sCDylZMmSlC5dmnnz\n5l3wPY4dO0ZiYiLXX389r7/+OitXrgSgU6dOvP3225nHnQ5a8QNpp2DBO/BmI+ziD1lVoQ+tjo9m\n2OFOPHZ9DFMebBuQgQF5oNMwxvQCugIlgI+stdN9XFKOaWn0nHn++edp0aIF5cqVo0WLFiQlJV30\n+N69ezNr1iwaNGhAtWrVaNWqFQChoaF88skn9O3bl7S0NJo1a8a9997rsTpP+7//+z/uvfdeTpw4\nQY0aNfjkk0/OeD4pKYmePXuSnJyMtZbXXnsNgLfeeov777+fmJgY0tLSaN++/XkH4SUXuVywdhLM\nGAEJfxFfvi1PJPZh1rbydI2pyDNdG1ChZGB3hj5ZGt0Y8zHQDThgrY3K8ngX4E0gGPjQWjsqy3Ol\ngdHW2gEXO7eWRpf8Rp/vXLJtLvw6DPYs51RYfd4Iup33d4ZT64piPNs9kra1y/q6whzx96XRPwXe\nAT47/YAxJhh4F7gW2AUsNsb8YK09PRI4NON5EZHcc2Q7THsa1v9IevHKTKr6NEO2NKBwSEGe6VaH\nO1pV98vVaL3FJ6FhrZ1rjAk/6+HmwGZr7VYAY8w3QE9jzDpgFDDVWrvsfOczxgwCBgGZl21ERHIk\nORHmvw5/vIcNCmZpzQe5d0sLjsQHc1uLajzcsQ5liob4uspc509jGpWBnVm+3gW0AB4EOgIljTG1\nrLXnXNS11o4BxoBzeep8J7fW+v3qkSLZFWg7b/qFtBRY+gn89jKciOdAjV48crAHC+JCaVe7LMO7\nN6DWFcV9XaXP+FNonJe19i3grZycIzQ0lPj4eMLCwhQcEjCstcTHx+uWXE+xFjb+4lyKOryFU1Va\n83rZO/hgbQmqlinMmNsbcG2D8vn+d4g/hcZuoGqWr6tkPJZjVapUYdeuXRw8eNATpxPxG6GhoVSp\nUsXXZeR9B9bDL4Nh62xcYbX5KeoNnlhZHmsNj3asxT+vrOH1PbrzCn8KjcVAbWNMBE5Y/AO4xRMn\nLliwIBEREZ44lYgEkpMJMOclWDQWChVjY5Oh3Le+EVuWpNAlsjxDu9WnSmn/3UXPF3wSGsaYr4Gr\ngLLGmF3AcGvtR8aYB4BpOLfcfmytjcvGObsD3WvVquWNkkUkkKSnwtJPYc4oOBHPsejbGXa0FxMX\nJFOjbEE+u7sR7euU83WVfskn8zS86XzzNEREMm2ZDVP/A4c24qrWhnFl72P44gIYDA92qMWAthEU\nKpD/LkX5+zwNEZHclWW+BaUjWHfl/7h/aXm2bjzBdVFX8Ey3BlQqVdjXVfo9hYaIBLaUEzD/Nfj9\nLQgKJqnNUww7cCWTpsVTPQw+uasZV9e9wtdV5hkKDREJTNbC2skwfSgk7sQVdSPflrqHkXMTSHUd\n4eEOtbnvqpq6KyqbAiY0NBAuIpn2xzm30G6bC+Wj2ND6vzyysAjrlsRzVd1yjOgRSfWwopc+j5xD\nA+EiEjiOH4LZLzh3RoWW5GTbIbywvzlfLNpDhRKhDOvegOuiKuT7CXrno4FwEck/0lNh0RiY8zKk\nHMM2G8i0cncxdNoeDh/fw4C2ETx6bR2KFdKvvJzST1BE8rbNM+GXIXBoA9TswO4Wwxg87xTz5m4j\nunJJPr2rGVGVS/q6yoCh0BCRvOnwVpg2FDb8BKUjSL3pK97bXZt3P9tCSHAQw7s34I5W4QQH6VKU\nJyk0RCRvSU6Eea/BwvcgqCB0GM4fV9zM01M2sfXQJrrGVGRYtwaUL6GFHL0hYEJDd0+JBLi0FFjy\nsbNk+cnD0LAf8S2HMHLuESb9tIJqZYrw6V3NuEpzLrxKd0+JiH+zFtZMgJnPQcJfEHEl6R1H8OWO\n0rw6bQOnUl3888oa3H91Lc25yAHdPSUied/+OPj5Cfjrd6gQDbdNZGmBxgybEEfcnj20qRXGcz2j\nqFmumK8rzTcUGiLif5L2w5wXYdlnEFoKur/Jwdo38/K0jYxf+gcVSoTydr/GdIupqDkXuUyhISL+\nIzXZGeCe919IS4bmg0hr+wSfr0ritdfmkpyazj+vrMFD19SmqOZc+IR+6iLie9bCuh/g1+FwZBvU\n7QqdnmdBQklGfLiWDfuTMvbnjqTWFboU5UsKDRHxrZ2LnEUFd/4J5erD7ZM4eEUbXvhpLd+vWE+V\n0oX54LamdI7U/tz+IGBCQ7fciuQx8Vtg5ghnJdpi5aH7W6Q3vIUvFu1m9OdzOJXq4qFravEv3RXl\nV3TLrYjkrpMJMPdV+PN/EBwCbR6CVg+wbH8qz3y/hrg9R2lbqywjekbqrqhcpFtuRcS/nJ6cN/cV\nOHEYGt8G1wwl3pTm5Snr+W7JLsqXKMQ7tzSma7TuivJXCg0R8S5rIW4izBiRMTmvPXQaiat8DF8v\n3sErv/zG8VNp/LN9DR7sUFsr0fo5/euIiPfsWeFshrTjDygfBbdOgFodWL37KEPfX8DKnQm0rFGG\n53tGUbt8cV9XK25QaIiI5x3ZDrNfglXfQpEw6P4mNL6dxGQXoyfH8cWffxFWtBBv3NyIno0q6VJU\nHuJWaBhjyrhxmMtam5DDekQkLzt2EOaNhsUfQVCwM8jd7jFcISX4bslOXpm2gYQTKdzZKpx/d6pD\nidCCvq5YssndTmNPxp+L/e9AMFAtxxWJSN5zKgn+eBcWvA2pJ51B7qsGQ4lKrNqVwDOTnUtRzcJL\nM6JHCxpUKuHriuUyuRsa66y1jS92gDFmuQfquWyapyHiA+lpsPwzmP0iHD8I9XtAh2FQtjYJJ1J4\nZdJqvl60g7LFdCkqULg1T8MYE2qtTc7pMblB8zREcoG1sOlX+PUZOLgeqrWCTiOhSizWWiYs282L\nP68j8WQq/VuH80jH2hTXpSi/5tF5Gu6EgT8Ehojkgn1rnGU/ts6GMjXgps+hfncwho37kxj6/RoW\nbTtMk2qleKF3NPUr6lJUILlkaBhjrgVuAt611q4wxgyy1o7xfmki4leOHYRZz8HyL6BQCej8EjQb\nCAVCOH4qjbdmbuKj+dsoFlqAl26I5ubYqgRpf+6A406ncTdwHzA04y6qRt4tSUT8SloKLBrjbLOa\negJa3Avtn4Aizk2V0+L28ewPcexNTObm2Ko8eV09yhQN8XHR4i3uhEZSxq20jxtjRgHNvFyTiPgD\na2HTdJj2NMRvglrXQpeXoGxtAHYnnGT45DhmrNtPvQrFeeeWxjSt7s7d+ZKXuRMaP53+i7V2sDHm\nQS/WIyL+YP9amPaUM24RVgtu+Q7qdAYgNd3FJ79v4/VfNwHw1PX1uKtNBAWDg3xZseSSS4aGtXby\nWV+/7b1yRMSnThyG2S84CwsWKgFdRkHsACjgXG5atuMIT01czfp9SXSsfwXP9oikSukiPi5aclO2\nlhExxsQCTwPVM15rAGutjfFCbSKSW9LTYMlHznyLU0nOAPdVQzLHLZKSU3l12gY+X/gX5YuHalOk\nfCy7a099CTwBrAZcni9HRHLd1jkwdTAcXAcRVzrdRfkGmU9Pj9vHsMlx7E9K5s5W4Tzeua5Wos3H\nsvsvf9Ba+4NXKskhzQgXyaYjf8H0p2HdFChVHW7+Eup1hYzuYW+iM9A9fa0z0P3B7U1pVLWUj4sW\nX8vWzn3GmA5AP2AmcOr049baiZ4v7fJoRrjIJaQch9/fdP6YIGj3GLR6AAqGApDusnz+x3ZGT99I\nmsvFIx3rMKCtBroDnbd27rsLqAcU5O/LUxbwm9AQkQuwFlaPhxnD4ehuiOoD1z4PJStnHrJu71EG\nT1zNyp0JtK9TjpE9o6gWpoFu+Vt2Q6OZtbauVyoREe/ZvczZDGnnn1CxIfT5CKq3ynw6OTWdt2Zu\nYszcrZQsXJA3/9GIHg21uKCcK7uhscAY08Bau9Yr1YiIZyXth5nPwYovoWg56PEONLoVgv6+1LRw\nazyDJ6xie/wJbmxahaevr09pzeiWC8huaLQEVhhjtuGMaeiWWxF/lJYCf34Av70CacnQ+kFn6Y/Q\nvxcPTEpOZdTU9Xz55w6qlSnClwNb0KZWWR8WLXlBdkOji1eqEBHPsBbW/wS/DoPDW6B2Z+j8IpQ9\n867C2RsO8PTE1ew7mszAthE81qkuhUOCfVS05CXZCg1r7V/eKkREcmj3Upj+DPz1O5StA7eOh9rX\nnnHIkeMpPP/jWiYu302d8sWYcGtrGlcr7aOCJS/SDB2RvC5hhzNusXocFCkLXV+DJndC8N//eVtr\n+Wn1XoZPjiPxZCoPdajN/VfXpFABdReSPQoNkbwq+SjMf93Zm9sYaPc4tHn4jHELgP1Hkxn6/Rp+\nXbufmCol+WJgC22MJJctu2tPDTvf49ba5zxTjohckisdln8Os0Y6+3LH3Ozsy12yypmHuSzfLN7J\nSz+vIyXdxZDr6jGgbQQFNElPciC7ncbxLH8PBboB6zxXjohc1NY58MtTcCAOqraEft9ClabnHLbt\n0HGGTFzFwq2HaVUjjJduiCa8bNHcr1cCTnYHwv+b9WtjzGhgmkcrukxae0oC2pHtzmZI63+EUtWg\n76fQoFfmOlGnpaW7GDtvG2/M2EhIgSBG3RDNzc2qapKeeExOxzSKAFUueVQusNZOAabExsbe4+ta\nRDzm1DGY/xoseAeCguGaZ85YJyqrNbsTeXLCKuL2HKVTg/I83yuK8iXOPU4kJ7I7prEaZ60pgGCg\nHKDxDBFPc7mcu6FmDIekvRB9E1w7AkpUOufQ5NR03pixibHztlKmaAjv39qE66Ir+qBoyQ+y22l0\ny/L3NGC/tTbNg/WIyLa5MH0o7F0JFRtB3/+Dai3Oe+iCLYd4auJqtsef4KbYKjx9fQNKFimYywVL\nfqLJfSL+4sA6+HU4bJoGJatC7zEQ3feMdaJOSzyRyos/r+PbJTupHqYlQCT3aLtXEV9L2udss7r8\ncwgpDh1HQIt7zztuYa3l59X7GP5DHEdOpPDPK2vwSIc6WgJEco22exXxlZQT8Mc7MP8NSE+B5v90\nFhUsGnbew/clOpP0ZqzbT1TlEnx6VzOiKpfM5aIlvwuY7V5F8gyXC1Z96yz9kbQH6nd3uouwmhc4\n3PL14h2M+nk9qS5N0hPfym5oDDfGfIgfb/cq4te2z4dpTzmD3JUaw40fQfXWFzx868FjDJ64mkXb\nDtO6pjNJr3qYJumJ72i7V5HccGAdzHweNvwEJSrDDWMh6sbzDnIDpKa7GDN3K2/O3ERogSBe6RND\n39gqmqQnPqftXkW86chfMOclWPkNFCoO1wyFlvdDyIX33V69y5mkt3bvUa6LqsCIHpFcoUl64ie0\n3auINyTth3mjYcknzkzu1g9A239DkTIXfMnJlHTemLGRsfO2UrZYIT64rSldoirkYtEil6btXkU8\n6WQCLHgLFr4Paaegye1w5ZPnncmd1YLNhxgyaTV/xZ+gX/OqDL6uPiULa5Ke+B9t9yriCSnH4c//\nwe9vQHKiM15x9VMXvCPqtKyT9MLDivDVPS1oXVOT9MR/ZTc0/mWtfTLrA8aYl4EnL3C8SGBLOwVL\n/w/mvgrHD0CdLnD101Dx0s331NV7GfZDHIePO5P0Hu1Yh9CCmqQn/i27oXEt5wbEded5TCSwpac5\ncy3mjILEHVC9Ldz8xQXXiMpq/9Fkhk1ew7S4/URWKsEn/TVJT/IOt0LDGHMf8C+gpjFmVZanigML\nvFGYiF9yuWDdDzD7BTi00Zlr0eNNqHH1OXtbnPtSy7dLdvLiz+tISXMx+Lp6DNQkPclj3O00vgKm\nAi8Bg7M8nmStPezxqkT80Y6F8MsQ2LMMytVzOot63S4ZFuDspDd4wir+3HaYljXK8NINMURoJz3J\ng9wKDWttIpBojHn+7JVujTFXWWvneKO47NDOfeI1R7Y7q8+u/R6KV4Re7zv7cgddevwhNd3Fh9pJ\nTwKIsdZe+qjTBxuzBvgceAVnj/BXgFhrbSvvlJd9sbGxdsmSJb4uQwJB8lFn17w/3nMCos3D0PpB\nCHGvQ8i6k16XyAo811OT9MR/GWOWWmtjL3VcdgfCWwAv44xjFMdZ9bZN9ssT8WOudGeZ8lkj4fhB\naNgPOgy75FyL087eSe+D25rQJUo76UlgyG5opAIngcI4ncY2a62WSJfAsXUOTHsa9q+Baq3glu+g\nchO3X75wazxDJq5m26Hj3Bxblaeur6+d9CSgZDc0FgOTgWZAWeADY0wfa21fj1cmkpt2LoY5L8KW\nWVCqurPFaoOebg1yAxxNTmXU1PV89ecOqpUpwlcDW9BaO+lJAMpuaAyw1p4eMNgL9DTG3O7hmkRy\nz64lzoKCm2dAkbLQaSQ0u+e8u+ZdyMx1+3l60hoOJCUzqL0zSU876Umgym5oLDXG3AbUsNY+Z4yp\nBmzwQl0i3rVraUZY/ApFwuDa56DZQLcHuQHij51ixJS1/LByD3XLF+eD25vSqGopLxYt4nvZDY33\ncPbRuAZ4DkgCJuBcrhLxf7uXOrO4N02HwmWg47NOZ1GomNunsNbyw8o9jJiylqTkVB7tWIf7rqpJ\nSAFN0pPAl+27p6y1TYwxywGstUeMMSFeqEvEs3YvywiLaU5YdBgOzQdlKywA9iaeZOikNcxcf4BG\nVUvxyo0x1Clf3EtFi/ifbN89ZYwJxtmtD2NMOf7ewU/E/+xZ7oTFxl+gcGnn1tnmg5wNkbLBWsu4\nJbt4/se1pLpcDO1an7vaRBAcpEl6kr9kNzTeAiYBVxhjXgBuBJ7xeFUiObV3pRMWG36G0FLOjnnN\n/wmhJbJ9qj0JJxk8cTVzNx6mro3DAAAUE0lEQVSkRUQZXrkxRvt0S76VrdCw1n5pjFkKdMDZgKmX\ntXadVyoTuRz7Vjthsf5HCC0JVw+FFpcXFtZavl28k5E/rcNlLc/1jOS2FtUJUnch+Vi2QsMY83LG\nfhrrz/OYiO/sWwO/jYJ1U6BQSbjqKWh5rxMcl2F3wkkGT1jFvE2HaFmjDK/0aUi1sAvv6y2SX2g/\nDcnb9sc5ncW6H6BQCbhyMLS8Dwpf3q2v1lq+WbyTFzK6i+d7RXFr82rqLkQyZHc/jRpn7adRAvjd\nG4WJXNSBdU5YrP0eQopD+/9Aq385g92XadeREwyesJr5mw/RumYYL/eJoWoZdRciWeV4Pw3giDGm\nlbX2D08XJ3KOfath3n8h7ntnIl67x6HV/VCkzGWf0lrL14uczZGstbzQO4pbmlfT8uUi55Gt/TSA\nfgDGmApAl4w/dYCFgEJDvGfnIpg72plnUagEtH3UWaY8B2EBZ3YXbWuVZVSfaKqUVnchciHuXp4K\nxlkCvQvObPAjwDTgWWvt+ou9VuSyWQvbfnPCYvs8Z1Le1UOh+T2XPWbx96nVXYhcDncvTy3CGbv4\nBXjeWnvSeyVJvudyOZPx5o12lv0oXhE6vwhN+2drbagLydpdtKkVxqgbNHYh4i53L0819XYhIrjS\nIW4SzHsNDsQ5S5R3ewMa3QIFCuX49Gd3FyN7RXFrC3UXItnh7uWpZdbai+5E484xIueVlgKrvoH5\nr8PhrVCuHvQeA1F9IDi7d4Wfn+6MEvEMd/+LrH/WrbZnM8DlzaKS/CvlBCz7DBa8BUd3Q8VGcPMX\nULcrBHlmxdis3YVL3YVIjrkbGvXcOCY9J4VIPnLsACz5BBaNgROHoHob6PEW1Ozg9k557sg6q1vd\nhYhnuDum8Ze3C5F8YM9y+PN/sGYCpKdA7U7Q9t9QvZVH38Zay3dLdvL8j5rVLeJpnrlgLHIh6anO\nelB//g92LoSQYs5dUM3/CWVrefzt9iaeZPCE1fy28SAta5Th1RsbqrsQ8SCFhnjH8XhY9iks/sgZ\nrygdDp1fgsa3XvYighdjrWXcUme/i7R0y4gekdzeUivSiniaQkM8a98a+PMDWD0O0pKhxlXQ9b/O\npaigYO+8ZWIyQyauYvaGgzQPL8OrfbXfhYi3uB0axpjSQHegN87SIduAycBka+0B75QneYIrHTZM\ndcJi+zwoUBga9nP2sbiivtfe1lrLxGW7GTEljpR0F8O7N+DOVuHqLkS8yN15GhOB0sBPwJPW2o3G\nmGpAT+ALY0yItfYqbxRojKkBPA2UtNbe6I33kMt08ggs/8K5CyphB5SsCtc+B41vz/GaUJdy4Ggy\nT01azYx1B4itXppX+zYkoqy6CxFvc7fTuNtam5D1AWvtDuBt4G1jTLYWAjLGfAx0Aw5Ya6OyPN4F\neBMIBj601o6y1m4FBhhjxmfnPcSLDm5wBrZXfg2pJ5xbZjuNdOZXeGgy3oVYa5m8Yg/Df4gjOTVd\ne3WL5DJ3b7lNyMnz5/Ep8A7w2ekHMhZFfBdno6ddwGJjzA/W2rXZPLd4Q1qKs8Lsko9hyywILgTR\nfaHFIKjYMFdKOJCUzNOT1vDr2v00qVaKV/s2pGa5Yrny3iLiuGRoGGOuBW4C3rXWrjDGDLLWjsnJ\nm1pr5xpjws96uDmwOaOzwBjzDc7lr0uGhjFmEDAIoFq1ajkpTc62bzUs/xJWfwcn4p3FA68ZCk3v\ngqJlc6UEay1TVu1l2OQ1nEhJ56nr6zGgbQ11FyI+4E6ncTdwHzDUGFMGaOSlWioDO7N8vQtoYYwJ\nA14AGhtjhlhrXzr7hRkhNgYgNjbWeqm+/ON4vHP304ovnNAIDoG610OjW6HmNV6/BJXVoWOneOb7\nNUxds4+GVUvx374NqXWFugsRX3Hnv/6kjMtPjxtjRgHNvFzTGay18cC9ufme+VJ6Gmye4QTFhl/A\nleqsBXXdqxB9o9cHts/np1V7eWbyGo4lp/Fkl3rc0y6CAsGeWZNKRC6PO6Hx0+m/WGsHG2Me9FIt\nu4GqWb6ukvGYeNOB9U5QrPwWjh+AImWh+SBnOfIKUZd+vRccPp7CM5PX8NOqvcRUKcnovg2pU764\nT2oRkTNdMjSstZMBjDFBwGBr7YteqmUxUNsYE4ETFv8AbvHSe+VvJxOc9Z9WfOlschRUAGp3dmZr\n1+4EwQV9Vtova/Yy9Ps1JJ5M5YnOdfln+xrqLkT8iNsXp621LmNMdyDHoWGM+Rq4CihrjNkFDLfW\nfmSMeQBnG9lg4GNrbVw2ztkd6F6rlufXMwoIrnTYOscJinU/QvopuCLS2REv+iYoVs6n5R05nsLw\nH+L4YeUeoiqX4IuBLahXoYRPaxKRcxlr3R83NsZ8AOzF2fLV5bWqciA2NtYuWbLE12X4j/gtTlCs\n/MZZA6pwaedW2Ua3OrfK+sG+EtPj9vHUpDUknEjhoQ61ue+qmhRUdyGSq4wxS621sZc6Lru3wZQB\nrgTuM8b8CawCVllrx11GjeItp5KcbVOXf+msLGuCoFZH6PyCcxeUB7ZO9YSEEymMmLKWSct3U79i\nCT67uzkNKqm7EPFn2QoNa+1NAMaYQkAkEI0zv0Kh4WsuF/w13wmKdT84M7XDakPHZyHmH1Cioq8r\nPMOs9fsZPGE1h4+n8HCH2tx/dS1CCqi7EPF3l3XDvbX2FLDMGFMUJzTEF1wu2LfKWSxw5VfO+k+F\nSkDMTdDoNqgS6xeXn7I6mpzK81PWMm7pLuqWL87H/ZsRVVk7BYvkFdkODWNMY5y7mm4C9uFsBXu/\nh+vKtnwzEH7sgLOMx+aZsHU2HD8IGKhxJVwzDOp3g4KFfV3lec3fdIj/jF/JvqPJ3H91TR7qUJtC\nBbyzXLqIeIe7q9zWAfrhhEUSzuWoq6y124wx27xYn9ustVOAKbGxsff4uhaPSjsFOxY6QbFlpjND\nG5z5FDWvgVodoMbVULy8b+u8iOOn0nhp6jq+WLiDGuWKMuG+1jSuVtrXZYnIZXC301iPM4/iRmvt\n6rOe07IdnmStc8fTlplON7F9PqQeh6CCUK0ldBjuhEWFGAjy/zGARdsO8/i4lew8coKBbSN4vHNd\nQguquxDJq9wNjRtwJttNN8bMAL4DfrHWpnqtsvwkORG2zXVCYstMZ2wCoEwNZ2Z2rQ4Q3hYK5Z1Z\n0cmp6bw6bQMf/76NqqWL8O2gVjSPyP2lSETEs9xdGv174PuMge+eOCvKfmiM+RnQPZLZ5UqHvStg\nc8Ylp52LwKZDSHGIaA9tHoaaHaBMhK8rvSzLdxzhsXEr2XrwOLe3rM7g6+pRtJB2FhYJBO6OaQRZ\na13W2uPAV8BXGdu/3ghUM8YYm51ZgvnR0b1OQGyZBVtmw8nDgHEm2LV9xAmJqs19uoRHTp1KS+fN\nGZv44LctVCgRyhcDWtC2du4sny4iucPd//371RhzCPge+MlaexQ4BRwG9gDL8d6S6W7xu7unUpNh\nx4KMS06z4EDGtiDFykOdLhkD2Ffl2p4U3ha3J5HHvlvJ+n1J9G1ahWe6N6BEaN4NQBE5P7eXETHG\nNMC5NHU9UBBnAHwa8IO1dpnXKswmny0jYi0c2vj3uMT23yHtpLMXRbVWTkjU7ADlI/1u7kROpKa7\neG/2Ft6etYnSRUMYdUM0Her7751cInJ+Hl9GJGPb1bXAS8aYwtbakzkpMCCcPAJbf8u402kWHN3l\nPB5WG5re6YREeBsIKerbOr1k4/4kHvtuJat3J9KzUSVG9IikVJEQX5clIl50uTPC82dguNKdpcRP\ndxO7l4J1QaGSUKM9tH/c6ShKBfaWs+kuy9h5W3lt+kaKhRbg/VubcF20fy1TIiLeoVtaLiVx198h\nsXWOc3ssBio3gXYZIVE5Nle3QPWlrQeP8fi4lSzbkUCXyAqM7B1F2WL+sQCiiHhf/vhNlx0pJ+Cv\nBX9Prju0wXm8eCWo39255FTjKp9sf+pLLpfl0wXbeWXaekKCg3jj5kb0bFQJE0DjMyJyaQqN0xJ3\nweQHnMBIPwUFQqF6a2hyh9NNlKsXUAPY2bHt0HH+M34li7cf4eq65RjVJ4byJUJ9XZaI+EDAhEaO\nb7ktUtYZ2G42EGpdA9Xb+O3Cf7kl3WX55PdtvDptA4UKBPHfvg25oUlldRci+Vi2du7LC7Rzn2ds\nOXiM/4xfxdK/jtCh3hW8eEO0uguRAOatnfskwKW7LB/P38bo6RsILRjM6zc3pFcjdRci4lBoSKbN\nB47xxPiVLN+RQMf65XmxdxRXqLsQkSwUGkJauosP52/jtV83UiQkmDf/0YgeDXVnlIicS6GRz23a\nn8Tj41excmcCnSPL83yvKK4oru5CRM5PoZFPpaW7GDNvK2/8uomihYJ5u19jusVUVHchIhel0MiH\nNuxL4onxK1m1K5HroirwfC/N6hYR9wRMaPjd0uh+KDXdxf9+28JbMzdTLLQA797ShK4xWjNKRNwX\nMKFhrZ0CTImNjb3H17X4o/X7jvL4uJWs2X2UrjEVea5HJGHqLkQkmwImNOT8UtNdvD/H2e+iZOGC\nWpFWRHJEoRHA1u45yhPjVxK35yjdGzr7XZQpqv0uROTyKTQCUEqai/fmbOadWZspVSSED25rSpeo\nCr4uS0QCgEIjwMTtSeTxcatYt/covRpVYnj3SEqruxARD1FoBIiUNBfvzN7Me7M3U7poCGNub0qn\nSHUXIuJZCo0AsGZ3Io+PW8n6fUnc0KQyw7o10F7dIuIVCo087FRaOu/M2sx7c7YQVjSEj+6MpUP9\n8r4uS0QCmEIjj1q9y+kuNux3uovh3SIpWaSgr8sSkQAXMKGRX2aEp6S5eHvWJt6bs4WyxUL4uH8s\n19RTdyEiuSNgQiM/zAjPOnbRp0kVhnVroO5CRHJVwIRGIMvaXWjsQkR8SaHh586+M0pjFyLiSwoN\nP6XuQkT8kULDD6m7EBF/pdDwIylpLt6ZtYl31V2IiJ9SaPgJdRcikhcoNHzs7O7iwzti6dhA3YWI\n+CeFhg+puxCRvEah4QNZV6Qto+5CRPIQhUYuO6O7aFyZ4d3VXYhI3hEwoeHva0+dvd+FugsRyYsC\nJjT8ee2ps7uLYd2134WI5E0BExr+SN2FiAQahYaXZN2rW92FiAQKhYaHqbsQkUCm0PCgrN1F78aV\nGa7uQkQCjELDA87uLsbeEcu16i5EJAApNHJI3YWI5CcKjcuUkubi3dmbeVfdhYjkIwqNy6DuQkTy\nK4VGNqi7EJH8TqHhJnUXIiIKjUtKSXPx3pzNvDNL3YWIiELjIrJ2F70aVeLZHpHqLkQkX1NonEdq\nujN28c6szZQqEsKY25vSKbKCr8sSEfE5hcZZ1u45yuPjVrJW3YWIyDkUGhnUXYiIXJpCI8PuIyd5\nf84WusVUVHchInIBARMaOd25L7xsUWb8+0qqlini2cJERAJIkK8L8BRr7RRr7aCSJUte9jkUGCIi\nFxcwoSEiIt6n0BAREbcpNERExG0KDRERcZtCQ0RE3KbQEBERtyk0RETEbcZa6+saPMoYcxBIABIv\ncljJizxfFjjk6bq87GLfjz+/V07Old3Xunu8O8dd6phA+3xB7n3G9Pny3eerurW23CWPstYG3B9g\nzOU+Dyzxdf2e/n799b1ycq7svtbd4905Lr99vjz9755b76PPl3f+BOrlqSk5fD6vyc3vx5PvlZNz\nZfe17h7vznH57fMFufc96fPl55+vgLs8lVPGmCXW2lhf1yGBSZ8v8abc+HwFaqeRE2N8XYAENH2+\nxJu8/vlSpyEiIm5TpyEiIm5TaIiIiNsUGiIi4jaFhoiIuE2hcRHGmKLGmP8zxow1xtzq63ok8Bhj\nahhjPjLGjPd1LRJ4jDG9Mn5/fWuM6eSJc+a70DDGfGyMOWCMWXPW412MMRuMMZuNMYMzHr4BGG+t\nvQfokevFSp6Unc+YtXartXaAbyqVvCibn6/vM35/3Qvc7In3z3ehAXwKdMn6gDEmGHgXuA5oAPQz\nxjQAqgA7Mw5Lz8UaJW/7FPc/YyLZ9SnZ/3wNzXg+x/JdaFhr5wKHz3q4ObA54//6UoBvgJ7ALpzg\ngHz4s5LLk83PmEi2ZOfzZRwvA1Ottcs88f76ReiozN8dBThhURmYCPQxxrxPYK4nJLnnvJ8xY0yY\nMeYDoLExZohvSpMAcKHfYQ8CHYEbjTH3euKNCnjiJIHKWnscuMvXdUjgstbG41xvFvE4a+1bwFue\nPKc6DcduoGqWr6tkPCbiKfqMiTfl2udLoeFYDNQ2xkQYY0KAfwA/+LgmCSz6jIk35drnK9+FhjHm\na+APoK4xZpcxZoC1Ng14AJgGrAO+s9bG+bJOybv0GRNv8vXnS6vcioiI2/JdpyEiIpdPoSEiIm5T\naIiIiNsUGiIi4jaFhoiIuE2hISIiblNoiIiI2xQaIiLiNoWGiIcYY6oYY87Z6MYYE26MOWmMWXGR\n1xY2xqwwxqQYY8p6t1KRy6fQEPGcDkCTCzy3xVrb6EIvtNaezHh+j1cqE/EQhYaIBxhj2gKv4exb\nsMIYU+MixxY1xvxkjFlpjFlzvu5ExF9pPw0RD7DWzjfGLAYet9auucThXYA91tquAMaYkl4vUMRD\n1GmIeE5dYL0bx60GrjXGvGyMaWetTfRyXSIeo9AQ8YCMwevEjCWqL8pauxFn7GM1MNIYM8zb9Yl4\nii5PiXhGOG4OYhtjKgGHrbVfGGMSgIHeLEzEkxQaIp6xHihrjFkDDLLWLrjIsdHAq8YYF5AK3Jcb\nBYp4gkJDxAOstceA5m4eOw1nhzWRPEdjGiLelw6UdGdyH1AQcOVaZSLZpO1eRUTEbeo0RETEbQoN\nERFxm0JDRETcptAQERG3KTRERMRtCg0REXGbQkNERNz2/5nliYoqU6eQAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots()\n", "ax.plot(tp.msd(traj, 1, 1)['msd'], label='random walk')\n", "ax.plot(tp.msd(noisy_traj, 1, 1)['msd'], label='random walk + random noise')\n", "ax.legend()\n", "ax.set(xscale='log', yscale='log')\n", "ax.set(xlabel=r'$t$ [s]', ylabel=r'$\\langle \\Delta r^2 \\rangle$ [\\textmu m$^2$]');" ] } ], "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.3" } }, "nbformat": 4, "nbformat_minor": 1 }