Gen_Data.ipynb 33.5 KB
 Maximilian Schanner committed May 06, 2020 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 ``````{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Imports\n", "import sys\n", "import os\n", "# relative import\n", "sys.path.append(os.path.abspath('') + '/../')\n", "\n", "import numpy as np\n", "import pandas as pd\n", "\n", "from matplotlib import pyplot as plt\n", "\n", "from scipy.stats import uniform, gamma\n", "\n", "import pyfield\n", "from corbass.utils import load, nez2dif" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# seed for reproducability\n", `````` Maximilian Schanner committed May 06, 2020 33 34 `````` "seed = 161\n", "np.random.seed(seed)" `````` Maximilian Schanner committed May 06, 2020 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 `````` ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Generate synthetic data for tests\n", "\n", "This notebook serves the purpose of generating synthetic data for testing the `CORBASS` algorithm. We therefore generate data from a given set of Gauss coefficients and add synthetic errors from the Fisher-von Mises and the gamma distribution." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# some basic parameters\n", "# the name of the output file\n", `````` Maximilian Schanner committed May 06, 2020 54 `````` "out = '../dat/synth_data_clean_complete.csv'\n", `````` Maximilian Schanner committed May 06, 2020 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 `````` "# switch for using locations and incompleteness structure from the example data file\n", "real_locs = False\n", "# the number of records to be generated, is only used if real_locs is False\n", "n_points = 412\n", "# the fraction of incomplete records, works as a switch if real_locs is False\n", "r_inc = 0.\n", "# switch for corrupting the data by noise\n", "noise = False\n", "# the error levels to be stored\n", "ddec = 4.5\n", "dinc = 4.5\n", "dint = 8250\n", "# the average concentration parameter from GEOMAGIA for the interval [750, today] is 650\n", "kappa = 650\n", "\n", "header = f\"# This file was produced using the notebook Gen_Data.ipynb with the following parameters:\\n\" \\\n", " + f\"# real_locs={real_locs}, n_points={n_points:d}, r_inc={r_inc:.2f}, noise={noise}, ddec={ddec:.2f}, \" \\\n", `````` Maximilian Schanner committed May 06, 2020 72 `````` " + f\"dinc={dinc:.2f}, dint={dint:d}, kappa={kappa:.1f}, seed={seed:d}\\n\"" `````` Maximilian Schanner committed May 06, 2020 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 `````` ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Sampling from the Fisher-von Mises distribution" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# The sampling process involves some rotations, thus we first define some convenience functions\n", "def angles(vec):\n", " return np.arctan2(vec[1], vec[0]), \\\n", " np.pi/2 - np.arctan2(vec[2], np.sqrt(vec[0]**2 + vec[1]**2))\n", "\n", "\n", "def rot_z(ang):\n", " return np.array([[np.cos(ang), np.sin(ang), 0],\n", " [-np.sin(ang), np.cos(ang), 0],\n", " [0, 0, 1]])\n", "\n", "\n", "def rot_y(ang):\n", " return np.array([[np.cos(ang), 0, np.sin(ang)],\n", " [0, 1, 0],\n", " [-np.sin(ang), 0, np.cos(ang)]])\n", "\n", "\n", "def rotator(vec):\n", " vec = np.asarray(vec)\n", " p, t = angles(vec)\n", " return np.dot(rot_y(t).T, rot_z(p))" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "# This cell contains the sampling procedure\n", "def sample_Fisher(n, mu=(0, 0, 1), kappa=20):\n", " \"\"\" Generate samples from the Fisher distribution\n", " \n", " Parameters:\n", " -----------\n", " n : int\n", " The number of samples to be generated\n", " mu : array-like of length 3, optional\n", " A vector pointing towards the center of the distribution. Its length is ignored.\n", " kappa : float, optional\n", " The concentration parameter.\n", " \n", " Returns:\n", " --------\n", " numpy array of shape (3, n) containing the sampled vectors\n", " \n", " Reference:\n", " ----------\n", " [1]: W. Jakob, \"Numerically stable sampling of the von Mises \n", " Fisher distribution on S^2 (and other tricks)\",\n", " http://www.mitsuba-renderer.org/~wenzel/files/vmf.pdf,\n", " 2015\n", " \"\"\"\n", " if kappa <= 0:\n", " raise ValueError(f\"The concentration parameter has to be positive, but kappa={kappa} was given.\\n\"\n", " f\"For kappa=0 use a uniform sampler on the sphere.\")\n", " trafo_mat = rotator(mu)\n", " \n", " # sample from the uniform circle, V in [1]\n", " angles = uniform.rvs(scale=2*np.pi, size=n)\n", " vs = np.array([np.cos(angles),\n", " np.sin(angles)])\n", " # sample W in [1] via inverse cdf sampling\n", " def inv_cdf(x):\n", " return 1 + np.log(x + (1-x)*np.exp(-2*kappa))/kappa\n", " \n", " unis = uniform.rvs(size=n)\n", " ws = inv_cdf(unis)\n", " ret = np.sqrt(1-ws**2)*vs\n", " res = np.einsum(\"i...,ij->j...\", np.array([ret[0], ret[1], ws]), trafo_mat)\n", " if n == 1:\n", " return res.flatten()\n", " else:\n", " return res" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We draw some samples and plot the result, to eye-check whether the procedure works:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { `````` Maximilian Schanner committed May 06, 2020 179 `````` "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlsAAAE9CAYAAAA8gKerAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deZCkd33n+ff3OfKozDq6DlW11HfrAIlDQFsWluSFMCCsYJDttRfYCduL/+hlY4iZ8Y43hh1iYcYTs7O7c+yGjxjQjhkbdgbLMSwYBpnLa4waEKglBOimT3Wru6vrzqqsPJ7jt39kVqm6VdVdraqnsjv784ooKo+nnudXT3UoP/yO78+cc4iIiIhINrxON0BERESkmylsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIaCTjcgC8PDw27Pnj2dboaIiIhcJ5544olJ59zIau91Zdjas2cPhw8f7nQzRERE5DphZifXeq+jw4hm9hkzO29mT6/x/jvMbM7Mnmp/fWKr2ygiIiKyEZ3u2fpT4I+Az17imEedc+/bmuaIiIiIbK6O9mw5574DTHeyDSIiIiJZuhZWI77dzH5sZn9lZnd0ujEiIiIiV6LTw4iX8ySw2zm3YGYPAF8CblntQDM7CBwE2LVr19a1UEREROQSruqeLedcxTm30H78CBCa2fAaxz7knDvgnDswMrLqyksRERGRLXdVhy0zGzMzaz++i1Z7pzrbKhEREZH16+gwopl9HngHMGxmp4FPAiGAc+5TwK8D/4OZxUAN+KBzznWouSIiIiJXrKNhyzn3ocu8/0e0SkOIiIiIXJOu6mFEERERkWudwpaIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhkSGFLREREJEMKWyIiIiIZUtgSERERyZDCloiIiEiGFLZEREREMqSwJSIiIpIhhS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREcmQwpaIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhkSGFLREREJEMKWyIiIiIZUtgSERERyZDCloiIiEiGFLZEREREMtTRsGVmnzGz82b29Brvm5n9gZkdMbOfmNlbt7qNIiIiIhvR6Z6tPwXee4n3fxm4pf11EPh3W9AmERERkU3T0bDlnPsOMH2JQx4EPutaHgMGzGz71rROREREZOM63bN1OTcBp1Y8P91+7VXM7KCZHTazwxMTE1vSOBEREZHLudrDlq3ymlvtQOfcQ865A865AyMjIxk3S0RERGR9rvawdRrYueL5DuBMh9oiIiIicsWu9rD1ZeC32qsS7wbmnHNnO90oERERkfUKOnlxM/s88A5g2MxOA58EQgDn3KeAR4AHgCPAIvDhzrRURERE5LXpaNhyzn3oMu874O9tUXNERERENt3VPowoIiIick1T2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREcmQwpaIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhkSGFLREREJEMKWyIiIiIZUtgSERERyZDCloiIiEiGFLZEREREMqSwJSIiIpIhhS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREclQR8OWmb3XzF4wsyNm9rFV3n+Hmc2Z2VPtr090op0iIiIir1XQqQubmQ/8MfBu4DTwuJl92Tn37EWHPuqce9+WN1BERERkE3QsbAF3AUecc8cAzOzPgQeBi8OWiIiIyKY4PrHAoaNTjM/VGe0vcO/+IfaOlDO9ZieHEW8CTq14frr92sXebmY/NrO/MrM7tqZpIiIi0m2OTyzw8OHTVOsxY/0FqvWYhw+f5vjEQqbX7WTYslVecxc9fxLY7Zx7M/CHwJfWPJnZQTM7bGaHJyYmNrGZIiIi0g0OHZ1ioBjSVwzxzOgrhgwUQw4dncr0up0MW6eBnSue7wDOrDzAOVdxzi20Hz8ChGY2vNrJnHMPOecOOOcOjIyMZNVmERERuUaNz9UpFy6cQVUuBIzP1TO9bifnbD0O3GJme4GXgQ8C/+3KA8xsDBh3zjkzu4tWOLxs/GzEKZ977OSWjseKiIjI1W20v8BCPaavGC6/tlCPGe0vZHrdjvVsOedi4KPA14HngL9wzj1jZh8xs4+0D/t14Gkz+zHwB8AHnXMXDzW+yky1ueXjsSIiInJ1u3f/ELO1iEotInWOSi1ithZx7/6hTK9r68gu15zdt73R/cvPfnX5eaUWUSoE/ObduzvYKhEREem0rFYjmtkTzrkDq73XyWHEzPjehXPvy4WAcxmPx4qIiMjVb+9Ied3harOCWVdu15OkF/bWbcV4rIiIiHSPzSwT0ZU9W0naGoctFwIW6jGztYj77xh9zefrRAE0ERER6ZyVZSKA5e+Hjk5dcQboyrC1rZSj1B46HO0vcP8do2vemMsFqaVkO1AMGWuvYnj48Gk+cGCHApeIiMg1Zj0dKMcnFvjWs+OAo78nx80jJYbLhQumJV18Hgty+bWu2ZVhKx9465oMv54gtZnJVkRERDpnrc/9e/YNcmKmxvhcncCDc5UG+cADg2aUcvjELAf2DJDzfUb7C6uex+vpX7UOKHRp2Bqv1Pk/v/kCBsRpq67Gnm3F5Ru5lGTXE6TG5+qMXTTfSxPuRURErj2rfe5PVxt8+tHj3HvzMGP9Bb7zswkqixGvHytzdHKRYi6gGHo8fabCvuEy998xunp+SJN4ret2ZdgC+MGxaRzwC/sHOTVV5QtPnOau3dvYNVxaTrIL9YjXbe+74OcuDlKdKoAmIiIim2u1DpSzlTqpc8uf8804JfDg8KlZioHPTK1JT+BTyPnLI1///tHjVOoR842YvkLIzSMlnEuTta7blWGrGaf09+Qw4OjkIgADxZBz8w32jJRbSXahwVOnZjk2UaWQay3KTB3kAo/bx3qXz3Xv/iEePnwaYNMm3IuIiMjWW60DZXKhwXAp335c5+xsjbNzNcLA584d/Qz6eWYWm+wZ7OHQ0Sn+5NBxHjs+xVBPjnIh4GfjC/zopRm8fGnbWtft2tIPhcAjH3pU2pVi+4oB8/UIgMn5Os+Nz5OmKadnqjx+YpqnTs1SWWxQWYw4V2ksL+3cO1LmAwd2LE+4LxUCTY4XERG5Bq1WQT7wPLb3FZhcqHP4xCxxmuKbhwecnFqkWo+oN2OeOVuhWo+Zq0cMlXK8NL3I0y/PMbPYYKEeYUGuuNZ1u7Jny/eMepxivDIPq7IYLT8+MlGlESU0EwgDj1IuoBGnvDRT5++8cYzR/p4L5m2tLIC2tPrgC0++rDIQIiIi15ClDpRDR6eWKxYcvHcv3z02zdNnKhRDD+egGHpEqWNyocFiM2HntiKebzTjhOfPzmNAPUmoNxNK+ZBiLgDn0rWu25VhKxd4zC02l+dsVRsJL00tcusNZVLneGmqyunZGsXAp5EkbOvJUQh85hsxU4sx+0dXnwD/6Avn+fSjx0lTx1A5RzNKeHimpp4uERGRa8RqFeR3DPbwT7/yDIYR+B5mCSO9BQIPpqtNnj83T0/OY6LSwOHI+z7NOCVNHYEHzSQFsFUvSJeGrYVGTCH02D3YQ5RAT87nzTf18cL5BX50aoapapPQ89hWCjk/75haaDJQDCnnfebr0fIE+JU1NHwPvntkknI+pL+Uox6nPD8+z+tGe1UGQkRE5Bq2d6TMu24fo1qPefzEFMcmUppxwvhixGIjxhlYBOWCEcUpcZKQpCkpUG0mFEMfnLu+JsiHvkcjSTk2tcj739jH0aka+dCnrxjw0nSV1DkwRz1KGegJGJ9rMFOLGB3IE/oes7WIO7b38q++/gLHJqssNmMacUoUJbxl9zbMrHVjaa1iyLUfi4iIyLXh4qKke7YVeeSZcxybqNKIYmbrMUnqKOdDSjmPc5UGaZpSygc049a8LkeKA4bKOcC5ta7VlWGrHiVUGwm4mP/1a8+T8z2asWOwHOKcoyf0qTYT6nFC4HuM9OVYbCbUGikHdvfyK3fexH/47nGeemmGUiGkEHhMzDdoxilPnJzhwO5BykuvLzS4a99Qp39lERERWYfjEwt86amXefTFSYb78tzYl+fF8/N8fnpxOS81EijnAupxynA5R5I6At/jbKVOIWh1sBRDn2rqSJPWCBnYmosOuzJseWYsNmLGK3VqzQTPg8D3OTfXwPeMkd4chZzP3qEeSoWQifk6t5by/KN337o8HPi9Y1OUCgG+Z0xXI0Lfw6UplVrM4yenGest0N8Tkgt87t2vsCUiInK1W6r8fmxigeG+PLVGzLeerzDWl2dqoUkzSegthEDKcG+R6WqTifkGqYO+vM9MLWKhEWFm5HyPIPAIfSN/mRGuriz94BzLKwhiB1ECSZoSpyn1KOZ8pUG1HvPi+Dz1KGHPUOmCoAVLRc2M+XqM7xmh79FMwAFRknJqZpGTk1UefNOY5muJiIhcA5YqvzeTlGLgM1ePCTzj5HSNOE0JfZ984JO6VhWDZpQwV4twLqUWJXgGgd9asVhtJCw2E+YWY8Zna2B2fU2Qb8QJc7WYpcFTB0SxwzPa4SuhNx9QCANwq4+y7h8pceT8Ao04JR94NOLWvLcwgFIuIHUwNlDgtLbtERERueqstuH0UgX5vkJIPU6pRwlRnOKcI/RbvVM9oU8h8Floxgz25GgmCbUoJU4dOd9jWzFgshqRAktBI3Vgnp9bqy3d2bPF8u+/LKUVtJben6tHvPHGXt5x2w3sGSpx6OjUBcf/zi/spZwPMaDWTGhECaFv3NTfw0hvgZu2FRntL/LUqdnsfyERERFZt6Xhwmo9Zqy/QLW9TZ/vtbbcu3mkRK0Zk6Qp09Umi82YZpIQpwln5+rM12OiJGVyoUEzoVXigVaJhzOVJs3kwpSxZpdWW1f2bEHrF19zWQCtocDvH5tmz0iZwVL+VXW17rvtBv4J8CffO87zZ+YJA5+hHp9iLiBKUm7s78Fc60qrpWcNLYqIiHTGK8OFCT88UaFSi8gFHtt78zTyjoFiyN6hHo6eXyBxjnLgUwo9ZmoxgedoxgkpkLbLlK5ZrbStlTeus2FEuHTQWurOm61FHJmocrvvv6qu1lJo+tMP/zzHJxb4/f/yDC+MLwCOXYNFAt9jZrHJ60bLPHz4NAPFkLH+Ai9NVvlfnj7HnqEebtvep+AlIiKyxcbn6oQ+PHFyjp6cT38xpN5M+OmZCr/7S7dwYqbGs2cr9BZDUpcwVY2Zq8V4BokHUXrpHHGlunIYES7dpZcCzQTqUcrPzlWYrUXs2VZctcvx+MQCe0fKfOJ9d/Bf3TrCjsEe4rRVWWOknOPk9CI/Pj3Ds+cqHD1f4fnxeQLPmKtHF5xDREREtsZof4Fnzs7Tk/MphD6GgRnD5TwnZmrcu3+IJIVizqhFjrxvpK4Vsmrxawxalwge12XP1spjXji3QBSnPH9ujtHePGfmau2Nq0PGevPL1eH3jpQ5eN++5Z6vwINzlQbT1YiR3jzNKOXRn02zvT9PfynHXO2VvRhVYV5ERGTr3Lt/iM//4CQOiJOU1DnixHFDb55vPTvO5Hyd4b48Tx6foREnxJvRk3WJE3Rt2Fovz+DMXJ2XZmrgHIPlPDf0FkiSlJmFJrXmK9X3V+6n9LnHTlIIAxabCY04pRj6pM4xW4so5cPloFUurL7PooiIiGSnEPrLE93na61t/PJB6+vRFyfZPVSg0ogAd9kJ7uvyWnu2zKyyjlOfdc7deuWtujrEDlyckrYT6fRik3qUEPoeY315ZnpaoenifRK/d2SKUj7AM5ivRwyWC5TyPnO1Jtt68tx+Yy/A8j6LIiIisjUOHZ3iwK5tPD8+z8RCg2IYkKSOc/MN3r93iKfPVnjuXJW+YsBcLSZNW+Wh0s2cqLXC5Xq2jjrn3nKpA8zsR5vYno4we6XWVhQ7qmmMb8Z0tclL04t8+D/8gGIY8KYdAwQePHZsmomFBvnAKBVyOBcRJUm7+KnP68ZaKxwrtYjZWsT9d4x29hcUERG5jozP1dk1XKJcCPjKT84SxQnNxJG6lCOTVW7sy/OTU7Ps2lbkufoCTXP4BumaW0lvzOXC1n+9jnOs55irWrxiTaejNUEuag+++l6rtMN0LebUdBU8Y6AQsnNbkbNzDW4u5BjszYODu/cNc8++QU7M1DjXXtF4/x2jmq8lIiKyBZZGoX768hwvjs9zU38eM5hZjMiHPtt6cjSjlOfmFtg92EN/KceOesLJ6UUADLepqxCXXDJsOeeOXe4E6znmWrRUp6tSj6k2W71WJyarpBhxmtJfDBkohuRCj7nFJmB84MAO9o6Uua/DbRcREek2l6tpuVTIdKAY8tZd/Xz1J+d47NgUrj05Pk0dhcCo1iMMuHWsl558yFt3buPc3CLfeHaciYUGUQa9W695gryZPeScO7iZjbmaOFqBK0nBN0iSlFoChdAIPKNSj5ivxxQCj7GBHl6/vVc9WCIiIhlYGaTG+gsstEsrLXVywCuFTPuKIZPzCY24tfPLYjOhJx/icKQOTs3WCDzjxFSVXYNFas2Es3MNKvWIpEPDiJfy6U1rxVVqqSsxStxy9dhG1NoMyAz68h7nFxoU8yHn5xvLNblERERk86wMUsCqpZWW9j0EODJRxfeMm7b1cG6uTl8hYHYxYnKhSZSk9OZ8+kshL55boFKPV93mbzO95qKmzrknNrMhV7P0osdLfxTzPBqx456bh9g9+Or9FUVERGTjxufqlAsX9g+VCwHjK0orjbZ7vAAq9YhSPqDeTCjlfc5XGlSbMc04xTOoJ46JSoNKPV7+XM/S5Uo/fOVSbXDOvX/TW3SV89srF3tCn8A3RvsLDJcLpM6pnpaIiEgGloLUUo8WvLq00r37h3j48GkAevMBi424NTE+MPKhR9JMcbTqbxVCn6mF5mX3PNwslxtG/Nft778GjAH/T/v5h4ATGbXpqpa6VndgPU5oJin7hktMLtTJtfdXFBERkc21Z1uRTz96nDR1DJVzbO8r4PneBaWV9o6U+cCBHRw6OkVfIWS2FnHP/kGeODlLM07oyQfkA0fge4SebVnQgsuvRvxbADP75865X1zx1lfM7DuZtuwqY0DOh9QZceqwFAo5Y64W8dfPnef2G/s4eN++TjdTRESkqxyfWOC7x6a5fXsvZyt1JhcazNYiDt6791WrEQ8dneL5MxXm6hGDxZBqM2VbT47A9xgu56jWY54/P08Sb2XUWv8E+REz27dU5sHM9gIj2TXr6uOAxLVqcPgehB7kQ5+zc3U8M/YOFTU5XkREZJOtnBy/Z7j1OVupRZyYqS2XWlparZgmKS/NLOKbMV+PKIYe1WaMb47zlQbNOKEvH3C+2dzS32G9E+R/F/i2mX3bzL4N/A3wDzd6cTN7r5m9YGZHzOxjq7xvZvYH7fd/YmZv3eg1NyJOW8OIeR8wj2ozIe979BV8nnhpjuMTC51snoiISNdZz+T4pUB2br5BTy4g8IzxSp0jE1VqzYRalDJfj6g3E6rNlMCD3GteInjl1tWz5Zz7mpndAryu/dLzzrnGRi5sZj7wx8C7gdPA42b2ZefcsysO+2XglvbXzwP/rv29I3Jea2uf1HkEvmFm1OKUGAgN/s03X2TPUGnVYmsiIiJy5dYzOX6p7EOlHuEDPz1TYb4WkZCS83yiNCVJHbnAI+d7ND2jGWe9BvEVV5LrbgFuA94MfMDMfmuD174LOOKcO+acawJ/Djx40TEPAp91LY8BA2a2fYPXfc2aaWslYpSmmNHetNIRxyn1OGW62mSsv0C1XWxNPV0iIiIbc+/+IWZrEZVaROrc8r7D9+4fWj5mKZB5Bscmq1QbMalLiZPWgjbPHEkKi81W6Yd4C4MWrDNsmdkngT9sf70T+D+AjZZ9uAk4teL56fZrV3rMlmq2C3K0QpfDM8PzoCcXMFzO45nR197KR3W3RERENmZplWGpEHBurk6pEFxQOR5eCWSNKMEZOOeIkla5Js+MOLXlY+cbyZauRIT1T5D/dVo9Wj9yzn3YzEaBf7/Ba9sqr10cNddzTOtAs4PAQQC/L7u5+0ZrorzvGb6DG/rzzNcSijmfm0dKy8eV2/8oREREZGP2jpQvOTVnKZA9eXKKNElJkrQVqFIAd0G4ipKt7dWC9YetmnMuNbPYzPqA88BG6xycBnaueL4DOPMajgHAOfcQ8BBAfvstmd3JpT0TPc+49YYynueRpHXetnOA4d5Xxo8vHk8WERGRbOWCgLG+Iphxerq25T1Ya1nvnK3DZjYA/N/AE8CTwA83eO3HgVvMbK+Z5YAPAl++6JgvA7/VXpV4NzDnnDu7wetuiNEq+5AmcGJqkZxvfOS+vQz2FqjUIibma3z7xfN845lzTM7XNW9LRERkCxw6OsXtY71MLzbJBz69xWDV4bFOuGzPlpkZ8C+dc7PAp8zsa0Cfc+4nG7mwcy42s48CXwd84DPOuWfM7CPt9z8FPAI8ABwBFoEPb+Sam8EBge9RyHm8aecAODg6VeOefYM88dIM33h2nNRBXz7gubPzjFeO8d//4j6tTBQREcnQ+FydXcMlhk7nWIwSfDO29QTM1mKcy37/w0u5bNhyzjkz+xLwtvbzE5t1cefcI7QC1crXPrXy2sDf26zrbRbnHMXQpxj4VOoRA8WQEzM1HNBXCOnvyVEIPOpxyonJKl966mV+9923dbrZIiIiXWtpReJof5FGnGJmVGoR9Til1uzsgOJ6hxEfM7Ofy7Ql15DEQW8xpB6n9BbC5eJqT52apa8npBj6mBnF0KevJ+SpU7OdbrKIiEhX27OtyHePTvLS9AJHzs9TWWzw8swijS3emmc16w1b7wS+b2ZH25Xcf2pmGxpGvBYtjf02E8dMtcGzZ+dYaES8NFVtT4Y37KJ+ytbzq2XUWEREpPss7Z/4+tFedm4rEfrGeKVJ4HuUcgGlnEfYwY/i9a5G/OVMW3ENmphvsGuwh4VazA9PzHDfzcPsGSzytz+bxDOjXPDpK4REiePt+wY73VwREZGudcH+iSNlzIzFZkKtmXDjQBHDaCYpZ2cXqUcpUZIS+B7NON2SFYvr3a7nZNYNudoZreJoiYOc33qhFiWkrslbdw5w+KUZmqmjvxjSjBPmazGVeswd2/t48M6O1mEVERHpakvb9Syp1CP6iyHVRkKcOELfCH3DOdgzXKKvGHBsospsLSLdgrpblxxGNLMnL3eC9Rxzrcr5tnyDWlvztEJXTy5goCfHtlKe/Tf0UotTnjo1S28uYKicI05bKxZvGihwyw2XLsQmIiIiG7M0OX5JXyGkGPiU8z6LzZhmnLDYjElxTC00eGlqEd8zhnrCS5x181yuZ+v1l5mbZUD/Jran4zwgF7SDlbHcvegZ+EubT0cpAz05GlFCIfCYWGgQxSnPjc/Tkwu4fXsf9ThlsRkzMd/s5K8jIiLS9e7dP8TDh08DrR1cCr5xZq5Ob8EncY5KIyKOUgZ6QiqLrdCVpI5GlGJkXxbicmHrdes4R7IZDbkaeAa9hYBm7PBJSBzkA1veKbyZpLRGED3majE7thWZa+/FNF+PGJ9vhbCRcp7eQkgjSpirR53+tURERLra0nY9h45O8dyZCmfnG9xz8yCLUcrkQoNRz2O4FHJmrkEzrrJQj4lTx1bt3HPJsHU9zdUyWr1ZcZLiGzRTw7PWcKBLExpRSiH0wIzAjIVGRJTkWGhEjPUXSB0kaUS9mfDSdJXRvgKGsa24NV2UIiIi17Ol/RM/99hJdg720Lfi87dSi/jPT5xm/0iJ056RDzyiZuvz3kH7cz+7tq13NWLXW5qb5RxEzhEnDt+HwLX2QfQB3/dwKRRyPr2FkP039LJ/pEQxDPD9CjMLDebqMQv1iLl6xNv3DrFzqHSpy4qIiMgmuniyPLSGFptJq5drptqkFqUXDB26jHu4FLbaEto3w1hemeBSiCwBZ4ShR+gZIwMFRsp53rprgDiFJG39EW8eKXG4GnFjf5HcUA+TC008z+Pe/UOd/LVERESuK0uT5Vf2bC3UY/oLHs+fmyfl1XO0si7/sN6iptcF86Dge8t/BEcr7fbkfHzPSFLY3l/gwO5tFMKA0f7C8h91uFzgwJ4BcmFrwvxgKccHDuzQSkQREZEtdO/+IWZrEZVaROoclVrEbC2imAsphh5p1t1Yq7hkz5aZzbP6JH2jtXVhXyat6pDUwVw9hvZKxKXXBntCGgn8wv5B3rRzGwv1mNlaxP13jAIsr4AYLOW53feZ7S8qaImIiHTAysnyL5ytMFOL6C+EnJpebNXb8jzMT2msMjs+q5WJl5sg35vBNa9aS/c957dud5S0erqK+ZCP3L2Legrn5uqM9he4/47R5TC19Edd7T0RERHZWkufwS/P1LhpWw/lQsC3nh+n3kwohD51g5SUJHHLQ4iB11oY18xgiaLmbK0iShyBb+T9VimIvmLI0anamr1VSysgRERE5OqwcgsfgNFyjhNTMbWoPYIFeB547YVwoe+RpI7AHA42tSzEdTVny1Z8XYoD0tQRBh7z9Sbz9YifnJrl33zzRY5PLGTfUBEREdmQ8bk65cIrfUp7Rnq5+YYSOAgDj5zvkfONXOhTDFtBy+HIhR6F0K48ILm1J4N1ZdgyILALn/sGod/qIgQIL/ObOweLzYRmAjnPGO7NMVVt8PDh0wpcIiIiV7mLt/C5eaRELgjYNVTi7r1DvGXnAH3FECOlHqeU8j7DpTyhZ4S+T+BtXkjqymFEMyMXeATtbsJyIcAM+ntCzs3WieIE328VNFv159v/44Bc4FGNEpqxY6S3wEAx5NDRKQ0bioiIXMVWbuFTj2KePTtPpdakvxByZGKBJHUUw4DeQshiI8b3fTxzeFbGq18AABJ+SURBVJHRXwgxYLH9+b++EUVbc+CsK8OWZ5A4B6mjnkJjoUkp7/O60TLztYhK4gAj9CBq562cb6TOkaSt21UIPQLzKOY8Jhaa3NCbcPuNvZQLAefm6h39/URERK5XxycWOHR0ivH2orR79w+tOZ/6Awd28KWnXubRFycZ7stz/xvGyAcB33ruHPnA49RMjXIYsHuoROB5NJOEk1OLTFeblPM+tShZsySD7xn5wHAOmklKqzrn6roybBVzPiPlHOfmG8tdgJ7BT16u8MYb+6jFKRPzdc7MJvQVPMDDDKqNmP6Cz7ZSnv0jZZ47V2GxkRD4xoE9AwyXC1RqEaMXVaYVERGR7B2fWODhw6cZKIaMtYcJHz58+pIL2IZ7C7znjrELipzmfJ/Bco6hcoFGnFIMfRyOepzwtl0DfO2n56g4RykfYCQ04hTPg0LgEaWtXWYCH3Zs6yHwPKpRzMm4UV2r3V05Z8s3Y7BcYHtvgULOIx94lPIh5ZzP2fkGb7yxj76eHOVCQE8ux8/tGeSdt42ye7AHzCjkfEr5gH1DJYqhx61jZQZL+eXCaKoKLyIisvVWrjD0zOgrhsvTey52fGKBzz12ki/96DTPnJljcv6VUamhco7JhQY3j5SoNWNqUUK9mZALPAbLBfYM95DzPZxrbdU3XA7J+R7NxOF7Rs6HQujTHkSjPx/gmvX5tdrdlWErdTC32GR6MSJNobcY0IhTGqmj3oyZWowZKuV5y84Bcr61Npt2ju39BTwzhoo55moRxXzAnbu28XO7Bzk3V6dUCFSsVEREpEMuXmEIrXnZ4xdN71nqAavWY24cKDJfjzl8cmY5cG3vK9CMUp4+U6EexZycqnJqusrtY73cs2+Q2EEpFzDWV6AQeMzVYkLPuHGgyNt2baMQhgRm1OMEh2PHYA9prTK9Vru7chgRoNpMSFxKT85nsJRnvh5TjxJSz5ivR/QVQyqLETePlskHHpV6hANuuaGMH3g4HK/f3suv3HmTwpWIiMhVYK19Dy+e3rOyB+yWG8ocPjGLZ8bPzi+QC3wqjZjt/QWaqaOYCxgo5Rkq5Xjwzps4dHSKHYM91Boxc/WY1Dk8z/B9o7cd9ALfKOdbc71ygUdvIeRSujJslQutNHp2ztFbCMkHHnHosdCIyNG6QWO9eV6aWuSu3dvYNVzipckqPzw5w117trFrqLS8Jc+S9U7IExERkWysXGFYLgSv2j5vyfhcnbF2AFvau/hn5xc4O1vnrn1DjPXlKQyVLghtlVq0/Dl/x/Zenjg5x439RWrNhCRNmKlG9IQ+s/WI20bLhIHPu28fW/5Zy5fW3MKwK4cR4yRlz3CJn987SOB7zDcievIBuwd76MkF9BVCdg6V+EfvuoWdwyXOzdU5U6lz155t7Bkuv2oceGV35Fh/gWp7Qp7qbYmIiGydpRWGpXZlgLWm91xcY2u4XOCO7f08+Jab+M27dxOnrDkcOdpfIB8EHNgzQK5dlNP3PN6yexsPvmUHQ6U8PYXwgt6sciHAPH/N7q2u7NmqNRNmqg2mqxFv29XPYpQyudAg8Dz+p/fs5b7bblg+9r7293/99ReWU/CSpTIPF5f8X/quelsiIiJbaz1b5F2uB+xSw5FLPztQDLlrzyBjvXl+eGKG4Z4c3z86ycmpKs7BO28dvuBnXZpErKEre7Z8zygXcuR8Y3y+we7BEg+88UZ+//13XBC0Vro4BcMrN369E/JERESk8y7XA3bv/iFmaxGVWkTqHJVaxMnpKpPzdb7w5MvkA6MexZybq7NzqMTf/bkdnJ1vMFltsHeoh75CyHPnFjhfqS1XKnCNamWt9nRlzxZAMfQZHSjSiBJ+7/7bgFeWga427+pSKfjQ0al1TcgTERGRq8OlesCWwtiho1Ocm6sTeK1t+ophcEEGWApon3vsJPfsH17OAZMLdZ4+U+FHp2Z51+1j3H/HKP9z3Gys1ZauDVsA5mBp2+nLFUK7+MaP9he4/47R5T/UeibkiYiIyLVhZRj73GMnKYTBmtOFxufqhD48drxCpdaqaHDH9l6iBH7z7t2XvVbXhq16lDBbi3j7vkGAdc27WisFXy6IiYiIyLVr5erFJSu35ws8+N7Rabb15OgvhjSilO8dnV7OGJfTlWErSR0pjn0jJR688ybg8jfyctYzIU9ERESuPZer39XaUbn1HXfR83XoygnyxZzPfbeMcPC+fcsB6eIJ8JMLdb7zswl++vIsn3vspMo4iIiIXKdWmzC/cnu+JIW79w0uF0HPBx537xskWXPr6Qt1Zc/WaF/hVWOoKyfAN+KY7x2dxmjdvOplNrIUERGR7nW56UKj7Rqbd+97ZW/kSi2iVFhfjOrKsLWalTfy+89O0V8MecONfQz3vjK0qLpZIiIi16dLTRdab+X6tVw3YQteuZFL87c8s+X3rmT+loiIiFw/NrpQriNhy8wGgYeBPcAJ4L9xzs2sctwJYB5IgNg5d2Azrr/ejSxFREREYGML5To1Qf5jwF87524B/rr9fC3vdM7duVlBCy4/EU5ERERks3QqbD0I/Fn78Z8Bv7KVF1/vRpYiIiIiG9WpOVujzrmzAM65s2a2+oaFrRIW3zAzB3zaOffQZjVAdbNERERkK2QWtszsW8DYKm99/ApOc49z7kw7jH3TzJ53zn1njesdBA4C7Nq164rbKyIiIpKFzMKWc+5da71nZuNmtr3dq7UdOL/GOc60v583sy8CdwGrhq12r9dDAAcOHFhvUVcRERGRTHVqztaXgd9uP/5t4C8vPsDMSmbWu/QYeA/w9Ja1UERERGQTdCps/W/Au83sZ8C7288xsxvN7JH2MaPAITP7MfBD4KvOua91pLUiIiIir1FHJsg756aAX1rl9TPAA+3Hx4A3b3HTRERERDZVV25ELSIiInK1UNgSERERyZDCloiIiEiGFLZEREREMqSwJSIiIpIhhS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREcmQwpaIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhkSGFLREREJEMKWyIiIiIZUtgSERERyZDCloiIiEiGFLZEREREMqSwJSIiIpIhhS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIY6ErbM7DfM7BkzS83swCWOe6+ZvWBmR8zsY1vZRhEREZHN0KmeraeBXwO+s9YBZuYDfwz8MnA78CEzu31rmiciIiKyOYJOXNQ59xyAmV3qsLuAI865Y+1j/xx4EHg28waKiIiIbJKrec7WTcCpFc9Pt19blZkdNLPDZnZ4YmIi88aJiIiIrEdmPVtm9i1gbJW3Pu6c+8v1nGKV19xaBzvnHgIeAjhw4MCax4mIiIhspczClnPuXRs8xWlg54rnO4AzGzyniIiIyJa6mocRHwduMbO9ZpYDPgh8ucNtEhEREbkinSr98Ktmdhp4O/BVM/t6+/UbzewRAOdcDHwU+DrwHPAXzrlnOtFeERERkdeqU6sRvwh8cZXXzwAPrHj+CPDIFjZNREREZFNdzcOIIiIiItc8hS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREcmQwpaIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhkSGFLREREJEMKWyIiIiIZUtgSERERyZDCloiIiEiGFLZEREREMqSwJSIiIpIhhS0RERGRDClsiYiIiGRIYUtEREQkQwpbIiIiIhlS2BIRERHJkMKWiIiISIYUtkREREQypLAlIiIikiGFLREREZEMdSRsmdlvmNkzZpaa2YFLHHfCzH5qZk+Z2eGtbKOIiIjIZgg6dN2ngV8DPr2OY9/pnJvMuD0iIiIimehI2HLOPQdgZp24vIiIiMiWudrnbDngG2b2hJkd7HRjRERERK5UZj1bZvYtYGyVtz7unPvLdZ7mHufcGTO7AfimmT3vnPvOGtc7CBwE2LVr12tqs4iIiMhmyyxsOefetQnnONP+ft7MvgjcBawatpxzDwEPARw4cMBt9NoiIiIim+GqHUY0s5KZ9S49Bt5Da2K9iIiIyDWjU6UfftXMTgNvB75qZl9vv36jmT3SPmwUOGRmPwZ+CHzVOfe1TrRXRERE5LXq1GrELwJfXOX1M8AD7cfHgDdvcdNERERENtVVO4woIiIi0g3Mue6bS25mE8DJK/iRYUCFUzef7mt2dG+zofuaHd3bbOi+ZudK7+1u59zIam90Zdi6UmZ22Dm35rZB8trovmZH9zYbuq/Z0b3Nhu5rdjbz3moYUURERCRDClsiIiIiGVLYanmo0w3oUrqv2dG9zYbua3Z0b7Oh+5qdTbu3mrMlIiIikiH1bImIiIhkSGGrzcz+uZn9xMyeMrNvmNmNnW5TNzCzf2Vmz7fv7RfNbKDTbeoWZvYbZvaMmaVmptVIG2Rm7zWzF8zsiJl9rNPt6RZm9hkzO29m2m5tE5nZTjP7GzN7rv3fgX/Q6TZ1AzMrmNkPzezH7fv6zzblvBpGbDGzPudcpf347wO3O+c+0uFmXfPM7D3A/+eci83sfwdwzv3jDjerK5jZ64EU+DTwe865wx1u0jXLzHzgReDdwGngceBDzrlnO9qwLmBmvwgsAJ91zr2h0+3pFma2HdjunHuyvY/wE8Cv6N/sxpiZASXn3IKZhcAh4B845x7byHnVs9W2FLTaSoBS6CZwzn3DORe3nz4G7Ohke7qJc+4559wLnW5Hl7gLOOKcO+acawJ/DjzY4TZ1Befcd4DpTrej2zjnzjrnnmw/ngeeA27qbKuufa5lof00bH9tOA8obK1gZv/CzE4Bfxf4RKfb04V+B/irTjdCZBU3AadWPD+NPrjkGmFme4C3AD/obEu6g5n5ZvYUcB74pnNuw/f1ugpbZvYtM3t6la8HAZxzH3fO7QT+I/DRzrb22nG5+9o+5uNATOveyjqt597KprBVXlPvtlz1zKwMfAH4hxeN0Mhr5JxLnHN30hqJucvMNjz8HWy8WdcO59y71nnofwK+Cnwyw+Z0jcvdVzP7beB9wC85TRK8Ilfwb1Y25jSwc8XzHcCZDrVFZF3ac4q+APxH59z/2+n2dBvn3KyZfRt4L7ChBR7XVc/WpZjZLSuevh94vlNt6SZm9l7gHwPvd84tdro9Imt4HLjFzPaaWQ74IPDlDrdJZE3tidx/AjznnPu3nW5PtzCzkaVV82ZWBN7FJuQBrUZsM7MvALfRWt11EviIc+7lzrbq2mdmR4A8MNV+6TGt8twcZvarwB8CI8As8JRz7v7OturaZWYPAP8X4AOfcc79iw43qSuY2eeBdwDDwDjwSefcn3S0UV3AzO4FHgV+SutzC+CfOOce6Vyrrn1m9ibgz2j9d8AD/sI59/sbPq/CloiIiEh2NIwoIiIikiGFLREREZEMKWyJiIiIZEhhS0RERCRDClsiIiIiGVLYEhEREcmQwpaIXPPMLDGzp8zsGTP7sZn9j2b2qv++mdk7zGzOzDZUi8jM/jsz+6NVXr/PzJ41sw1VmxaR7qKwJSLdoOacu9M5dwfwbuAB1t5u61Hn3AMXv2hm/kYb4Zx7tH1tEZFlClsi0lWcc+eBg8BH21uarKnd0/U3ZvafaFXixsy+ZGZPtHvJDq449sNm9qKZ/S1wT5a/g4h0l+tqI2oRuT445461hxFvoLVFzKXcBbzBOXe8/fx3nHPT7X3RHm9v5ZUD/hnwNmAO+BvgR9m0XkS6jcKWiHSrS/ZqrfDDFUEL4O+3950E2AncAowB33bOTQCY2cPArZvWUhHpagpbItJ1zGwfkADn13F4dcXPvQN4F/B259yimX0bKLTf1kayIvKaaM6WiHQVMxsBPgX8kXPuSgNSPzDTDlqvA+5uv/4D4B1mNmRmIfAbm9diEel26tkSkW5QNLOngBCIgc8B//Y1nOdrwEfM7CfAC8BjAM65s2b2T4HvA2eBJwEfwMzeDxxwzn1io7+EiHQnu/L/4ycicm1qDxP+nnPufRleYw/wX5xzb8jqGiJybdEwoohcT5rAGzZa1HQtZnYf8BVgMovzi8i1ST1bIiIiIhlSz5aIiIhIhhS2RERERDKksCUiIiKSIYUtERERkQwpbImIiIhk6P8HQHkS8F5HMBIAAAAASUVORK5CYII=\n", `````` Maximilian Schanner committed May 06, 2020 180 181 182 183 184 185 186 187 188 189 190 `````` "text/plain": [ "