{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Linear vs Circular, Normalized vs non-normalized correlations" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# we should test linear vs circular to understand the differences\n", "# the main difference in the new version was a bug as \n", "# the zero-padding for linear case does not work as long as you leave values on the \n", "# borders not zero. \n", "# see new normalize_intensity function that should take care of uneven illumination\n", "# both between the frames and at different regions of the image A or B\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "from skimage.registration import phase_cross_correlation\n", "from skimage import feature\n", "# from skimage.feature.phase_cross_correlation import _upsampled_dft\n", "from scipy.ndimage import fourier_shift\n", "from skimage.exposure import match_histograms" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "from openpiv.tools import imread\n", "from pylab import *" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "from numpy.fft import rfft2, rfftn" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "from openpiv.pyprocess import find_subpixel_peak_position\n", "from openpiv.pyprocess import normalize_intensity, fft_correlate_images" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "a = imread('../test11/A001_1.tif')\n", "# b = imread('../data/PIVChallenge2001_A/A001_2.tif')\n", "a = a[:32,:32].copy()\n", "a[16:18,16:18] = 255\n", "# b = b[:32,:32]\n", "\n", "# should be in the order of y,x:\n", "# so it's about 5 pixels upwards and about 3 pixels to the right\n", "shift = (-12.035, -10.92)\n", "# The shift corresponds to the pixel offset relative to the reference image\n", "b = fourier_shift(np.fft.fftn(a), shift)\n", "b = np.fft.ifftn(b).real\n", "b = match_histograms(b,a).astype('uint8')\n", "# b = b + np.linspace(10,85,32)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAATEAAAD5CAYAAABPqQIFAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAArJklEQVR4nO2de3Rd1XHwf6O3JdmSLRm/8AOMIdgEDDaEFkgIaRqHkEJaQqBZ5LlwaGAt+BZNA6RtaBNWk6Yh5WsCiRMIJKU8AoRQSgIO0AD9gsE2xjY2D5vYxrb8kGzLlqz3ne+Pc9R1Je05upKur3Rgfl5n+d49Z5+9z0Nz99kze0ZUFcdxnLRSNNodcBzHGQmuxBzHSTWuxBzHSTWuxBzHSTWuxBzHSTWuxBzHSTUlI6ksIkuAW4Fi4Ceq+q2k/Yurq7SkblJQVtRh1yvqCZdnihP6lrFlmqC6F0zZa8pebQ33XdvsjhR12m1R222K5lfuN2Xr9042ZVphnHhGzDqlB21Z13hTRFGpcWMAaQ5fk0ypfTzsbiTez0zCUyxGFyXBsyjpniU9Oz3V9kGl0z45LQ3Xk+6EC2I01XVgHz2trQkVB+cjH6zSpn32vc1m1dqOJ1R1yUjaGynDVmIiUgz8APgwsB14SUQeVdUNZmN1k5j21WuCsvGbbUVQ1hy+Y5019r1KehAzZbbsxa/cZspOeuHTwfKetTVmncqd9oOtF+yz+7H4flM2/wdfNmUdJ7YFyzOt9q2e9rR97XedY/e/cnqLKSv/zYRgedvkpD9mU0TJYVvWXm/3sbQ53F5xl328qp22xuyqtPu/70z7oBXb7IeufVq4XmmTfc+KDaW49Ue3mHVypWlfDy8+MSunfYunvVk/4gZHyEheJ88ANqnqW6raCdwHXJifbjmOM1ookMnx31hgJEpsBvB21vftcZnjOClGUbq0J6dtMERkpog8IyIbRORVEbkmLr9JRHaIyJp4Oz+rzg0isklEXheRjwzWxojmxHJBRJYCSwGKJ9Ue6eYcx8kDeRxldQPXqepqERkPrBKR5bHse6r6L9k7i8h84FJgATAd+K2IHK9qa8yRjMR2ADOzvh8dl/VBVZep6mJVXVxcXT2C5hzHKQSK0qO5bYMeS7VBVVfHnw8BG0l+Y7sQuE9VO1T1D8Amoqkrk5EosZeAeSJyjIiUEWnPR0dwPMdxxggZNKcNqBeRlVnbUuuYIjIHOBVYERddLSJrReROEZkYlw15mmrYr5Oq2i0iVwNPELlY3KmqrybVKT0EU58LW1Uypfbw9cSrw4dd8/P3mnU6au1+tM2yrUiLbvorUya14b6XJVi69p9sn9fEx8IuGwAnrrYtkO0J/S/dVhEWzGo363TU2I9B3WrbGlf160pTtu2j4fOuaLAtoRW2dwslbfavft2rtqvK1j8Ll094zTaFNp08PA+Fj5/8iil75rXTTVnt2nBffvLX/2rWufiZ8POhZSOPSqNAj+XDMZBGVV082E4iUg08BFyrqgdF5HbgG3Fz3wC+C3xhOP0d0ZyYqj4OPD6SYziOM/bI5K7EBkVESokU2D2q+jCAqu7Okv8YeCz+mtM0VTbuse84Th8U6FLNaRsMERHgDmCjqt6SVT4ta7dPAOvjz48Cl4pIuYgcA8wDXkxq44hbJx3HSReKDuV1cjDOAi4H1onImrjsRuAyEVlIpDO3AF8CUNVXReQBYAORZfOqJMskuBJzHKc/Cj150mGq+jzhBWXmNJSq3gzcnGsbrsQcx+lD5LGfHlyJOY7TD6EnaTX+GKOgSqx7HOw/cei2hP9+5T3B8voEs3t3VcLi8MO2mT9TatfrqA23N26vXafkoH2+SYuhkxY8T5lpR7hobjgqWK47DdcLoNxYYD8YTQtsN4W61eFjNp5tr8wvX2Uvkt53nu0i8tJ5PzVlJy4LuyIcfK/dj6P/034+ti+xxyjPPGC7USRhuY987vZrzToTz24Klu8tyS36RBLRxL4rMcdxUkrkJ+ZKzHGcFJPxkZjjOGnFR2KO46QaRehJkR+8KzHHcQbgr5MG0gNlzWHZX35heVgA/PvPPxwsL20dnjdLSav9K5O0cNwiKZzxK1/5oSmbe/+VpmzR6W+aslUvzTNlNaeFrVanHLXTrPNs6QJTVr3NvlYfvOQlU/bcT8OWuj8s+YlZ58ZTTzZlv/zV2aZs/kZ7sXzFgXB5WbNtCW08xRQxbpttuWydZVsGk8Kvt8wKK4wkC3XZAxOD5bJv5H/SitCpCQksxhg+EnMcpw+Rs6u/TjqOk2J8Yt9xnNSiKvQk5aYbY7gScxxnABkfiTmOk1aiif30qIb09NRxnILgE/sJSMbOzH3Hf/6JWW/6q+GY8odm2N2v3Gu7PYy//G1Ttmm3ndC4ZF04W1NbvT30PvVm2/wfzpEdsa75BFNWbbipALS0heP2r95bZ9aZnHCtZl35uil77IXTTFnJjPCi5je6Ws06jzxou1F01tt9LGu2/+AWXr4uWP6eql1mnd/smm/Ktr883ZQVH7b7kUnIbm5liV92w61mnaX/dE2wPF9TWT3uJ+Y4Tlpxj33HcVJPxq2TjuOklWgBuCsxx3FSiiJ0pWjZUXrUreM4BUEVerQop20wRGSmiDwjIhtE5FURuSYu/46IvBZnAP+liNTG5XNEpE1E1sSbvfg4xkdijuP0Q/Lp7NoNXKeqq0VkPLBKRJYDy4EbVLVbRL4N3AB8Na6zWVUX5trAiJSYiGwBDgE9QPdg6cyLO5WaLWF3iXGN9vB1z8KwfbrifeGIDQD7u+1Ta//5bFNWfdE+ux5hF4vWc1rMOuWrw3UADp/cZsqWnvK8KfvJr213lJLD4YevbbIdR79lUYcp27d9himb9eukKCJh2V8c+IpZo2OafbzKBvtXv73ePreN3w9H6FhfdpJZp7PG/gOuTchHUH2ZHSlk93/b19HKcfB3Wy4acp2ikYfYj+bE8jSxr6oNQEP8+ZCIbARmqOqTWbu9AFw83DbyMRL7oKo25uE4juOMEYYwsV8vIiuzvi9T1WWhHUVkDnAqsKKf6AvA/VnfjxGRl4GDwN+q6nNJHfDXScdx+qDIUIIiNg72BgYgItXAQ8C1qnowq/xrRK+c98RFDcAsVW0SkUXAIyKyILtOf0aqxBR4UkQU+JGlgR3HSQ9Ryrb8jW9EpJRIgd2jqg9nlX8OuAD4kKoqgKp2AB3x51Uishk4HljZ/7i9jLSnZ6vqDhE5ClguIq+p6rP9TmApsBSgfFztCJtzHOfIk7/kuSIiwB3ARlW9Jat8CfA3wAdU9XBW+WRgn6r2iMixwDzgraQ2RqTEVHVH/P8eEfklcAbwbL99lgHLAMbXHj28LK2O4xQMJa8e+2cBlwPrRGRNXHYj8H+BcqLBD8ALqnol8H7gH0Wki8g6dKWq2tY2RqDERKQKKIotDlXAnwL/ONzjOY4zdsjXSExVn4fgwR439n+I6NUzZ0YyEpsC/DLWoiXAf6jqb5IqdE6ErX8RHoxNWDd0zd9yuNyUTXqs0pT12DkiaF9hR3voqgn3vfo5240ik9BWzQQ7E8Td94WTowDUbbVdEVqnWw+f/VBOmnTIlHU9MtmU7T7dPmbFXlNkkhSNYuoK2x2l6Tr7OrY1hO/nobm2L0JFg+3u01Vpn3Pn47YbBfbjSFt9+LzbHpxj1uk4Llze/Tu7nVxRlXfH2klVfQtIyAvjOE4aiSb207PsyF0sHMfph8fYdxwnxUQT+x4U0XGcFOOheBzHSS1D9NgfdVyJOY4zAE8UYlDULlS/FvY5ODzN9oOd+vuwS0HLAdu1oafMPl57QmKP8gOmiLoN4X7s/pRt/reSiwB8bOZGUzZtrp0N5LZffMyUWUknWmfZbhl7X7eTo0xPSCJS3mxfx72LwrLKBrvOMR/YYh9vmx15pH3FOFPWYZx33Wr7j7S40z7nJBeLJIo67XpFXeF71l1l1xm314hiEQ4SMyRUoSvjSsxxnJQSvU66EnMcJ8Xky2O/ELgScxynD+5i4ThOyvHXScdxUk4eY+wfcQqqxFSgx0jnPmGTXa9lengdV+sM2wKZFNNtxu9sE05ndcLiX8NaVPNklVmn6TR7ofGDry80ZeNW2FbNqefvMGX7jEXIE9fav6xF3aaIlukJse3tteFELyUD+cW13zFrfOqf7fj7B9+bFM9/6BGekuLoVxkWXoDKBGvtOX//e1P2yINnm7LaTeH2Dpxo92PKC+HyvMTYV+jK+NpJx3FSiju7Oo6Tevx10nGc1OLWScdxUo9bJx3HSS2qQrcrMcdx0oy/ThoUd0FVQ9hsPPnyrWa93feGF/92TbVdJaY9YZ9a0wLDzwOo2plkyg+TSbiKSa4NmVLbjaKkzTavNy63Y7mXtQ7dXF+7ceiLkwE2Lr3dlM3/wZeD5Rfc89dmnZqEc560zu5j9zhbZl3HdnvNuxnzHqC82X4+7n35DFNWYT9yputO9Tb7vA4cyRj75E+JichM4GdEOTmUKEP4rSIyiSjr9xxgC3CJqu6PU7zdCpwPHAY+p6qrk9pIz5jRcZyCkVHJacuBbuA6VZ0PnAlcJSLzgeuBp1R1HvBU/B3go0S5JucR5au1fyljXIk5jtOHXj+xfCgxVW3oHUmp6iFgIzADuBC4O97tbuCi+POFwM804gWgVkSmJbXhc2KO4wzgSPiJicgc4FRgBTBFVRti0S6i102IFNzbWdW2x2UNGLgScxynD6rQnXtQxHoRWZn1fZmqLuu/k4hUEyXFvVZVD8b5auP2VEVk6GvHYlyJOY4zgCFM7Deq6uKkHUSklEiB3aOqD8fFu0Vkmqo2xK+Le+LyHcDMrOpHx2UmPifmOE4f8jknFlsb7wA2quotWaJHgc/Gnz8L/Cqr/DMScSbQnPXaGWTQkZiI3AlcAOxR1ZPisqB5dLBjKbY7guVGATDlsrD7RSYhzfvBT9kx6jOrakzZ7nPsMADjtoXt5JpgPi85bMuqzt1jyvatsUNEHP10hynbclH4Alc02FEJuu0gHIxLiNow9/4rTVml5f3Saj/4PRc1mbKu5+tMWdkB+02k7YKDwfKkZ6AonAZiUKZMPWDKDm4+ypQd84U3guUbHz3BrFO7KXxfStrNKkNC8+cndhZwObBORNbEZTcC3wIeEJEvAluBS2LZ40TuFZuIXCw+P1gDubxO3gV8n8jXo5de8+i3ROT6+PtXcziW4zgpIF8T+6r6PJgH+1BgfwWuGkobg75OquqzwL5+xZZ51HGclKOaVz+xI85wJ/Yt86jjOKlH6Hk3pWwbzDwqIkuJPG8prZ440uYcxykAeZwTO+IMV93u7vWi7WceHYCqLlPVxaq6uKQiYQbZcZwxQe/aybS8Tg5XiVnmUcdx0o5G82K5bGOBXFws7gXOJfLM3Q58Hds8mkhRt51owUoGArDlt3OC5d0JiUKK19om9GojkkaE7S9hpY5vm2z/In38sudN2S+Wn2XKpr5kuza019l9HL85/Lt0aK7tOvK+xRtM2cbvLzBl1duG/hvYMteOPFKc4EbRUWvfs6SIFNVPhZ+DpGgUB2fb51VWlhBN479sN4ruc1pM2bYfzguW1+23r1VRR7j/0p0fzfKOCk+tqpcZogHmUcdx0o++2yb2Hcd55zFWXhVzwZWY4zgDSJN10pWY4zh9iCbtXYk5jpNixor7RC64EnMcZwA+JzYMxjXaJu9yIyDFnjPtK/3Gn//IlF3ylm1YXfdbO3LA3M+Eow2s2jLLrPPE7bYbRWmCa0ZPgil//4m25ajzuLZg+Vvn/dSsc/7r55uyItvKT6ftxUJnTfh+zn7UrtM6xb6fRV329UiKtNFVabdnkRSVpLjT7mOpkaQFoGWb7ehtJWPpKRuGhTAPRkVFyLh10nGcNJOigZgrMcdx+uET+47jpJ4UDcVciTmOMwAfiTmOk1oUyGRciTmOk1YU8JGYgdiuA0luA0c/E3YbODh7nFln0U1/ZcoypQk3yM7PwZrtM4Ll5RvtfhQlRBVIMuU3nmL3ccpLdkSK7bPC0UCSknpkarpNWWlCP2o22efWfFy4XnmTneTk0Az7cSzqNEU0H2c/OxuX3hYsP+eqL5l1uiuH50ZxeLLdjwmbTBEHjP53nBh+7gGmPxTOZqKSr9j4eTlMQUiPM4jjOIVDc9wGQUTuFJE9IrI+q+x+EVkTb1t6syCJyBwRacuS/TCXrvrrpOM4/ZB8TuzfRb9saar6qf9tSeS7QLY7+2ZVXTiUBlyJOY4zkDy9TqrqsyIyJySLE+teApw3kjb8ddJxnL4oaEZy2kbIOcBuVX0zq+wYEXlZRH4nIufkchAfiTmOEyBnBVUvIiuzvi9T1WU51r0MuDfrewMwS1WbRGQR8IiILFDVcBr3mIIqsUwxdNSEL461YBhgy8cqguXVW+22imyDG1Mu22LKHj/hcVNmWTyTLJBJ1L9iWxmTFoA3zbfzEcx+KHwdK3Yeso936gRT1nycKUq01B37i3B7rXOqzTqTXm83ZYdmlpuypMXhJy77crB8Sqe9sr1ubf7dC5KeR8sKWbwt/NwDHJgX7mPP/xtSt2xyf6QbVXXxUA8vIiXAnwOL/rdJ1Q6gI/68SkQ2A8cDK4MHifHXScdxBpIn62QCfwK8pqrbewtEZLKIFMefjwXmAW8NdiBXYo7j9KXX2TWXbRDibGm/B04Qke1xhjSAS+n7KgnwfmBt7HLxIHClqu4brA2fE3McZwD5cna1sqWp6ucCZQ8BDw21DVdijuMMxNdOOo6TZiRFy45ciTmO05eRT9oXlEGVmIjcCVwA7FHVk+Kym4ArgL3xbjeqqu2bEKPF0F0VHqZ+94Kfm/Webj4xWP7cT08367TX28Pht389x5Qd88YVpmxac/jOJsVdPzjHdofY9167j5NXmSLqNtiuGe0Tw+2N22rXySQ8BfWv2OfWOt22C+0+c3ywPOm8mubbC+knbLVdcA4muIHM+F3YlaJpgb36vuxA0sJ2+5yTXGYaT7Gfg0xnWFbWmZCDweh+flYL5TZpP1bIxTp5F7AkUP49VV0Yb4MqMMdxUsSRd7HIG4OOxJLWPjmO8w7FHviOOUbiJ3a1iKyNQ21MzFuPHMcZXfLoJ1YIhqvEbgfmAguJ1jt919pRRJaKyEoRWdlzuHWYzTmOU0hEc9vGAsNSYqq6W1V7VDUD/Bg4I2HfZaq6WFUXF1faCUQdxxlDpGhObFhKTESmZX39BLDe2tdxHOdIkouLxb3AuUQhN7YDXwfOFZGFRLp4C2AHLM8iUwYts8Izhjfe9RmznhmLPiEefs0me2YyKfpC3av2MUtbwvHhd5xtuwaM22u3NWGT/RvSU2bXa0uI5V61M3zeLfNqzDo1b9kB7J/+2R2mLAkr4sfuc+xwDtOeto+3d5E9/1LRaMt2fCD8iE9IyA9QuTfJnSPBreT0hPuywxRRcjgcLz/J1cOKilFsB+cYEmPlVTEXcrFOhtY+De/Jdhxn7KP4siPHcVLOO2kk5jjOu4931Ouk4zjvQlyJOY6TalyJOY6TVsaSI2suFFSJFbdB7UYjUUitXa9lbtgF4A9LfmLWOfZh2+tjxlN2W53VdrQBK0JE6WH7eEkJIroTfH9r19tJM7Z/0HbpKGkNX99Jr4fdQ44UlptC5W9sN4TiTjsKREWjHXVixvJmU9ZZF0620TXefvSTErF019vuKEnP49z7rzRllqtN22TbQjjl3LDPRtFzefKxcOuk4zhpxkdijuOkG1dijuOklpTNiXnKNsdxBpKnBeBxqK49IrI+q+wmEdkhImvi7fws2Q0isklEXheRj+TSVVdijuMMQDK5bTlwFzlGhhaR+UT5KBfEdW7rTaabhCsxx3GOGKr6LDBoAtyYC4H7VLVDVf8AbCIhzFcvBZ0TkwyUHg6PQa1ygLIDYfP6/M1fNuvUNdjH66pKiGCwyzZR7z+hPFjesqjNrPPAWT8yZUv/6RpT1j4pHNkAkhNSbP94WDbpdbMKPaX2b9mSj33arjfB7mPLSeEf0Oqddt+LOuyf9smrbReRplMnmLKq3WEfl92fsu/ZpMcqTVntJvtPZm6z7UYxbvYhU3bOOW8Eyzded5JZZ/lVvwqWn1FxwKwzJHKfE6sXkZVZ35ep6rIc6l0tIp8BVgLXqep+YAbwQtY+2+OyRHwk5jhOX3KM6hpP/jf2Bj2Nt1wUWM6RoXPBlZjjOAM5gpFdEyJD7wBmZu16dFyWiCsxx3EGcgSVWEJk6EeBS0WkXESOAeYBLw52PPcTcxynD0LOlsfBjzWEyNCq+qqIPABsALqBq1TVnkSNcSXmOE5f8ujsOtTI0Kp6M3DzUNooqBLTIuiqDC8snfDJnWa9ppenB8trX7Ov9EvfvN2U/fH/sa1ImXL7DTtjGOOW/dHPzDoXP3mVKasttRfZ7vojux+1r5kiShrD9YoP2ta4rum2NW7rhXZs/mnP2xbD+vXh9jZ/0rZoVm0LW3/Bzh0AULHfljUtCFu233j/j80687aF8wMA1K21n7mauftNWfEjk0zZ/5QsDpaXT7Xbsp7h197+V7POkEiRx76PxBzHGYgrMcdx0kya1k66EnMcZyCuxBzHSS2aP+tkIXAl5jjOQHwk5jhOmnlHzYmJyEzgZ8AUIv28TFVvFZFJwP3AHCKHtUviRZwmWgKdtWG3ggOP2Os8a9vCV7T5OLut0//WNpNnbK8BGs6zZZed8Xyw/Pqbl5p1JtTYbhQHj7PH7NXbbBeLllmmCOkOX6sDN9vuEEV3j7cPmMChWba7xPht4Vj0E9cmLRKx/3J6yuzrmCSre3XoMefv++Stpuy639kuM12P2W4Uje+zky2UNIf/DI992E7ecGi27RaTF1KkxHJZdtRNtMp8PnAmcFUc9+d64ClVnQc8FX93HCft5LrkaIwoukGVmKo2qOrq+PMhYCNReIwLgbvj3e4GLjpCfXQcp4AIQ4piMeoMaU5MROYApwIrgCmq2hCLdhG9bjqO8w5grCioXMg5ioWIVAMPAdeq6sFsmaqag0sRWSoiK0VkZc/h1hF11nGcAvFOep0EEJFSIgV2j6o+HBfv7g2pEf+/J1RXVZf1BkwrrkzIFus4ztjhnaTERESIVp1vVNVbskSPAp+NP38WCMfLdRwnXQwtsuuok8uc2FnA5cA6EVkTl90IfAt4QES+CGwFLhnsQNIDJcYbZc0W2xS+4wPhbiZFczjwHls2YZMtmzjtoClb/m9nBcv3nddu1pnyn3Zkhrp19jnvep9dr+aM4KAXgJ4HJwfLmzuPMutU1Jsijn7Gjn7ROtXuoxUNZNLr9rXad0KFKetMcFUpabX/mnZeHnYtOfNv7Egm/GWjKeqYE46KAcnP8IRP2vfsmQXh3/8P/9fnzTrtnw57M2VW2q4cQ2KMKKhcGFSJqerzRAaLEB/Kb3ccxxkL+LIjx3FSzVh5VcwFV2KO4/RlDE3a54IrMcdxBuJKzHGctNLrsZ8WPGWb4zgDkIzmtA16HJE7RWSPiKzPKvuOiLwmImtF5JciUhuXzxGRNhFZE28/zKWvY2Yktv9423Q957/CZv6iDjub06QNdltd1XZbe5+vM2VV+8PtJblRTHjDTl+fRN0G+9Y0YbtLbPjmbcHyYx/+klnn+5f/xJR95+JLTdnBa8KRKgAm3B+OjFHaYrsAFCV4B7Qssl09irfZrhkLjw7nXt3eNc+s02W4qQCMa7XNdodm2Pes8rappmzJlk8Hy9vn2FFCup4fFyzXljz8Sed3Tuwu4PtEkXB6WQ7coKrdIvJt4Abgq7Fss6ouHEoDPhJzHGcA+XJ2VdVngX39yp5U1d6frBeIMn0PG1dijuMMJPdlR/W9a6PjzQ6uF+YLwK+zvh8jIi+LyO9E5JxcDjBmXicdxxk7DGFiv1FVw4kzB2tD5GtE8QrviYsagFmq2iQii4BHRGRB/4AT/fGRmOM4AznCC8BF5HPABcCn4yg4qGqHqjbFn1cBm4HjBzuWj8Qcx+nLEc52JCJLgL8BPqCqh7PKJwP7VLVHRI4F5gFvDXY8V2KO4/Qhn35iInIvcC7R3Nl24OtE1shyYHkUJIcXVPVK4P3AP4pIF5ABrlTVfcEDZ1FQJZYphXbDel1i50SgfVLY1Lz947aLxZSn7VPrqrQjIiRR3Bn+eWqfaLd18Hg7CUdxp/2ktE+03/RrN9k/k2euuThYnpSg4++eu8KUFc+x+3j0N1pMWWddOKJD40lh1wCAjO1RQNFO243i4o/+jyl74vZw5JFMvf0MtM6wz7n+FVNE22T7mPsW2ffs+J8WB8sTn4E3w89+sR0kZGhofrSYql4WKL7D2PchoriFQ8JHYo7jDCBNHvuuxBzH6YsvAHccJ+14PDHHcVKNKzHHcdKLkreJ/UJQUCVWchjqXwlbVRpPCVtoACobwot/J6yrNus0nWzfhKm/t39muqvsfliLl5uNOO6DMf5pO/tTebPd/64q2wpW/B/hgPlTVuw26xw4zV5QfnCOfT16ymzLq9XHltn2eSXlTDg8zZZZuQ8A2qaF+3HUGjsefvVO2ypYtcW2yIJ9PWo32dfx0OzKYHlSnojqreHjZezYBkPCJ/Ydx0k3rsQcx0kraQuK6ErMcZy+aG4BD8cKrsQcxxlIenSYKzHHcQbir5OO46QXBd5Jr5MiMpMoPvYUotNbpqq3ishNwBXA3njXG1X18eSDQU9Z2OQ9+1fNZrXdf1wTLB+313aVSFqMe3iybUJ/5SvhGPUA514RXig976t22vtdHxle5N2kxeEVTbZ7wK73heP91/y7HdEk8z7bxaJqp32Ny1rsBfgdNeFHq26tfV6Vu+zzqnvZXtncPj3sogBQ3Bl2Rdh9uu2LUJOwwH7v6RNM2fgdCUkCEugpCz+PdWvtOhX7w20lPTdDIj06LKeRWDdwnaquFpHxwCoRWR7Lvqeq/3Lkuuc4zmjwjnqdVNUGorCxqOohEdkIzDjSHXMcZ/RIk3VySOGpRWQOcCqwIi66Os4dd6eITMx35xzHGQVyDU09RvRczkpMRKqJApZdGwfuvx2YCywkGql916i3tDcTSldH68h77DjOESVydtWctrFATkpMREqJFNg9qvowgKruVtUeVc0APwbOCNVV1WWqulhVF5eW22sFHccZQ2Ry3MYAgyoxiYJg3wFsVNVbssqzl+R+Aljfv67jOOkkTSOxXKyTZwGXA+tEZE1cdiNwmYgsJHoz3gJ8abADSUZNs3xnnR1DXT7SFBY8YE/DacKZ7Vtkm8Ln/+DLdsWF4eIJE207R5LZvbN66JE7ADLldr2pK8IRNeTUBWaduud2mLKDp003ZUnUrw1PHZQ07DfraLOdXrDnPbNNWdMC211i8urw9Zi0wu7H25+YasoqGu0/3PImO5pJcZMd/aJ9dvg5PjAv7C4DULMpnJRCuvOgWMbQfFcu5GKdfJ7oNbk/yT5hjuOklPytnRSRO4nyS+5R1ZPisknA/cAcogHQJaq6P37ruxU4HzgMfE5VVw/WhifPdRxnIKq5bYNzF7CkX9n1wFOqOg94Kv4O8FGiXJPzgKVExsNBcSXmOE5f4uS5uWyDHkr1WaB/7sgLgbvjz3cDF2WV/0wjXgBq+829B3El5jjOQPI3EgsxJXaiB9hFtKQRIif6t7P2204OjvW+ANxxnIHkrp/qRWRl1vdlqros52ZUVWRki5xciTmOMwDJ5OwE1qiqi4d4+N0iMk1VG+LXxT1x+Q5gZtZ+R8dliRRUiUmPmsk2Wqfa5uQyw5Vi0oo9wXKAnrIppqzu5bB5GqBnQpkpK90ZjrRx6KTJZp2q9btM2c4v2REuGj9hiih+bZwpm/WbsCk/KfpCzVu2e0tSYozWOXaiFssNJDPednjOTLNdZjrq7Odj1v1vmzLLRaQ9ITnK1BW2e0tPqT0Dc3iafV+qOuyIH+11YReRlvd0mnXenBVuq+PbeZghUo60I+ujwGeBb8X//yqr/GoRuQ94H9Cc9dpp4iMxx3H6IOTPkVVE7gXOJXrt3A58nUh5PSAiXwS2ApfEuz9O5F6xicjF4vO5tOFKzHGcgeRJianqZYboQ4F9FbhqqG24EnMcZyBjZElRLrgScxynL0d+TiyvuBJzHGcAQ7BOjjquxBzH6ceIHFkLTkGVWHdFEftOCJvza96yzcmHZoXdHg4fX2e3Nc5OFLLnH+zIEvpEOCkJwMGLwmb+E26zXT26E9wGjrt1s92PyZNMWdtsO6HGodnhpBlTfmW3tfvCuaasyLj2APUv2pEgLLrrbDeErmr7cdx/vB2ponWK7aqy77xwgpHyjXY/al5NcLGYOryYeNJlu1jUrg4/Px01tptQpZEkp9EOBJI7iisxx3FSTnreJl2JOY4zkLES8DAXXIk5jjMQV2KO46QWVehJz/ukKzHHcQbiIzHHcVKNK7EwmXI4NCcs666yoxRU7QwPbVun2N1/+Wu3mbJzr7jClLVPtIfRf/fRh4Pl/7bxL8w6bZNtV49Zj9vRL4pabDN/El1V4fba3zszWA7JyUySXBsy1bb7RfFrW8P9mPces87BOQkJUH5vR9NIov7FsOvO5r+0I3c0nmG7xSS57hy1yu6jdNguRA1/Gk5MUt5sP4sV+8LHK+rOw2ugAinKAO4jMcdx+qGgPifmOE5aUXxi33GclONzYo7jpBpXYo7jpJd32AJwEakAngXK4/0fVNWvi8gxwH1AHbAKuFxVbRMMUFzZzaSFe4OyPW/Zi7khbLUqSmjt1Ju/bMomdtjp5rsq7UvyvR9dHCyvSrAiJS1sbzrVjntf/6K9YLiow26vdVrYejZpg71ovKPctgyPMxYagx1HH6Br8XHB8p6yBOveatsiW9JgLzbXMtuC2lMXzgPQPcE+r/oXw7kUALTUPueDx483ZZnjwgvzAerXh8+7KCEuf7NxvJ61+Yqxn545sVzOuAM4T1VPARYCS0TkTODbwPdU9ThgP/DFI9ZLx3EKy5HNO5lXBlVicTbeXgeY0nhT4Dzgwbg8O4uv4zipJl52lMs2BshpTkxEioleGY8DfgBsBg6oaq+XZE6Zeh3HSQEKmic/MRE5Abg/q+hY4O+BWuAKoHd+6UZVfXw4beSkxFS1B1goIrXALwHb7bofIrIUWApQdpQ9B+Q4zhgiTx77qvo60TRU72BoB5EO+TzRdNS/jLSNIc0CquoB4Bngj4BaEelVgmamXlVdpqqLVXVxSY09uek4zhjiyMyJfQjYrKrhNWnDZFAlJiKT4xEYIjIO+DCwkUiZ9ZrrsrP4Oo6TZlQj62Qu29C4FLg36/vVIrJWRO4UEXvB6iDk8jo5Dbg7HgoWAQ+o6mMisgG4T0S+CbwM3DHYgVShsztsop6YYBqe9HrYBH30tzeZdZ5/9iRTVtRlL1zOlNouANZCdGvRNUBFkyli8nO7TVnXdDvWf3GX/fBUNYR/HbctCbsaAEx73nY5qdxlt3Vgrr2Iet97w/0o6rKv1e7ppojiEvsZL37NPreJG8P9r1tt9yPJjaJnQsKi9057ZDJ+q+0+smdRuP/Tntxl1qkzXE5KWhK9nHIn91FWvYiszPq+TFWX9d9JRMqAPwNuiItuB75BZCT8BvBd4AvD6eqgSkxV1wKnBsrfAs4YTqOO44xlFO2xfdT60aiqi3PY76PAalXdDdD7P4CI/Bh4bMjdjMmDZ5zjOO8oekPx5LLlzmVkvUqKyLQs2SeA9cPtri87chxnIHkMxSMiVURz6V/KKv5nEVlIpDK39JMNCVdijuP0QQHNY1BEVW0lWp6YXXZ5vo7vSsxxnL6oB0V0HCflDGFif9QRLeAiThHZC/Q6utUDjQVr3Mb70RfvR1/S1o/Zqmonb8gBEflN3F4uNKrqkpG0N1IKqsT6NCyyMkfTrPfD++H9GCP9GIu4i4XjOKnGlZjjOKlmNJXYgKUJo4T3oy/ej754P8Y4ozYn5jiOkw/8ddJxnFQzKkpMRJaIyOsisklErh+NPsT92CIi60RkTb+V+Ee63TtFZI+IrM8qmyQiy0Xkzfj/YYcmGWE/bhKRHfE1WSMi5xegHzNF5BkR2SAir4rINXF5Qa9JQj8Kek1EpEJEXhSRV+J+/ENcfoyIrIj/bu6PI0M4qlrQjSh10WaiMLVlwCvA/EL3I+7LFqB+FNp9P3AasD6r7J+B6+PP1wPfHqV+3AT8dYGvxzTgtPjzeOANYH6hr0lCPwp6TQABquPPpcAK4EzgAeDSuPyHwF8V8j6N1W00RmJnAJtU9S2NUrzdB1w4Cv0YNVT1WWBfv+ILiRKuQIESrxj9KDiq2qCqq+PPh4iCbs6gwNckoR8FRSM8OU+OjIYSmwG8nfV9NJOMKPCkiKyKcwGMJlNUtSH+vAuYMop9yUvEzeEgInOI4tetYBSvSb9+QIGviYgUi8gaYA+wHE/OY/Jun9g/W1VPIwrYdpWIvH+0OwTRLzGRgh0NbgfmEiV3aCCKuFkQRKQaeAi4VlUPZssKeU0C/Sj4NVHVHlVdSJS/4gyGkJzn3cZoKLEdwMys72aSkSONqu6I/99DlIFlNCPV7u4NFBf/v2c0OqGqu+M/oAzwYwp0TUSklEhx3KOqD8fFBb8moX6M1jWJ2z7AEJPzvNsYDSX2EjAvtrSUESUPeLTQnRCRKhEZ3/sZ+FNGEF0yDzxKlHAFRjHxSj4jbg6hTSHK0bBRVW/JEhX0mlj9KPQ18eQ8Q2Q0rAnA+USWn83A10apD8cSWUZfAV4tZD+IwvQ2AF1EcxtfJAoa9xTwJvBbYNIo9ePnwDpgLZESmVaAfpxN9Kq4FlgTb+cX+pok9KOg1wQ4mSj5zloihfn3Wc/si8Am4BdAeaGe2bG8uce+4zip5t0+se84TspxJeY4TqpxJeY4TqpxJeY4TqpxJeY4TqpxJeY4TqpxJeY4TqpxJeY4Tqr5/3r5wYkE9WcGAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "imshow(a)\n", "colorbar()" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAATEAAAD5CAYAAABPqQIFAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAsHklEQVR4nO2de3xdZZX3fyv3a3NveqctbcGC0kKpIIIoMla8AI4vgg6COhZ9QWEGdRBmhMHXeUUEZUbFNw4d1JcBKhepiEJBoaAUeqH0CvZO06ZJk+baXM/Jmj/2jp+T5Fk7O8nJSXb5ffvZn57zrPPsZ51n77Oyn2etZz2iqiCEkKiSNt4KEELIaKARI4REGhoxQkikoREjhEQaGjFCSKShESOERJqM0VQWkWUA7gGQDuA/VfW7QZ/PzM7X7PxSpyyebdfrzXGHgZxUWGvWqe4pNmXdtTmmLFCPTHd5VosdptKbLrYsoPfTYwF6pNsyibvLg75XevfwzwcAsXxbltFhCHoD2uoNCPexuxGaFiA0RLFcu8rJZfZ9dbCn0JS1NeaZsuyiLlPW2eO+sdKb7WcM617saT6KWPuxgA4Zmg+9P18bjgZc+AQ2bO56WlWXjaa90TJiIyYi6QB+DOBCANUA1onIKlXdbtXJzi/FOy+83ilrmmf/Mtvnu39lqz74Q7POTQc+Zsr2/3CBKWs+0dajY4r7FzhztX3BO4sDzldu36T5tfavvbPEvkezm9yGoHmu3VbBQdt4ZDfbetSeaZ+zfLP7nOlddlsZHXZbvZn2d45nBfyhMOrVn2bXeebvvm/Kbqn5oCl76fHFpmz+st2mbMfhSmd58W/tvxJt093671txt1knLA1H43j16VmhPps+dWf5qBscJaMZTi4FsEtV96hqN4CHAFycHLUIIeOFAugN+W8iMBojNh3AgYT31X4ZISTCKBQ9Gg91DIWIzBSRP4rIdhHZJiLX++W3ichBEdnkHxcl1PmmiOwSkTdF5ENDtTGqObEwiMhyAMsBICuveKybI4QkgSQ+ZcUA3KiqG0WkEMAGEVnty36gqv3G7iKyEMDlAE4BMA3AsyKyQNW2mKN5EjsIYGbC+xl+WT9UtUpVl6jqkszsglE0RwhJBQpFXMMdQ55LtUZVN/qvWwHsQPCI7WIAD6lql6ruBbAL3tSVyWiM2DoA80VkjohkwbOeq0ZxPkLIBKEXGuoAUC4i6xOO5dY5RWQ2gMUAXvGLrhORzSKyQkRK/LJhT1ONeDipqjERuQ7A0/BCLFao6ragOrE84Mhit92c+rIdU3BgutvDd8GDX7fbmmSP1zOW2La7dJv91yWe5a6X2dxp1ukuCIiHCKB9sq3jn7/xQ1N24T98xVmuAVc6q80eOvTk2168SXvsc7bOcOvfsdiKvQCmPJ5lyoK8mg2n2n08c/UxZ3nrDNvzd9E/f82UWaENANA7zZZtOWALc7a54z16M+zvLGOYfEYBxBG6gXpVXTLUh0SkAMCjAG5Q1RYRuRfAt/3mvg3gLgCfH4m+o5oTU9WnADw1mnMQQiYeveGN2JCISCY8A/aAqj4GAKpamyD/GYAn/behpqkSYcQ+IaQfCqBHNdQxFCIiAO4DsENV704on5rwsUsBbPVfrwJwuYhki8gcAPMBvBrUxph7Jwkh0UKhwxlODsU5AK4EsEVENvllNwO4QkQWwbOZ+wBcAwCquk1EVgLYDs+zeW2QZxKgESOEDESBeJJsmKq+BPfiL3MaSlW/A+A7YdugESOE9MOL2I8ONGKEkAEI4kEr7icYKTViGceAynVuG99wiq3KzKfd4Rf5u+rNOo2LykxZ03z7AqX12H+DZj/W6Cw/Nm+SWaeg2g6/AOxsGrFcW8czf3yDKSvvcPdVyRv290qL2WOHWLbt+wnKfnHsne7vnbHf/s7NcwIWtjcG6Jhry3Zd4W6veHtA5pGAX4Ve3GDKtp7+kClb8m/XmbLsJve1CVr0XvqG+zq/1Tn6caA3sU8jRgiJKF6cGI0YISTC9PJJjBASVfgkRgiJNApBPEJx8DRihJBBcDhpoOlAV5Hbwlds6jHrtU92q5mXbi/8DUpZXLzT9tR1Fdv13rymyFle+We7TsN5dt71or22Ho0nmyJM+5O9WL5pvruvMgP2AWiaZ98GZdvtYOmm+Xb/a6+7TzIC0r8HeQWDFl7//UeeNWV13e6c+M9vfbdZJ5Zn6yi/t7MxL9j1v01ZZYN9rTONBfjNc+wOaVzo7vuedaM3PgpBt44sccF4wCcxQkg/vGBXDicJIRGGE/uEkMiiKogrn8QIIRGml09ihJCo4k3sR8c0REdTQkhK4MT+UBie/qPvsH3ohQfcbv72OfbC6+Kd7absRw/+xJQtyLRzr5/5z192lh9+nx2GMP1ZO7ThyCL7Rsmptx/nDy+13d8F1e72Cg7ZYRktJ9q3weGzbB0L99nfbdKmbGd5RoddJ7M9II/+O+3++P+/vNCUxY3bKj4jQPe9tizviH2tJ+03Rcg9bN+P1R9038e5dbYeRX9xl9cE5RsYBnHGiRFCogoj9gkhkaeX3klCSFTxFoDTiBFCIopC0BOhZUfRMbeEkJSgCsQ1LdQxFCIyU0T+KCLbRWSbiFzvl98pIm/4O4A/LiLFfvlsEekQkU3+8dOh2uCTGCFkAJLMYNcYgBtVdaOIFALYICKrAawG8E1VjYnIHQC+CeCf/Dq7VXVR2AZGZcREZB+AVgBxALGhtjNP7+xFyZvubeUPnl9g1muvcFv8su22Pzmt23aFH4jZoRknZtjZNLJa3dkGZj3pLAYAZBg57wEgq8kdhgAAM1a78/kDQKzYzlPfU+C+pEcW2SEs8Vw7w8Irn7zLlJ351A2mrPh1d3vdk+wfR+UVdoxCz0MnmLK43Y04Nst9H6z5uP29zl19gyk74RFb/54C+8mkdY4dupPeYYpMrGwryRgFKpC0ZUeqWgOgxn/dKiI7AExX1WcSPrYWwCdH2kYynsTer6r2jh2EkMgxjIn9chFZn/C+SlWrXB8UkdkAFgN4ZYDo8wAeTng/R0ReA9AC4J9V9cUgBTicJIT0QyHDSYpYP9QIDABEpADAowBuUNWWhPJb4A05H/CLagDMUtUGETkDwK9F5JTEOgMZrRFTAM+IiAL4f5YFJoREB2/LtuQ934hIJjwD9oCqPpZQfjWAjwK4QFUVAFS1C0CX/3qDiOwGsADA+oHn7WO0mr5XVQ+KyGQAq0XkDVVdM+ALLAewHABystyZUQkhE4nkbZ4rIgLgPgA7VPXuhPJlAL4B4H2q2p5QXgHgqKrGRWQugPkA9gS1MSojpqoH/f/rRORxAEsBrBnwmSoAVQAwqWD66Hf2JISMKYqkRuyfA+BKAFtEZJNfdjOAfweQDe/hBwDWquqXAJwH4HYR6QHQC+BLqno0qIERGzERyQeQ5nsc8gH8DYDbR3o+QsjEIVlPYqr6EuA82VPG5x+FN/QMzWiexCoBPO5b0QwA/62qvw+q0Jsp6JjiDg8o22qHIrz1MfcDXPFO+6+FZNm+5js/ebkp+37cDs3oPNvdXl6dHaKgAVkFpqyzfesNi4tNWX6NHQZS/Wl3Pxb/0b7UkzfafX9W/EZThlK73iXXPO8sf/nqxWadH39lpSn7RPo3TJmV5QQAMtrd1+wju+3zlR6zBwydpbas4FC3KWuaa8eBtJzmrvf0BfeYdT51x9ed5WLfiqFRlbfH2klV3QPgtCTqQgiZAHgT+9FZdsQQC0LIAJhjnxASYbyJfSZFJIREGKbiIYRElmFG7I87NGKEkEFwoxCDngJBzXvcXo/SrXa9DCOhQ2abOyMGAMRz7a924EI7o0BQuEHFevfyrd2fsrNiZDXarvWivbY/PLvZlrWcYGek2PC+/3CWL2r/qlknt8H2RE1aYGfTSHui1JT9fvV5znJZYFbB1df+oynLLLdDG7r+3o6FPLrPreO0P9h6/P4eO7Th/d/6B7tiAPEc+8nm1vescpYve+ErZp3SLqM/khBOrgr09NKIEUIiijecpBEjhESYZEXspwIaMUJIPxhiQQiJOBxOEkIiThJz7I85KTVi6V1AwT535zTPs90qWc3uOt3fbjbrdFdNMWUZATnN2ytsT11WS5azvGS7fb7edPt7xbPsG6U735aV7OwyZbt63PqvW/ZDs84nfmt73LIeLjFlde8JWPie5v7eJ93X7iwHgLYTbK9xR4XdH9mPl5uyyZ1uPepPs8/3gVsCPJABSwr3fML2Gn/8Pa+asv+4+2+d5aUBi7mbFrq/VzwwBUM4PO8k104SQiIKg10JIZGHw0lCSGShd5IQEnnonSSERBZVQYxGjBASZTicNNA0IJbnluUvbrDrPVXmLK9/Zrrd2BxbNO35VlPWOtd288dy3G7nor12Iv3WmfYC8GPT7L92BQds/3rd9XaMyPVfdy8aPnq5vVg+fZ59G1RsssM5MsvtXP8bznFvQfqxJ+xFzTWX2jnqc7blmrKSnXb/7/24u/8r19r9W/tu+wdcttmWFc2yQ35evudMU9Z4gbuPJ692h/QAwNxH3fdAfePok+wnc05MRGYC+AW8PTkU3g7h94hIKbxdv2cD2AfgMlVt9Ld4uwfARQDaAVytqhuD2ojOMyMhJGX0qoQ6QhADcKOqLgRwFoBrRWQhgJsAPKeq8wE8578HgA/D22tyPrz9au8dqgEaMUJIP/rixJJhxFS1pu9JSlVbAewAMB3AxQB+7n/s5wAu8V9fDOAX6rEWQLGITA1qg3NihJBBjEWcmIjMBrAYwCsAKlW1xhcdhjfcBDwDdyChWrVfVgMDGjFCSD9UgVj4pIjlIrI+4X2Vqg6aEBWRAnib4t6gqi3+frV+e6oiMuJ0jjRihJBBDGNiv15VlwR9QEQy4RmwB1T1Mb+4VkSmqmqNP1ys88sPApiZUH2GX2bCOTFCSD+SOSfmexvvA7BDVe9OEK0CcJX/+ioATySUf1Y8zgLQnDDsdDLkk5iIrADwUQB1qnqqX+Z0jw51rvRuoGive8v5zkZ3GAUArL39R87yk391rVnnFx//iSm7/fmrTNnhi2w3f9nzbnd92Wbbxa9pOaas7R12+EJXse1en1dku/KPTCp2lhc/UmDWgdpu+bYZth4VK+3b57J/+Yy7qfm2GoWv2mEUaT32aOPIu+x6cx53Z81I77T3Upj0F/c9CgCxYvt61j1v7zmQecw+56yV7tCd7CMBYTFH25zlErPbGQ6avDixcwBcCWCLiGzyy24G8F0AK0XkCwD2A7jMlz0FL7xiF7wQi88N1UCY4eT9AH4EL9ajjz736HdF5Cb//T+FOBchJAIka2JfVV8CzJNd4Pi8ArCfThwMOZxU1TUABm4nY7lHCSERRzWpcWJjzkgn9i33KCEk8gjib6ct24Zyj4rIcniRt8jKKx5tc4SQFJDEObExZ6TmtrYvinaAe3QQqlqlqktUdUlmdsDkMiFkQtC3djIqw8mRGjHLPUoIiTrqzYuFOSYCYUIsHgRwPrzI3GoAt8J2jwai6UBXkdtu5jXYruHFr3zWWV683f5LcMvT15iy2Czbdp/8XTtSpHuyu7x+UaFZp6fQ1jFnr53h4qKL15qyP9x3linLMPaxbzjV1qNikx1i0Vlq99WRc+x6ubXuPukos2+5wgP2PdA03964Im7vz4GeIrdw1+dsPcr+bJ8wLSCCoTfg15R5LCCMZZq7YmexnVElu9kdVtJbl5z49eMqPbWqXmGIBrlHCSHRR99uE/uEkOOPiTJUDAONGCFkEFHyTtKIEUL64U3a04gRQiLMRAmfCAONGCFkEJwTM5A4kN3idjV3F9jekLxVk5zlFS8cMus0nz7FlB2rDHLX25kIrF2smgMyM5Rus++Gjgq73p/uXmrKemaaIhQecscA5NfaeuRWuzMiAEB6lx0+Mmm/fc26yt3ZL5oWmFWQW2tfl7It9qYkRxbZIRHZDe5MIRUv2uELR0+1+2rGc3aMRUe5rUdmq61/vRExU1BhZ7GI/bbIWR7PHP0TlELQS+8kISTKROhBjEaMEDIATuwTQiJPhB7FaMQIIYPgkxghJLIogN5eGjFCSFRRAHwSC8AYaxf9xXYnH7jQ7eYv2bvfrJN+qp1sNrfBziiQV2tvFNI8x511YvJ6+3w5R+0NKcrWtZiyYyeWmLKMTnvCov5U9yWdtC8gi8LUYlNWcNDWPyMWMHFieOhLdth1Cg7ZG6e0T7Y3LDnh10dMWc0H3HEsufV2f6R32D/gY1Ptn8xz199pyj5y9Gum7KSvbHSWN3zmdLOOGJEeI9+9sT9RihOLTjAIISR1aMhjCERkhYjUicjWhLKHRWSTf+zr2wVJRGaLSEeC7KdhVOVwkhAyAEnmxP79GLBbmqp+6q8tidwFIHEPwt2qumg4DdCIEUIGk7Rhqa4Rkdkumb+x7mUAPjCaNjicJIT0RwHtlVDHKDkXQK2q7kwomyMir4nICyJybpiT8EmMEOIgtIEqF5H1Ce+rVLUqZN0rADyY8L4GwCxVbRCRMwD8WkROUVXbA4YUG7G07l7k1XQ6Zb1Z9uLfaX9y15EzTjHrFGysNmXt75xuynryA3Kvv+7uy/R6u4+1zfa66nQjaT+A5rm2HiVv2h7UOf91wFl+8BOzzToZ7fbYIaeu3ZSlNduyeLnbo3zkXfaOV0W7bU9odrN9f/SU2Yu5ey5odpZnPuFOKgAAaT32D7h+qa3jhd//uimL25caRy93eyELDtlt9RgJE8R2ug6P8MPJelVdMtzTi0gGgE8AOOOvTap2AejyX28Qkd0AFgBY7zyJD4eThJDBJMk7GcAHAbyhqn992hCRChFJ91/PBTAfwJ6hTkQjRgjpT1+wa5hjCPzd0l4GcJKIVPs7pAHA5eg/lASA8wBs9kMuHgHwJVU9OlQbnBMjhAwiWcGu1m5pqnq1o+xRAI8Otw0aMULIYLh2khASZZK1fCkV0IgRQvoz+kn7lDKkERORFQA+CqBOVU/1y24D8EUAfStvb1bVp4Y6V3clsO+r7t5Jf8N2k8/6vTsHfNM77PzvOVPyTFnu/lZT1j7Hdr3Hc9051KXI1h0Bso6Ztv7TH7MXt7cttkNEGs6f5SyfvMEO9QgKb2mfaeufF+DO7yzPcZ/vBDtsYOeVuaZszmN2ve5ie3F43hPuc1Y8b4fgpJ1r9+/sJ+x7p6fE1j97f4MpO7bQHX+Rt+kts86+z53oLI+/aFYZBuEm7ScKYbyT9wNY5ij/gaou8o8hDRghJEKMfYhF0hjySSxo7RMh5DglWUGzKWA0cWLXichmP9WGnfyKEBItkhgnlgpGasTuBXAigEXw1jvdZX1QRJaLyHoRWR9vsedlCCETB9Fwx0RgREZMVWtVNa6qvQB+BsDc6VVVq1R1iaouSZ8UMAFOCJk4RGhObERGTESmJry9FMBW67OEEDKWhAmxeBDA+fBSblQDuBXA+SKyCJ4t3gfgmjCNZVfHMf+mJqes80R7C/i2mW7Xddmze806DR+cY59vqj2FN/nPtiscMXdi81iFHSrRU2h/r6a5tqy93B0qAQBH3muHG+RUu/8ulbxu56/vLbfDUcRuCtJjJHoHULClxlleUT7DrJN3xG6s5QQ7jKLwLTurR3q3uz9aTp9m1oll23M9rXfY/dj5SLEpO/rJKabspJ+6M23EZtmpL05YsctZfqje1m84TJShYhjCeCdda5/uGwNdCCETAQWXHRFCIs7x9CRGCHn7cVwNJwkhb0NoxAghkYZGjBASVSZSIGsYUmvEBNAMd8aE3gzbG9KT53aTt58206xTUG273Y+c5s6wAACxEjvcIGObO6QjNqvYrHP0JDuMYtrzjaZMM+3MEhUv2BuT7LnaHcJQv9QOK+kNuAtKt3eYsrRjtqz2QrceaT12W1mN9jU7/G5bya7ibFNWcHD4iwB/feudpuyyb3zNlBU22l/u/9z0S1N2y5bPO8u7i+zfxNQ/ue8PbbPvt2FB7yQhJMrwSYwQEm1oxAghkSVic2Lcso0QMpgkLQD3U3XVicjWhLLbROSgiGzyj4sSZN8UkV0i8qaIfCiMqjRihJBBSG+4IwT3I2RmaBFZCG8/ylP8Oj/p20w3CBoxQsiYoaprAAy5Aa7PxQAeUtUuVd0LYBcC0nz1kdo5sbQ0aI7bHZ7VbLvX6053ZzAo3mVnPYjl2gY8v3ZkuXc7l8xzlncX2d1YsbnTlKU120kiNdM+Z+8kOwykp8D93co3uDMlAIC025kPWk+tMGXHZkw1ZaU73OEXae12GMKx2QWmLG+R/TtoPFhkytJ63P0Y9BRxye1fN2VBYRTtlXZ4w813ucMoACCvwa1MyZt2Wy1z3GFC8e1Jei4JPydWLiLrE95XqWpViHrXichnAawHcKOqNgKYDmBtwmeq/bJA+CRGCOlPyKyu/uR/fV/SU/8IY8BCZ4YOA40YIWQwY5jZNSAz9EEAiRHsM/yyQGjECCGDGUMjFpAZehWAy0UkW0TmAJgP4NWhzsc4MUJIPwShPY9Dn2sYmaFVdZuIrASwHUAMwLWqaqcP9qERI4T0J4nBrsPNDK2q3wHwneG0kVIjphlp6DHyuYvavTZpn/vPwoELbS/d1Jdsj1tere0JbZpnb0Xf+A53edAFj9mp3JGdZ3vV4jttT93kjfafydKt7oW7mmHPHPSW2W0FMenNVlNWe7b7u015yfYyFv7Z3jOh/tOVpqxoh30bT33B3V7hT46Ydda/5vZCA0BazPZAdpXYi6YLqu0Hiu4C97XJarLPV/pKnbM841jACvvhEKGIfT6JEUIGQyNGCIkyUVo7SSNGCBkMjRghJLJo8ryTqYBGjBAyGD6JEUKizHE1JyYiMwH8AkAlPPtcpar3iEgpgIcBzIYXsHaZv4jTZP6cI/jdAz9zyt5V9RWzXucsIyQize7prm32V6tbYocb5CywF0oX/dYdNhDkWk/vsPP5ZzfZsq5SU4TG+fbi9nQjsuToqZPsEwZQ/wE7VCV7rx0i0lXiHo+0zbZz/We22F+6u9oe32QFJGs5YuwtUPdLW4/igHUsQfsR/OYfvmfKPr39s6as7cUpzvJYrjvxAQC0fM7d9523DJm5JhwRMmJhlh3F4K0yXwjgLADX+nl/bgLwnKrOB/Cc/54QEnXCLjmaIIZuSCOmqjWqutF/3QpgB7z0GBcD+Ln/sZ8DuGSMdCSEpBDBsLJYjDvDmhMTkdkAFgN4BUClqtb4osPwhpuEkOOAiWKgwhA6i4WIFAB4FMANqtpv40NVNR8uRWS5iKwXkfVHGoZcy0kImQgcT8NJABCRTHgG7AFVfcwvru1LqeH/71zMpapVfQnTKsqSNOlICBlbjicjJiICb9X5DlW9O0G0CsBV/uurADyRfPUIISlneJldx50wc2LnALgSwBYR2eSX3QzguwBWisgXAOwHcNlQJ9p+eDKWfM8dStF5ir36Pnev29XcdbI7jzsAdBXZ9rl0mylCwTN2Fou3jA2kMuxU+civtq+0BjyYTllr7x+w/2P2d8vf6z5p83y7rSA9ls7bZ8pe332yKZtspLL73rfvNet8fu3VpqxylR2OktFp99WUr+12lu944iSzjgTMerSeaWdA+eS/2Ln5O8vsMJxsYxuGlhPte6f4D+7+ONya8hz7486QRkxVX4LnsHBxQXLVIYRMBLjsiBASaSbKUDEMNGKEkP5MoEn7MNCIEUIGQyNGCIkqfRH7UYFbthFCBiG9GuoY8jwiK0SkTkS2JpTdKSJviMhmEXlcRIr98tki0iEim/zjp2F0Te2TmAC9xsL8yhl2AoymwxXO8rTqILe73cFZbbbr5fDSbFNWttldr/Zc28Wf0W5vLNE+3dbjxW/92JSdUXWDKWtb4A5VmfmU7eI/9F77b9mWJ+0wiqx2U4RYjru9f/y3L5t1Cj9m3wONJ9mbwkBt/as3uTf9KGuw749Yrt1XGUfs69mTb9frLrbby25010vvsM93bIZb1murF57kzondD+BH8DLh9LEawDdVNSYidwD4JoB/8mW7VXXRcBrgkxghZBDJCnZV1TUAjg4oe0ZV+/7yr4W30/eIoREjhAwm/LKj8r610f6xfJgtfR7A7xLezxGR10TkBRE5N8wJOLFPCBnEMCb261V1yYjaELkFXr7CB/yiGgCzVLVBRM4A8GsROWVgwomB8EmMEDKYMV4ALiJXA/gogM/4WXCgql2q2uC/3gBgN4AFQ52LT2KEkP6M8W5HIrIMwDcAvE9V2xPKKwAcVdW4iMwFMB/AnqHORyNGCOlHMuPERORBAOfDmzurBnArPG9kNoDVXpIcrFXVLwE4D8DtItIDoBfAl1T1qPPECaTUiGUU9WDyh6qdsn1bppn10he400ScfcI+s85LuQtNWcl2O23DBX+7zpS9+F9nOsv/sOz7Zp0n204xZff+90dM2bv//QZTltNmipDV7PaxH1lk18mpt2XtU+0/yZN227MR7VPdIQBBGSIKf1VsytKz7V/VBde9bMqeXPkeZ3l2s61Iercd2oC37O8ctyN+ELcjd1C63Z0Z4+H/vMesc/ajNzrLNVm/aE2OFVPVKxzF9xmffRRe3sJhwScxQsggohSxTyNGCOkPF4ATQqIO84kRQiINjRghJLookjaxnwpSasR6WjNx6EX3MqmCJrted2OBs/z1l08160yrtr1PpV/dZ8p+s2GRKSvId5e3BySpv/dB2wMZKwjIv99pe8hO/Ig7bzwAnFbs9v6uqXMvhAaA6tdsz3BmS8Ci5kJThJI33f1/1x0/Met88Sfu/RcAIK/O7qvHf3e2KZuxrstZ3jzXyEQAoPCgvaA/68paU9bUYbsns14psevNc+ty0W1fM+sUG7dcrb3txLDgxD4hJNrQiBFCokrUkiLSiBFC+qPhEh5OFGjECCGDiY4NoxEjhAyGw0lCSHRRAMfTcFJEZsLLj10J7+tVqeo9InIbgC8COOJ/9GZVfSroXCUlbbj0kpecsl+tPsesZ7n5O8vtjm49w51rHgDqDlWasllP2iEF6R1ud/1n7nYvxgWA7ll21GDeIXsxcdB3O1w1x5TVxdyy9kq7rbKGgMjGz9irwxs3uPc+AICMDrf+t+y51KyTf9jWI7feDpnJarW/W90SY+X10mazTk6eHafQe7997/T+r1ZTFhCFg5az3e1l78g163TMdIeBxH+XJOMTHRsW6kksBuBGVd0oIoUANojIal/2A1W1UzgQQiLJcTWcVNUaeGljoaqtIrIDwPSxVowQMn5EyTs5rPTUIjIbwGIAr/hF1/l7x60QETskmRASHcKmpp4gdi60ERORAngJy27wE/ffC+BEAIvgPandZdRb3rcTSnuje06JEDJx8IJdNdQxEQhlxEQkE54Be0BVHwMAVa1V1biq9gL4GYClrrqqWqWqS1R1SV5JQHpLQsjEoTfkMQEY0oiJlwT7PgA7VPXuhPKpCR+7FMDWgXUJIdEkSk9iYbyT5wC4EsAWEdnkl90M4AoRWQRvZLwPwDVDnaixsQCPP/5ep6xyu+1Cb5nl9k9L3A6HKCqz3d3xVeWmrP7UgC3sO93dld5pX8y0gHzt0/5ou/kP3Wr/mTvWZE8/dlS662U1mlUQy7P/lmU+OdmUZbmTiwAAWmYZffWg7RPqmGn3Ve377Psjf7cdv5DR7i5v63TvRQAAkx6203P0Zts6ytoiUxafZN8jk152h1Kkd9l1spe6t2I8kpWEx6MJNN8VhjDeyZfgDZMHEhgTRgiJKslbOykiK+DtL1mnqqf6ZaUAHgYwG94D0GWq2uiP+u4BcBGAdgBXq+rGodrg5rmEkMGohjuG5n4AywaU3QTgOVWdD+A5/z0AfBjeXpPzASyH5zwcEhoxQkh//M1zwxxDnkp1DYCBe0deDODn/uufA7gkofwX6rEWQPGAuXcnNGKEkMEk70nMRaUfRA8Ah+EtaQS8IPoDCZ+rRojAei4AJ4QMJrx9KheR9Qnvq1S1KnQzqioyukVONGKEkEFIb2gvZ72qLhnm6WtFZKqq1vjDxTq//CCAmQmfm+GXBZJSIzajrAH/98pfOGXfeORKs17RTrehbjzFNuDd2+wwisp6+wIVVtuyg+e6uyv/kO12n3vWflNWd2CWKUt7zj5n22xbx7LX3fXSu+w6sVy7rXhAfHJWky0TIyIilm+3lVdrX8/MNjskon2aXa/kTff3zjhmZ4josiMlEA8IsUjvtOtNe93OqnLgKndGiuwteWadpZXu3/bBjG5bibAoxjqQdRWAqwB81///iYTy60TkIQDvBtCcMOw04ZMYIaQfguQFsorIgwDOhzfsrAZwKzzjtVJEvgBgP4DL/I8/BS+8Yhe8EIvPhWmDRowQMpgkGTFVvcIQXeD4rAK4drht0IgRQgYzQZYUhYFGjBDSn7GfE0sqNGKEkEEMwzs57tCIEUIGMKpA1pSTUiN26HAZvn2nO5Si+GJ7Q4r2ZveGFOWv2R0dCwgN6CoK2lgiIKNDm9u9vvLGO806n/re101Z+/QA/fPtv4RBmTHaK92ywrfstvIPu138ADDvX7ebshf+8C5TNm2N+5zVH7QzTkx52daxY7Jdr3CvKUJHmft6HpsRFFZi6zHtJTtUIp4VsPFLqf1Tq/iNW1Z/mq3Hi3tOdJa3diUhZ5+CRowQEnGiM5qkESOEDGaiJDwMA40YIWQwNGKEkMiiCsSjM56kESOEDIZPYoSQSEMj5kYApBne/K5n3WEUAJDT4u7QhoBNPQoO2LL8w/amE+uu/KEpe/eKf3SWX/pTO4wiJxaQYWFHQIhFQGaJNNvLj65Sd7kVagAAaXFbj2dfX2jKcjoCMm1Md99aRX8xq+DoybaOJW/a1yzvS4fsc66c4SzvrLTDSia/bIdzNC6ws2lkNdv9mNVqD886S9zfu+QNswp63nJnuEhvS0KeUwUQoR3A+SRGCBmAAso5MUJIVFFwYp8QEnE4J0YIiTQ0YoSQ6HKcLQAXkRwAawBk+59/RFVvFZE5AB4CUAZgA4ArVTUwwXcsB2hc6O6cuY8Z+80D2PO37nzoubW2JyYoN3xvhu1VW3q/2wMJAFnGAnANcAj1fLjJlHWsKzFlBdUBi6EvaTZlva8WO8vT7a8cSOlUu62Ot+x9DEquqHaWV6+Z6SwHgMpXbbdry2z7VpV7p5myE766yy14YJ5Zp/HD9r2YuTnflLXOtq9ZTp3t8bT2bwzaAyjrPQ1uwdO21zU0CiBCqXjC+GO7AHxAVU8DsAjAMhE5C8AdAH6gqvMANAL4wphpSQhJLWO772RSGdKI+bvxtvlvM/1DAXwAwCN+eeIuvoSQSOMvOwpzTABCzYmJSDq8IeM8AD8GsBtAk6r2PbuG2qmXEBIBFNAkxYmJyEkAHk4omgvgWwCKAXwRwBG//GZVfWokbYQyYqoaB7BIRIoBPA7g5LANiMhyAMsBIL3EngMihEwgkhSxr6pvwpuG6nsYOgjPhnwO3nTU90fbxrDWKKhqE4A/AjgbQLGI9BlBc6deVa1S1SWquiQ9354UJYRMIMZmTuwCALtV1d5RegQMacREpMJ/AoOI5AK4EMAOeMbsk/7HEnfxJYREGVXPOxnmGB6XA3gw4f11IrJZRFaIyIiHaaJDWFMReRe8ift0eEZvpareLiJz4YVYlAJ4DcDfqWpX0LkKSmboovdf75T1BsQAtE1zu6fbZtud+DfnbTJl63+02JT15Nt6qOElb5ln65FfPbIwkJ4C+7oU7bTrZbe4dWmZbbv4s4/abaUHLDYPCkUo/a07LCa/JuCEAQTlr2+dac+KdE9yX8+eSfZ3lnhA8oCA0JeuMrte0OJw6547/9PrzDq/+8MSZ3n1PT9A14GA7AchKEov17PzPxbqs0+33r8fQOIGGVWqWjXwcyKSBeAQgFNUtVZEKv16CuDbAKaq6udHou+Qc2KquhnAoF+9qu4BsHQkjRJCJjIKjdtZQwZQr6pui9qfDwPYqKq1AND3PwCIyM8APDlsNX2SkLeDEHJc0ZeKJ8wRniuQMJQUkakJsksBbB2pulx2RAgZTBJT8YhIPry59GsSir8nIovgmcx9A2TDgkaMENIPBaBJTIqoqsfgLU9MLHNvQDsCaMQIIf1RJkUkhEScYUzsjztDhlgktTGRIwD6At3K0d81O15Qj/5Qj/5ETY8TVNXesCIEIvJ7v70w1KvqstG0N1pSasT6NSyyPqRrlnpQD+oxQfSYiDDEghASaWjECCGRZjyN2KClCeME9egP9egP9ZjgjNucGCGEJAMOJwkhkWZcjJiILBORN0Vkl4jcNB46+HrsE5EtIrJJRNansN0VIlInIlsTykpFZLWI7PT/H/MMkoYet4nIQb9PNonIRSnQY6aI/FFEtovINhG53i9PaZ8E6JHSPhGRHBF5VURe9/X4V798joi84v9uHvYzQxBVTekBL6XPbnhparMAvA5gYar18HXZB6B8HNo9D8DpALYmlH0PwE3+65sA3DFOetwG4Gsp7o+pAE73XxcC+AuAhanukwA9UtonAARAgf86E8ArAM4CsBLA5X75TwF8OZXXaaIe4/EkthTALlXdo94Wbw8BuHgc9Bg3VHUNgKMDii+Gl7cNSNHGK4YeKUdVa1R1o/+6FV7SzelIcZ8E6JFS1IOb84RkPIzYdAAHEt6P5yYjCuAZEdng7wUwnlSqao3/+jCAynHUJSkZN0eCiMyGl7/uFYxjnwzQA0hxn4hIuohsAlAHYDW4OY/J231i/72qejq8hG3Xish5460Q4P0lhmdgx4N7AZwIb3OHGgB3paphESkA8CiAG1S1JVGWyj5x6JHyPlHVuKougrd/xVIMY3OetxvjYcQOAkjcBtrcZGSsUdWD/v918HZgGc9MtbV9ieL8/+vGQwlVrfV/QL0AfoYU9YmIZMIzHA+o6mN+ccr7xKXHePWJ33YThrk5z9uN8TBi6wDM9z0tWfA2D1iVaiVEJF9ECvteA/gbjCK7ZBJYBW/DFWAcN15JZsbNYbQpAO4DsENV704QpbRPLD1S3SfcnGeYjIc3AcBF8Dw/uwHcMk46zIXnGX0dwLZU6gEvTW8NgB54cxtfgJc07jkAOwE8C6B0nPT4JYAtADbDMyJTU6DHe+ENFTcD2OQfF6W6TwL0SGmfAHgXvM13NsMzmN9KuGdfBbALwK8AZKfqnp3IByP2CSGR5u0+sU8IiTg0YoSQSEMjRgiJNDRihJBIQyNGCIk0NGKEkEhDI0YIiTQ0YoSQSPM/YGb8ur640s0AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "imshow(b)\n", "colorbar()" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "# like moving window with the 0th dimension the IW no.\n", "a = a[np.newaxis,:,:]\n", "b = b[np.newaxis,:,:]" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "# %%timeit \n", "c1 = fft_correlate_images(a,b,'circular',normalized_correlation=False)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "# %%timeit\n", "c2 = fft_correlate_images(a,b,'linear',normalized_correlation=False)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "# %%timeit\n", "c3 = fft_correlate_images(a,b,'circular',normalized_correlation=True)" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "# %%timeit\n", "c4 = fft_correlate_images(a,b,'linear',normalized_correlation=True)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-12.02031891 -10.95825492] 0.05293600647141261\n", "[ 0.03563358 -0.00205613] 22.98857744116232\n", "[-12.0226627 -10.95144842] 0.043785720871595046\n", "[-12.01144001 -10.94500473] 0.048564728262190826\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAzkAAACyCAYAAABhjQGfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABSqUlEQVR4nO29e/wdVXnv//4Q+CZcwkWjqIQEDkS5qdxBrcpBsUDPAfsTK7QKeKiWAhWr9iXWvrDV9hyt1R49opYKclFBC1TTgmKKWLQKJiK3EC7hFr5cjAgJQSBfEp7fHzM7TnZmZs/svea2v8/79dqv73fPrJm1Zu/9zFqfeZ71LJkZjuM4juM4juM448JmTTfAcRzHcRzHcRwnJC5yHMdxHMdxHMcZK1zkOI7jOI7jOI4zVrjIcRzHcRzHcRxnrHCR4ziO4ziO4zjOWOEix3Ecx3Ecx3GcscJFjtMoks6XtFLSbQXK/qOkm+LXXZJW1dBEx2kNbi+O4ziOUwz5OjlOk0h6A/AUcJGZ7VPiuD8D9jOz/1VZ4xynZbi9OI7jOE4x3JPjNIqZXQc8ntwmaTdJ35P0c0k/krRHyqEnAJfU0kjHaQluL47jOI5TjM2bboDjpHAucKqZ3S3pEOCLwOG9nZLmA7sCP2iofY7TJtxeHMdxHKcPFzlOq5C0DfBa4F8k9TbP7Ct2PHCZma2vs22O0zbcXhzHcRwnHRc5TtvYDFhlZvvmlDkeOL2e5jhOq3F7cRzHcZwUfE6O0yrM7EngPklvB1DEq3v74/kGOwA/baiJjtMa3F4cx3EcJ52RRI6kIyXdKWm5pLNCNcqZPki6hGgA9gpJk5JOAf4IOEXSzcBS4NjEIccDl1rJtICSdpZ0raTbJS2VdGZKmT0k/VTSWkkfKnNswTa4vTgj0TJ7kaTPx7/nWyTtP/yVpbbB7cXpBIN+q1np3CXtG/c5S2MbekdVbXCctlDAXj4Q9z23SLomnlea3L9t3P99YWBdw6aQljQDuAs4ApgEFgMnmNntQ53QcSpE0kuBl5rZjZJmAz8H3pr8vUp6MTAfeCvwhJn9Q9FjC9Tv9uJ0hoL2cjTwZ8DRwCHA58zskED1u704naDsbzWZzl3SywGLk4a8jMjO9jSzVVW2wXGaoshvVdJ/B24ws6cl/SlwmJm9I7H/c8CLgMfN7Iy8+kbx5BwMLDeze81sCriUjZ8gOk5rMLNHzOzG+P81wDJgp74yK81sMfBc2WML4PbidIaCv/ljidbrMTO7Htg+FkchcHtxukLZ3+qGdO5mdpeZ3R3//zCwkmjwVnUbHKcpBv5WzexaM3s6fns9MLe3T9IBwI7A94tUNorI2Ql4MPF+kvIDP8epHUm7APsBN9R4rNuL00lyfvNV/qbdXpyuUPi3mpfOXdLBwARwT5VtcJyGKftbPQX4LoCkzYDPAB/KKb8RlWdXk/Re4L0Am82YOGCrrYd5SOFMB5568qHHzCzzB/K6w2bZqsefT913+63PLQWeTWw618zO7S8Xp9y9HHh/PGm7MKMcW6KODfYyg80P2HrzHaqoBps1Ucl501g/S4MLNciMZ4cL2Q2Fnp0qfcyT636VaysArz9slj2RYS9La7CXqqnLVvLot6NhvsvpQN79psxnNuznXZe9FCA1nXvsAb0YOMnM0hsxInn2Mkp/UOQ7KHN+PTs1sHyvTxn23p3XJzXdH1RN7/vK+oz17FQhe/mdjPFYQFtB0juBA4E3xptOA64ys0mp2LhiFJHzELBz4v3ceNtGxBd3LsDs7ebaAa/5sxGqdMaZ/7z6rAfy9q96/Hm+8e87pu7bd/7ks2Z2YN7xkrYgGrB93cyuKNO2UY6NKW0v223xYnvtnLcPUVU2U3vMHVwoIKt361+ypb1sd8/aRuufuGOycNnvPfrFXFsBeOLx57n8yjmp+/aY90gIeyn0mx6SgefexFZeeFygqosztefOm2ybWPZgSkkn+Vlt9BltXf4cG44veOz3fvmlKu2ljB1sks5d0rbAlcBH47DPYShvLwP6lqJ9xaD7VhV9zurdZm50v+5/P+jYIjTdH1RF7/vK+l4m7pgs1L+sevx5vnXlpjpo73kPD+pbCtmLpDcDHwXeaGa9L+M1wOslnQZsA0xIesrMMhNtjCJyFgMLJO0aN/B44A9HOJ/jVIYi2X8esMzMPlvXsQkatxcXOPn02ttU5za1x9xSQqdKCv7mFwJnSLqUKPHAajN7JFATgtvL1J471yJANhmIV1xPEbLaUldbQ9XRUgFZ6LeqlHTukiaAfyWa23ZZ1W0ow8QdkwP7jKbuV/0Cp/d30L27TJ/UdH9QBS3pXwb+ViXtB/wTcKSZrextN7M/SpQ5GTgwT+DACCLHzNZJOgO4GpgBnG9mS4c9n+NUzOuAdwG3Srop3vaXwDwAM/uypJcAS4BtgeclvR/YC3hV2rFmdlXRypu2lzoFTtfETT9lngqOMQPtBbiKKLPacuBp4N2hKg9tL73BfF1CJ1ln2foyPR4jtiV5rn6BVOfnMm5k/VYlfRxYYmYL46Jp6dz/AHgD8MJ40AZwspndFKINQ19UTJbQ6R8s95fJ8hZk9Q3D3m/7z1fFvTutzW3oH7oqwgray6eJPDX/EoelrTCzY4apb6Q5OfEgr/BAz3Gawsx+DOQGcZrZoySyeCQYeGzBNjRiLy5wytNUB9L7rpp+4lbQXoy+0JvAbRiL/qWMgKhSfHRByHShjWmk/VbN7Oy+93+dctzXgK9V1YYQJAVL2n0prX8J2ecU8SglyRI6Ifumrj8Iy/s8p/aYC49WW/8gezGzNxc4xwXABYPKjbQYqOM47cYFzmg0dU11hxaOOxPLHtzw6hpTe+6cOfen/5VHcn+bP4es63WaJevBy8Qdk4UfygwjDLLuhXn35izvS+81Km0QOCE/y3Gm8uxqjuM0w7gJnDXzN3YszH6gniw4TXp1mvboOKMzrBcnbX8ZITMo7G1i2YNDhdTltXMU8eQhdN2mqMelrvtonrelt73rD+Zc6AzGPTkj8MQr6kvD6zhlGHeBk7WtSproEKdbh9Q1QoaSFS1XtGwRD09ZcVNEiIXEPTrtZWqPuZlzcepmzXyl9gerd5uZe9+uWnD16k97OfXgImdIegLHhY7TNuoaGNdxs87qvIruD01TQsfFTnsJJSKaokgby4iNuoWOh7c1S5vuT1l9QRNCZ1BfUbYvqUpAjnu0gIscxxkj6hQ4VVNGvNQpdnyeTvuwWROVD3TzBtNdEDN14h6d8SYr2UDT96hhQ5jLCp1QfUDR84y7EKmSzoucJ14x0Yg3ZYc7pza8km1xnKYYF4EzimCpU+h4yMH40xM2yUF2UaFTRvg0lRQh7fqaoOn6neI0LWR6JIXJ7AcsV+AUETFNCZ0iVPmZt+X7rIrOJh7oFxTJ90nhUTdPvGKi0fqd6ck4CJxQAqV3njoSE9SdStSTEeQTarL6sEkA+rf1n6do4oA6vFJZ26ZTyurpTqg5Nf39Ql33xO3uWRu0Typ7vhBJaYr2IeMuRqqik56cQR6TJrw77sVxmsIFTn3nTKNuj453dvmMKhDKJAEou7+o16RKcVB0jktbvDxF6GWK60p7m8ZmTQy1vk3a4p5p978674mDBEJZATKMYGnrgqFOB0VOv5hYs9vzrNnt+YZa4zjN0nWBU/Vcmrrm6oxb+Jqk8yWtlHRbxv7tJP2bpJslLZX07rrbWAXJAfIo32mIgbZ7QYrTtLCRdKSkOyUtl3RWRpk/kHR7bC/f6Nu3raRJSV+op8XZZPUpg/qa/nttG4ROntDIs28XKNUyyF4kfSC2lVskXSNpfmLf9yStkvTvRerqlMhJCpyTf/4D9n7+LrbZdTXb7LqaNbs9z0EP3s3JP/9Bavmq8RA1p27qEDhVDt7rzIo2rl6dCrkAODJn/+nA7Wb2auAw4DOSWuHOHlYcJAfKx/7mv9hr/XIee/1zTP73WazebSb7PnEP71hxXeYxeduGbV+RkLiyFF04NO9VtJ4ybRqWYdsXCkkzgHOAo4C9gBMk7dVXZgHwEeB1ZrY38P6+03wCuI6a0LPlxitpfU3vXnfiLT9gr/XLWTtvirXzplgzXxzw8HJOvOUHtQudIot99vdpWX1cqIVDBzHdBFURewF+ARxoZq8CLgP+PrHv08C7itbXiTk5aWJl8X5z+cJFF/HNVx/MPQe9mCdumcmnvn8FH37LiYXPF1qYuNBx6qIugVMFda9v019v1XN1mlo8NCRmdp2kXfKKALMlCdgGeBxYV0fbsiizEGaSflGyereZLFk/j3+8/CI+sc/RfHv+ARxy3wP81Z2X8revOL7UuUKQXLSzSN1ViJCQx+ads2mPzJAcDCw3s3sBJF0KHAvcnijzHuAcM3sCwMxW9nZIOgDYEfgecGBdjS7KoL5myd6RrXx8n9/jF/vOY6dlq/k///kNPvLGaCxW97xFSL/3Fk3p3H9s//sy/eKgOT5d7iNGYKC9mNm1ifLXA+9M7LtG0mFFK2u9yEkTOGt2e57bdn0R33z1wbznw9fx2Du3ZpuLpzjzI3/AtbN3Y/Y9Gx+flQHNkwRMHyTtDFxE1JkYcK6Zfa6vzB7AV4H9gY+a2T8k9h0JfA6YAXzFzD5ZV9v7cYETpg11iJ0x7sS+ACwEHgZmA+8ws8bihvXsFGwd/d8/UM4a/KcNqFfvNpM188R18+fz3Tfsw9998Dvs93sP8j++cytnHn8it8/YfaTvNIQXJo3+7G9dDXXrqNDZCUh+4JPAIX1lXg4g6b+I+pC/NrPvSdoM+AzRIO7NNbR1AxN3TOb2JXn7ev3Dmvniunnz+fg+v8f//rtvc+M75vHKSx7iT//sndy49YIN99gm74Vl+7JBbe3tS553UHmfs7MRRewlySnAd4etrHUip2iI2Vvm3clrX72cx27dmpd8fg2L3rsXi1+1C9wXiaDZ9xSLxHOhM21YB3zQzG6UNBv4uaRFZpZ82vY48D7grckDE+7VI4gMcrGkhX3H1kJXBU4bxE0/dYidJjv3J9ZvxWVP7p+x98o5kpYkNpxrZueWOP3vAjcBhwO7AYsk/cjMnhyutdVTdPC8dv4Uu8z9FQ+8ZA7PnDjBSZ+7ngvfeSjXHTafOT8a/vxVeU7qyJRWJ00JnYrtZXNgAVFo51zgOkmvJBI3V5nZZOQUrZcsoVO2n/nFvvN4+sQJXv+55fzo1N156IjtWDs5BUw0KnSG7cuKtDUpdgbN7emSoCmaYe+J9Vvxr0/ul7Ln4VFtZQOS3knk3XzjMMdDC0VO0XVnvr/iFey2eCXv+dp1PPq+2bzm4ns4aN79XDv71YUFTn99zvhiZo8Aj8T/r5G0jOiJQtJFuhJYKen3+g4vEo5QOVULnHH23uRRtdhpafjaY2Y2SmjMu4FPmpkByyXdB+wB/CxI64YkxOB45gMT3M+LmL/ix2x50RQXvvNQ/sd3buWazffj9hm7b1Q2KSKyBudNCo2uenXKCJ2aBF2evTwEJBs7N96WZBK4wcyeA+6TdBeR6HkN8HpJpxGFfU5IesrMUpMXtJX9blrBVhdN8aNTd+eVlzzETtut5pdbv2yje2pXBE7y+CLze4rub9n9P5OpPebCoyOdYlDfUsRekPRm4KPAG81s6A+vdSInSU+A9Iud2fdsxt53/4p3XPQz/vkzb+BLO76eg158P3//f67gw2/ZhsU7L9jkHFnvnbGh8NODeK7BfsANBc9d1r0aHBc41bNmvsbWq1MBK4A3AT+StCPwCuDeZpuUT5FwtWiyMey0y4McdfltfPRjx/LtOQdwzeb78fmvX8gn9jqBm3bYLfV8Wdv6KTNfqAiDBEEbhM4w15y8rpaH6S0GFkjalWiwdjzwh31lvg2cAHxV0hyi8LV7zeyPegUknUw02br1AqcXgjX7AeOQex/g7Muv5C8/9tZoTs52qznn/0Vzcn7+st0Hn6yjhA6Dm0YMtBdJ+wH/BByZnL82DK0WOT12uHNqE6Fz0C8mOePEE7ltxxfx1H3bce3sV/Pht2zD3isf3CByXNC0g1CJHn69fhsufuK1GXu/VejJtKRtgMuB97c5tCZJlQKna2vfrJ0X/YZmrqgmkdc09epsgqRLiEJr5kiaBD4GbAFgZl8mygZ1gaRbAQEfNrPHGmouUN6LkycOfue6+/jzt53IDb9ZwJxlxkP378wn9jqBV6x5iJt22C1IBrfe+6LnGmUQX/bYYT7Loucp43kpG6bXBGa2TtIZwNVE823ON7Olkj4OLDGzhfG+t0i6HVgP/IWZ/bq5Vo9OT+gcuHQFf/62E7luznyYhF9u/TI+8sYT2fuxFRtEzrD3u2E9IaH6taz7ddb5k31eWh/iQqewvXyayLP5L3EY5wozOwZA0o+Ioga2ifumU8zs6qz6OiFyYFOvzgUHHA7AU/f9dv7N4p0XuMBpMU3Pf5K0BZHA+bqZXVHi0ELu1SroosCpUtz0v69S7Exnr46ZnTBg/8PAW2pqzshkDZazhM43570B7ofZM2zD93TTDrsFFTjDHFuH0BmmnaMcE9L70pQ3x8yuAq7q23Z24n8DPhC/ss5xAVHq9toI0b9c9KpoLDZzxW8fEP38ZbuP7MXp75+K3jOr6NeSYqeIwOm9rzq5TVcpYC+ZSTjM7PVl6uqMyOnR79VJm3/jAsfpJ051ex6wzMw+W/LwIuEIwbFZ1S070pXkAv3iJm2/e3WcQQw7+E1+N02HfDnFaEN43nQhOejPuk/WdX+rY02eEHW0MdPcOPdBnVoMtMcOd06lCpms7U47aPi7eR3RAlKHS7opfh0t6VRJpwJIekns/vwA8FfxCtTbmtk6oOdeXQZ8y8yWNnUho9IFgdNbWC502WGoel7RGC0g2iiDFoLM8jaMGjZVlCY8QMPUH2qx0SYW5uwxtefOrQlnG3fSFs3MW5izaPauMrThHtp2r03eZ9SGz68qOufJSeKCphukfU91h66Z2Y+J5g7klXmUKBQtbd8m7tUu0hWBM8pxVXh26ghfg/F+olYXeXNtskKk0o6pY8HLUHW0YX2ZvJBAaCattXt16iP0vas/NCzr/G0aoM9+wGpbhy00bQ+hHpZOenKcbtMLNyy6JpIThrYLnFAemao8O2vmy706HWHQwDZPEFTtgeidu82D76oEXhPX3LT4ayt1rLmWRVFvzqgpnPPo3c/7XyFom8Ap+jmNY//jIsdxpgFdEDihqSqEzYVONygrdIYZhA8bFlX3YL/NgsqZXkztMbeUwOqad6GrAmdc6XS4mtNNPMywXtoscKqcS9M7f1Xha1BtUoKude5dZNRwpjrEQ5H1YpL7Qs/hGXS+pgRUG0L0uszEHZONenPKUOe9sOrQ5KLX0t9v52V2SzvvdBc3PVzkOM4YE/pG13bvTV49XZur40LH6VFGSIQe/CfrHmXeUpmFPYdZMNRxypCcP1MXwwqcrG1l9hdhHPsbFzmOM4ZMZ+9NXr1d8+q40OkWTUyuT6Oo92eY8w5DUoikec/SFkkt2o48j1PT30OXmQ7phusUOqMInNDtSKtjnL7XJD4nx3HGDBc4+fV3ba5OE2EHks6XtFLSbTllDotTsS+V9J91tq8IRZ/wh/IE9A/kRz3XqOdoMoXzIKrwvvRfb53XLulISXdKWi7prJT9J0v6VWL5gj9O7Jsn6fuSlkm6XdIutTU8g+kU6pT2cKpoEoJkquw0kZC3L406BE7vb5OipoC9vEHSjZLWSTqub9/fx33OMkmfj9dAzMRFjuOMES5wiuFCZyAXAEdm7ZS0PfBF4Bgz2xt4ez3NKk6XJ/cXEShdDtOqOntdzQJnBnAOcBSwF3CCpL1Sin7TzPaNX19JbL8I+LSZ7QkcDKysvNHkZ1frwlP9kPfELC/8ILGTNm+mrLBpkibaWtBeVgAnA9/oO/a1RGsevgrYBzgIeGNefS5yHGdMqGL+zTgKnB5VppqugjqFjpldBzyeU+QPgSvMbEVcvpaBWVmKLEpZVbrkKukJnK4InUEelrZ6nApyMLDczO41syngUuDYIgfGg7vNzWwRgJk9ZWZPV9fUCJs1OGw3yzPRBsreC4ukiM4LNy4jdIaljvt7Szx0A+3FzO43s1uA5/uONWAWMAHMBLYAfplXmc/JcZwxoK0JBkKKiF3m/mqj9/dPvijIeauYq9PVBeFK8HJgC0k/BGYDnzOzi5ptUj51TFLvymC9ifkrg87flc8uhZ2AZOMngUNSyr1N0huAu4A/N7MHiexolaQrgF2B/wDOMrP1Fbe5ED1R05LBMTCcwClK3hydvCQzvphzKYrayyaY2U8lXQs8QrS4+xfMbFneMS5yHKfjTEeBk9wWQuxUlYEtdPa1MokInlw3i0WP7pGx98o5kpYkNpxrZueWaMrmwAHAm4AtgZ9Kut7M7ipxjtoZh2xco7S/69deJRXby78Bl5jZWkl/AlwIHE5kR68H9iMK0fkmUZjOeWXbXyWhBu/JvirtnGnhX1n7ipDWl41yTx50bNlEMU2Ix7w2bmjPDwef58l1s7j60T1T9oxsK5lI2h3YE+jFWi6S9Hoz+1HWMR6uNiZsedtDTTeh1UjaWdK18cTOpZLOTCmjeCLbckm3SNo/sa/UZLe6GHeBs8vcX6UKnP4yoehC+Fqg7/wxMzsw8SrbCU0CV5vZb8zsMeA64NUhGlY1HfYYOM2RZy8PAUn1ODfetgEz+7WZ9UaWXyF6QACRHd0Uh+6sA74N7E8Hmbhjsukm1EKI+/nq3Wa2yjsWmEF9y0B7yeH3gevjsM6ngO8Cr8k7wEXOGOFCJ5d1wAfNbC/gUOD0lMluRwEL4td7gS/BcJPdukgbBU6ZsqHEThVzdepej6EGvgP8jqTNJW1FFG6QGzbghCE5z6isYHOBVwmLgQWSdpU0ARwPLEwWkPTSxNtj+K2tLAa2l9RzRx8O3F5xeyuhzMKiRbw4WdvKMPsB28TzMsiLU8TLM8r9vM3ipqZwu4H2ksMK4I1xv7MF0Tgst99xkTNmuNBJx8weMbMb4//XEBnGTn3FjgUusojriTqflzLEZLc6CHmzbJPAGUWwtNmr0yWhI+kS4KfAKyRNSjpF0qmSTgWI46C/B9wC/Az4ipllppt22oMLnbDEHpgzgKuJ+pVvmdlSSR+XdExc7H1xJMDNwPuIQtKI5958CLhG0q1E8wz+ue5raJqqB/49sdPkHMmueG+qFjpF7EXSQZImibJ2/pOkpfHhlwH3ALcCNwM3m9m/5dXnc3LGkC1ve4hn9ukfvzs94nUI9gNu6NuVNiFup2Emu1XNOAucUOdo41yd0HN0qsLMTihQ5tPAp2toTqWkLU457vQvoDndrj80ZnYVcFXftrMT/38E+EjGsYuIogTGnnGamJ91L0+b89IFcVMnBexlMb+dd5Mssx74kzJ1ucgZgZ7XpI2Cok6hU6aupKepbPuenJrF91e8Imt3oclukrYBLgfeb2ZPFql3mMluVeICp/j52piBrStCp6v4wL0445CIwXFCk5dlrQguatrDwHC1tJWvJb1A0iJJd8d/d6i2me0jOVjf8raHWhcmVofASV537/+szyFtX+DPbOBE6jiG83Lg62Z2Rco5sibEFZ7sVrW9jKPACTmfJu3coQgZvtal0LUqqcJeisxV6e2ve+HINjLdr79LjOt4LMvDM06enzaRN44YN4FWZE7OBWy68vVZwDVmtgC4Jn4/7Rk00K+Snqh5Zp+dKhc4g64x+Tm0RQDG2dDOA5aZ2Wczii0EToyzrB0KrDazRyg32e0CKrKXcRU4VeNCp9VcQEP9iw/unQ5yAWM6HtvunrUbiRoXOPUzbgIHCoSrmdl18RyGJMcCh8X/X0iUVfvDIRs2iDJhT/2D7Lq8HGlUWXcT19Wrsw1CZgCvA94F3CrppnjbXwLzAMzsy0QxokcDy4GngXfH5S4jynpzK1ESgu9lTXZrq70kmU4CJ1mXh661jy7Yi+O0helgL20RN4NC1pq+b4deqHUcBQ4MPydnx/gJN8CjwI7DnKTsvJG8UKis86Qdk7atjfNXqqbMnKJBIuaZfXZqtdAxsx8TJQ3IK2PA6SnbS09262Nkewl1A5qOAqe/zrYlJJjuQieFIP2L40wT3F4qYtS5OVXQLwKHFTvJBAnjKnAgQOIBMzNJmT20pPcSrTnCzFnbZ87LyBtoFxk8p4mHMoPuUSbEl6VpoVP2Oyj6OSaPL3JMW8RenZSyl5nbA+MlcJoQN2ltaJtXx4VOOnn2krSVWZttU2u7HKeNFLWXXt8ynen1h1n33baJmx698UAIgTLO4qbHsOvk/LK3uFX8d2VWQTM7tzcZfOb67MFACC9AfzKAUc4zaM7JqAyaoF/VXJbQ52wiLK+DDGUvE1tsHawB4yJwjnjJHRzxkjtGPk8b5+m0tVNtgEL2spGtbLZlrQ10nBZR3l4C9i11EHownrzXrpmv3PdptOGBVFfW3WmaYT05C4GTgE/Gf78TojFtC3dK83CEbGPSo5MnevrbkHe+0IxyzrwQtmkmgoa2lxA3sXESOMn/Fz26x0jnc49Oa6mkf3GcMWWs7aVKgVNkez959+fkOfw+3g6KpJDeZOVrImM6QtLdwJvj96XoZQELPditwgNSpfgq2t5BZaq45rQkA2W/s7Sy4yxwQtrL+lmji5NxEDhZ3ptQHp1QXp1x8uikparNKHeQpHWSjhuynkr6F2d8aeO6PpKOlHSnpOWSMrObSXqbJJN0YPx+C0kXSrpV0jJJqQuGJo4PYi96doqJOybLXKKTQ5r3p+n7eK9NZdtRR7sH2YukN0i6Ma1vkbRe0k3xa+GguopkV8ta+fpNg44dJ9rgZaprLk+RxAxlsqolPTrjLHCgXfYyLgKnyP62eHXGyKNzAfAF4KKsApJmAJ8Cvj9sJW2yF6f9tFTgzADOAY4AJoHFkhaa2e195WYDZwI3JDa/HZhpZq+UtBVwu6RLzOz+tLrcXoYnRH/o3pnRKWgvK4CTgQ+lnOIZM9u3aH0jJx4YljaIhlFpIn1yv9AJHT5XlrLpu5vMbOeUpymBU9ZL06bwtZApppsiI1VtP39GtLjuQdW3yHFau7bRwcByM7sXQNKlRGmdb+8r9wmihwJ/kdhmwNaSNge2BKaAJytvMTC1x9w6qqmFfgGTFCN1eCZ69dUZrjYoccKwmeFqEHID7aUn8iU9P2plwyYeaISqQtyGbUvy/yba1XahWLR9bVkwdJwIcWPvisAZ9bgkbQpdazrcIQ9JOwG/D3yp6bY4TsPsBCTV12S8bQOS9gd2NrMr+469DPgN0Ft0+h/M7PEK2wqAzerOQ5hB6+ak3SeHDdUalZ5AqEIoJK+pP1FCXnta6H0aaC8DmCVpiaTrJb11UOHGPDllyQqXgvoH+3lipg7vTpuuNwsXLc0xHQVO8vhx8uiMErY2NbV53nXMkbQk8f5cMzu3xOn/L/BhM3teaq8Yc7pDMhStCY9NVfYiaTPgs0ThN/0cDKwHXgbsAPxI0n/0nnJ3hWRygDoX82zDg6D+e3RVAqdtTE1tzorJOWm7Ru1bBjHfzB6S9N+AH0i61czuySrcCZEzaJBd50KURQf8bV8csyqm4zW3iekscPrPM4rYCbVwaNNCJ4fHzOzAEY4/ELg0FjhzgKMlrTOzb4donDO96J9rM7Xnzm0LTcuzl4eA5AXMjbf1mA3sA/wwtpeXAAslHQP8IfA9M3sOWCnpv4hsq/UiJyvrWXKRyS4xzALPvf52UOhY1nH9DHufD9lHBPjuBvUtg+wlFzN7KP57r6QfAvsBmSKnE+FqdWcWczbGP99u0FWBE2rtm7TzjkqI8LVQWdfahJntama7mNkuRCE3p7nAcYalX9C0TOAMYjGwQNKukiaA44nSOgNgZqvNbE7CXq4HjjGzJUQhaocDSNoaOBQIfzMMSNXrs/TOn3ylEdq7MXPFRJAHUnkpqqsMoQt13hrW3sm1lzwk7SBpZvz/HOB1bDr3bSM6IXIgfaBd91yOMmFb4yYMun49knaWdK2k2yUtlXRmShlJ+nyc1vCWOI66t2+epO/HaT5vLzApu1a6LHCqZFyETt3hCmmpaiWdKunUWhviTBsmlj244dUlzGwdcAZwNbAM+JaZLZX08dhbk8c5wDaSlhIN/r5qZrdU2+LyDBIcaeWHraet9HtK8u7J/fNnRr1/j+ohahNF7CVemmCSKPvgP8X2AbAnsETSzcC1wCf7sxj204lwtR69zGJNDLjbkOygaYqkgW5xmN464INmdmOcyvPnkhb1GchRwIL4dQjRpOpD4n0XAX9nZoskbQOMnPUjFC5witUzavha06FrdaaVzklVm1b25Aqb4kxjWhi2loqZXQVc1bft7IyyhyX+f4poINdK6hQdZepqKqFA1RS9x6dlc0ueI1mmjQyyFzNbTBTG1n/cT4BXlqmrM56cHi0dQE8ruvgdmNkjZnZj/P8aoicI/WrtWOAii7ge2F7SSyXtBWxuZovi458ys6frbH+VjLPACVlniIVDu+bRcZymaePaOONOiJC0Msc3IXBGEQFV3ofz5uv0Z0vLu4Ymssu1kc6JHKcddFHo9IhDzfZj40XZIDu14cuBVZKukPQLSZ+OF7RqnFFvYtNF4LSh7h7jOEfHcRxnGJry4JQ5V8jQszL1FSFP6LTZm1MXnQpX6xItDtsKRpULkw5BobSFcajZ5cD7zazoomubA68nEkYrgG8SpQM9b6QWj0jXBM6oAuO4bW/c8P9lT+6fU7JYO4YNX2s6dM2fzjnThS6Eqo0jvQxbIbw5Wdm62jz/pijJPrSKBaDLhK/lLYg6nXGRUyFVCZ1hhEVV6/fUKWzWr53BU/dtl7V7YEpcSVsQCZyvm9kVKUWyUhtuDtyUWKH320RZcBoVOaPQVJroYUkKnOT7UcVOk0LHcRynzYQQO2WPzRrYl3m4M0w66LL096EhlgsYhZ7QcXGzMY2Gq02Hyfwhr/GZfXZKXRQ1q47evuT+6fCZp6FogYLzgGVm9tmMYguBE+Msa4cCq83sEaKsN9tL6o1qD2dA2sKqGeVpfpdC1I7b9sZNBE7//lEYxbvU9PwcxwnN1J47b3g504dBQmS7e9Z2ZpHP5H21qnts1nmrqK/MZ+ECZ1MaEzldGmyP2tYQ11pkQdRQ5xpTXge8Czhc0k3x6+i+tLhXES3Cthz4Z+A0ADNbD3wIuEbSrYDi/Y0wnQROyHJZuNBxHGe6U8Tjkid0qhRCRfu8tPvp2nlTQe+zg87V9Xt6FxdyzaOxcLX++RxtJVQbq56jk7WOUFr7x32uUBpm9mMicZJXxoDTM/YtAl5VQdM6QVsFTn/5YcPXPHTNcaI5MO7FmX6EGNiGmmNTZj2aMqSJj7LhZUUETJMha6MybgIHGg5Xa/tgO7QIG/Z8eccNWhC17Z+xU44mvDhtFzihjm3Ko9P1J3/OeNHFBTmd+hingXDPy1PE29OUwMman+QpoovhKaQzqMrL1Dtvcr5M3ryaPIoKmGS5UUTPKG11RqepMLVhaELghDiHC50ISedLWinptoz9fyTpFkm3SvqJpFfX3UbHaQuSjpR0p6Tlks7KKfc2SSbpwMS2j8TH3Snpd+tpcTeoexCfJXja5MHp/0y6KHQG2YukN0i6UdI6Sccltu8r6aeSlsb9zzsG1eUip4+ig/hRxULR7f3bep6bsvUPc0yyDS5smqVL83CGEQqDEgwMc75hcaEDwAXAkTn77wPeaGavBD4BbJKu3akWDytrB/GaaecARwF7ASfEC0j3l5sNnElifba43PHA3kT29sU2rMEWYjHQUSnb54W+f/bON+i8M1dM5Aqc/gU8y5A8Ls9zU3bNnyzq+M4L2ssKomU6vtG3/WngRDPr2cv/lbR9Xn0uchIUHcj3xELdoWBNhZ55yFt36YrAqYIuCp22YGbXAY/n7P+JmT0Rv72eKN26UyMeUtYaDgaWm9m9ZjYFXAocm1LuE8CngGcT244FLjWztWZ2H1HSm4OrbnAeyYFuT+z0ttUlfrronchilLCy0GFpRc5Tw/c70F7M7H4zuwV4vm/7XWZ2d/z/w8BKIHdCbOMppEOEbYVqSxH6B/whPSpZyQOc6U0TN/y2CZwXfOkptvzJxrHgW/5kLS/40lNBzp/GqIuXDkPLvDlFOQX4btONcJyG2AlIKs7JeNsGJO0P7GxmV5Y9ti007dlpgiL3467cs5PjiIZFZJDfvKSDgQngnrxyrUwh3bTgKUuIuTGDzuFixylLnYkGqvbgPPPqLXjZaas2CJ0tf7KWl522imdevUXQevoZVui0IWxNU9oQStH/AuZIWpJ4vXeoOqT/TiRyPhyk0Y7TEJoSMx+YSH0xgr1I2gz4LPDBqtoeijYImWEG4F0RGlWS97lVIWqy+hcC9S25dUsvBS4G3m1mz+eVbSSF9DBrugwa5KfNXQndniJZzMpcWy/FcxcETsiU38nztOHa2sywN6dxEjgAz7x2Jg9/cXtedtoqVr1rK7a/+Gke/uL2PPPaYp3ycdveOHR66WFpeWrpx8zswMHFspH0KuArwFFm9uswzXKcVpJnLw8ByQlSc+NtPWYD+wA/jNak5iXAQknHFDi2NpoWOFl9Xa8vCz2xv3e+EAJp7bypodpXVbrsIudcM19VLR46qG8Z6TcvaVvgSuCjZnb9oPJjMSdnlMF36EU0s8LRyoSoNUUT3rOueOuc5nnmtTNZ9a6tmPO5p1j1rq0KC5wu0vYnk5LmAVcA7zKzu5puj+M0yGJggaRdJU0QJRJY2NtpZqvNbI6Z7WJmuxDNYTvGzJbE5Y6XNFPSrsAC4Gd1X0CTAmeUOSdZqZ8HCY4617LJamNFAqML5NpLHnH5fwUuMrPLihzT2GKgIUkKhdCiJeuYLnhfhiXrGrt8TdONmSsmhhoo3z/5otLenEWP7lHam3PZk/uX9uZs+ZO1bH/x0zx25jZsf/HTPP2aicJCp24vDjCSF6fpBeUkXQIcRhR6MAl8DNgCwMy+DJwNvJAoGxTAulE9Q47TRcxsnaQzgKuBGcD5ZrZU0seBJWaWOYCLy30LuB1YB5xuZutraXjNDCtksvqyIlnP0sol762hHiY1fb/OE0xtE1NF7EXSQURiZgfgf0r6mzij2h8AbwBeKOnk+JQnm9lNWfU1InKGCX0qEtbVO3f/cSHqT2tPWn1lj0+2KSSjtq93jjqETV31dJHZD1jTkwQro4zQ6c3B6YWoPf2aiY3eD6pnWBY9usfQxw5L0x0mgJmdMGD/HwN/XFNzHKfVmNlVwFV9287OKHtY3/u/A/6ussYNoAvZ0soKk2T4WFLshBI4w9yjB4W0DfsZtU3EFGGQvZjZYlIydprZ14CvlalrLMLV0sgLvWoyRCqrXW0N26qrXW29/i4z7GB5GA/EsGKgqADZ8ubnNhI0vTk6W978XJDzpzGKwGnxXBwnIFN77rzhNep5HKdtJFMYj5IGOQRZYV9FaULg5NUT4rPtosCpm7EUOdN9wDyMl2zcPzNJO0u6VtLt8Wq5Z6aUkaTPx6vw3hKn/Uzu31bSpKQv1Nfy7tzIqhQ6j//pNpt4bJ557Uwe/9NthqqzSrocpjadCC1Mss43SAT19rnQceomz4tTRQTBsEIl5DF1e3CqoivjgqZpTOR4eNLGhPg8yoiVJtN0N/TdrwM+aGZ7AYcCp6essnsU0cTPBcB7gS/17f8EcF3VDQ1Jnd4cqN6jU9c5p2uY2nQi9GKaw56vd5wv7um0idkP2IZXKBJphmshZAKXsu0eZh5RUVzgFKdRT06ZwW5VZeuml32t/zUKw4ibIu2rkrq/IzN7xMxujP9fAyxj0wWojiXK2mFxasLt43zsSDoA2BH4fo3N3sAoN7XpKHQ8TM2pmqQoGSRQRt3vOE0yitgZpe8aNUStd45RCSFwQuECpxyNh6u1WZCUoSmvSJXZ5KoWO01995J2AfYDbujblboSb7yY22eAD9XSwGlMCKHTRYHjXpxuMrHswVICJtQcnrYzHa7RqZa2pNEvOwdnULtHua461tUZN1qRQnpQtrM2C6H+dteZKayIaCm7OOmw5xv2mstkupuxFmbfk6nL50haknh/rpmd219I0jbA5cD7zezJgs08DbjKzCbjdLmNMEqmtTpTSsNwaaVD0LUQNWe8SBvYTyx7MHX71J47F/LgJI9twuPTq79o3cn2Fr1Gp37qXBunbN81ihCYN/exDf+vmJyz4XxlHySFSjBQpnxanVmem97nmfa5urfntzTuyemRN0hu46T4OkK6ukQLvqPHzOzAxCtN4GxBJHC+bmZXpJwjayXe1wBnSLof+AfgREmfDH4FY0bdYWtNrIXTw704Th4+0G8O9yhtyvpZ9T6sq9Pj0BM2aZQRI0VTVQ8jyIre8wcJnCzcw/NbWuHJKUJZD0mItXDSzpm3rW1enLpp83o3ilww5wHLzOyzGcUWEomZS4FDgNVm9gjwR4nznAwcaGZnVdzkVLrkzRmFsouFjipwfB6OUzWjJiZoiibrL+tFctrHMP3VoL4qa6HPHmlCp66sbmmkiZph++Oi9D736e7VqdWT8/yWW+TuHzRAbmpgX9RrU8UAP5kFrUxGtKraEqJMQ7wOeBdwuKSb4tfRkk6VdGpc5irgXmA58M9EYWoO9SchgOLCpcsCp61eHEnnS1op6baM/bnp1uvGZk00/sS+6frHjf6wtzYj6UhJd8b2sMkDsLifuTXud37cy+wp6QhJP4/3/VzS4XW0d32LbjtlB/rJe2bd2dqGIa99w7S/rGhso1engL3MlPTNeP8N8TxqJE1I+mpsLzdLOmxQXbV7cnqD4GEH4WW8BSG8OaHEQrIdWeeswvPUYtFRK2b2YyDX2s3MgNMHlLkAuCBYw2qm6qdHaYwyP6esR6csTc7DaXnnfAHwBeCijP3JdOuHEKVbP6SWluXgc0C6T9sFTT+SZgDnAEcQJatZLGmhmd2eKPYNM/tyXP4Y4LPAkcBjwP80s4cl7QNczaZZPyuhjif9gwbYvb6o6JyZrDL927uWtCCr3Jr5GjtPTEF7OQV4wsx2l3Q88CngHcB7AMzslZJeDHxX0kFm9nxWfa2Zk9Oj6OT3QYP3ENnOQqV2TktOUIf4GGZR0CbqdcrR1E1vFK9FVR6dJhMNjHOYmpldBzyeUyQz3brjVEGLxc/BwHIzu9fMpoBLiexjA31JbrYGLN7+CzN7ON6+FNhSUn0ZASpizXyN7EEY5SFQz0PS5IOklj/EapKB9hK/vzD+/zLgTfGUg72AHwCY2UpgFXBgXmUDRU7WSvGSXiBpkaS74787FL/GbEZJiRxigcu61omBTcVOyDqH/Qxc6IxG3fYyDF28+aaJmS4LnFG+g5Y82UtNt172JF2wl6KkeZCmo1epzPpBg47vZxihU8N3UMgWJJ0u6R7g74H3pZznbcCNZrY2rZIu2EpZcVNWjLTFQ1Mnyc9zWOHYspC1IvayoYyZrQNWAy8EbgaOkbS5pF2BA9g4WdQmFPHkZK0UfxZwjZktAK6J3xci9GA+lNemKZFRldAZlrqFThuuOSDB7SWN6ebNgY1FzXQVOGWYMbXxquV9K5jPkbQk8XpvLY3alKD20rSo6NVfZP2ctlCFh2TU62/is5sxBbNXWOqLAPZiZueY2W7Ah4G/Su6TtDdRSM6f5JwieN8Sqh8J4bkZJ0Lf49v42Wb1L1Tbt5xPJIqWAP8X+AmwPu+AgSInZ6X4pDvpQuCtRVrYxgFt6DYNIxJc6DR/zSEIbS9VMcpNuGmh02Sq6KYJKG4HplwfQFa69VK02V6GXdCyrgF6iAU3e8dXHQo2TFuzPseGwtby7KWsLVxK4vcsaS7wr8CJZnZP1kGhbSWkwKmT0GmgQ9J0iFxLGNS3FLGXDWUkbQ5sB/zazNaZ2Z+b2b5mdiywPXBXXmNKzcnRxivF7xin1wV4FNhx0PGbPfNc7v66Q5zyvDdNh1uNy6B/OjOqvQyiybClJoXOKEwXL04NLCRaL0qSDuW36daHZlR7qUpchB5U9wb8oc4b8lxV03Ghk8ViYIGkXSVNAMcT2ccGJC1IvP094O54+/bAlcBZZvZfRSsc1VZmBBr7N+VhGEa8zJv72EYLhRahzP266nt7G705QzLQXuL3J8X/Hwf8wMxM0laStoYoMyGwri9hwSYUzq6mvpXilVj5Pa48dcQVu6reCzBrxuzUczctKLLob1fdwqPp7GhtXvem7YSwl4mtqw2tbiLTWo9RMq6NUucoNC1w6hS1ki4BDiMKPZgEPgZsARBniboKOJoo3frTwLtHrK+0vSRtZebM7UepPpcqPTMh1oHpHTvMuSaWPdj6jHT915fc1gbMbJ2kM4gyo80AzjezpZI+Diwxs976a28GngOe4LcDuDOA3YGzJZ0db3tLPKk6lTb1LaOs21Y3/WvnZAme/nJF+8miGeKmOwXt5TzgYknLiRLgHB8f/mLgaknPE3l73jWovkIiR+krxf9S0kvN7JE4q06qUcauqnMBZm83dxPja3Ltmyyy2lR00N+0OAlJ8pqHua7pKJJC2cvWL9p54Kh21E5mFKEz6iKhdQqdrgucujGzEwbsH5huvSjD2kvSVrbb4sW5tjKsAAhNlgdiGKExyMtR9Hx1CIYQYqpNwqYfM7uKSPgnt52d+P/MjOP+FvjbovXU2bcUpffwpb8fqloAFREVoR7kdVnotCRRzUYUsJdngbenHHc/8IoydRXJrpa1UnzSnXQS8J0iFZZd1DKNZAa00NnQiqSmLkKo9rRBJPTPFyqzTtF0I7S9tJ0upFFuMjwuFG3sqEJQh710ZVHJUdqWduyw56sq/K3NIqULtL1vSZmAPtJ9q6ioGESa6CgbtlaGtmSAS34P05kic3JSV4oHPgkcIelu4M3x+8rJGzjXlfq5KE0LnSqzxbXts24RtdvLqDeyJp88VS1AQpx/unlxaqYye8karLdZ6ISm7LX2C8Lp9Fl1gFaNxYoyLgPtMvfxpoROv8h0CoSrDVgp/k1hm5NPGQ9CW8LFQrWl7HnqEiAudDamTfZShnEMWxsXgTPOHVZX7cVx6qbLttL03J26w9bqYpz7hlCUyq7WNHUM8qvwUHRJCKS1dRy8NlkLqfWVkaTPS1ou6RZJ+8fb95X00/i4WyS9o/4ryKfpm92oYiC0R6cNAicETX+vTliqCtnKWojUQ8ScrhLaAx7qfEXOU7X33r01xSmcXa2L5Hk/BiURCO0NGnS+JkVEf91dFzQZ9BZSu1HSbODnkhb1pR88ClgQvw4BvhT/fZpoDYO7Jb0sPvZqM1tV7yXk02QSAmiPR6ctAsfD1Jw00rKFVXH+tp3LccrSu4dm9UvD3mNXTM7ZaF5Of0a1Uajqvu+iZjjGWuSMSpvC3oahSNtHFTS9Ogadp2mBF68j8Ej8/xpJvYXUkiLnWOCiOHPU9ZK2j7PW3JU4z8OSVgIvAlZV3vCa6brQGSeB451acWzWRGtTDGfRyzgW8nxNESIV9ih1d+H7doaninCzQcIm7x6e1Z6QAsfv/2HoVLgalE85PUyq6KLHl2XQ4qMhMs+VbU8dFMlYV6eY7FtILclOQLK3nIy3JY89GJgAMlembpI23BibCl1zgeP06MqEeR+chyH0AqtOeMrc09Im0M9cMTFtvON+/w/H2HtyigywBw32kx6dEMIg7Rxp7Uy2LXTq6lAenKJ1hRAxM541drgz82nOHElLEu/PjdcF2Ij+hdTK1B+vQXAxcJKZPV/m2C4R4qlZ3R4dFzhOj7Q1WdI8DU16H5K0aUHOpHep6TWFulC3U5wyodRr5mujNNTJY3v907D32kEhcMkyg86TPEcoAeb3/7B0zpNTBW0IScvz8BQ9PnmOKsPDhvm8avAaPWZmByZeaQInbSG1JA8ByUeBc+NtSNoWuBL4qJldH7754Qhxk2zDE7OiwqUtAmeckHSkpDvjJBxnpeyfFyfy+EWcjOPoJtq5oT3PTm0y0b4/fC1rXxto00DdExaUp4C9fCBOenOLpGskze/bv62kSUlfqK/VxVkzXxteddaZRYj+qecZSns51VLAXmZK+ma8/4Y4+qa371WJRFC3SpqVV9dYi5wyg/GyA/c6Q6wGCZYy3puqBE7R0L+mkhrkLKSWZCFwYpxl7VBgdbyK9ATwr0TzdS6rqcmdJ4RwGCRg2iRwxsWLI2kGcA5RIo69gBMk7dVX7K+Ab5nZfsDxwBfrbeXw9Auctgkep1sUtJdfAAea2auAy4C/79v/CeC6qts6DP1iYxihM6w46h2Xd19cO29qw6tuevf8ur04TYjOUBS0l1OAJ8xsd+AfgU/Fx24OfA041cz2Bg4Dnsurb6xFTlmKiIl+cRM6A1uZsmXLj0qoa21I6KQupCbpVEmnxmWuAu4FlgP/DJwWb/8D4A3AyYlj9625/aVoizenSqHjAqcyDgaWm9m9ZjYFXEqUlCOJAdvG/28HPFxj+xynTQy0FzO71syejt9eTxQlAICkA4Adge/X1N6RqXJwXUas9Jdtwxo2VaV3zhI2HRQ6RfqXY4EL4/8vA94UP6h+C3CLmd0MYGa/NrP1eZWN/ZycrhE61XRbFwWtO3PdgIXUemUMOD1l+9eInh50ihALsLVhfg5sOkdn3AROy0hLwHFIX5m/Br4v6c+ArYlWWm8VWSFXobOaOdOeIvaS5BTguwCSNgM+A7yTFtpQFlU+kOna/TRke5NzkcaYIvayoYyZrZO0Gngh8HLAJF1NlOH2UjPr94puhIucPupePyerjl5b0raPC2WTKjjN0CahE4q2zcGpomOb8ayx3T1rs3YXStSRwwnABWb2GUmvAS6WtE9bEnIMmlNS5ZyTtiQRcMpRsb0AIOmdwIHAG+NNpwFXmdlk9KC6ffQ/LBvmXhXigVvaOaqa/J9Fr64q68kK0Uv73JsURTn2EsRWMtgc+B3gIKL1C6+R9HMzuybvgNYQehDf9XVuRml/FYIorT1ZmeLqDqVz0gnRuUB7hE4IQgqcDoepPWZmB2bsy0zAkeAU4EgAM/tpPPlzDrAydEMdpwWMai9IejPwUeCNZtYbHb4GeL2k04BtgAlJT5nZJpOxm6SpgXSRekOtoTNILCX3j5LdDYr1y0UETBXfS47QL0qerUAxe+mVmYzn4WwH/JrI63OdmT0GIOkqYH8gU+S0Zk5O2YFu3WvKJAk5KB90Db25N1ULgWE+y7xU2E18L+PGjGfD3MBC3QjbMkenLfV3WOAMYjGwQNKucdKN44mSciRZAbwJQNKewCygeQXbEjwcblox0F4k7Qf8E3CMmW14EGBmf2Rm88xsF+BDRMltWiVwQhHiXpd3jlHux8PM7RlWWCVTYw+ig/NtilCkf1kInBT/fxzwg3g6wdXAKyVtFYufN7Lxgu6b0AqRU+UAvqpzT0fvQ1K0lFnrxwVP87jQaZ/AaStmtg44g6hDWUaURW2ppI9LOiYu9kHgPZJuBi4BTo47IceZVhS0l08TeWr+JU5c0z+oq5W8B2jb3bM2xNP8kSnTZ7Xpfhw661mdQqeO772gvZwHvFDScuADwFnxsU8AnyUSSjcBN5rZlXn1tSpcrQzJ+RxFFvPslS163jJtCEXZMK8m6V8ctchn26Xrawvb3bOW1bvNDHKu6Ry61sYkAy314gBgZlcRZRtMbjs78f/tRBkLnT58Ps70o4C9DEwqYGYXABeEblsWgwa0IfueHsl7XlZfVCREqwkPx6jzf/qvq21emjqFbQF7eRZ4e8axpRJBNe7JGXXQW2Xa5enEoM8lb32cKtcjcuq9+RSlSx4dFziO4zjlqbLv6aVa7n+VPTYE/ff2Qff6UQVO20h+z6GFbdM0KnKaEhzJeS7JNrRBAA0rAqoUD6HP7UKnPKE6m5A321BCp0qx4wLHcRxnfBlW7AxzzMwVE0H6grq8OEVC58ZZ4ECD4WptEBQ92tSWsvQvTFpFCF0VpJ63u19DLYQKHwjp7g+V2aaK8DUXOI4zviSTO3iIYHVUEbZWBWX6tTrvw03d88uIGxhPgQMNiJwuC4qqCfHZlBE6/fNq0o7ztWzGk3EWOm1NMOACx3HC4cLG6adov5ZWptePtSmJwSj0Pou0fme6CByoWeRs9sxzMLvOGttPFaKviNAp6wHK2+8CqD7amIgA2iN0XOA4juOEpWi/0xs8d3XQHFrgNH3fL1J/V7+ronQ2u1poupr5K2vB0KzryUsgMOz1J49zwVM900HoAKXETuh5PW0WOG1MROE4SXqhZO5tmT4k70t596guDaqL9I/D9KFNhNalfSdd+i6GpfHsam1gXAfm/WvUDLrOEJ9DW4WipJ0lXSvpdklLJZ2ZUkaSPi9puaRbJO2f2HeSpLvj10n9x9ZNyIFu25IR9CgqXEJ7b1zgOE4YfFHU6UGb7kuhFh0ts2hnVR6bQeftrWlU5NXPdBA44CIHaO/APDSjhKSFrKch1gEfNLO9gEOB0yXt1VfmKGBB/Hov8CUASS8APgYcAhwMfEzSDnU1vA66KnTaGp4GLnAcxxkfsu4/43ZfKnPfTpZNy/JWRRa1POFShNW7zZw2Agdc5GygroF5FV6jvLanpcnOSzAwrpjZI2Z2Y/z/GqKVdvu/jGOBiyziemB7SS8FfhdYZGaPxyvuLgKOrLH5qYRembrNQqdfzIROPd12gVM3ko6UdGfs1Twro8wfJDyj36i7jU57mVj24IbXdGCQvUh6g6QbJa2TdFzfvnmSvi9pWWxPu9TW8JL0D7DL9j9VDa77UyWHTiudtr3o+bPSOGe1MW8NoP50z2UFS1vETQF7mSnpm/H+G3o2IelgSTfFr5sl/f6gunxOTo30h42FFBZFzpWcd1O1qMmaK1QhcyQtSbw/18zOTSsYG8x+wA19u3YCkr3yZLwta3sraPMcHSDIPB2obvHQLgicOp+WSpoBnAMcQfRbXyxpoZndniizAPgI8Doze0LSi2troOO0iCL2AqwATgY+lHKKi4C/M7NFkrYBnq+4yUHIy9DVluxdIfuz5LmGucdnZToreq689WxW7zZzYB/RIoFTxF5OAZ4ws90lHQ98CngHcBtwoJmtix9A3yzp38xsXVZ9LnJqokxygLrqr5rQQmezZ57L+7weM7MDB50j7kQuB95vZk8Ga1zDtFXoQLiEBFXgAieVg4HlZnYvgKRLibycyU7oPcA5sWcTM1tZdyOdjfHJ/o0x0F7M7P5430YCJg6Z3tzMFsXlnqqpzUFJG3TXRbKvShMSo/Rn/aJm1HC0LKEziCJ9QJbQaYu4SVCkfzkW+Ov4/8uAL0iSmT2dKDMLGPhherhaDeQN9Mc16UEbkbQFkcD5upldkVLkISA5U3ZuvC1re6toa+gahBcTIeiCwGmIIp7LlwMvl/Rfkq6X1Hj4phPhk/1rZxRP/8uBVZKukPQLSZ+On3Q7ASkavtYLLesPMatibk0/RTLSDQpPa6GgSaOIvWwoE3tpVgMvBJB0iKSlwK3AqXleHBhDT04VoWCjMsijUYdHp02fRxNIEnAesMzMPptRbCFwRvxk4RBgtZk9Iulq4H8nkg28hShUZ6wZZ49OVwTOKMJVz04xccdk1u7C4Z0ZbE6UoOMwItF/naRXmtmqYdrqhMM9OcOhZ6fyPrtR7SWLzYHXE4VPrwC+SRTWdl6Ac08LRg0jC02vDVnzcLIYJFDKrFNUBzn9S1W2AoCZ3QDsLWlP4EJJ3zWzZ7PKj5XIKbvAZZvoWns7yOuAdwG3Srop3vaXwDwAM/sycBVwNLAceBp4d7zvcUmfABbHx33czB6vr+nFCRm2Bvk37GEIPU9n2PpD0kaBU4C88M4instJ4AYzew64T9JdRKJnMU4juLiplFHtJYtJ4KZE6M63ibJ/dkrkhO53ylKHuBkUCpe2ry7y+ooGvptBUweK2EuvzKSkzYHtgF8nC5jZMklPAfsAS8hgbERO1pwX6I4Xw4VOdZjZj4HckbqZGXB6xr7zgfMraFpwqripjYNXxwVOYRYDCyTtStTZHA/8YV+ZbwMnAF+VNIco7ObeOhvp5JMMW3MBVClF7CXv2O0lvcjMfgUcTs6AzWmWLKHTJB1M4V3EXhYCJwE/BY4DfmBmFh/zYJx4YD6wB3B/XmVjMSenjkUu66JLbQ2Ji7uwVHHj6/I8HRc4xYljnM8AriZKtf4tM1sq6eOSjomLXQ38WtLtwLXAX5jZr9PP6NRN/7wcn6dTHUXsRdJBkiaBtwP/FM8pwMzWE2Vcu0bSrUQP4v65iesYlabvW6Oydt5UoQdvdczPKUrRz7xN303B/uU84IWSlgMfAHpppn+HKKPaTcC/AqeZ2WN59Y2FJ2fQnJcmB9BtFi1t8Ry1oQ3jSFc8OlBd+FpVQmpcBU4PM7uKKHwzue3sxP9G1Pl8oOamOQWYWPage3JqpIC9LCYKy0k7dhHwqkobWAPD9jV1hVOFTCXdRdI+56b6mwL28izRA4H+4y4GLi5T11iIHGjnQHkYgVN3SumqhU7y3GmfRxu/t3GiC0IHqglfc4HjTGd6QscFjuMMpkh/0TaB0+vbu+jRqYuxETl1khysZw3SuyBw6qZ3fV2bK9V1ppvQqTIMzgWO0yVc4DhtpytenLYJnCR5C7FOd8ZiTs6oPLPPTgNTPCdfg44dJCLSBvdNCpyQomfQZwnRtbrAqZcuzNGB0QWKCxxnujC1584+18ZpBaOEqnWB0Pf+/nV4QtKRtXJqwz05A6jK61FmkF+H56No2FpV4W3P7LNTC5fXdAZRlUcHys3T6aK4ge508o7jOKFp+4C8int/mxIXTAcGenIkzZL0M0k3S1oq6W/i7btKukHScknflNS+Jc1TSPPIjOpZ6D82pAgI7WUZdn/yM8vz1gzzWbY5OUNZ2mwvVQ2oqxIBRYWLC5zu0mZ7aTsTyx70cLRpRJttpe33stkP2MB7ea9M8uV0nyKenLXA4Wb2lKQtgB9L+i5RVp1/NLNLJX0ZOAX4Ut6Jnt9yi5EbPAzDDqKLHJc1oB+U8W2Y+gd5UUKIhTLnaEt2tpYRzF6qoKpMNqEXDe0xyKvjAqfztNpeHKdFuK2MSPKenrWIZ5XU4cUZ1Ldk9f/j2icN9ORYxFPx2y3ilxEtWnVZvP1C4K1FKgwx8C87EA+9r8egAf4oSQnKXmNT3pDQQq7rhLaXKqjyZlaXV2fmigkXOGNAF+zFcdpA221lu3vWbvIKea7Q981x9NYMK3B6+9oePjgMhRIPSJoRL76zElgE3AOsihf1AZgEUkerkt4raYmkJc9N/QYYblDef0yRcxQVE/1hWCEETrJcyPCtsp9BHbShDW0ilL1MPfebytrYVaFTtbgBFzh1M6y9bGQrzz9TW3sdpym60LeMSt49chSxM8733iKfyzgKmCIUSjwQr8q7r6TtiVYZ3aNoBWZ2LnAuwOzt5m40eiiShSuEkClC1QP1kNfSRlGRlVa7SLrtOpB0PvA/gJVmtk/K/h2A84HdgGeB/2Vmt8X7/hz4Y6KnZrcC744Xq0ollL1st8WLK33MVOUibFUkJKiaqp/qdamTlXQk8DlgBvAVM/tkRrm3ET1FPsjMlgxT17D2UqetOE4eg+xF0kzgIuAA4NfAO8zs/jjs7CvA/kTjsYvM7P9k1ROqb9l29txK7SVkvxLqXG0Y5CfDukP0N2X6lCL9fV191LD2Eu/7CFFI5nrgfWZ2dV5dpVJIm9kq4FrgNcD2knoiaS4V5MZq42B+FHqD/H4R0BaPTChaek0XAEfm7P9L4CYzexVwIpEBImkn4H3AgbE4mgEcX6TCEPYyccdkkWKtpEuhAC5wfoukGcA5wFHAXsAJkvZKKTcbOBO4IUS9dfcvjhOCgvZyCvCEme0O/CPwqXj724GZZvZKogHdn0jaZVCdbbaVUcVE7/hxDZ+CMP1NyD6linDALEaxl7jc8cDeROO5L8bny6RIdrUXxU8NkLQlcASwjMjAjouLnQR8p8D1jQ1Z681kvXpkeTmqpg1r0zQpeszsOuDxnCJ7AT+Iy94B7CJpx3jf5sCWcUeyFfBw1km6Zi9V39jaLnTqiMvuksCJORhYbmb3mtkUcClwbEq5TxB1PplezUF0zV768bVqHIrZy7FE82Ug8ny+SZKIogO2jvuWLYEp4Mm0StpuKyFFSVPiprd+Td1RCHWJjDrmOhVgFHs5FrjUzNaa2X3A8vh8mRTx5LwUuFbSLcBiYJGZ/TvwYeADkpYDLwTOK3R5Y0rZULS6BU4baZmn52bg/wOQdDAwH5hrZg8B/wCsAB4BVpvZ93POE9xeqvbmTFehU0e7OihwIIrpT+Ym3iTOX9L+wM5mduWIdXW2f0mKm6JCpyeKpqMwGuPrHmgvyTLx/JnVRL/ry4DfEPUtK4B/MLOsh3GttZVx8Lj0C5s6BE+/yCgqOjrar/QYxV6KHLsRMqtvACJpDXBnbRWOzhzgsaYbUZAutRXS2zvfzF6UdYCk78XHpTGLjZ8onxvHICeP3wX494w5OdsShajtRzTvZg/gPcADwOXAO4BVwL8Al5nZ17LaGYqO2cs4/P7aTH97c20FhrcXSccBR5rZH8fv3wUcYmZnxO83I/J6nhzPK/gh8KFh5+SEoGO2At3//bWd1thLvO22uMxk/P4e4BDgFcBpwMnADsCPgKPM7N6C1zkUbi+V0/X2jmIvuWOxEe3lr4Hre+MvSecB3zWzXnbBTSiUeCAgd5rZgTXXOTSSlnSlvV1qKwzXXjPLm1MzEmb2JPBugNgteh9wL/C7wH1m9qt43xXAa4HKRQ4dspfp8Ptrkprt5SEg+ci9P85/NrAP8MPIVHgJsFDSMQ0Knc7YCkyP31+TtMxekmUm49C07YgmVP8h8D0zew5YKem/gAOJ+p4qcXupkOnQ3obspcixG1Eq8YDjjCuSttdvV4r+Y+C6WPisAA6VtFUsft5EFAftOOPKYmCBopXUJ4gmei7s7TSz1WY2x8x2MbNdgOuBJgWO4zRJrr3ELCSaLwPR/JkfWBRGs4JonRskbQ0cCtxRS6sdpxlGsZeFwPGSZkraFVgA/Cyvsro9OY7TCJIuAQ4D5kiaBD5GtJgaZvZlYE/gQkkGLCXK7oGZ3SDpMuBGYB3wC+I0nI4zjpjZOklnAFcTZRM838yWSvo4sMTM+jskx5m2FLSX84CL43kzj/PbDJ3nAF+VtBQQ8FUzu6X+q3CcehjFXuJy3wJuJxqPnR6nVc+k7jk57+2fJ9FmutTeLrUVutfeJujSZ9SltoK3d9zo2ufj7a2WrrW3brr2+Xh7q6Vr7S1DrSLHcRzHcRzHcRynanxOjuM4juM4juM4Y0UtIkfSkZLulLRc0ll11FkGSedLWhmnrette4GkRZLujv/u0GQbk0jaWdK1km6XtFTSmfH2VrZZ0ixJP5N0c9zev4m37yrphvh38c3ExP9pjdtLWNxexhu3l3C4rYw/bi/hcHtpP5WLHEkziCbXHUW0qvwJkvaqut6SXAD0p8M7C7jGzBYA18Tv28I64INmthdRNpbT48+0rW1eCxxuZq8G9gWOlHQo0Wrp/2hmuwNPEE/2n864vVSC28uY4vYSHLeVMcbtJThuLy2nDk/OwcByM7vXzKaAS4Fja6i3MGZ2HVEGhyTHAhfG/18IvLXONuVhZo+Y2Y3x/2uIUhrvREvbbBFPxW+3iF9GlDqzt4hTa9rbMG4vgXF7GWvcXgLitjL2uL0ExO2l/dQhcnYCHky8n4y3tZ0dzeyR+P9HgR2bbEwWknYB9gNuoMVtljRD0k3ASmARcA+wyszWxUW68ruoGreXCnF7GTvcXirCbWUscXupCLeXduKJBwoQL0LUujR0krYBLgfeHy9cuYG2tdnM1pvZvkQr1B4M7NFsi5yqaNtvr4fbi9NG2vbbA7cVp7207fcHbi9tpg6R8xCwc+L93Hhb2/mlpJcCxH9XNtyejZC0BZFRfd3Mrog3t7rNAGa2CrgWeA2wvaTegrRd+V1UjdtLBbi9jC1uL4FxWxlr3F4C4/bSbuoQOYuBBXH2hgmilUu7sGL2QuCk+P+TgO802JaNkCSiFWGXmdlnE7ta2WZJL5K0ffz/lsARRLGr1wLHxcVa096GcXsJjNvLWOP2EhC3lbHH7SUgbi8dwMwqfwFHA3cRxf59tI46S7bvEuAR4DmieMRTgBcSZcW4G/gP4AVNtzPR3t8hcn/eAtwUv45ua5uBVwG/iNt7G3B2vP2/AT8DlgP/Asxsuq1teLm9BG+v28sYv9xegrbVbWXMX24vQdvq9tLyl+ILdBzHcRzHcRzHGQs88YDjOI7jOI7jOGOFixzHcRzHcRzHccYKFzmO4ziO4ziO44wVLnIcx3Ecx3EcxxkrXOQ4juM4juM4jjNWuMhxHMdxHMdxHGescJHjOI7jOI7jOM5Y4SLHcRzHcRzHcZyx4v8H1GHZOAkqN2sAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig,ax = plt.subplots(1,4,figsize=(14,2.5))\n", "\n", "counter = 0\n", "for c in [c1,c2,c3,c4]:\n", " s = ax[counter].contourf(c[0,:,:])\n", " ax[counter].invert_yaxis()\n", " colorbar(s, ax=ax[counter])\n", "\n", " default_peak_position = np.floor(np.array(c[0,:,:].shape)/2)\n", " i = np.array(find_subpixel_peak_position(c[0,:,:]))\n", " ax[counter].plot(i[1],i[0],'rx')\n", " print(np.array(i - default_peak_position),\n", " np.sum(np.abs(np.array(i - default_peak_position)-np.array(shift))))\n", " counter += 1" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcwAAAChCAYAAABKxo+2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA02klEQVR4nO2deZgV1bX23326oRmbOSioOEXFAUecjYjiQIwSruKAnxo1osZoNA5R891o4sX74b1XjF4ThxtnjXFK1IQgRpwV0EscwSkOiIAyNzP0qe+PKvSs364+XRi7jU/W73l49O2qU7Vr166qc+rda62QJIkcx3Ecx6lO6atugOM4juN8HfAHpuM4juMUwB+YjuM4jlMAf2A6juM4TgH8gek4juM4BfAHpuM4juMU4B/2gRlC6B1CeCqE0BBC+M+vuj1fhBDCr0MI//erbofz9xFSbg4hLAghTM7+dnoIYU4IYUkIoUcL7tvHkPO1I4TwRAjhlC/42Y2y66rmy27X30uLPjBDCO+HEJZnBz87hHBLCKFTwY+fKmmupPokSX7cgs1sMZIkOS1Jkl981e1wqhNCODGE8GoIYVk2Tn8VQuhascrekoZI2iBJkl1DCG0k/ZekA5Mk6ZQkybwvuN9BIYSPqq3jY+jrRQjh2BDCi9k9b1YIYVwIYe+vul3/yGTPiQPW6iRJPsyuq8avsl15tMYvzO8kSdJJ0g6SdpR0UcHP9ZP0RvIFMiuEEGrX9TPOPychhB9L+n+SzpfURdLuSsfehBBC22y1fpLeT5JkaaZ7S2on6fVWbq7zD0wI4VxJYyWNVjpGNpJ0naTDc9b92tyj+Evv69T2L50kSVrsn6T3JR1QocdI+mOF3l3Sc5IWSnpZ0qDs77dIWi1plaQlkg5Q+nD/iaR3Jc2T9DtJ3bP1N5aUSDpZ0oeSnsr+fpKkaZIWSBovqV/FvhNJp0l6O9v/f0sKFcu/n322QdIbknbK/t5H0v2SPpX0nqSzqhz/LZIuz/5/kKSPJF0g6RNJsyQNkzRU0luS5ku6uOKzu0p6PmvbLEnXSmpbsfxASW9KWqT0onxS0ikVy5s8dv/3WR/VZ+NrBP7eKTu/J2VjaoWkxmzduyUtzcbPEkmPSwqSrsrO62JJr0raNttWnaT/yMblHEm/ltReUkdJyyWVs+0skdTHx9DX85/SL1tLJB3ZxPJLJd0n6Y5sjJyS3Useys7bO5K+j3P3YrbuHEn/lf29XbaNedl5nSKpd5V2NXUf6y/piWwbr0s6DGPuV5L+lI31A5Teyy+U9IqklZJq1cT9O9vGE2vHkqTNsutkntK3hndK6potuz27BpZn/XeBPr+f12brVOunS5U+C27LjvF1Sbu02Hlu4UH0vrIHpqQNlN5Irs5036wDhyp9GA7JdK+Kk3Z5xbbOlvRCtp06SddLujtbtraDb1N6I2qv9FvdO9nAqJX0U0nPVWwvkfSIpK5Kvwl+KungbNmRkmZKGqj0Zri50l8ZJUkvSfpXSW0lbSrpb5IOauL4PzsGpTe7Ndln2ygdyJ9KuktSZ0nbZINmk2z9nbMBWZsd3zRJP8qW9VR6IQ3Plp+t9AvG2gFa9dj932fn5+DsnNTmLLu1YnydKOmZimVrx9vaC/qgbFx0zcZLf0nrZ8uuUnqxd8/O88OSrqgYEx8100YfQ1+Df9XGUrb80qx/hym9j7SX9JTSLyrtlL6B+1TS4Gz95yX9n+z/O0naPfv/UdkY6iCpJjvH9U3ss6n7WJvs3F6s9D42WOnDZsuKMbdI0l5ZW9spvZf/VdKGWdubu38/UTGWNs+W10nqlR332Ip2vi/7w2pj2eurWj9dqvQL7dCsP66Q9EKLnecWHkTvK/3W0JB1wF/0+TeLCyXdjvXHSzqh4qRVPjCnSdq/Qq+fDcC1N4NE0qYVy8dJOrlClyQtU/YtOVt/74rlv5P0k4p2nJ1zPLtJ+hB/u0jSzU0c/2fHoPRmt1xSTaY7Z23YrWL9lyQNa2JbP5L0YPb/x0t6vmJZkDSjYoBWPXb/91m/HCdpdhPL/l3ShOz/T1T1B+Zgpb/wdpdUwnlZKmmzir/tIem9ijGxrg9MH0P/gP8kjWxqLGXLL1X25ivTGyp9a9G54m9XSLol+/+nJF0mqSe2c5LSX3UDCrSpqfvYPpJmY6zeLenSijF3Gz7zvqSTKnRz9+8nVPG2AusNkzQV2859YBbop0slPVaxbGtJy1vqPLeGhzksSZLOSi/2rZR+s5XSbzpHhhAWrv2ndHLF+k1sp5+kByvWnaa0I3tXrDMD619dsf58pTeFvhXrzK74/2VKv8lJ6Ul6t4k29EGbL0YbqjEv+dzIXp79d07F8uVr2xBC2CKE8Eg2CWWxUl9kbd/1qTzWJB0plZNHihy7k74e6tmEJ7N+trxZkiR5XOnrzv+W9EkI4YYQQr3Sb9MdJL1UcS7+nP39i+Jj6B+TeWp6LK2l8v7UR9L8JEkaKv72gT7v35MlbSFpeghhSgjh0Ozvtyt9MP02hPBxCGFMCKFNCGGfbKLRkhDCWm+9qftYH0kzkiQpN7FvtjXvb4Xv31nEw29DCDOzcXiHPh+HzdFcP0nxfbxdS/msrRZWkiTJk0q/ufxH9qcZSr+hdK341zFJkn9vYhMzJB2C9dslSTKzcjdYfxTWb58kyXMFmjtD6Xv3vL+/h212TpJkaIFtriu/kjRd0jeTJKlX+mAO2bJZSl9NS0rDHiq1/r5j/2fieaV+zPDKP2YzuQ9R+kakEEmS/DJJkp2VfsPdQukkorlKH2DbVJyHLkk6CU6y47Ul8DHUeqwdS8OqrFN5vj+W1D2E0LnibxspfYWqJEneTpLkGEnfUDop7b4QQsckSVYnSXJZkiRbS9pT0qGSjk+S5OkknVnaKUmSbbLtNXUf+1jShiGEyvv/Z/vOaWve39bl/j06++x22Tg8Tp+Pw6b2VdnWJvuptWntOMyxkoaEELZX+i3jOyGEg0IINSGEdtk0+w2a+OyvJf1bCKGfJIUQeoUQDq+yr19LuiiEsE22fpcQwpEF23mTpPNCCDtnMXibZ/udLKkhhHBhCKF91u5tQwgDC253Xeis1GNaEkLYStLpFcv+KGm7EMKw7JvUDyStV7H87zn2fxqSJFmk9LXXNSGEg7Nv6hsrfT3/kdJv880SQhgYQtgtCzdZqtRTKWff4G+UdFUI4RvZun1DCAdlH50jqUcIocuXemCf42OolcjG0r9K+u+sTztk4+mQEMKYnPVnKH21ekV27xug9FflHZIUQjguhNArG0MLs4+VQwj7hRC2y2auLlZqS5W5/Yym7mOTlP4SuyBr4yBJ35H023U45HW5f3dWas0tCiH0VfplspI5SueDRDTXT61Nqz4wkyT5VOnEnH/NOuJwpd96P1X6jeX8Km26WunkiUdDCA1KJwDtVmVfDyr9Zvbb7DXAa0p/NRRp572S/k3pZIoGSb9XOiO3Uek3uh2UzpCdq3RQtsQN7zxJx2b7v1HSPRXtm6vU0B+j9FXQ1kpn1K3Mln/hY/9nI0mSMUrH4H8ovQFNUjoW90+SZGXBzdQrPUcLlL4umifpymzZhUonWLyQnYvHJG2Z7Xu6Uu/ob9lrrT5fykF9jo+hViRJkv+UdK7SCVJr72lnKr1/5HGMUr/uY0kPSvpZkiSPZcsOlvR6CGGJ0nvf0UmSLFf6peY+pWN1mtKZzblf7Krcx1YpfUAeovQedp3SX6nT1+FY1+X+fZmknZROJPqjpAew/ApJP82ugfNyPl+tn1qVkBmlzteY7NXKR5JGJkky8atuj/P1w8eQ4zTPP2xqPKc62auQriGEOn3uTb3wFTfL+RrhY8hx1g1/YH592UPpDLi5Sl+vDMte2ThOUXwMOc464K9kHcdxHKcA/gvTcRzHcQrgD0zHcRzHKUDVbAijR48272vffvtts3zRokVGd+lioytWroxn5dfV1Rl98803G33bbbcZPXXqVKNnzrTxqt/97neNPuaYY4w+88wzozZss802Ri9dutToCRMmGL3vvvsavcEGNtTokUceMfob3/iG0W3atInawH3yMwsXLjR61apVRn/88cdGd+jQwei99trL6Pfffz9qQ9++NmnLp59+WnWfV111VVArMGTIEDPuNt98c7N8q622Mvrkk082evLkydE2r7jiCqM328zGc/fpYyM6Hn30UaPr6+uN7tXLJuuZM2eO0d26dYvawHPKNnz0ka30xetr1113Nfrll182mtfbihUrojbU1tpLnmOzbdu2Ru+www5GX3rppUa/+OKLRt91111GH354HCrN4/zDH/5gNK+v6667rsXH3YG7XlbVmyrX4VZZtquXVqyJPlOzyF7jSckexsp+3Y1uO3eZ3UCw66/q3r5aE9VmcXy/5T7L7e1xhNU2hLMErbLVSRtbojKsthW4Skty2oC+K6FfGnvaa6vc1q6f1NhjCOj72gW231b3jCtIhjX2OGoabDsbO9vn0mPP/jR3zPkvTMdxHMcpgD8wHcdxHKcAVV/JPvnkk0bz9c3559sMR9dff73R3bvbVw6StPHGGxs9YsQIo7t27Wr06tWrjeYrogcffNBovo7j/qT49STX4atlvnb605/+ZPTQoTaVLF+X8nWdFL8ue+ONN4weNmyY0a+/bmsVsx8IX9/x85L0v//7v0b/8pe/NPrOO++suo+WYqeddjL6ueds+tINN9zQ6NGjRxvN8ydJAwfa7IWvvfaa0bQKGhoajO7YsWOVFsev1H/zm99E6xxyiE2Uw9ejS5YsqbrPd9+1ebT56nKXXXYx+r777ovawNe2AwYMMPovf7Hpc7mPo446ymjeE/ia/29/+1vUhldffdVo9kO53FSmt9YjqbW/JVZ3sm2sWclXmfbVpCSt3Mje/2rw2nZNB/t6s7GPfZVYs8ruo80ie+5Wd7FjdlVXqyWpzWJrqyR4zdvYAY8AvMKtbbD331Xd7Cv7dnPs69CkXWw/ldvbvyV45VrGa95GvDYuNdp+aGxr1y+1s20Kjc2Pn9U9rYVVs2x1E2ta/Bem4ziO4xTAH5iO4ziOUwB/YDqO4zhOAap6mFtvvXXVD0+ZMsVo+nIMS5DicAr6F5yO/8knnxjNqfb0jrj9fv36RW1gaAo190HYL/Q4GYYyadKkaBtDhgwx+vTTTzf6nHPOMZqe5e9//3ujv/3tbxvNcJ08eL4uvvhio/fbb79mt9ES0IdmGNC4ceOM/vDDD43OG7dbbLGF0fTGH3roIaM5buhLk/bt7ZT/H/7wh9E67O9p06YZvWaN9bgY4sHriWPivffeM3qTTTaJ2rBgwQKjGY5Ej5Led02N9Y+OOOKIqp8fOXJk1AaG1/A4p08vXDTjS4N+YM1yey7q5tlzt7rerr+iVxzyUbfAeo5rOlovr2a59T3Ldfb3S7mN1St7tLM7gB+5ulP8+6exnT1fdfNtm2qXsA32frymM9q8onqYyequaKOkNe3tOkLARt08tAl+In3XpLPd3sredgzTt5Vi7zWqwFkwcMl/YTqO4zhOAfyB6TiO4zgF8Aem4ziO4xSgqodJr4G+2rXXXms0/cM8GGuWF6tZCb0jenOXXXaZ0UwJ9+yzz0bbZIwi06Ixtoxxm9/73veMvuaaa4w+++yzjabHKcUxjkzhN3/+fKN5LuhfNRdbyrSGUuzT8fyxjUw72FIwjpIp4JiWjqnW6IlJ0ltvvWU0432ZXo/jZp999jGa8YX0TZ9++umoDUzxR6/8m9/8ptH0D7mc8Yo8f3kpGfm3o48+2ujly211L8ZMMjaU3u+bb75pNMdxXjs5VrfddtvoMy3NqnqkDISnFdohNnCV9f4YUylJn+xs4yo7fGI/02YJ09BZ2VhnG1GDrHO1S6zXF9bEbVjZHakPF9nfSGs6VfcXE+g18A/L6KjQGGcYZLq9xna2DYy7rFlh+4k+LNP9Me3dqi7x9d8If7h2uf0M0+01hf/CdBzHcZwC+APTcRzHcQrgD0zHcRzHKUBVD5N5J+lpMf8pvUCWi5Li/LP0KP/6178a3aNHD6OPPfbYphus2OujlqRly2z+Q5ZyYrt/9KMfGT1x4kSjmduU3tQtt9wStYE+HNs5atQooxmfypg7xsTyGPP8LMafXnfddUafccYZ0WdagySxfgLLec2YMcNo5jvNy1/6zDPPGE2fdPjw4UYzBpVjZMcddzSa8cK77bZb1Aa2e9asWUbTo9xjjz2M5jl/6aWXjKb/yDEgxbljOQ+B/i/HNvuW1wqPMe8eMGjQIKPZd19FHGbdQusHruwaXy+VJCXEK7aLA/mmXmKvp22vttdTr5etj7ayq/XqIo8T0Hcrt41//7RdbP1A5sTlNpgjty1iSRs72H7hPktrYi+wDM+x7SIb48ptNNbZfqhtsHGVZeSSXdOx+jFJ8XHRmy2K/8J0HMdxnAL4A9NxHMdxCuAPTMdxHMcpgD8wHcdxHKcAVSf9MKn1vHnzjGbwPJMGMMG1JN10001GH3nkkVUbOHjwYKOZmJuTL9iG3/3ud9E2GTjNBO5MbMCkACwgTTixiUV1pXjCDdtw5ZVXGs3JGJy0wolGJ5xwgtETJkyI2sDPsG+bS0LfUrzzzjtGM2kAJ0xxkkjeuOM5r6+vN5oFox9++GGjmTidScyfeOIJozt37hy1gUH/vH6YCP3CCy80mhN0WACcx3TiiSdGbbj11luNZuFsJhE466yzjOb1e8899xjdqZMN1uekISmenESYVKI1aMTEE04cKa1Kqi9fE//22O0ntqBCtwYkLkBxZgbsM4kAJyLVrKyBjicJMSED1wmYH8MkAZyQU7PMTtgJmNCzqkv8SGHfqRbJ1HGcjR1RULqdTeheu8z2U1TMe1VczHvpBjY5fhsknVdOwoU8/Bem4ziO4xTAH5iO4ziOUwB/YDqO4zhOAap6mPQXmVSA0CdikWMp9uKY5Jr+xeOPP240PReuf8MNNxhN/1GKEyyQnj17Gn3vvfcaPXDgQKPnzJlTdXt5yROYaOCCCy4wmh4ng9qZ8IH+GPuey6U4KcQdd9xRdZ+tBT1KJiE//PDDjabnmZekgUW9mbh8l112MZrB9IRJA5gwfPbs2dFnjjvuOKMfe+wxo6dOnWp0Y6P1WXi+6JWPGDHC6KuuuipqA88xEw0w2QELYfMaZ0INXgt5Bdx32mkno+m1sph6a0Avr2YlkgLAd6OHWbs09s3q5tmAewbYUzNJQFJj92nvfDEsOC3lJCJYiDZ1wCMAniSLMbOYM4s9N7aP28DjEKYY8LjLSMbeiONisvWalSiCXVskgQPXqZ6oYi3+C9NxHMdxCuAPTMdxHMcpgD8wHcdxHKcAVT1MekeMJ6T/yGKxeYV86Q/yM9znyJEjjaaPSp9t++23N/rGG2+M2kB/kLFnTCrPBNKMT6Snye0PGTIkasMHH3xgNPuSbWIR6lWrrLfAGD7CRPqStGaNjamiz9mcj9dS7L333ka/8sorRnPM0Hdrh7gtSdp///2NpgfJcUX/njHITITOMbHRRhtFbTjppJOMZrsZ+9m/f3+j6aX/+c9/NnrTTTc1eujQoVEb6IPSL3z33XeNZvzu+PHjjT7qqKOifVSSdy7Yrttvv91oxry2Bm3n2djQNZ2tYxgQj7gGXl27T+Lri3GUkUcJr21FT9tXdVifBaOFIgVK4vttCf5e5IvOtce9upttQyMKZ7N4M/3DdrPjhP8qwYNcDb8X8wlWoh/ou7JN9IJrVsXxqG0Xoi+XI4F7u7j4dh7+C9NxHMdxCuAPTMdxHMcpgD8wHcdxHKcAVT1MepBvvfWW0euvv77RzBuZFw/HnJv0zaiZt5XxcvSe6AUOGzYsagNjFLkNens8zt69e0fbrGTs2LFG58WV0WttDuaOZRwmfVbmos0rYv32228bvc022xj93HPPrVMbvyxYpHjzzTc3mnGto0ePNpo5V6XYL2Ss55IlS4x+7bXXjC7Bh+G1QQ+ZXqAkLV682GgW6H7ooYeMprfHnLps080332x0Xrwxj5u5lukfMt/tEUccYXRNjfV+vvvd7xp9wAEHRG245JJLjGZM7HbbbRd9pqVZBe+u7SLrm5U+tXplb5szt3ZRnB83LLN9We5qPxMXVrb7YLwhPctynb19t5kf+4dJnb0HN7ZFseVV9l5Hn5TFnGuQpzWKX51rx7gklTvbe3pAfDHbWDfX9ltorJ4rll5wuU3zfmTNEsaXNhflmu270FqO4ziO80+OPzAdx3EcpwD+wHQcx3GcAlT1MOmTnXLKKUbTF6Nnklc3cuutt666j6iBtbaJ9Bu5nB7MH//4x2ibzIdJ746xoPQ0maPzzDPPNJqeWl5c2csvv2w049X4GcavMp6O/hi93Dw/mX4x8+5+Vblk6QeyPukf/vAHox955BGjN9lkk2ibTz/9tNHMccocxYyTZQwyxzpzsp5zzjlRG44++mijmVuWXi3HCGumHn/88UbTm8+rNcuYVnqthx12mNH0sffYYw+jeT0yFpTXgiQdeOCBRjOe9IEHHjD6sssui7bxZcOcqYz1Cyutt8ccquX28fVVove2yN53wkrETSJPa8L4xRU2lrDU1u4zLI191KSNPY7SYvh78I/DSrvPWhx3uYNtc1iO2NByXFeyhHaxTWxDzRJ7v20OHne0/bzP4LhKDe5hOo7jOM6Xhj8wHcdxHKcA/sB0HMdxnAKEhPkIKxgxYkTTCxXXsWOtvYMPPjj6zOTJk41mzlT6aPTe6D8yTytzrrKephTHitFH5T4Y49hc7ljG6O21115RGyZOnGg0YwsPPfRQo9lvzbXxhBNOMPrJJ5+M2sAYVcbk0ccbP348smO2DIMHDzbjjjVP77vvPqMvv/xyo0eNGhVtc9CgQUb/4Ac/MJrjhPmEt9hiC6MZX8gxwPqbknTTTTcZzdhQXovt27c3mj4142jpy9JPlKTVq+HFYQ4APU7mIKaf37lzZ6MZt8kxJMV5kJmnl8c5c+bMFh93B295oel8xlAmHe25YExkaLDXoySpzt4HEnqOy5F/tga/X+BhJrXwH+HVhcU5beA26TFiedIeOXTh9SXt4Luuht+Y9zyBN0sd9TX3QdjX9Ta+Nayw4ysX+KZJvZ3P8efXR+eOOf+F6TiO4zgF8Aem4ziO4xTAH5iO4ziOU4CqHuZ+++1nFjKnavfu3Y1mXB+9Cin2NBjPxhqMjD2jb8rYM8Yf5tXk5D4Zb0jPkXlZmY+TMK6TOXglqVu3bkbTq2U+Vfo69H7pwzKuc/DgwVEb2Hf0UelfTZkypVU8zOHDh5txx/qXrCvJOEzGRErx2HznnXeM5vnYbLPNjKZXV19fb/Sxxx5r9NVXXx21gbG13MfUqVONZk3Nbbfd1uhbb73VaHqYDQ0NURs41nl9MK6SuWOfeuopo5kHltcGc/JK8X2E/cJ2P/300y0+7g7pc6YZc0knGwsc6NWRVavjv9VU9xyTBpu/OLB2KHKu0tPk9iI/UVJ5se3LUg87zqNtcp+MDWWsZ3eb4zrk9cOa6vUvo312gF9Mv5HeLfbJ3LWSVFqC+QLN+MHuYTqO4zjO34E/MB3HcRynAP7AdBzHcZwCNJ90r4Lm4qX23HNPoy+++OJoG4899pjRzPV66qmnGs1YszFjxhhNj5OeDHPXSnFMHetVMgcuY9cYF0hNb5BtzNsm4yiZT5XQ86RfzDhNxrfmtSEv3+xXAcdVp042zoq5RidNmmR0//79o23yHPBYA7wajhF6cdwe853Sf5Rin5l+/ezZs42mh8xzSN91yy23rLo9SXr11VeN3nHHHY1esGCB0axfO2DAAKMnTJhgNHMes86oFM8RYF8yV3BrwLkczE9aXrDQ6FIve5xJV+txS4o8RW4z8iwxBhvn23NRQr7jAK+vvJ7NLy1J6ma9di1Avcq2zVzzyF+r1ajZSU+TMZeStDrH16ykHfK40tNE7GjSBevDhy0tisdc5JsusffHZGWB2E35L0zHcRzHKYQ/MB3HcRynAP7AdBzHcZwC+APTcRzHcQpQddJPcxNoOCGAkxpGjBgR7xDJnnv37m00JzY0V8yZkzeY/JlJy6U46J8TNjiJh9vkRCT2EydOMFGCFE/oYEA5i1Tvt99+RjMQn5N8GECelzyB+2TCBSZjby0YDM8ix/vvv7/RnCTCQsxSnPyAk6wY9D9lyhSjmdSBxZ05ySevv7fbbjuj33vvPaPXW289o5nQnWO/BoHxvFbmzJkTtYGTcjjpjhNyONGPCeOZhOCaa64xeuTIkVEbLrroIqNZaKCRkz5aAU6gSRDwH7pg8gwnCS2Li8RHk1UwgSbBRL1SvZ04VLOenQyVcKJKF7t+UhNPuCmh0HWUBGCevVeFjnayYXmpPd8RmCwTOuVMVmSSgOXoKyY2aCZhQ/h4rtFJHzumWWhbkpIV6AdcO9EErCbwX5iO4ziOUwB/YDqO4zhOAfyB6TiO4zgFqOphduliE+vS2/v2t79tNBOG03OR4uTpQ4cONZq+Z942KmGyaPqJ9OUk6d577zWaXg+9Jnp99AsZSL/BBhsYTc9Tkvr27Ws0j/Occ84xmkkGeNzNkdeP7Bv6xzzu1oL+Is8pE+wzycC0adOibfJYmBiC/h89y+XLbYA2PeJNNtnE6DwP8+677zb6kEMOMZrjir4rkwDweuS4zSvezKQQ3AaLvt9///1GcxzuvvvuRp922mlG5yXt4HHPmjXLaPrzrUGC5NzqZu99mmvPRXmR7cfI45RiDxNjqIyk86XuSIxOb49JAeD1lRpyfFR+pr316gKSBrAfomTtICqszaLYUuT3JkzYjmQK9HajpPJMvj7TJgBRp3j8BBbKxv0w1BbL4eO/MB3HcRynAP7AdBzHcZwC+APTcRzHcQpQ9cUtvQT6j4y7ZAwWY76kODH5G2+8YfSiRYuMZlJxxlBuv/32RjN+jt6UFPs2jHdjrCj9QsZVMuk825Dno9IvpA/HGFfGwNJjo0fHNjOWVIr7ljqv+HZrcPbZZxvNc/joo48azWPPK/DNc0yfmQW3mfCdmvGH48aNM3qrrbaK2rDrrrsazXPEccd98BhYwJ0+al4CePrtLIRNz/KFF14wmjGSjONkoe28BPAsnP3mm28afdRRR0WfaXHoaa1E3F57eHU4V6FA4YJQ29VqxDgmiHENLMDQ3X6eMZV5xZsbu9nzUVqK+FH6pPAXGWfJwsthsfXd80jo99IXZfJ19EuAn0zfNPnbh/bzOfM1Aq5fFu+O2tAE/gvTcRzHcQrgD0zHcRzHKYA/MB3HcRynAIGFUyu57LLLzELGxzVXcDgv/pCxZozr4nLmtmTRXfo+KxDbxFhSKfZe6d0x3o2eJGHsKL0qFtmVpIEDBxrN42L8GuP6eC7YT+zXGTNmRG2gL0dflLGbDzzwQE512C+f9dZbz4y7DTfc0CynB/3SSy8ZTU9Zij3KF1980Wh6c9wH+4r5ahnfuO+++0Zt+PBD67UwB25z+/zoo4+MZl7XPfbYw+i5c23OTSkeVxwXnHfA64lFq4877jij77rrLqN32mmnqA033nij0Tw39D0XLFjQ4uPukL4/NGOOBaHL7ayfX1qCYtB58YcozpyUcBjzFtptMCYSeVmTjrafSosRO8oiyVIUy5m0RxFqthsxjlHeV+TcjfaZk8eVsZ9Rnl7OlcD6jMNM6lD8HbGjSWPcD1Hf4nxyG+NmXpM75vwXpuM4juMUwB+YjuM4jlMAf2A6juM4TgGqxmEy/o0eC+PfGD83evToaJv08xjryX3Q06Svxs8zj+igQYOiNnCbeV7rukA/kbFp9BOl2LOcOHGi0c35WYztZOwhzwVrfkpx/Bv9rSFDhkSfaQ1Y55GaeWAZw8p4RUkaMGCA0YwXbK42JffBfLb00idNmhS14ZhjjjF6xx13rLoNnh/6ifvss4/RHGeMaZZif725mGRy2223GX3EEUcYTR987Nix0TZYS5TzDJ5//vmqbWgJmstfWkJeWPqTeTlUI68OOU7L/WycrZZV9xPDcsREIj+q6JEqjrMMDaipyfqVyD0b6qzXlzBGkjlYca3mtTOqPYnjjLxeepaMkSU5/RCdX2wzyonb1KYLreU4juM4/+T4A9NxHMdxCuAPTMdxHMcpQFUPk3Fe1PQr6AUy3lGKPUfWoqQvSs/lyCOPNPqMM84wmp5MXt3IYcOGGc1200tim+gn9u/f32jW91u5MvY3pk6dajRrGdLjHDlypNHMyctj4LlhnlFJOvDAA41mX1577bXRZ1oD5rBtaLBeEL0/1mFlXJ8UjzPCc0iP87XXXjOasbocM3k1HbnNV155xWh6lBx3zC3bnMfJfLlSHP/L2FD6oIxPZT5pjjP64Oeff37Uhnnz5hnNGON+/fpFn2ltkln2+gvr2+MKjA1Ercv8jcIHXWPjBZN29OqwD3qaiHkM9FkVx10ms228cIk5VFm7Ep6lGOOIWPwy4lclqbTMjsukA+Iy4ZOWFiHPK/LXRrGikdcb1wWNPEzUN1XvOP90Hv4L03Ecx3EK4A9Mx3EcxymAPzAdx3EcpwBVPUzGu9FHY45V1rZk7T1JOuigg4yuq6urqg8//HCjb7nlFqPpHTGnKtskxXGT9KN69+5t9KGHHmr0M888YzRjP+mzsvalFPu7CxYsMJr1L0855ZSq26RXy7jBe++9N2oD+47+E8/V+PHjo220BPShCY/9gAMOMJrnR4rP+aabbmo041rLyJFJP5Ae5fTp041mflsp9iwZX5qX97iSd99912j6kYz9zPPv6XMy1vq5554z+vrrrzea9TEZF81cs8x/K8VzABgvmtfulibpa8d+ab71zYXYv6QT4va6xecurIBvxhhHeJBretnapLXwMOkFslZlHoH77GlrqEYxkmwza1difeZtDXltYr5Zeo7wgxlHmSDPq7rZfiojx25NTn70QG8W3nxe7GYe/gvTcRzHcQrgD0zHcRzHKYA/MB3HcRynAFXrYZ566qlmIb2+5uLG8nK00jdjnCVzvzIfLX2ejTfe2Gj6W3l+GOP8dtllF6Nvuukmo+knfvDBB0azpmBezCOh90MYh0lfjrUOGVvap08fo9nPUhxD9+CDDxpN3278+PGtUg+zXLYBZeeee65ZznNeKtnvfQ8//HC0Tfp9jLtkbUp6cfTz6Tn36NHD6LzYW14vjEckzOVLD7N7d+tHcZ95uWRZ75LXC2Ne6bOyrxm3yTkGjz76aNSGxx9/3Gj2A2uJnnbaaS1fD3Oz82w9zAZ7rwudre9aRu7RQE9MkhgXSW8OfUtPUwttzCTbwFqXUZxnzt/KvbpiG/AgcRxhoR0PCfzDaJ85saCRhwlflZrxpqUFth+ivK+I0yy3j+s013yCa6HKc0+Sxn041uthOo7jOM4XxR+YjuM4jlMAf2A6juM4TgH8gek4juM4BaiauICTfJiomRN0mCSAk2WkONFy3759jeYkH058YNJxJsXebbfdjN5zzz2jNtx9991Gr0JiXk7geOKJJ4xmgD+TLXCyEydDSfHkJRa1ZpJrJk844YQTjN5oo42M5rnKawOTQrDwNdvQWnzve98zmsWfOfmF4ywvAQDHCQuZM0lDc2Ofk9cY0H/YYYdFbWC7mVyCE2w4MYnjioWWf/Ob3xjNiUxSPLGLE5F++tOfGn3PPfcYzWTsu+++u9FMbsJrKW+fHJucKHTaaadF2/iyiYoxd+9q5OoedsJN7Xx7vUYTdqQ4UTmTiCMxOos9C0UEokTqmIAVZuYUGMA4Lc2zE2jKPTCJBxN0EujQYJMIMJF6bhHrZpLIRwkdmNgeRa7DUiS6x/Zqcop5NzfJp9nlGf4L03Ecx3EK4A9Mx3EcxymAPzAdx3EcpwBVPUxCr4GeF71AFriV4gDwiRMnGs2AcHqWLBZMmEw6DyaFp2/KorcMpD7nnHOMvvPOO42mH5mXpICJBBhAzuB8BojTh+U+mbhg7NixURvOO+88o5mogOfmmGOOibbREjCxBP3ZmTNnGs3EEfTIpNiTZLIKJvymP8/+3XDDDY2+6667jP7xj38ctYGeJBMsfOtb3zL6vvvuM5oF3Fko+8wzzzR6+PDhURu22moro5m4ZPLkyUYzqQD3ycLcTEqflziE1ziPi+e7NUhQSJk+Wu1iFEFmsedV9vqUmi++zMToSVsE3C+0wfbJgoX28xvYguL0KyVJvPesQSKKRXZcR0kFOtr7tZgonb4sl+cQFYxusG1IFlmfVTX0fpG4AG0ufxonBCn1sl56lHSe/nET+C9Mx3EcxymAPzAdx3EcpwD+wHQcx3GcAlT1MFlImT4O/QkmAGesmxR7mIReHZNgs/AvfVOun5cAnsWU6V8xDpDe3i9+8Quj2Q/0qhinKcU+HItOjxkzxmjGs9FD4zEwVjQv+frNN99sdJG+aw2+853vGD116lSjGXe5+eabG92eHofiOMzjjz/e6GuuucboUaNGGU2f9LrrrjP6jDPOMJqeqRQnu99xxx2NZsJ9xijzfHA5vXL6qlLsD1900UVGc0ywzbx+6b3T46SW4oLr9MqvvPLK6DMtDT1I+omRR1lr+zGK45Qi7y20sduMYhiZpHx9ew2HuQurtinpnlPEehGSpzOucjn8fvqg6IekDv3CotbcnqSks/W1mUw9isPkmIHn2VwcZ6m+c9QGFtsOPC56sU3gvzAdx3EcpwD+wHQcx3GcAvgD03Ecx3EKUNXDpC+2ww47GE3Ppbn1JenZZ581mnFb9EmHDh1arYlR3kkWzeX2JenWW281mu3eZpttjGZMJH0g7oM+K3PsSs3n0GUcJds4a9Ysoz/66KNoH5Xk+ZH0ZukNNhfz2lKwwDDj9ubPn2/0008/bTT9Rkm64IILjOY47NTJ5gqlH/irX/2q6vqM/fyXf/mXqA29evWq2s4pU6YYzbFcX2/zfk6fPt1oxgNffPHFURt4zTLWk21gMfSPP/7Y6JNPPtloXr/sZym+Pph/lnl6W4Okvb2mwwrrDwbECuoTW8A9dOsab7MeBZ+R4zRAl+uRM3WVHYP03SK/MM9HRXwpP5Mw321nxKPCJw2C14f4xWR1HI+aF6NaDXqeWoxi3ozLxDyWyBvOawPijxl32xT+C9NxHMdxCuAPTMdxHMcpgD8wHcdxHKcAVc0CenWPPPKI0QMHDjT6008/NZqxgFLs4zDua++99zb6kksuMZp5XQnjMul3SXHNTObEpfdK34Y+Desv0m9kvk1JevXVV41mzCuPg31NT5IxfWxzHm0QF8Z2MwavtTjuuOOMZj5henscU+wrKfYLWSuSHuTZZ59tNH1sesZsw9VXXx21gbUkec7nzrW+GOOBJ02aVLXNrMk5e/bsqA3MDdvQYOP0dt55Z6NZT5Z+/Q033GA08w3n+ZH77LOP0cz//FXkkg3LEaeHuD3GL6psPbBkGWo0Ko67pG+WLLZ9T580wItLmCOZuWhnWX9akkJn+KhodxnzA0qM5YRfmMy199PQDeujjqgkJcjtGuD/03sN9FWRzzaKw6R3m1PbsrGn9f9rFiyN1imC/8J0HMdxnAL4A9NxHMdxCuAPTMdxHMcpQFUPs7ncscz5yTgwemRSnNOUNfzoH3Kf9OboqTDOjD6sFHtBzEXKWDTGJxJ6MPSi2I+S1K9fv6r7YP5U+ouMT+U+GbtI30+S1l/f1tNjvdP33nsv+kxrcPvttxvdnMfFmMk87/z+++83ettttzWafiD1ypU2Zo6xnz/72c+MZv9LcX/TX2ee1gcffNBoxn4yn/Att9xSdbkkDR482Oi33nrLaNamXL7cenO8xpkn+ZlnnjGa/r4U9x39evqmrQLzm9JXa2/9xDWf2PkaUU5WSaUuNqdpgJ+bwCeN6GLPd+SJkk5xzDnrnYaONs9yDY6THmWZcZiMBeUx1MX3W6EWJWM3y2hTaQnqiJL59l4WGGtKj1NSabEdx5EPSo+6CfwXpuM4juMUwB+YjuM4jlMAf2A6juM4TgHWycOkN0Q/gz5bnn/I+EP6MozlpK/G+MMTTzzRaNaVpP8lxflnGcNI744xd8zByuOkF9WzZ8+oDaw9ybhAeprMb0uflfFuPFf0iSSpe/fuVffB2MPW4txzzzWaeVzpu2233XZGv/LKK9E2Oa5Yj/Rb3/qW0TznjB388MMPjR4/frzRzJcqxb40fVHW8aQ/yLhYere8lvLif1mnkzmN33zzTaN5XPQk2UZ6v3nngvGlL7/8ctU2tQoJPCzUu6TvVtPFxvUpxL89Auuywi8s0e9jfct5C+1yeKIJfVPE+ab7RN5VeI6N69t7U81sxEzyHl6q/hsrLI7rwDZXzzLMsftMcF0kOK5SVzsGE8QSK6f+MPuu3MmuU1rhuWQdx3Ec50vDH5iO4ziOUwB/YDqO4zhOAQLjdCoZMWKEWUhvj3FejIHMywnJ2M2HH37Y6C233NLoG2+80Wh6KIzLZNxmnh/C2M3JkycbPWTIEKPpFdGLImxDXk1OerOMk2Qb6MUy/y19UvpZ77//ftSG1ahdx/g3+n4TJ060yStbCI67I444wixnvUz6jXneOdliiy2Mpl9/+eWXG3366acbzXHF/qY/KcX9SR908WIbCzhmzBijec6ffPJJow888ECj82pRMuaYsaBlxKPxGm8u7pKx1+3axbUJCecQMEZ51KhRLT7uDu5xqhlz5SV2PNR8A/MQapFjdUV8vkMzOVBZg7E8294/S/XWdyPlpTZeMeT0dUA7xVjQdvZaKXdCnOZCeJKMV0Q+2yguU4pzuzK/MPuyAfvE+lHcJX3VNbGXyxqZ5Q52G/Rux828JnfM+S9Mx3EcxymAPzAdx3EcpwD+wHQcx3GcAvgD03Ecx3EKUDVxwSabbGI0J69wwg2D5e+6665omywQzQkyxx57rNGc7MJ9MLCayaLz4HFwIgSPixN0mJjggw8+MJoB/3kJADjRiHBi0YgRI4zmRAlOrmAwOCf4SHGwPiedfFUFpDnJaty4cUZzgtSAAQOMzpvgxGLKP//5z43+yU9+YvRFF11kNIs/c1xyTOT1NxNF9Ohhk1JzHxwjHNuceDR9+nSj+/fvH7WB65x11llG33HHHUavt956RnMiUXOT11igWoqTdrzxxhtG8/yNGjUq2saXTUDi8poOSDqAiSUJJ9xggo+keILMQiR4xz5KSFIeTRJaas9/6GMLkOeyzCY3KCOhewnbrFnUTGFlJkLISTofwUk6IGlvlwdOoGLSeU7yYWHuujhJffKxTeJRw/PF890E/gvTcRzHcQrgD0zHcRzHKYA/MB3HcRynAFU9zF122cVoJj4nDNamH5m3Dr2f0aNHG01PhMkTmHi7uQTXktSnTx+jmXCawdkMOGdgPL0iLmdSeiku+MyAbyZooOZxMvk6k3/Tj5akadOmGc0E7XmfaQ1Y3JmFlIcPH270rFmzjM4Llmf/XnLJJUYzKXkJPgkTbrzwwgtGz5gxw2gmW5Ck559/3ujAoG9Aj5Le3gknnGD0X/7yF6NZKF2K/fSrr77a6AsuuMDo//mf/zGa/cBxyHGXlwCeyRGo99133+gzLU3ko62KPWizHInVo4LEkkKjPa5AHxTJ2IXizPQshSLWUZHrRTmJz5EUoLTQesrRcbN4M/zHyFddhaQDHXISVVRJjiNJYgKHntbjDkyGgGMKK1HkekWcPCH0tPMHuM883zMP/4XpOI7jOAXwB6bjOI7jFMAfmI7jOI5TgKrJ1++++26zkF4SPc2xY8caPWHChGibLP5L723ePJsElzGT9IroX/Xq1SvaJ6E3x2TPTChNmJB6++23N5oxeHkJ4OmpMa6PsWnsh+YKcTO+lYnxpbjdTz31lNE8jjlz5rRK8vVBgwaZcUc/kf1J34yJ6yXpyiuvNJr9WYP4MiZnZ6J0+qSdO9tE2Yx3lOKYRPqJHPuMR6V/yKLjeX4h4fXGsVxfX19V029kv7AIAOcgSHHff//73zeaSeYvv/zylk++vuWF9kbYtpin9Rl5Sccxppjwu9zenouaJXa+xZp6O8YCc5gvQCzosjgmMmHydVxL9GrpaTYb81iL7a2Ivd/AZOjYZ7mrvXYaOyMx+nJsE23I22cEP7McRarr7Zj886v5Y85/YTqO4zhOAfyB6TiO4zgF8Aem4ziO4xSgqofpOI7jOE6K/8J0HMdxnAL4A9NxHMdxCuAPTMdxHMcpgD8wHcdxHKcA/sB0HMdxnAL4A9NxHMdxCvD/AZdUdFxh4wwgAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "Detected pixel offset (y, x): [12. 11.]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAS8AAAChCAYAAACWEznEAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAiJ0lEQVR4nO2defxVc/7HX59IQoUsZZlEFFqUkFGjkGhQ0jLSQ/ZCI8ZazG+mmSbzq/lNksmSsStLZMnka2tSCmkaQraRFlMRUchW5/fHPbjv5+d0b5i+5ny9n49HD16de8/53HM+99M9r/NeQpIkchzHyRvVfugBOI7jfBd88XIcJ5f44uU4Ti7xxctxnFzii5fjOLnEFy/HcXLJf+3iFULYPoTwZAhhVQjh/37o8XwXQgjXhBB+/UOPw/l+hAI3hhBWhBCeTf/uzBDCshDCRyGEuhvw2D6H1kHYkHFeIYS3JG0vaY2kjyQ9LGlAkiQfrcd7fy2ppaTjEg9GczYgIYSTJJ0vaTdJKyVNlDQoSZIP0u3tJI2X1DhJko9DCNXT17VJkuT573Hc9pJuS5Jkp+8z/h8rlfHL6+gkSbaQtI8Ki9Gg9XxfA0kvf5eFK4Sw8bd9j/PjJIRwvqT/lXShpDqS2qgw9x4NIWySvqyBpLeSJPk41dtL2lTSS5U8XKeYJEk22B9Jb0k6rEgPl/RQkW4jaYakDyQ9L6l9+vc3SfpC0ucq/GI7TIWF9hJJ/5L0nqS7JG2dvn4XSYmkUyUtlPRk+venSJonaYWkCkkNio6dSOov6fX0+H9R+ks03X56+t5Vkl6W1Cr9+x0k3SPpXUnzJZ1T4vPfJGlo+v/tJS2WdJGkdyQtkdRVUmdJr0l6X9LgovfuL2lmOrYlkq6StEnR9sMlvSrpQ0ljJE2VdFrR9nV+dv/z9Tmqnc6vnvj7LdLre0o6pz7VN3cP4yV9nM6fjyQ9ISlIGple15WS5kpqmu6rhqQ/pfNymaRrJNWUtLmk1ZLWpvv5SNIOPoe+xfXbwJPjLaWLl6Sd0os6KtU7qrAIdVZhYeqY6m150VI9UNLT6X5qSLpW0vh02y7pZLolnRQ1JXWR9IakPSVtLOkySTOK9pdImiRpS0k/SSfrEem2HpLelrRfOjEbqfCvbzVJsyX9j6RNJO0q6U1Jndbx+TnxvkzfW12FxfFdSeMk1ZK0dzqZG6av31eFxX3j9PPNk3Ruum0bFb4k3dLtA1VY7E9Lt5f87P7n6+tzRHpNNs7YdnPR/DpJ0vSibV/Nt41T3SmdF1um82VPSfXTbSMlPSBp6/Q6Pyjp8qI5sbjMGH0OrevcbODJ8ZYK/6KsSi/245K2TLddLOlWvL5CUl9etFTPk3Roka6fnuyvLkwiadei7ZMlnVqkq0n6ROm/Hunr2xZtv0vSJUXjGJjxeQ6QtBB/N0jSjes58VZL2ijVtdIxHFD0+tmSuq5jX+dKmpj+/4mSZhZtC5IWFU28kp/d/3x9XvpIWrqObX+U9Gj6/yep9OJ1iAq/fNpIqobr8rGk3Yr+7kBJ84vmxLddvHwOpX8qw/PqmiRJrfTEN1FhxZcKv2R6hBA++OqPpLYqLEpZNJA0sei181T4Kb990WsW4fWjil7/vgoXaMei1ywt+v9PVLhdkKSdVbg9zRrDDhjzYIyhFO8lSbIm/f/V6X+XFW1f/dUYQgh7hBAmhRCWhhBWShqmb87dDsWfNSnMrsUYZ7nP7kjLJW2zDo+0frq9LEmSPKHCLdlfJL0TQrguhFBb0raSNpM0u+haPJz+/XfF51BKpYVKJEkyVYV/Rf6U/tUiFX55bVn0Z/MkSf64jl0sknQkXr9pkiRvFx8Gr++H19dMkmTGegx3kQpPnrL+fj72WStJks7rsc9vy9WSXpG0e5IktVVYJEO6bYkKt8+SCo/yi7W+32f/MTFT0mcq3Dp9TQhhC0lHqnCnsF4kSXJlkiT7StpL0h4qPABYrsJisnfRdaiTFB5gSXa+bgiq9Byq7DivKyR1DCG0kHSbpKNDCJ1CCBuFEDYNIbQPIazrsfE1kv4QQmggSSGEbUMIXUoc6xpJg0IIe6evrxNC6LGe47xe0gUhhH3TGJ9G6XGflbQqhHBxCKFmOu6mIYT91nO/34ZaKngSH4UQmkg6s2jbQ5KahRC6pr8azpZUr2j79/nsPxqSJPlQ0hBJo0MIR4QQqocQdlHBQlgs6db12U8IYb8QwgFpCMXHKhj8a5MkWStprKSRIYTt0tfuGELolL51maS6IYQ6/9EP9g1Veg5V6uKVJMm7Kpjq/5MkySIVTMHBKpiOi1T412pdYxqlgvH5SAhhlQrm/QEljjVRhUfgd6Q/mV9U4V/T9Rnn3ZL+oIIRukrSfSo82Vwj6SgVwj7mq/Av6/UqPGL/T3OBpN7p8cdKurNofMtVeKgwXIWHHHtJek6FXxHf67P/2EiSZLgKc/BPKnzRn1FhLh6aJMln67mb2ipcoxWSFqhwTUak2y5Wwfh+Or0Wj0lqnB77FRWeXr6Z3p7t8B/5UN9QpefQBg1SdSqHEEI1FX4pnJAkyZQfejxO/sjjHPqvTQ9ySpPebm8ZQqihb7yMp3/gYTk5Iu9zyBev/HKgCk9El0s6WoWnuqtLv8VxDLmeQ37b6DhOLvFfXo7j5BJfvBzHySUlqy8MGzbM3FO+/vrrZvuHH35odJ06NmLgs8/iJ801atQw+sYbbzT6lltuMXrOnDlGv/3220Yfe+yxRh9//PFGDxgwIBrD3nvvbfTHH39s9KOPPmr0wQcfbPROO9lQtEmTJhm93XbbGV29evVoDDwm3/PBBx8Y/fnnnxv973//2+jNNtvM6IMOOsjot956KxrDjjvaYOl333235DFHjhwZVAl07NjRzLtGjRqZ7U2aNDH61FNPNfrZZ5+N9nn55ZcbvdtuNgZ5hx1slMIjjzxidO3atY3edlsbJL9s2TKjt9pqq2gMvKYcw+LFi43m92v//fc3+vnnbTUeft8+/fTTaAwbb2y/8pybm2yyidH77LOP0b/97W+Nfu6554weN26c0V26xKGY/Jz333+/0fx+jRkzJnPe+S8vx3FyiS9ejuPkkpK3jVOnTjWaPzEvvPBCo6+99lqjt95662ifu+yyi9E9e/Y0essttzT6iy++MJo/YydOnGg0bxl4PCm+heJrePvLn8Z/+9vfjO7c2aY28paOtxRS/JP+5ZdfNrpr165Gv/SSrXvH80B4i8H3S9I//vEPo6+88kqjb7/99pLH2FC0atXK6BkzbDrdzjvvbPSwYcOM5vWTpP32sxlcL774otG0M1atWmX05ptvXmLE8W3/DTfcEL3myCNtgDpv4T76yBYY5jH/9S9bK4C3V61btzZ6woQJ0Rh4a9m8eXOjH3/cpnPyGL169TKaawKtiDfffDMaw9y5c43meVi7dm30niz8l5fjOLnEFy/HcXKJL16O4+SSkp7XXnvtVfLNs2bNMpo+Dh+1S3GIAO93+Yj5nXfeMZqPj+k1cP8NGjSIxsBwC2oeg/C80BNjaMUzzzwT7aNjx45Gn3nmmUafd955RtPjuu+++4z++c9/bjRDULLg9Ro8eLDRHTp0KLuPDQF9S4a2TJ482eiFCxcanTVv99hjD6PppT7wwANGc97QxyQ1a9Y0+pe//GX0Gp7vefPmGf3ll18azbAFfp84J+bPn290w4YNozGsWLHCaIbY0NOiV7rRRhsZ3b1795LvP+GEE6IxMGSEn/OVV16J3pOF//JyHCeX+OLlOE4u8cXLcZxcUtLz4r0pfZirrrrKaPpNWTCWJSsWrBh6DfRyhgwZYjTTYp566qlon4yBYmoIY1cYF3byyScbPXr0aKMHDhxoND0xKY6hYhrT+++/bzSvBf2OcrFrTO2SYl+H149jZOrVhoJxWkyDYWoO003ooUjSa6+9ZjTjCZlixHnTrl07oxm/RJ9t2rRp0RiY5kRvdffddzeafhO3Mx6K1y8rLY1/94tf/MLo1attRRzGZDH2jF7hq6++ajTncdY4OVebNm0avScL/+XlOE4u8cXLcZxc4ouX4zi5pGQl1d69e5uNzCmjL0DviCVWpDgfkp7WP//5T6Oz7plLsf32tv8rvSIpjg1jLBnHfe655xrNOCR6F4yXuvTSS6Mx0LfhOJnbyPg3xvQw5o7wM2Ydc8yYMUafddZZRt91112VUhKnS5cuZt7RF120aJHRzNfLuubTp083mr5at26mdaPatGljNEvktGzZ0mjGT9FTk+JxMz7tZz/7mdGMqeI1p7dHv+qTTz6JxsBcxuXLbV9d7nPTTTc1ml5f1ne83Pb27dsbzblZrZr9TTVhwgQvieM4TtXBFy/HcXKJL16O4+QSX7wcx8klJYNUaSi+9957RjNQkgGiTH6VpOuvv97oHj16lBzgIYccYjTNchqnHMNdd90V7bNc/XgGsTIAlMUICR9KsECbFCeDcwwjRowwmkYqDWc+NOjbt6/RrMuf9Z5yDyIqizfeeMNoBojyYQcTebPmHa85a9Kz+OCDDz5oNJOqadD//e9/N7pWrVrRGBjgye8Pk6QvvvhioxkUzmKS/EwnnXRSNIabb77ZaBZhZMDoOeecYzS/v3feeafRW2yxhdF8yCfFDxZI1sOOLPyXl+M4ucQXL8dxcokvXo7j5JKSQaqjRo0yGxlAWu4engXzpNi7YWAr73cZFMd79MMPP9xoBs5m9SvkMfk5fvrTnxp99913G00Piw02eN+f1byBHslFF11kNIMimZRbrjgePZ6sMZQLAGZw6FVXXVUpQart27c3827NmjVmO3sB0nfJ6lfIv2NSM5tXMHCSBQVYOJBJ7kuXLo3GQA/qscceM5oJ5pzbdevWNZpBrWxmQ09Mkvr162f0kiVLSh6D3h3nFc8LvwtZxUBZrJPnhYU5x40b50GqjuNUHXzxchwnl/ji5ThOLikZ58ViaPR66C/RQ8kqCrfNNtuUfA+PyQL+9N3oBbVo0cLosWPHRmNgjBVjW9gAlsmljH9iQ1Pun802JGnBggVG81xyTCxoyGYMjBEiTF6WYt+GvlhWMndl0LZtW6NfeOEFozlnmIDMZGJJOvTQQ42ePXu20ZxXjD9kjCOTpDknfvKTn0RjOOWUU4zmuBlbtueeexpN7/Xhhx82etdddzWazZCl2NOiv8TGtowPrKioMJoNN0jWteC4br31VqOzPMss/JeX4zi5xBcvx3FyiS9ejuPkkpJxXgMGDDAb6YGwYSzzmLIaADC+iT4L30PfjPE4zAmkd8SiflIcf8bYFXpBvG9nwUN+BuY20leQpMaNG0d/VwrGqzHOix4YcyNvuummaJ9sysH4thkzZhhdUVFRKXFeLIJZr149s52eSKdOnYxmDqAUN2BlET56WIxPYoE8zhF+F+hpSrFfxHjAI444wmj6mMzx5JjYnIbxjJK02267Gc0cTZ5bNuulJ8aCia1atTL6sMMOi8bw5z//2WjG3LHJ8JAhQzzOy3GcqoMvXo7j5BJfvBzHySUlPa9u3bqZjaeddprZTh+FtbSy6l7RSyiXo1cufoleA+/ZGVcmxQ0aeI/NWDPCGmFsGEsPLKvhK6GvVu5zMKds5cqVRjO2LCvHk74Ma4Txc1RWbuNpp51m5h19k/vvv99oekchxMNkE1jOQ84JxsUxxpF+EptrMD5Rir8vffr0MZpNQthst3///kafeOKJRjPXkdulOGbu6aefNvqYY44xmr7ngQceaDTj49jE9o477ojGwJxNxqtxXs6dO9c9L8dxqg6+eDmOk0t88XIcJ5eU9Lx69uy57o2K6/DcdtttRjNuRYrrpDOHj94M474Yk0VvhzmATzzxRDSGY4891mj6bjwGcxvL5TIyNo1eoCRNmTLFaHouRx11lNE8b+XGyBr2U6dOjcbAGDh6efQeKivO65BDDjHzjjXeJkyYYPTQoUONZs0qKW50evbZZxvNecL8VsY70YfjHGDNKimu/87cRn4Xa9asaTR9TXqp9PGY6yjF9f3pGdMTo/fH+DbW+2JNv6x4N8YkMm+Un/Ptt992z8txnKqDL16O4+QSX7wcx8klJT2vDh06mI3169c321njnHFBvLeV4ntgxsswf5K5i/TZmKfI3MasmmI8JmuC0aNiniD7BhLGjb322mvRa7baaiuj6e2xdj99AHqF9O0YN8b+l1J87srFOs2aNatSPC/GF7J+F+tiTZo0yWjGXEnx3GRvSF4P5gDS22GPxN69exs9atSoaAyM3eMx5syZYzRrgjVt2tRo9mCk58VelFI81/n9YNxW9+7djX7yySeNZl4ivxsvvvhiNAauIzwvHPe0adPc83Icp+rgi5fjOLnEFy/HcXJJyRr2pFw8ButBDR48ONoHe9U99NBDRp9xxhlGM5Zl+PDhRtMT4z38XnvtFY2BMTust8WcTMbGMO6Iml4Sx5i1T8ZplcvppEdGf5FxYFm5jRxDVv21HwLOK+Yu3nvvvUY/88wzRrP2uxRfA35W5kNyjtC74f6GDBlidFYNe/qS9HfZ65GeI68hfTrWiMvqHTl37lyjW7ZsaTTre7H+XvPmzY1mjftyObhS7CnzXDJ3dV34Ly/HcXKJL16O4+QSX7wcx8klvng5jpNLShr25cxvmnk0JHv27BkfEImgbGZBU5KFAZkYSuOViaFMaJbiAE+arTTguU8+ROB5oumZ1QiBZiyDB0eOHGl0hw4djGbQJQ16BgtmBcrymAyuZaJ2ZcHARzZ9YANZGryNGjWK9slAVz4gYYDnrFmzjGYALwsF0qDPOt/NmjUzmk0/2GiEyd6c+2x+we/KsmXLojHQUOcDM5rpfEjHZHIGnI4ePdrorKKMgwYNMppFCNasWRO9Jwv/5eU4Ti7xxctxnFzii5fjOLnkWxUjpCfCIEcmEzPhUoqbeHTu3JnHNJr3+YT+Ev0n+jpSnHxKb4DeBPfB9/Nz77TTTkbTL5Hihg4MSuV2ejD05cp5XlnnsVyTDp7b6667rlISs5s1a2bmHZurUtOHoU8jxdeQ3isDRnnNVq9ebTTPL72iLM+LyeBHHnmk0byG9LAY8Mk5wAYcTKKW4oRynrvzzz/f6Hvuucdo+sWcpxwzCypK0pIlS4zmuaef+9hjj3lituM4VQdfvBzHySW+eDmOk0tKxnnx3pOFABnXxRiPrIavTFpmo4MPP/zQaCYc8567RYsWRtMbymogS6+A8TT0Q+hxMW6LCekcA71CKfaXmAjMGDrG2NGTKefLMVZNis8tdVYhx8pg4MCBRvMaPvLII0bzs2cVi+Q1pi/J4o1MBqdmfNPkyZONbtKkSTSG/fff32heI847HoOfgcVAGzZsaHRWcjiT2umB0eNiU1rGYDFOjEUbs5LDWYTx1VdfNbpXr17Re7LwX16O4+QSX7wcx8klvng5jpNLSsZ5DRkyxGxkPli54nXMAZTiWBbGinE742cYE0KfgLFlderUicZAr45eD+Np6GERxqbR22DBNiluXFuu8QjjhngteJ54XrOaUtDHoY/G2LB77723UuK86tWrZ+bdzjvvbLbTs5w9e7bR9CCl2NNic2J6OTwGzxXzJ1lAMavp7MKFC41mTma5Yy5evNhoxmgdeOCBRi9fvjwaA+cV5wV9an6fWACxT58+Ro8bN87oVq1aRWMYO3as0bw29MlWrFjhcV6O41QdfPFyHCeX+OLlOE4uKRnnxfga3pMzvobxOcOGDYv2Sf+HsWTlcvbow/D98+bNM7p9+/bRGLjPLG/u20D/ibEv9J+k2OOaMmWK0eX8D8aOMbaJ14I1y6Q4voZ+CHMfKwvWqaJmHihj5BgPJcWNIxiPVK62Fo/B5sb0XtkURJKOP/54o9n8gvvg9aH/1K5dO6M5zxgzKcV+bLmYR3LLLbcYzaa09E2vuOKKaB+shUZfeubMmSXH8BX+y8txnFzii5fjOLnEFy/HcXJJyTivo446at0bFXsP9EwYTyXFMVWspUUfjff1PXr0MPqss84ymv5TFl27djX6pptuMpreA+/j6T+xySnrFWU13iSMZaInxrxS5ojS72CsE/PepPhc81zOmTOHb6mUOK+GDRuaecd5xrpV06dPN5pxQ5JUt25do+nt8D30xOgflcuHrVGjRjQG1nNnTBQ9LV4f7pNj4PuZvynFTZ0Ze8bmuox/4/eXsWdcA+hJS/F3lPXXWCNs0qRJHuflOE7VwRcvx3FyiS9ejuPkkpJxXvQaeM/OnD/W5mLtIEnq1KmT0byPp+7SpYvR5fwp5vhxTFJ8z01Pir0k2VeOHgtjy+jLsXaXFPuBrBfP+l2s/c99sjYU45LuvvvuaAw8d/QreK0qKiqifWwIGHtE+NkPO+wwo3l9pPia77rrrkYzbo69I+kf0ct55ZVXjGa+pSS98MILRjN+LSsPtxj2r2R+LGPLOCek2BdjLOeMGTOMvvbaa41mfS/GXTL3kZ6YFHvG9BOzxp2F//JyHCeX+OLlOE4u8cXLcZxcUjLO64wzzjAb6Q2Vi0vJyhmkz8JcK+YiMj+SvgBr2tMPyfJPWJu9devWRl9//fVG039asGCB0ayJlBVTRcr1o2ScF30c1mpi7Brr7PM8S3H82sSJE42mz1NRUVEpcV5r16418+5Xv/qV2c5rzj6ODz74YLRP+kOMzWNtLXo39H/pUTKO7LPPPovGwO9LVn/JYphbSs+LNex5zKzcRtbr4vdl1apVRtOX47lmfCE9afYbkKQnnnjCaJ4H1kLr37+/x3k5jlN18MXLcZxc4ouX4zi5xBcvx3FySckgVRr0TOKkuc6AUBrdUtxIgkmYNOhpWjIhmYmkBxxwgNFZzRjGjx9vNJNNab6yASyDORlYywcVfJAhxQ8eWCCRidgMlO3bt6/RbDDKa5U1BgYAs4gix1BZnHzyyUazkCCNa86zrGBPzhMWxWRAbrm5zwdPDN485phjojFw3AwkpjnOhwqcVyzad8MNNxjNhxBS/FCGDxEuu+wyo++8806jmajdpk0boxnIzu9S1jE5N2ny9+/fP9qH5L+8HMfJKb54OY6TS3zxchwnl5T0vAjvTemR0DtiAwEpDvZj4wkG/9HjevPNN0uOkYmmWTBhnD7bhRdeaDSD5s477zyjb7/9dqPLNQ2R4qBRBgsyEJPBgPTteEwGqWY1QrjggguMZlAqrw0bSGwoGERMP49F+BgkTE9Fij0sBiYzGZh+Ls8vi0eysOD5558fjYEeFoNpWWRxwoQJRrPJLAsoDhgwwOhu3bpFY2jSpInRDFJ/9tlnjWYAKY/JYqLlGupI8Xecn4vXe134Ly/HcXKJL16O4+QSX7wcx8klJT0vFuXjfT/vZ5kczFgaKfa8CL0dJsiyiBx9Nr4+KzmchfnodzDOiF7Q73//e6N5HuhtZDVjoG/DAobDhw83mvEy9Fz4GRiLlpWYfeONNxq9PueuMjj66KONZiMQxnU1atTI6Jo1a0b7ZJzXiSeeaPTo0aON7tevn9H01caMGWM0m5dkNV1hIjybzjIZnzGQvB7cTm+VPpwU+4mDBg0ymnOCY+b3l14tPbGsZigs3klvdcSIEdF7svBfXo7j5BJfvBzHySW+eDmOk0tKel70UfbZZx+jeY9e7vWS9NRTTxnNuBD6ap07dy41xCgPigXYsppe3nzzzUZz3GzGwJgr+gY8Bn25rCYg5XI6GafFMbKxbVajg2Ky/Ct6efSSysXUbShYrI5xQe+//77R06ZNM5r+lCRddNFFRnMesikE/aOrr7665OsZW3bcccdFY9h2221LjpMNbTiXa9eubTSbfjDecPDgwdEY+J1lLBnHwMKabBxz6qmnGs3vL8+zFH8/mA/JvNF14b+8HMfJJb54OY6TS3zxchwnl5S8uaS3M2nSJKPZ1ICNPRlrJMX3/Ywradu2rdGXXnqp0cwzJIz7oj8ixTW/mKNJr473+byvZ/0o+lPM/5KkuXPnGs2YOn4Onmt6WIwZ4pizqF69utEcN2N8Kos+ffoYzfxWekGcUzxXUuwvsdYVPauBAwcaTd+THiPHMGrUqGgMrIXFa86mKow3ZFNZjpk1xZYuXRqNgbmKbLix7777Gs16ePR3r7vuOqOZ/5rlX7Vr185o5iN7bqPjOFUaX7wcx8klvng5jpNLSjadbdu2rdnIppZsBso4k6w4L+YBMjeKcDt9GN6D04ugbyfF9+Hc58qVK41m/BPzLxmHRO+CHpkkNWjQwOiDDjrIaObzMa6L8W88JmupZTUgrV+/vtH0zRij89JLL1VK09mBAweaeVfOE2HuaIcOHaJ93nPPPUazYSv9I2oekzXhfvOb3xj95JNPRmPg+aYfy3nFecjYsnJNoOmpSXFNPvaUYG0tzkt+bp77evXqGZ3VT4B5vfR3WYesd+/e3nTWcZyqgy9ejuPkEl+8HMfJJSU9r5YtW5qNjNsqFyeUVUOK8U28j2esGO/RmSfI3oJ8fVZNe+ZDMpdx/vz5RtMTY04gt5ergSSVr51Fn43eHXPU6OPxvGblPjLHjD4avbrBgwdXiue1YMECM++YV0ifplmzZkYznkqSZsyYYTR9FtaP5/nq2bOn0fR3ef55bqXYx6SHzDpk9JP4XWX+ZTlPTIo9YubYMpaM85Sfk2MMwU4RXispjl97/vnnje7Vq5fRBx98sHtejuNUHXzxchwnl/ji5ThOLinpefXs2dNspBfEOur0YbJylFifnL3rGjdubPTYsWONZtwIc/gYR8Z7eim+72evuo4dOxpNb4HeBeEYsmqK0ZtjHBbHQO+O+ZjbbLON0fQ/suLpGFfEPDb6j1OmTKkUz4vzrnv37mY7632tT2wf2WOPPYymvzd06FCjzzzzTKM5r3i+6WdJ8flcuHCh0YzrYh8DXvOpU6caffjhhxudVUuL9f4Za7Z27Vqj+R2nxzV9+nSj6YNn1bAnjC9kDGS/fv3c83Icp+rgi5fjOLnEFy/HcXKJL16O4+SSksUIWXyfxjPNcgaxZjW9ZLFBmtu9e/c2mkY1j1FRUWH06tWro2MSfg6amPxcNNcZdLpgwQKjWbiOWoofEhA+FGCQJE1OGqMM/KM5L8WBmTSMf6hihHxAMnnyZKP5cKN58+ZGZz2cYGG+3/3ud0ZfcsklRrMZKwNfOS85J7LON5PB2UiYx+Ac4dzmQwM25GDz46zXnHPOOUbfdtttRjPRmg8Byj14YrFDKQ58ZeA5rx8bAH+F//JyHCeX+OLlOE4u8cXLcZxcUtLzat26tdFZjUuLYWAe/aus19ArYPNV3kMzUJZJueWSX6W4oSsL9zEQj8GFDIKkt8DtTFiX4uKBDO5jMC41P2e5xGD6l5I0b948o5msnPWeyqBp06ZGT5w40WgWq2OhxqzASJ5fNnZhEnO1avbfdQZXs9kFG3wwsFaSZs6caTSTmAk9LXpBffv2Nfrxxx83mknWUuy/slEIm/P+9a9/NZrngfOQ8y6r+QwDYanLNdn5Cv/l5ThOLvHFy3GcXOKLl+M4uaRkYvb48ePNRnoP9MCuuOIKox999NFonyxORq+GTTEZk0VvgX4Hm39mQS+HiaBMNiVMVm3RooXRjPHJSg6nB8O4Ica+8DyUK+rI+DkmzUvxuNk0gp9j2bJllZKY3b59ezPv6D/xfNJnYVK7JI0YMcJons+NNtrIaCZuM4mavlqtWrWMZjyVFMc80X/i3Ge8G/0mNs7N8pcIv2+cy7Vr1y6p6U/xvLBAQFbjF577008/3WgmoA8dOtQTsx3HqTr44uU4Ti7xxctxnFxS0vNyHMf5b8V/eTmOk0t88XIcJ5f44uU4Ti7xxctxnFzii5fjOLnEFy/HcXLJ/wOKWh4h4u0UKgAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "Detected subpixel offset (y, x): [-12.031 -10.93 ]\n" ] } ], "source": [ "image = a[0,:,:]\n", "offset_image = b[0,:,:]\n", "# pixel precision first\n", "shift, error, diffphase = phase_cross_correlation(image, offset_image)\n", "\n", "fig = plt.figure(figsize=(8, 3))\n", "ax1 = plt.subplot(1, 3, 1)\n", "ax2 = plt.subplot(1, 3, 2, sharex=ax1, sharey=ax1)\n", "ax3 = plt.subplot(1, 3, 3)\n", "\n", "ax1.imshow(image, cmap='gray')\n", "ax1.set_axis_off()\n", "ax1.set_title('Reference image')\n", "\n", "ax2.imshow(offset_image.real, cmap='gray')\n", "ax2.set_axis_off()\n", "ax2.set_title('Offset image')\n", "\n", "# Show the output of a cross-correlation to show what the algorithm is\n", "# doing behind the scenes\n", "image_product = np.fft.fft2(image).conj() * np.fft.fft2(offset_image)\n", "cc_image = np.fft.fftshift(np.fft.ifft2(image_product))\n", "ax3.imshow(cc_image.real)\n", "ax3.set_axis_off()\n", "ax3.set_title(\"Cross-correlation\")\n", "\n", "plt.show()\n", "\n", "print(\"Detected pixel offset (y, x): {}\".format(shift))\n", "\n", "# subpixel precision\n", "shift, error, diffphase = phase_cross_correlation(offset_image, image, upsample_factor=1000)\n", "\n", "fig = plt.figure(figsize=(8, 3))\n", "ax1 = plt.subplot(1, 3, 1)\n", "ax2 = plt.subplot(1, 3, 2, sharex=ax1, sharey=ax1)\n", "# ax3 = plt.subplot(1, 3, 3)\n", "\n", "ax1.imshow(image, cmap='gray')\n", "ax1.set_axis_off()\n", "ax1.set_title('Reference image')\n", "\n", "ax2.imshow(offset_image.real, cmap='gray')\n", "ax2.set_axis_off()\n", "ax2.set_title('Offset image')\n", "\n", "# Calculate the upsampled DFT, again to show what the algorithm is doing\n", "# behind the scenes. Constants correspond to calculated values in routine.\n", "# See source code for details.\n", "# cc_image = _upsampled_dft(image_product, 150, 100, (shift*100)+75).conj()\n", "# ax3.imshow(cc_image.real)\n", "# ax3.set_axis_off()\n", "# ax3.set_title(\"Supersampled XC sub-area\")\n", "\n", "\n", "plt.show()\n", "\n", "print(\"Detected subpixel offset (y, x): {}\".format(shift))" ] }, { "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 [conda env:test_robustness]", "language": "python", "name": "conda-env-test_robustness-py" }, "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }