{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Generating a number range with normal spacing density\n", "\n", "We want to generate N numbers between x_min and x_max with both limits included. This is trivial to do with a uniform distribution, for example with `np.linspace(x_min, x_max, N)` in NumPy.\n", "\n", "But sometimes it can be helpful to have a range of numbers that are more concentrated around a specific area. For example if they are used as input for monte-carlo simulations. Ideally their density would be proportional to a [normal distribution](http://en.wikipedia.org/wiki/Normal_distribution). That shall be the goal of this notebook! Source is [on github](https://github.com/s9w/articles).\n", "\n", "We'll assume a normal definition with a standard deviation of $\\sigma$ and a center of $\\mu$:\n", "\n", "$$f(x) = \\frac{1}{\\sigma \\sqrt{2\\pi}}\\exp\\left(-\\frac{(x-\\mu)^2}{2\\sigma^2}\\right)$$\n", "\n", "Let's have a look:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib as mpl\n", "import matplotlib.pyplot as plt\n", "from scipy.stats import norm\n", "%matplotlib inline\n", "mpl.rcParams['lines.linewidth'] = 2\n", "mpl.rcParams['axes.grid'] = True\n", "mpl.rcParams['legend.labelspacing'] = 0.0" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAEACAYAAACuzv3DAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXd8VMX6h59JgSQQCF2lhV4VbKhYiGKnWa6KHStgb/dn\nuV6JBeu9igqIHQFFroiCVEUJIlKkBJDeIWggIZQQ0vf9/THZlCXJ7mbP2ZLMw+d8smfPnJnvmT28\nO/ued95RIoLBYDAYQpewQAswGAwGg28YQ24wGAwhjjHkBoPBEOIYQ24wGAwhjjHkBoPBEOIYQ24w\nGAwhjltDrpS6Uim1SSm1VSn1dDnHGyul5iqlkpVSfyqlhtii1GAwGAzloiqLI1dKhQObgUuBfcAf\nwM0isrFUmUSgtog8q5RqXFS+mYgU2CncYDAYDBp3I/JewDYR2SUi+cDXwCCXMn8D9Ype1wMOGiNu\nMBgM/iPCzfHmwN5S+ynAOS5lPgZ+UUr9BcQCN1onz2AwGAzucDci92T+/nNAsoicAvQExiilYn1W\nZjAYDAaPcDci3we0LLXfEj0qL01vYCSAiGxXSu0EOgErShdSSpmkLgaDwVAFRERVdtzdiHwF0EEp\nFa+UqgXcBMxwKbMJ/TAUpVQztBHfUYGYoN9GjBgRcA1Gp71bXkEeXUZ3gUSIez2OmEtjiH4lGhLh\nyXlPBlxfKPWl0Wn/5gmVGnLRDy0fAuYBG4ApIrJRKTVUKTW0qNirwFlKqTXAfOD/RCTDo9aDkF27\ndgVagkcYnVVn9PLRbEzfSPuG7Ul9MpUbmt/Ab3f/hkLx7rJ32Zi20X0lASAY+7I8jE7/4zaOXETm\niEgnEWkvIq8VvfehiHxY9DpdRAaISA8ROVVEvrJbtMFQVVKPpTIiaQQAo64YRe2I2gCccfIZ3HfG\nfRQ4Cnhk7iMej4QMhmDAzOx0YciQIYGW4BFGZ9V4Zv4zZOZl0q9DP/p17AeUaBzZdyQNohowf8d8\nvtv0XQBVlk+w9WVFGJ3+p9IJQZY2pJSYUY4hkCzZu4Ten/WmVngt1j+wnvYN259QZuwfY3lw9oO0\nqt+KjQ9uJCYyJgBKDYYSlFKIjw87axxJSUmBluARRqf3vLjwRQCeOu+pMka8tMahZw6lR7Me7Dmy\nh89Xf+5viRWilDJbDdiqijHkhhrB4ZzD/LzzZ8JUGE+c90SF5cLDwnnyvCcBgs69EujoCbMFPjql\nIoxrxVAjmLxuMrdMu4U+rfuQNCSp0rIZ2Rk0faspSinS/plGXFScf0RWQtHP60DLMNhIRZ+xca0Y\nDEXM2KKnPwzq5Joq6EQaRjfkwtYXUuAoYM7WOXZLMxh8xhhyF4LJp1sZRqfn5BXmMXvrbAAGdhp4\nwvHyNA7sqMtN3zzdVm0GgxUYQ26o9izctZCjuUfp1qQb7Rq28+icQZ31yH3OtjnkFebZKc9g8Blj\nyF1ISEgItASPMDo9Z8bmyt0q5Wls26At3Zt252juURbuWminPIPBZ4whN1RrRKTYPVKeW6UyjHvF\n/4wePZqzzjqLqKgo7rrrrkrLZmRkcO2111K3bl3i4+OZPHlyjdVnDLkLweDT9QSj0zOSU5PZe3Qv\nJ9c9mbObn11umYo0Ot0rMzbPMBEjfqJ58+b8+9//5u6773Zb9sEHHyQqKooDBw7w5ZdfMnz4cDZs\n2FAj9RlDbqjWOEfTAzoOIEx5d7ufdcpZnFz3ZPYe3cvq1NV2yDO4cO211zJo0CAaNWpUabmsrCym\nTZvGyy+/TExMDOeffz6DBg1i4sSJNVKfMeQuBINP1xOMTs9w+scrc6tUpDFMhTGg44Ay9QQrSlmz\n2UX//v1p0KBBudvAgSd+Nu5+AW3ZsoWIiAjaty+ZodujRw/Wr19fLfW5wxhyQ7Xl78y/WZ26mpjI\nGPq27VulOpzulVlbZ1kprVpSWFjIBRdcULx/zz33sG3bNgBmzpzJoUOHyt1mzDjxS9LddPVjx45R\nr169Mu/FxsaSmZlZ4Tlbtmzhxhtv5OKLLyY2NpYBAwYwbty4oNHnC8aQuxBon66nGJ3uWZKyBIDe\nLXsTFRFVYbnKNF7Y6kLCVBjJqckczz9utUTLELFm84UlS5bQunXrIj3CkiVLyoxIvbueysXUrVuX\no0ePlnnvyJEjxMaWv8pkRkYGw4YNY8KECSxYsIC+ffsyadIkhg0bFhT6fMUYckO1ZWnKUgDObX5u\nleuIrR1L96bdKXAUsOrvVVZJq5bMnTuXK664AoDVq1dz6qmnFh+76qqriI2NLXfr16/fCXW5G/F2\n7NiRgoKC4hE/wJo1a+jevXu55ceMGVP88BEgNzeXmJiSzJaB1ucr7tbsrHEE2qfrKUane4oNeYvK\nDbk7jec2P5e1+9eyNGUpF7S6oNKyNZl58+YxePBgAGbNmkXfvn2ZMWMGAwcOZM4cz1IdFBYWkp+f\nT0FBAYWFheTm5hIREUF4eHiZcnXq1OG6667jhRde4JNPPmHVqlX88MMPLFmypNx6MzMz6dq1KwDr\n16+nW7duREZGFh8PtD5fMSNyQ7UkvzCfFX/p9b/PaXGOT3Wd1/I8oOSLwXAiaWlp7NmzhxkzZjB7\n9mxiYmJIS0srM+r1BGeUxxtvvMGkSZOIjo5m5MiRAFx99dW8/vrrxWXHjh1LdnY2TZs25bbbbmPc\nuHF06dKl3HqHDx/Ojz/+yLfffsv8+fPL1BMM+nzFZD90ISkpKSRGu0Zn5az6exVnfnQm7Ru2Z+vD\nWyst607jpvRNdBnTheaxzUl5IsVipZ4R7NkPJ02axMaNG4uNmsF7bM1+qJS6Uim1SSm1VSn1dDnH\nn1JKrS7a1imlCpRSgc/7aajReOpW8YSOjToSFxXHvsx9pBwNjCEPdpYvX851110XaBk1lkpH5Eqp\ncGAzcCmwD/gDuFlEyl1mXCnVH3hMRC4t51hIjMgN1YM7vruDiWsnMvqq0TzY60Gf67ty0pXM2z6P\nb274hn90/YcFCr0j2EfkBt+xc0TeC9gmIrtEJB/4GqgsofMtgP0JDwwGN1g5Ii9dj/GTG4IRd4a8\nObC31H5K0XsnoJSKAa4AvrVGWmAw8dnWEgidB48fZGvGVqIiojit2Wluy3ui0RhyQzDjzpB781tu\nAPCbiBz2QY/B4DPL9i0DdK6UyPBIN6U9o1fzXgCs/HulyU9uCDrcGfJ9QMtS+y3Ro/LyGIwbt8qQ\nIUNITEwkMTGRUaNGlRkJJSUlBcW+M3ohWPRUtO98L1j0BFN/Lk1ZCjvhlPRTPCqfkJDgtv61y9bS\nMqMlOQU5rN2/NiCfd00gWNPEVqXNhIQEoqOjiycWeRp6mJSUxJAhQ4rtpUe4WdU5AtgOxAO1gGSg\nSznl6gMHgehK6hKDwR9cNuEyIRGZun6qpfXe+d2dQiLy/rL3La3XE2rK/59p06bJ999/L8OHD5ch\nQ4ZUWnbw4MEyePBgycrKkt9++03q168v69evt1WfN20mJCTIp59+6nHdFX3GRe9XaqsrHZGLSAHw\nEDAP2ABMEZGNSqmhSqmhpYpeA8wTkWzPvj6Cl1AZBRmd5eMQR7FrxdMHnZ5qNH5y+wnWNLFVbVP8\nFGnkdoq+iMwB5ri896HL/hfAF9ZKMxi8Z1P6Jo7mHqVFvRY0r1fuc/kqE8yGXL1oTQ5aGWGP4enf\nvz+LFy8u99iFF154QoZBdwawojSxVR04eKKvKm0+++yzPPPMM3Tq1ImRI0fSp0+fKulzh8m14kIo\nzJYEo7MiqhJ26KnG7k27ExMZw/ZD20nLSqNJnSZVkVhtKSwspE+fPvz222+ATmP77LPP0r59e2bO\nnOlVXXalsX3++edJS0tjxYoVJCQk0K9fP4YNG+aRPm/bfOONN+jWrRu1atVi8uTJDBgwgOTkZNq2\nbeu2LW8xhtxQrXBmKDz7lPKXdfOFiLAIzjj5DH7b8xurU1dzebvLLW+jqtg1kvaGUEhjO3v2bKKi\norjmmmv44osvqF+/vseavG2zV69exa/vuOMOJk+ezOzZs3nooYc8btNTTNIsF4zv2Vr8rXPt/rUA\n9GjWw+NzvNHorNfZjqGE6p7G1t+pab3BjMgN1QYRKTawnkwEqgrOeo0hP5HqnsbWmzaPHDnC0qVL\n6dOnDxEREUyZMoVFixbx/vvve9QPXuMurMWqjRoSPmUIHLsP7xYSkcZvNhaHw2FLG0v2LhESkR4f\n9LCl/ooI9v8/Bw4ckKZNm8rIkSNl1qxZ8p///EdeeeUV+emnn7yqZ8SIEaKUKrO9+OKLIiJy1VVX\nyWuvvVZcNiMjQ6655hqpU6eOtG7dWiZPnlxhvTt27JBRo0bJ1KlTZdSoUZKfn1+l66yszdL60tLS\n5Oyzz5bY2FiJi4uT8847T+bPn19p3RV9xngQfmjS2BqqDTO3zGTA5AFc0uYSfr7jZ1vaOJZ3jNjX\nYokMiyTruSzLZo66I9iTZpk0tr5jaxrbmobxPVuLP3UWu1WaeudW8UZj3Vp1adegHfmOfDYf3OxV\nO9UZk8Y2sBhDbqg22O0fd2L85Cfy3nvvceaZZwZaRo3FGHIXTHy2tfhTZ1UNubcajSE3BBvGkBuq\nBTkFOWw+uJkwFUbXJl1tbcsYckOwYQy5C8b3bC3+0rkhbQMOcdCxUUeiI6O9OtdbjcaQG4INY8gN\n1QJ/+ccB2jZoS0xkDPsy93Hw+EHb26tJeJMmNpjT2NatW7fMpKKIiAgeeeQR23QZQ+6C8T1bi790\nVjViBbzXGKbCOLWpnrW47sA6r9szVIxz9uWBAwf48ssvGT58OBs2bPC5rL/1HTt2jMzMTDIzM0lN\nTSU6Opobb7zRNl3GkBuqBf4ckZdux7hXrMObNLGhksYWYOrUqTRr1owLLrjANm3GkLtgfM/W4g+d\nIsKa/WuAqhnyqmh0jsiDxpArZc1mE/3796dBgwblbgMHDgQqTk27fv36E+rzpmwg9JXmiy++4I47\n7qiSLk8xuVYMIc/+rP2kH0+nXu16tKrfyi9tmhH5ifiaxtabNLGhkMYWYPfu3fz66698/vnnbuv3\nBWPIXTC+Z2vxh87SbhV3WfPKoyoaT22mR+R/HviTQkch4WHhbs6wmSCYvu9rGltv0sSGQhpbgIkT\nJ3LhhRcW94tdGNeKIeTx5UFnVWkY3ZAW9VqQXZDN9kPb/dZuMONrGltv0sSGShrbCRMmcOedd1bc\naVbhLquWVRtBnr3NyYIFCwItwSOMzhJun3a7kIiM+2Nclc6vqsarv7xaSES+Wf9Nlc73hlD4/3PW\nWWfJunXrRETkpZdekg8//FCmT5/uVR2DBw+Wm2++WbKysmTRokVSv3592bBhg89l//nPfxYf+/PP\nP+XJJ5/0SldV2ly8eLHUqVNHjh075lHdFX3G+Lr4MoBS6kql1Cal1Fal1NMVlElQSq1WSv2plEqy\n8ovGYHCHvyNWnDh/AaxJXePXdoORtLQ09uzZw4wZM5g9ezYxMTGkpaWVGfV6wtixY8nOzqZp06bc\ndtttjBs3ji5dugBw9dVX8/rrr3tU1pXhw4fz448/8u233zJ//vwy9dilb8KECVx//fXUqVOnSm15\nQ6VpbJVS4cBm4FJgH/AHcLOIbCxVJg5YDFwhIilKqcYikl5OXVJZWwZDVSh0FFLn1TrkFuZy5Jkj\n1Ktdz/1JFvHl2i+57bvbuK7LdXx747e2tmXS2FZ/7Exj2wvYJiK7RCQf+BoY5FLmFuBbEUkBKM+I\nGwx2sfPwTnILc2lRr4VfjThQnNNlY9pGNyWrPyaNbWBxZ8ibA3tL7acUvVeaDkBDpdQCpdQKpdTt\nVgr0NyY+21rs1uk0ol0al/+T2hOqqrFT404oFFsztpJfmF/l9qsDJo1tYHEXfujJb7lI4AygLxAD\nLFFKLRWRra4FhwwZQnx8PABxcXH07NmzOPTL+Z8p0PtOgkVPRfvJyclBpSdQ/bkhQk+Pjv0rlqSk\nJL9fX+u41uw6vIuvfviK1nGtbWvPUHNISkpi/PjxAMX20h3ufOTnAokicmXR/rOAQ0TeKFXmaSBa\nRBKL9j8B5orIVJe6jI/cYDlDvh/CF2u+YFy/cQw9a6jf2+/3VT9mb53Ntzd+y3Vd7HMtBLuP3OA7\ndvrIVwAdlFLxSqlawE3ADJcy04ELlFLhSqkY4BzA3sw1BkMRG9L0rdalSdVdK77gdOk4dRgMgaBS\nQy4iBcBDwDy0cZ4iIhuVUkOVUkOLymwC5gJrgWXAxyISsne18T1bi506RYRN6ZuAwPjIS7e7Md08\n8DQEDrdT9EVkDjDH5b0PXfb/A/zHWmkGQ+Xsy9xHZl4mjWMa06ROk4BoMJErhmCgUh+5pQ0ZH7nB\nYn7c/iNXTLqCC1tdyK93/RoQDYdzDtPgjQZER0Rz7LljhCl7sl5UJYeMIfSoqo/cJM0yhCxWhB76\nSlxUHCfVPYnUY6nsPrybNg3a2NKOGQQZKsMkzXLB+J4t4sgR2LePpAULbGvC6Zf2dbFlX/uy2L1i\ns5886D/zIoxO/2MMucE6HA74+We4+WZo1gxatIAbb9T748ZBdralzQU6YsWJiVwxBBrjIzdYw7p1\ncP31sLVoHphSEBsLpfM3n302TJ8OJ59sSZNN3mpC+vF09j6+lxb1WlhSZ1UY+8dYHpz9IHf3vJtP\nB30aMB2G6okVceQGg3vWroVLLtFGvFUrSEyEnTvh0CH480/44ANo0wb++APOOUeX95G0rDTSj6cT\nWyuW5rGuWSP8S/GIPN2MyA2BwRhyF0LFbxY0Op1GPD0drr4aNm+GESOgdWsICyMpLQ2GDYOlS+G8\n82DvXjj/fJg716dmnf7ozo07+xzR4WtfOl07G9M22vpQMmg+czcYnf7HGHJD1Vm3Thvxgwe1Ef/2\nWyhageUEmjaFX37R/vJjx+CGG2DHjio37YxY8fVBpxU0q9OMBlENOJJ7hNRjqYGWY6iBGEPuQqgk\nKQq4zpwc+Mc/3BrxMjqjouDLL/V5x47BHXdAYWGVmi9+0GlB6KGvfamUKh6V2/nAM+CfuYcYnf7H\nGHJD1Rg5ErZsgc6dKx+Ju6KUjmA5+WRYvBjeeMP9OeXgdK0EOmLFiZmqbwgkxpC7ECp+s4Dq/PNP\ncC5p9fHHlRrxcnU2agRFaToZMQJWrvRaglUx5GBNX/pjqr65N60lVHR6gjHkBu8oLIT77oOCAv0Q\n84ILqlbP5ZfDww/rem67zasY86O5R0k5mkLt8Nq0ibNnJqW3mMgVQyAxceQG7xg9WhvgU06BDRug\nfv2q15WdDWeeCRs3wttvw+OPe3Ta8n3LOeeTczi16amsHe57KKMV7D68m/h342lWpxmpT5kHngbr\nMHHkBmtJT4fnntOvR4/2zYgDREfDm2/q16+/DllZHp1WnLo2SPzjAC3rtyQ6Ipr9Wfs5lH0o0HIM\nNQxjyF0IFb9ZQHSOGgWZmXDFFXDttR6d4lZnv37QqxccOABjxnhUp9OQd27U2aPy7rCiL8NUGJ0a\ndwJg88HNPtdXHubetJZQ0ekJxpAbPOPIET0KB3jhBevqVQpeekm/fvNN/UXhhmJD3tgaQ24VTj1O\nfQaDvzCG3IVQiS31u84xY7QxT0iA3r09Ps0jnZdfrmd7HjwI777rtrjVhtyqvnT+QrDLkJt701pC\nRacnGENucE9WFrzzjn79r39ZX79S8PLL+vV//wuHD1dYtMBRwLaMbQB0bNTRei0+4Pxiscu1YjBU\nhDHkLoSK38yvOj/6SD/o7NUL+vb16lSPdV58sd4OH650VL7z0E7yHfm0rNeSOrXqeKXFZ41ucPrI\n7RqRm3vTWkJFpye4NeRKqSuVUpuUUluVUk+XczxBKXVEKbW6aHveHqmGgJCbC/8pWo71+ef16Nku\n/v1v/XfcOMjLK7dIsPrHoeQXwraMbeQX5gdYjaEmUakhV0qFA6OBK4GuwM1KqfJivhaKyOlF2ys2\n6PQboeI385vOSZPgr7/gtNOgf3+vT/dKZ0ICdOsGqanw/fflFrHDkFvVlzGRMbSu35oCRwE7DlU9\nIVhFmHvTWkJFpye4G5H3AraJyC4RyQe+BgaVU86sDFtdGTdO/33qKXtH46Drf+AB/Xrs2HKLBPOI\nHEzkiiEwuDPkzYG9pfZTit4rjQC9lVJrlFKzlVKBzyvqA6HiN/OLzlWrYMUKaNBAZyysAl7rvO02\nqFsXFi6E9etPOOx8kGilIbeyL+184GnuTWsJFZ2eEOHmuCdz6lcBLUXkuFLqKuB7oNxwgiFDhhAf\nHw9AXFwcPXv2LP554+zUQO87CRY9Fe0nJyfb397bb5MAcMcdJC1bVqX6nHjV/u23k/TBB/Cvf5FQ\n5GJxHneOdA9uOEjS7qSg+Tyc+50a6QeeCxYsoFd+L0vrT05ODvj1Vaf9YO3PpKQkxhcllXPaS3dU\nmmtFKXUukCgiVxbtPws4RKTC3KNKqZ3AmSKS4fK+ybUSShw7pvOpZGbqbIfduvmv7T//hFNP1SPz\nv/7Sa38C6cfTafJWE+rWqsvRZ476vDKQHSzYuYBLJlzCeS3O4/d7fg+0HEM1wIpcKyuADkqpeKVU\nLeAmYIZLI81U0f8opVQv9JdDxolVGUKKKVO0ET//fP8acYDu3eGii/SXyaRJxW+X9o8HoxGHsj5y\nM3Ax+ItKDbmIFAAPAfOADcAUEdmolBqqlBpaVOwfwDqlVDIwChhsp2C7cXUJBCu26/zoI/33/vt9\nqqbKOks/9CwyiHY96LSyL0+qexL1atfjUM4h0o6nWVYvmHvTakJFpye485EjInOAOS7vfVjq9RjA\ns2xHhtAgORmWL4e4OL22ZiC49lpo0kS7WVatgjPPZHN60YNOi5Jl2YFSis6NO7N833I2p2+maZ2m\ngZZkqAGYmZ0uOB8+BDu26vz4Y/339tt1qlkfqLLOWrX0Qs1Q7F7ZdNCeEbnVfWlXCKK5N60lVHR6\ngjHkhrLk5sJXX+nX990XWC233ab/Tp4MBQXFhtE5FT5YcUaumFhyg78whtyFUPGb2aZz7lyd76RH\nDx054iM+6TzrLOjYEfbvJ2/eHHYc2kGYCqN9w/Y+6yqN1X1ZPCI/aK0hr/H3psWEik5PMIbcUJbJ\nk/XfW24JrA7QMz1vvx2A459/iEMctIlrQ1RExYs9BwNmdqfB35g1Ow0lZGZCs2Z6Lc3du6FVq0Ar\ngh07oF07CqJrE/dYLgnd+zHzlpmBVlUpuQW51Hm1DoKQ9VxW0H/xGIIbs2anwTumT9dG/IILgsOI\nA7RtC+efT0R2LtdsCt4cK6WpHVGbtg3a4hBHce50g8FOjCF3IVT8ZrbodD7ktNCtYonOooeet68t\neZBoJXb0pfMLZ2PaRsvqrNH3pg2Eik5PMIbcoElLgx9/hIiIwMWOV8SNN5Ifrrh0B5zmaBJoNR5h\n/OQGf2IMuQuhEltquc6pU6GwUK+f2bixZdVaoVMaNGBexzDCBbolnZgR0Vfs+MztiFypsfemTYSK\nTk8whtygcbpVnJNwgoi/Mv/iy66FANT9YV6A1XiGGZEb/Ikx5C6Eit/MUp179sBvv+lZnIPKWzek\n6lihc2P6RmZ2hNwIpXX+/bfvwkphp498U/omHOKwpM4aeW/aSKjo9ARjyA3w7bf6b//+xSljg4lN\n6Zs4Vhs2nNlKJ9CaNi3QktzSMLohTes05Xj+cVKOpgRajqGaYwy5C6HiN7NUp9OQX3+9dXUWYYVO\nZ+RH6hXn6zemTvW5ztLY9Zl3aayXt7UqcqVG3ps2Eio6PcEY8prO33/D779D7dpw9dWBVlMuzgeG\n4YOugchI+PVX2L8/wKrcY/zkBn9hDLkLoeI3s0zn999rd8Xll9viVrFCp9MQdmhzltbpcGjdFmHX\nZ+4ckVtlyGvcvWkzoaLTE4whr+nY6FaxgiM5R/gr8y+iIqJoVb9VySLQ33wTWGEeUDwpKN26SUEG\nQ3mYXCs1mYMHdW4VpeDAAWjQINCKTmD5vuWc88k5nNbsNNYMWwMZGVqziHYLNQneCUK7D+8m/t14\nmtVpRupTqYGWYwhRTK4VQ+VMn64nAV1ySVAacShxSzjdFDRsCJdeqnVPnx5AZe5pWb8lMZEx7M/a\nz6HsQ4GWY6jGGEPuQqj4zSzR6XSrXHed73VVgK86nREfZZJlOd0rFkWv2PWZh6kwSxeZqFH3ph8I\nFZ2e4NaQK6WuVEptUkptVUo9XUm5s5VSBUop+6yCwTqOHIH587Vb5ZprAq2mQpwRK8UjctCTlsLC\n4Jdf9HUEMV2aWPvA02Aoj0oNuVIqHBgNXAl0BW5WSnWpoNwbwFygUl9OsBMqsaU+65w1C/Ly4MIL\ntc/ZJnzVWe6IvHFjnWo3P1+vaOQjdn7mzoWirXjgWWPuTT8RKjo9wd2IvBewTUR2iUg+8DVQ3hzu\nh4GpQJrF+gx28d13+q+NbhVfyS/MZ/uh7SgUHRt1LHvQmUrAwjBEOzCx5AZ/4M6QNwf2ltpPKXqv\nGKVUc7Rx/6DorZAOTQkVv5lPOnNyYM4c/dpmt4ovOrcf2k6Bo4D4uHiiI6PLHnQa8tmz9S8LL8jP\nh9Wr4bPPYOxYeOyxJMaOhf/9D3bt0gExVmGla6VG3Jt+JFR0ekKEm+Oe3NKjgGdERJRSikpcK0OG\nDCE+Ph6AuLg4evbsWfzzxtmpgd53Eix6KtpPTk6u+vm//EJSVha0a0dC69a26nVSlfMX7V4E6FHt\nCcf37oU2bUjYuROSkkiqVavS+n74IYnZs2HdugRWrYLsbKe+BKfC4v2mTaF9+ySuuAKeey6BiIiq\nX/95F5xHmApj26pt/Pjzj1ze9/Iq90dycnLQ3H/VYT9Y+zMpKYnx48cDFNtLt4hIhRtwLjC31P6z\nwNMuZXYAO4u2TGA/MLCcusQQJNx/vwiIvPBCoJVUyqu/viokIk/MfaL8As8/r69j+PAK61i3TuS+\n+0Sio3VR59a+vcjNN+tThw8XGTZM5OqrRRo1KluuZUuRV18VSU+v+nW0f6+9kIj8uf/PqldiqLEU\n2c7KbXWlB/WIfTsQD9QCkoEulZT/HLiugmN+umxDpRQWipx0kv7oV64MtJpKuX3a7UIi8tGKj8ov\nsGKFvo5jV71FAAAgAElEQVTmzfV1leLIEZGHHxZRqsQoX365yLffVm6UHQ6RbdtERo0S6dCh5NyG\nDUU++0wf95b+X/UXEpFv1n/j/cmGGo8nhrxSH7mIFAAPAfOADcAUEdmolBqqlBrq2Zg/tHB1CQQr\nVdb5xx+QmgotW8Lpp1uqqTx86c/iyUBNTgiU0pxxBrRoAfv2wcqVxW9/9x107Qrvv6+jFIcPh40b\nYd48/Wy3UaOKNSoF7drBo4/Cpk06KCYhQU8ovftuuPhi2LzZu+uwKudKtb83/Uyo6PQEt3HkIjJH\nRDqJSHsRea3ovQ9F5MNyyt4lIsGfLLomM2OG/jtwoLZaQYqIFBu+MqGHpVGqTPRKbi7ce6821vv2\nQa9esGKFfqDZuYIqKiMsDK64QoerT5yoox4XLoQePfS+p5icKwbbcTdkt2rDuFaCg27dtK/gxx8D\nraRS9hzeIyQijd9sXHnBn34SAcnr1E1699aXFhUl8v77IgUF1mpKTxe5884Sd8uTT4rk57s/7/c9\nvwuJSM9xPa0VZKgR4KtrxVDN2L4d1q+HevWgT59Aq6mU9Wl6keVuTbpVXrBPHwpi44jcvJ79v2+j\nRQtYvBgeegjCw63V1KgRjB8PH3wAERHw3/9Cv35w+HDl53Vt0hXQrpVCR6G1ogwGTK6VEwgVv1mV\ndDqTTF19NRSF69lNVftz/QHPDPni5ZFMy74KgEfif2DFCu069wZvNQ4bBj//rF0tP/6oc45lZFRc\nvn5UfVrUa0FOQQ47Du3wTpwPOgOF0el/jCGvSTgNucULLNtB8Yi8acWGfNEi7cOeVjAQgIdaz7Az\n20AZLrpI+9/bt9eTi/r21VmBK8L5heS8LoPBSowhd8EZoB/seK3z4EG9An1kJFx1lS2ayqOq/enO\ntbJwob6MrCyof+OVSEQEYb8tqnxobLHG1q0hKQk6dIDkZD0yT6sgSUWxIT9QdUNebe/NABEqOj3B\nGPKawuzZeom0hASoXz/QaipFRNiQtgEof0S+dKn2DmVlwR13wNiv4lB9+ugc5c7UA36ieXNtzDt1\ngrVrdar08hIyOq/DjMgNdmAMuQuh4jfzWmfpsEM/UpX+3HNkD8fyjtG0TlMaxzQuc2z7dhgwAI4f\n10b8s8+KHmo6r8t5nTZrLM0pp8CCBdCxozbmN9yg87mUxgrXSrW9NwNEqOj0BGPIawK5uSXpXgcM\nCKwWD6jIrZKert0p6enaN/7JJ6UiU5zXNWeO10m0rODkk3XTTZvCTz/B0KFlk2+VjlwpcBT4XZ+h\nemMMuQuh4jfzSmdSEhw7Bqedph27fqQq/VlexEpOjk7UuHWrnpDzzTfa3V9MmzZw6qmQmakd6DZr\nLI+2bWHmTIiOhs8/h5dfLjkWWzuWVvVbkVeYx/aM7VWqv1remwEkVHR6gjHkNYEAuVWqimvEigjc\nd5+OD2/RQq+JERtbzok+uFes4uyz4euv9azQESNgypSSYyZyxWAXxpC7ECp+M491igTUkFelP11d\nK2PHwqRJUKeOfmbbvHkFJ5Y25F4kFbf6Mx84EN5+W7++5x49Bwt8j1ypdvdmgAkVnZ5gDHl1JzkZ\nUlK0E/fMMwOtxi0OcRQv79ataTd+/x0ee0wf+/RT7T2pkLPOgpNOgj179FPHAPLII3DrrTqy5rrr\n4OhRE7lisA9jyF0IFb+Zxzp/+EH/HTBA/973M972554je8jKz+KkuieRf7QhN9wABQXamN90k5uT\nw8JKHno6Jz/ZoNETlIIPP9RfPFu2wF13QdfGvhnyandvBphQ0ekJxpBXd0LNP17kdujauBuDB8Nf\nf+l1lt9808MKgsBP7qROHZg2TYftT5sGP03W6Ww3p28mvzDfzdkGg+cYQ+5CqPjNPNKZkqLzdMfE\n6GmHAcDb/nSOVrN3dyMpCZo102tplolQqYy+ffX1rlypr98Gjd7Qvj1MmKBfJz5Xl1Oi48l35LMt\nY5vXdVWrezMICBWdnmAMeXXG6Va57DIdExcCOA350h+0G2LiRO3e95joaLhcr4vJzJkWq6saAwdq\nn3lBARzeavzkBusxhtyFUPGbeaTT6V4IYJIsb/tz7d/awMn+bjz5pP4O8hqne8VDP7k/PvM33tD+\n8uO7qh65Uq3uzSAgVHR6gjHk1ZXMTL20jVLQv3+g1XiEQxysS9URK6ee3JWRI6tYUf/++rp/+UX3\nQxAQFQWTJ0PkYT3Dc+4qMyI3WIcx5C6Eit/Mrc65c/VU9d69oUkTv2gqD2/68/2JuygMO446djLf\nTGhA7dpVbLRJE33deXl6oU4LNfpCt27w+G16RL5813pPXfjFVJt7M0gIFZ2e4NaQK6WuVEptUkpt\nVUo9Xc7xQUqpNUqp1UqplUqpwDxVM5QlhHKPg14P+vn39Ci1U8NudOrkY4VBFL1Smn8P05Erjrgt\n3H1fnjfzlgyGClFSyZ2klAoHNgOXAvuAP4CbRWRjqTJ1RCSr6PWpwHci0r6cuqSytgwWkp+vszcd\nPqyXfO/YMdCKKkVEf9/8cPhV6PsvHj3nMUZd+Y5vlW7erFdcbtgQ9u/Xa7MFCW3eac+uo9th7Fo+\nHnkq994baEWGYEYphYhUulK6uxF5L2CbiOwSkXzga6DMEM9pxIuoC6RXRazBQn77TRvxTp2C3oiD\njkz54QeIbJkMwOkn9fS9Uue1Z2ToJC1BxBnNe+gXJyXzxBOwe3dg9RhCH3eGvDmwt9R+StF7ZVBK\nXaOU2gjMAR6xTp7/CRW/WaU6g8it4q4/U1J0aB5Aw67akPe0wpCDx+4Vf3/mPZvp62t/QTKZmTof\ni8Ph/rxqcW8GEaGi0xPc/d70yBciIt8D3yulLgQmAuV6OIcMGUJ8fDwAcXFx9OzZszgEyNmpgd53\nEix6KtpPTk4u/3ifPjB9OkkALVuSEODrcVLecRF4660EjhyBs3vP4Y+tW4lsF0mXJl2sab9VK339\n06eTVBTJEgyfX8+TesJOiG6SROPGeiHnp55KYuDAys9PTk4OCv3VZT9Y+zMpKYnx48cDFNtLd7jz\nkZ8LJIrIlUX7zwIOEXmjknO2A71E5KDL+8ZH7g/WrtUJu5s21fPbi1deCD4mTYLbb4e4OPjilyUM\nmtGbHs16kDws2ZoGCgv1bKK0NN0vlWbc8h97juyh9ajWNIpuxJjWaQwerIiNhQ0bdJpeg6E0VvjI\nVwAdlFLxSqlawE1Amd+pSql2SilV9PoMAFcjbvAjTjdC//5BbcTT0kqyGv73v7CvwGK3Cujrd7pX\nvv/eunp9pGW9ljSIasDB7IP0vnIfgwbpcPfhw73KvmswFFOpIReRAuAhYB6wAZgiIhuVUkOVUkOL\nil0PrFNKrQbeBQbbKdhuXF0CwUqFOp0GKwj841CxzkcfhYMHdWqUu+6C5FQbDDnoZYUAvvvOa412\noZQqvs61+9cwdqxOrDVzZtmFKFwJ+XszyAgVnZ7gNo5cROaISCcRaS8irxW996GIfFj0+k0R6S4i\np4vIhSLyh92iDRWwe3dJkqwqzW33D7Nm6VmOMTHw0Ud6Euaa/WsAGwz5pZdC3bqwejXs2mVt3T7g\nvM7k1GROOQXeeku//8gj+gvOYPAGM7PTBefDh2CnXJ3O0fjVVwdNkixXnUePwrBh+vXLL+t1Lgsd\nhazdrxeC6NGsh7UCoqL0is1QoXslEJ+58zqT9+tfIvfeCwkJ2uX0+OPlnxPS92YQEio6PcEY8urE\ntGn677XXBlZHJTz/vA45PPts7V4B2JqxleyCbFrVb0WD6AbWN+qBe8XfOEfka1L1LxGl4OOP9ffO\nxInw00+BVGcINYwhdyFU/GYn6DxwABYt0om7+/ULiKbyKK1z+XIYPVo/g/z445Jnsbb5x53066f7\n5bff9JC3Eo3+okuTLkSGRbItYxuZuTqxV/v28MIL+vjw4ZCdXfackL03g5RQ0ekJxpBXF5wLDl96\nqX5yFmQUFMD992uJTzyhIySdOEellrtVnNSvrxfWcDhKcrQHmFrhtejapCuCsO7AuuL3n3oKuneH\n7dvhlVcCKNAQUhhD7kKo+M1O0Ol0q1x3nd+1VIZT56hRsGYNtG4NI0aULeP0E9s2IodK3SuB+sxd\n3Sugfzg4HwC/+Sb8+WdJ+ZC9N4OUUNHpCcaQVweOHIH58/Xiw0G4NueuXSXG+4MP9FqWpbHdtQI6\nHFMp7XwOkhzlpSNXSnPeefqBsPNXjCfT9w01G2PIXQgVv1kZnbNn64yHF1ygZ3QGEQsWJPHgg3D8\nONx0U0kAiZP9x/aTeiyVerXrER8Xb5+Qk0+Gc8+F3FyYM6fMoUB95sWGfP+JM1lfe01LXrJEP0+A\nEL03g5hQ0ekJxpBXB4LUrQLw66/6e6Z+fe1eccUZP35as9MIUzbfjv/4h/77zTf2tuMhzmcC6/av\no9BRWOZY/frw7rv69TPP6Ey8BkOFiIhfNt2UwXKOHxeJiREBkV27Aq2mDEeOiJxyipY2dmz5Zd74\n7Q0hEXlo1kP2C9q9W4uJjhY5dsz+9jyg1TuthERkY9rGE445HCJXXaUl33JLAMQZgoIi21mpfTUj\n8lBnzhzttzjzTP0kMYh4/nmdt+ucc2Do0PLL+MU/7qRVKy0mO1v/TAgCKvKTg3bpjxmjY8u/+srE\nlhsqxhhyF0LFb1as83//039vuilgWspjxQodMx4WlsSHH+rnsOWx6u9VgJ8MOcANN+i/pdwrgfzM\nnbnJV/61stzjbdqUxJYPGZJETo6/lFWdkPs/VA0whjyUycoqiYu+8cbAailFQYEegYtot3SPCsLD\nD2UfYvPBzdQOr82pzfyUYtbpJ581S/+SCTC9mvcCYPlfyyss8+ST0LWr/nXz6qv+UmYIKdz5Xqza\nMD5y65kyRTtQzzkn0ErKMGqUltWypUhmZsXl5m2bJyQi531ynv/Eiej+ApFvvvFvu+Vw4NgBIRGJ\nGRkj+YX5FZZbtEhLjowU2XiiO91QjcH4yKs5QehWSUnRvnHQrpW6dSsuuyxlGQDnND/HD8pK4XSv\nOPsvgDSp04S2DdpyPP846w+sr7DcBRfoJeHy83WMuclbbiiNMeQuhIrfLGn2bO0egBJ3QRDw6KNw\n7JieSDlwYOX9uWxfkSFv4WdD7uJeCfRn7vwic/ZHRQwYoJeGW7gQJkzwh7KqEej+9JRQ0ekJxpCH\nKkuWQE4OnH8+tGwZaDWAdtdPm6ZH4e+9V3lZESkx5P4ekbduDb16aR95EESvFBvylMoNef368Pbb\n+vWTT0J6ut3KDKGCMeQuhEr+hYR1RYmWguQhZ1YWPPSQfv3yyyXfLRX1587DO0k/nk6TmCb2zuis\nCGe/TZkS8M/c+YvE3Yg8ISGB227T+b8OHoT/+z9/qPOeQPenp4SKTk8whjwUOXpUx48rFTRulREj\nYM8eOP30EoNeGcX+8RbnULTkq3+58UbdfzNn6lw1AaTnST2JDItkQ9oGjuYerbSsUjpfTa1a8Pnn\nUI28AwYfMIbchZDwm02fTlJeHlx0EZxySqDVsGoVvPOOjhX/6COIiCg5VlF/Bsyt4qRlS+jTB3Jy\nSApwTF9URBQ9T+qJIKz4a0WF5Zx92bEjPPecfm/oUIIutjwk/g8ROjo9wSNDrpS6Uim1SSm1VSn1\ndDnHb1VKrVFKrVVKLVZKnWa9VEMxX36p/wZBtErpDH2PPAJnneXZeQE35AC33qr/BsGUSU/95E6e\neQY6d4YtW3SCLUMNx118IhAObAPigUggGejiUuY8oH7R6yuBpeXU46+wy+rNvn0iYWE6oDg9PdBq\n5J13PIsZL01uQa7Ufrm2kIgcyj5kr8DKOHRIpFYtEaVEUlICp0NEJq6ZKCQigyYP8vicX38tiS1f\nv95GcYaAgkVx5L2AbSKyS0Tyga+BQS5fBktExOloXAa08OnbxVAxX32lh7/9+0OjRgGVsmdPScz4\nmDGVx4yXZk3qGnILc+ncuDNxUXH2CXRHXBwMGKCDsidPDpwOyoYgiodB4hdeCPfdp2PLhw41ectr\nMp4Y8ubA3lL7KUXvVcQ9QOBjuqpI0PvNJk4EIOmMMwIqQwQeeEBHq/zjH9oelkd5/RkUbhUnt95K\nEpS4qwJE+4btaRjdkNRjqew9urfcMuX15RtvQLNmejlSZ97yQBP0/4eKCBWdnhDhvggezyFTSl0M\n3A2cX97xIUOGEB8fD0BcXBw9e/YsDgFydmqg950Ei54y+9u2kbB2LTRsSHJUFCQlBUzPCy8kMWsW\n1K+fwLvvetefy/Ytg53QoFmDco/7df/qq6FOHZKSk+Hzz0m4666A6Fm4cCHtjrQjo1YGy1KWsWP1\njhPKJycnl3v+u+/C4MFJPPEE9O+fQPPmQXK/Bvl+Rf0Z6P2kpCTGjx8PUGwv3eLO9wKcC8wttf8s\n8HQ55U5D+9LbV1CP7b6kas8TT2in6AMPBFRGWppI48ZayiefeH9+h/c6CInIyr9WWi+uKtx/v76Y\nZ58NqIwRC0YIiciT85706jyHQ2TgQH0JAwbofUP1AYt85CuADkqpeKVULeAmYEbpAkqpVsA04DYR\n2ebZV4jBKwoKtH8c4PbbAyrlscf0rMJLLoG77/bu3IzsDLZmbCUqIopTm/op46E7nNErX34ZUEez\n09W0NGWpV+cpBWPHQr16enZtEKSQMfgZt4ZcRAqAh4B5wAZgiohsVEoNVUo5lwt4AWgAfKCUWq2U\nqjgnZ5Dj6hIIGubPh9RU6NABzjknYDpnzdL2Ljpa+2TdzeVx1blo9yJAp2+NDI+0SaV3JBUUQHy8\nfnr7888B03Fui3NRKP746w+y87NPOF7ZZ968OfznP/r1ww/rmZ+BImj/D7kQKjo9waM4chGZIyKd\nRKS9iLxW9N6HIvJh0et7RaSRiJxetPWyU3SNpOghJ3fc4d562sTRozB8uH79yivQtq33dfyy8xcA\nLom/xEJlPhIWVvLT4pNPAiajQXQDTj/5dPIK8/h97+9en3/vvZCQAGlp+leToQbhzvdi1YbxkVed\njAyRqCjtBN25M2Ay7r1XS+jVS6SgoGp1dB/bXUhEFu5aaK04X9m7tyQ+/8CBgMl4ct6TQiLy3Pzn\nqnT+1q16SVIQmT7dYnGGgIDJR15NmDBBz8O+7DLtAggAc+fqwWrt2jB+PISHe1/HgawD/HngT6Ij\nooMj9LA0LVrAVVfpoOwA5oi9pI3+pfLLrl+qdH779iUzPYcODayLxeA/jCF3Iej8ZiIwbpx+PWxY\n8dv+1Hn4sP7ZDjqzYZcunp9bWmfSLv36glYXUDuitnUCfaRY43336b+ffBKwlRsubHUh4SqcP/b9\nQWZuZpljnn7mDz+sJwulpuq0Cf4m6P4PVUCo6PQEY8iDnUWLYNMmOOmkimfd2Mzjj8O+fXDuufDE\nE1Wvx+kfvzj+YouUWczVV+t+3rQJFi8OiITY2rH0at6LQilk0Z5FVaojLExnRoyJ0YFO06ZZLNIQ\ndBhD7oIzQD9o+PBD/feeeyCyJMrDXzpnztSulKioqrlUSussftDZJogedFJKY2QkFE0ICuQ0SecX\nnbO/nHjzmbdrp2d9gv4hd+CAVercE3T/hyogVHR6gjHkwUx6OkydqqNUnD/7/cj+/SXBHCNHQqdO\nVa8r5WgKWzO2ElsrljNPOdMagXZwzz367zffaJ9SACj2k++smp/cyQMPwMUX6yiW++4z63xWZ4wh\ndyGo/Gbjx0Nenn4I17p1mUN26xTRfvG0NG0MqhrO5tS5YOcCAC5qfRERYZ5khvAfZfqyXTs90yk7\nuyTk08/0btmbWuG1SE5NJiM7o/h9bz/zsDB9C9WvDzNm+O9HRlD9H6qEUNHpCcaQBysOh16lAXT4\ngZ/56CPtVomLgy++0EbBF5xRGMHmVimXBx7Qf997LyAzPaMjozmvxXkIwsJdC32qq1UrvaIQ6Gcd\nW7ZYINAQfLiLT7Rqw8SRe8e8eToYuHlzkfx8vza9aVNJLPKUKb7X53A4pNU7rYREZNVfq3yv0G7y\n80Vat9YdMGNGQCS8mPSikIg8NOshS+q79VZ9OWefLZKXZ0mVBj+BiSMPYf77X/13+PCya6fZTG6u\nTj2Sna1TulixtvOOQzvYc2QPDaIa0OOkHr5XaDcRESVxe++8ExAJvsaTuzJ6tB6d//GHXl/VUL0w\nhtyFoPCbrV0LP/6o48ecc+JdsEvn//0frFwJbdrA++/7Xl9SUhILdmn/+MVtLiZMBd8tV25f3nOP\nXiljwQJITva7pl7NexETGcOGtA2kHksFfPvM4+K0yz8sDF5/Xd9edhEU/4c8IFR0ekLw/a8ylIzG\n77kHGjb0W7PffafdwpGRMGWKfkhmBbO2zgKgb5u+1lToD+rXLwnZefddvzdfK7wWF7W+CIA5W+dY\nUudFF0Fion6Qfdtt8PffllRrCAKU+CkmSSkl/morpElJ0cNhhwO2bdOv/cCuXXD66Tri7p13rEu6\ndDz/OI3fbEx2QTZ7H99Li3ohtArg9u0622RkJOzerScL+ZFxK8YxfNZwBnQcwIybZ7g/wQMKC+GK\nK3SSx4sv1utOVyXdgsF/KKUQkUoz5ZkRebDx3ns69/g//uE3I56XB4MHayM+cCA8+qh1dc/bNo/s\ngmx6Ne8VWkYcdCjiwIG6g5yhH35kUKdBKBQ/bv+RY3nHLKkzPBwmTdLLwy1YoFMuGEIfY8hdCKjf\n7OjRkpmcTz1VaVErdT7+OCxbBi1b6qndVmbJHTt1LADXdr7WukotptK+fPxx/Xf0aMjMrLicDZwc\nezLntTyP3MJc5mydY9lnftJJOqe8UvDSSzDb4hV2Q8X3HCo6PcEY8mDi44+1Mb/oIjj7bL80+dln\nenWZWrX0JFIrXfL5hfks2bsECG5DXikXXQTnnw8ZGdY8/fUSZ799t+k7S+vt21cbcRG45RbtxTOE\nLsZHHixkZemVGg4c0Ot19e9ve5PLl+sseXl58Omn3i/b5o6ftv/E5ZMup2uTrqx/YL21lfuT+fN1\nCuGGDWHnTr2mmp/YnrGd9u+3p17teqT9M41a4bUsq9vhgOuvh++/h27dYOlSHahjCC6MjzyUeO89\nbcR79YJ+/Wxv7sAB/Z84L08nVbLaiEPJKDJkR+NO+vaFCy7Qo/LRo/3adLuG7Ti16akczT3qc+4V\nV8LC9Kzdzp1h/XqdL8yMtUITjwy5UupKpdQmpdRWpdTT5RzvrJRaopTKUUo9ab1M/xEQv9nhw/Dm\nm/r1yJEeOal90ZmTA9deqwNkeve2J7rOIQ6+3/Q97Ax+Q+62L5XScXugF8Y8etRuSWW4rst1AIz+\nn/VfIvXq6bDT2FjtWnvpJd/rDBXfc6jo9AS3hlwpFQ6MBq4EugI3K6VclxY4CDwM/MdyhTWBt9/W\nxjwhQY/+bMThgCFD4Pff9aI4U6dq/7jVLEtZxt/H/qZp3aaccfIZ1jfgby65RPuhDh3yu6/c+UW4\neO9iCh2FltffubPOWx4Wpr+vJk2yvAmD3bibww+cB8wttf8M8EwFZUcAT1ZwzNZ8BCHLgQMidevq\nRBiLF9ve3L/+pZuKjRVZs8a+dv754z+FROTROY/a14i/+fln3XkNGogcPuy3Zh0Oh7QZ1UZIRBbt\nXmRbO++9py+vVi2RhUG2pGpNBotyrTQH9pbaTyl6z2AFb7wBx47p1Wl697a1qc8/156b8HCdbvu0\n0+xpR0T4duO3QPC7Vbzi4ouhTx89KrfCB+EhSqnifpy6Yapt7Tz8sE4xk5enXW8mU2Lo4Ikhr1GP\nP/zqN9u5E8aM0a9fecWrU73V+cMPJWtTjBmjZ/fZxYJdC9hxaAfNY5tTsKPAvoYswuO+VEq7wZTS\nD6c3bbJVV2luPvVm2AkT104kpyDHtnbefluvKJiRoe+Rffu8ryNUfM+hotMTPEmrtw9oWWq/JXpU\n7jVDhgwhvmgV+Li4OHr27Fm83JKzUwO978Qv7T37LAk5OXDrrSQdOQJJSR6fn1yUyMmT8klJcP31\nSRQWwrPPJjB0qL3XN27FONgJl/W8jPCwcPv6L1D799xD0iefwB13kLBsGShle/vHthyjxfEWpGSn\nMHXDVFpktLCtva++gl69kti4ES67LIFff4U//7T3+gKxn5ycHFR6nPtJSUmMHz8eoNheusWd7wVt\n7LcD8UAtIBnoUkHZRIyP3DOmT9cOyXr1RP7+27Zmli8vccEPHy7icNjWlIiI/J35t0S8FCHhL4ZL\nypEUexsLFPv3i9Svrzv1hx/81uzHKz8WEpHzPz3f9rYOHhTp3l1f4llniRw5YnuThgrACh+5iBQA\nDwHzgA3AFBHZqJQaqpQaCqCUOkkptRd4HHheKbVHKWWmFlTE8eMl+a5fftm2ZEzr1ulV4o4dg5tv\n1iHQVk6/L4/PVn9GgaOAAZ0G0LxeNX2U0rRpSVLvxx/XSdz9wODug4mtFcvivYtZt3+drW01bKhT\n3bZtCytWaHdLVpatTRp8wZ2lt2ojREbkCxYssL8RZ+hIz55VXv3Hnc5Vq0QaNdLN9Ovnn1VhCgoL\npPU7rYVEZO7WuR7pDAaqpDEvT6RzZ93BI0darqk8FixYIA/MfEBIRB6c9aBf2tyxQ+SUU/RlXnih\nyNGj7s8Jhc9cJHR0YlYICkI2b4a33tKvx4yxZfWf5ct12PPBg3qS6NSpOhOr3czbPo/dR3bTJq4N\nl7W7zP4GA0lkZEk8eWKiXgzEDww9S6/fOmHNBMsyIlZGmzY6S2Lz5rBoEVx+uZ7yYAgy3Fl6qzZC\nZERuK7m52uEIInfdZUsTv/2mY8RB5NprdZP+YsBXA4RE5PVFr/uv0UAzbJju7NNOE8nJ8UuTvT/t\nLSQiH6/82C/tiYhs316yjOmZZ4qkp/ut6RoPZkQeZLzwgnY4tm6t47ws5rvv4NJLdbbVm27Sq/zU\nsmHWZnnsObKHWVtnERkWyV2n3+WfRoOBt97SecvXroUXX/RLk8POHAbABys+cA6SbKdtW1i4UP9d\nuULQWR0AABFUSURBVFInhNy50y9NGzzAGHIXnGFAlvPLLzqfSliYTgYdF+dTda46339fJ8HKydHx\n4pMm+ced4uTlhS/jEAc3dLuBpnWaVqgzGPFJY926OvOUUnpy1++/W6bLFafOG7rdQJOYJqz6exUz\nt8y0rT1XWreGX3/VE8k2b4Zzz9Xjkop0BjuhotMTjCH3BwcP6iXpReDf/9bDGYsoLNRrUDzyiK7+\nlVf02hQ2uN4rZFP6Jj5L/oxwFc6IPjVwifbzz9erVjsc+nPOyLC1uaiIKP514b8AePbnZ23Jv1IR\nTl/5pZfqDJp9+ujJZoYA4873YtVGTfWRFxSI9O+vnYu9e1c5SqU80tNFLrtMVx0RITJhgmVVe8V1\nU64TEpH7Z9wfGAHBQE6OyOmn6w+jb1/bw4Ry8nMkflS8kIiMXz3e1rbKIzdX5M479eUqJZKYKFJY\n6HcZNQI88JEbQ243jz4qxYmWdu60rNqVK0sePjVpIhKoSKqle5cKiUj0K9Gy7+i+wIgIFnbvFmna\nVH8oDz1ke3MTkicIiUird1pJdn627e254nDoyEulSsJcMzL8LqPaYwx5FbA0tvTdd3UXR0aKJCVZ\nUqXDITJ2rEhk5AIBkV69RPbssaTqKmhxSJ/P+wiJyLPzny23TCjE6lqqcfFinT4QRD780Lp65USd\nBYUFcurYU4VE5O3f37a0LW+YO1ekYUN9yW3biowZsyBgWrwhFO5NERO1ElhmzIDHHtOvP/9cOxN9\nJDVVrwD3wAOQnw/3368fPrVs6f5cO5i7bS4Ldy+kQVQD/u/8/wuMiGCjd++SBbQffBDmzLGtqfCw\ncF7r+xoAryx6hSM5R2xrqzKuuEJHspx+OuzYobMoJibqe9TgJ9xZeqs2QmREbgnz54vExOghyksv\n+VydwyEydapI48YlXpr//c8CnT5wOPtwcY7stxa/FVgxwcg//ynFyb1nzbKtGYfDIRd9fpGQiNz1\nvT1zEzwlO1vkqadKXC29eols2BBQSdUCjGslAMyYIVK7tu7ae+7xOUvVjh3a96hjUvTDzZQA56Jy\nOBwyeOpgIRE5fdzpkpPvn4kwIYXDof3kTmM+c6ZtTf25/0+JeiVKSES+WvuVbe14yoIFIi1blngV\nn3tOJCsr0KpCF2PIq4BPfrOvvhIJD9fd+sADPj3Gz87WD5Kio6U4SeKYMSVVBtK/9+mqT4VEpM7I\nOrI5fXOlZUPBD2mbRodD5OGHS4z5d9/5VF1lOsf9MU5IRGJfjZVtB7f51I6vLFiwQA4dErnvvpIB\nSHy8HuPYnX3TG0Lh3hQxPnL/IQL//S/ceqsO7H76aZ1qMMz77i0s1C71Dh3gX/+C7Gy45RY9AeOB\nB6pUpaVsTNvIw3MeBmBsv7F0bNQxsIKCGaX0ytbOZXeuu04H+jscljd1/5n3c32X68nMy+Tmb28m\nrzDP8ja8IS4OPvoIFi/WE4h27YKBA/WytEuXBlRa9cSdpbdqI0RG5F5z6JDINdeUDD2qmAmvoED7\nvbt1K6mqRw/tbg8WMo5nFEdJ3DbttkDLCR1c4/Suu86zNIJeknE8Q1q906o4O6IjSIa/+fki77xT\nEtnizAO0alWglYUGGNeKzaxYoeOtQC80UIWfztnZIuPGibRvX/Zn6KRJwTXBIj0rXc748AwhEenw\nXgc5mmO9Iar2zJpVsiBFly4iS5ZY3sTiPYsl4qUIIRF5ZPYjQWPMRfR61c89V+IudD7zmT8/uFwu\nwYYx5FXAI7/Z4cMijz1W4g8//XSRbd75Jbds0YENTZqU3NRt2+oYcU+S6PnTv3fg2AE57YPThESk\n3bvtZPfh3R6fGwp+SL9q3LxZG3HnlMiHH/Z4dO6pzumbpkutl2sJicjwmcOl0OHfEYE7nfv2iTzx\nhEidOiX3/qmniowerX/g+otQuDdFjI/cehwOnY2qc2cYNUrfg48+qhMltWvn9vRDh+Czz3Su8I4d\ndeK8tDQdf/v119oPPnw41K7th2vxkD1H9nDxFxezdv9aOjXqxMIhC2lVv1WgZYUuHTvqoOtnntEP\nPN5/H7p1gwkToMCahaoHdhrI9MHTqR1emw9WfMC9M+4lt8A/qxh5wimn6EdKe/boBbKaNtWrWT30\nkD52550wd66JQ/cKd5beqo0QGZGXS26uyGefiXTqVDKE6N1bZPVqt6f+9Zc+tX9/HYrlPD0mRuTu\nu0WWLg3On5UOh0M+WfmJxL4aKyQiXcd0lb8z7VtbtEayerVO7u28KTp00AlzLMrH89P2nyT6lWgh\nEek2ppus2LfCknqtJidHZMoUkUsuKekK0Ctc3X+/9kjV5PBFjGvFR7Zv1xN6mjcvubtatRL5/PMK\nHdiZmSI//qh9gWecUfbGDAsTufRSkY8/1t6ZYGXP4T1y1aSrhESEROTar6+VtKy0QMuqnuTni4wf\nL9KuXcmN0rKlyL//rScR+MiylGXS8f2OQiIS/mK4/PuXfwckL4unbN0q8uKLJd4n5xYVJXLVVSJv\nv63zDBUUBFqp/7DEkANXApuArcDTFZR5r+j4GuD0Csr46bJ9wOGQBRMnirz3nsj555e9k7p1E5k4\nsUxWu4ICkfXrRb74Qrs6e/UqcZs7t+hoPaFn7FiR1FTrpNrh31ubulbu/O5OiXwpUkhEGrzeQL5c\n+6VPD8xCwQ8ZFBqdBr1Dh7I3UEKCtl6bN1dZZ1Zeljw25zFRiUpIRJq+1VReSnrJti9nK/rT4RBZ\nu1ZkxIiSRbVKb/XqiVxxhf6+mz5dT5Lz9jYNis/dA3w25EA4sA2IByKBZKCLS5mrgdlFr88BllZQ\nl98u3GPy8kSSk0U+/VTPXoiPl3dK3y0xMSK33iqZU+fKiuWFMmWKjiK79Va9bnJU1Ik3WHi4yNln\nizz5pMjs2SLHj9sj/Z133rGknl2HdsmY5WOk7xd9i0fgYS+GyU3f3GRJNkOrdNpJUGksLNRTI2+7\n7YQb7J1GjfRs4Y8/Flm3zuth6cJdC+X0cacXf87Rr0TLnd/dKd+s/0aO5Byx7BLs6M/UVP09d9dd\nJYFirlvDhnqB6GHDdL66mTNFNm6sOHggqD73SvDEkLtbfqAXsE1EdgEopb4GBgEbS5UZCHxRZKmX\nKaXilFLNRGS/l+56ezh+HPbtg5QU2LsXtm9HtmzBsWkLauMGwnJzyhQ/EBHN8lYDSarbn2/yrmHL\nD3U5+mXF1bdqBWefDWeeCWedpVdNiY21+ZqAw16ugCsiHMo5xKb0TSSnJpOcmszve39nfdr64jIx\nkTHc3fNuHjv3Mdo1dP/w1g6dgSCoNIaF6VkzCQn6Qejs2TBrFsydy+GDB+HTT/UGEBUFnTpBly76\nAXx8vL4hW7aEk0+GOnXKVH1R64tYef9KFuxawH9+/w9zts3hizVf8MWaL4gMi6R3y96cdcpZ9Dyp\nJz2a9aB9w/ZER0Z7fQl29GezZvoh6J136v2UFFiyRD83XrFC/83I0IteLFpU/vktW+rtlFP0/vLl\nh2nTBho1goYN9d/69XW3hhruDHlzYG+p/RT0qNtdmRbACYb8+/98Do6iL1CHNi44HFAoiEOg0IEU\nOlCFhXq/oBAKClH5BZBfiMrPJywvH/LyCc/LJTw3l/CcHCLycqmVfZxaOVlE5WRRJ/sIdXOOEJt7\nhJiC7BMuSqF/agBsC2/BisgurIjsyoLaZ7EqeyGvZQ+DbIDNUAeiGumVUVq00FubNtC2HbRto1f6\nKs2WTCATBKmwU/WXbNFrl3IigiBl/jrEUbwVSiEFjgK2HNzC95u+J68wj7zCPHIKcjief5zj+cfJ\nzM0kIzuDQzmHOJh9kJSjKew9spes/KwTtMTWiuXydpfTr0M/BnUeRMPohhXqNviZuDg9rfeWW3RE\ny7Bh0L27nhq5dCns3g1r1uitPKKjdUiI00LFxaHq1+eSOnW4JKY76aodfx7dxtrDm9mcuYvsFQvZ\nH76QWWHwfRgUhEGd6Fga1mlCXN1G1ImqR93a9agbXY/aEVHUjoymdmQUkRG1iAyvRXhYBBHhEaSs\nX8of00YTpsIICwtHKUWYCtMzXQHl/FfJNGVVVLZop+QlJa+7x0L3BLgzQZuUQ4dhXwqk7IPUv2H/\ngf9v715DpCrjOI5/f866aqXdjLIUCjKooDIqpIJWohCJzBdRQnSD6EVRREQXXxRE2OVFQREV3YhI\ng0rRNGy7SEJkBJqWu2wX1krLC5V0L/XXi+fsNu3O7IyxznMm/h94mDkzz5nz27PMM+fyPOfA9m2w\ncyfs3QN7+2FzP2wu5t/6ywa2fvDisGV3dKRVd8AE6BwH48elxr2z858ytjPVGzs2lUoFOipQ6UiP\nYyrpN7lSgTFKr0vptYHHMQKUHgdWjwb+XBXTGhavpkYNef3W6N+GLq7mfJfcdm2THzd6/qjAlonw\nzSTYMgm+PBT6Dk+ldzL8OOEb0m9Pd5phCTDv6X99xu/AF0UBYEdRcg41XgOLJi/ap1kmdk7k+MOO\nH9zimjFlBjOnzqSzsv/u0Nzf37/fPnu0tENGOjro3737n0sjA+zaBT09qfT1pf58A2XbtnR9h82b\nU6lhMtBVlPqKLRO+bDrqU8CZr6xqun4uVwPP/7Ykd4xRoeqtw2FvSjOBe2zPLqbvBPbafqCqzhPA\natuLi+le4Lyhh1YkNfujEEIIoYrtEbfNG22RfwRMl3QssBW4DJg/pM4y4EZgcdHw/1jr+HijICGE\nEP6bERty27sl3QisIh1WfsZ2j6Tri/eftL1S0hxJnwO/ANfs99QhhBAGjXhoJYQQQvm19Forku6V\n9LGk9ZLelpTpbpP1SXpIUk+R8zVJB+fOVIukSyV9KmmPpNNz5xlK0mxJvZI+k3R77jy1SHpW0jZJ\nG3NnGYmkaZLeLf7fn0i6KXemWiSNl7S2+H5vkrQwd6Z6JFUkrZO0PHeWeiT1S9pQ5PxwpLqtvmjW\ng7ZPtX0asBS4u8XLb8abwMm2TwX6gDsz56lnIzAPeC93kKEkVYDHSKOCTwLmSzoxb6qaniNlLLu/\ngFtsnwzMBG4o4/q0/Tswq/h+nwLMknRu5lj13AxsovmeeTkY6LI9w/ZZI1VsaUNu+6eqyYOAna1c\nfjNsd9seuIXLWlKf+NKx3Wu7L3eOOgYHktn+CxgYSFYqttcAP+TO0Yjt72yvL57/TBqQd3TeVLXZ\n/rV42kk6r/Z9xjg1SZpKGpH+NMO7TpdNU/lafhlbSfdJ+gq4Cri/1cvfR9cCK3OHaEO1BokdkynL\n/0rRg2wGaSOjdCSNkbSeNCDwXdubcmeq4WHgNmD077k3ugy8JekjSdeNVLFR98N9JqkbOKrGW3fZ\nXm57AbBA0h2kFdryXi6NMhZ1FgB/2n6ppeGqNJOzpMq8u9q2JB0EvALcXGyZl06xN3tacW5plaQu\n26szxxok6SJgu+11krpy52ngHNvfSjoC6JbUW+xFDjPqDbntC5qs+hKZtnYbZZR0NWnX6/yWBKpj\nH9Zl2WwBqk9kTyNtlYf/SNJY4FXgRdtLc+dpxPYuSSuAM4DVmeNUOxu4WNIcYDwwSdILtq/MnGsY\n298WjzskLSEdsqzZkLe618r0qsm5wLpWLr8ZkmaTdrvmFidv2kHZjvMNDiST1EkaSLYsc6a2pXTh\nkWeATbYfyZ2nHkmTJR1SPJ8AXEDJvuO277I9zfZxwOXAO2VsxCUdIGli8fxA4EJSB4eaWn2MfKGk\njcUxtC7g1hYvvxmPkk7Edhfdfh7PHagWSfMkfU3qxbBC0hu5Mw2wvZs02ncVqWfAy7Z7Rp6r9SQt\nAt4HTpD0taSyDmY7B7iC1AtkXVHK2NtmCvBO8f1eCyy3/XbmTI2U9TDgkcCaqnX5uu0361WOAUEh\nhNDm4ubLIYTQ5qIhDyGENhcNeQghtLloyEMIoc1FQx5CCG0uGvIQQmhz0ZCHEEKbi4Y8hBDa3N/i\nHS6byWjnnQAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x = np.linspace(-3, 5, 100)\n", "parameters = [(1.0, 1.0), (1.0, 0.5), (0.0, 0.7)]\n", "for mu, sigma in parameters:\n", " plt.plot(x, norm.pdf(x, mu, sigma), label=\"$\\\\mu$={} $\\\\sigma$={}\".format(mu, sigma))\n", "plt.legend();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The key idea for our range is that the points divide a normal distribution into segments of equal areas. So first we need to look at the integrals. Normal distributions are always normalized, so $\\displaystyle\\int_{-\\infty}^\\infty f(x) dx=1$. The integral up to a certain point is given by the cumulative distribution function $\\displaystyle\\text{cdf}(x)=\\int_{-\\infty}^x f(x) dx$. For the normal distrubution, that is\n", "\n", "$$\\text{cdf}(x) = \\frac 1 2 \\left(1+\\text{erf}\\left(\\frac{x-\\mu}{\\sigma\\sqrt 2}\\right)\\right)$$\n", "\n", "Let's have a look:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAEACAYAAACuzv3DAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmczdX/wPHXuWNmDMOMfYvsa5ZKhDBfVLZSUmgdSVQk\nWiyVrVWlhKgsqV+hVEKEyCBE9n3JvosxYzBjzNzz++PMRpi55n7mc5f38/G4j+Yz9zOfc5w+855z\n35+zKK01QgghvJfD7goIIYTIHgnkQgjh5SSQCyGEl5NALoQQXk4CuRBCeDkJ5EII4eUyDeRKqUlK\nqRNKqc3XOWeUUmq3UmqjUupW91ZRCCHE9WSlR/4V0PJabyqlWgMVtdaVgGeBcW6qmxBCiCzINJBr\nrZcBZ65zyv3A1ynnrgLClVLF3FM9IYQQmXFHjrwUcCjD8WHgJjdcVwghRBa462GnuuJY5v0LIUQO\nyeWGaxwBSmc4vinle5dRSklwF0KIG6C1vrKzfBl3BPJZQE9gmlLqTiBGa33iGpVxQ3Heb8iQIQwZ\nMsTuangEaYt0ntIWcRfj2HxyM9v+3caOUzvYfmo7u0/vZn/Mfi45L2X683kC81AkTxEK5SlEwZCC\nFMhdgAK5CxCWO4yw4DDCcoeRLygfoUGhhAaFkjcoL3kC85AnMA8huULInSs3I98fyZAhQwgKCEKp\n68Ywn5eVf3+mgVwpNRVoChRWSh0CBgOBAFrrL7TWc5VSrZVS/wDngS7ZqrUf2L9/v91V8BjSFuns\naIu4i3GsObqGvw7/xZpja9h4fCN7zuy55vklQktQrkA5yoSVoUz+MpQOK03JfCUpEVqC4qHFKRZa\njDyBebJdr2OHjxGcKzjb1/EXmQZyrXXnLJzT0z3VEUJYKTo+mqUHlhK1P4qo/VFsPrkZp3Zedk5Q\nQBDVi1TnlqK3ULVQVaoVqUblQpUpX6C8W4K0cD93pFaEiyIjI+2ugseQtkhnRVs4tZN1x9YxZ9cc\n5v4zl7+P/I3OMBYhlyMXt5W4jfql6lOvVD1uLX4rVQtXJTAg0O11cYXcF65ROZW3VkppyZELYb1k\nZzIrDq1g+rbp/LT9J47GHU17LyggiAY3NSCibAQRZSOoX6o+IYEhNtZWZEYplenDTgnkNoiKiiIi\nIsLuangEaYt02WkLf38g6CuuFiOzEsgltSKEF7uYdJHvt34PyKgwb5edP8bSIxfCCx2LO8a4NeP4\nYu0XnDx/EoZIIPd2KT3va31feuRC+Ir9MfsZ/udwJm2YRGJyIgC1i9VmIxttrpmwk6xHboOoqCi7\nq+AxpC3SXa8t9sfsp8vMLlQaXYnP137OpeRLtK/WnqWRS1nffX3OVVJ4JAnkQniw6PhoXlnwClXG\nVGHyhslorXmi1hNsfX4rPz3yE41vbiwPOt0gIiKCiRMnZuncqKgoSpdOX5XklltuYenSpVZVLUsk\ntWIDGaWRTtoiXca2SHImMXrVaIYtHUZMQgwAj9V8jKERQ6lQsIJNNbxxZcuWJT4+nn379pEnj5lU\nNGHCBL777jsWL15sc+1MHvpG/yBu2bIl03P2799P+fLlSUpKwuFwf/9ZeuRCeJjlB5dz+5e303dB\nX2ISYmherjlrn13Lt+2/9cognsrpdPLpp59m+zpaa699sGtVvSWQ20DywumkLdL9uuBXus7syl1f\n3cWmE5soG16W2Z1n8/sTv3Nbidvsrl62KKV45ZVX+Oijj4iNjb3qOStWrOCOO+4gPDycevXqsXLl\nyrT3IiIieOONN2jUqBGhoaHs3bsXh8PBuHHjqFSpEvnz52fQoEHs2bOHBg0aEB4eTqdOnbh0ySzy\nFRMTQ9u2bSlatCgFCxbkvvvu48iR/yzSelXx8fFERkZSsGBBatSowd9//33Z+2XLluWPP/4AYPXq\n1dStW5ewsDCKFy/OK6+8AkCTJk0ACA8PJ1++fKxatcq1BsyEBHIhPMD8f+bT5ZcuTNowiaCAIN5o\n/AZbn99K28pt3ZIDV8p9rxtVt25dIiIi+Oijj/7zXnR0NG3atOGll14iOjqavn370qZNG86cSd+c\n7Ntvv2XChAnExcVRpkwZABYsWMD69ev566+/GD58ON26dWPq1KkcPHiQzZs3M3XqVMB8GujatSsH\nDx7k4MGDhISE0LNn1paIGjp0KPv27WPv3r3Mnz+fr7/++rL/Jxm/7t27N3369CE2Npa9e/fy8MMP\nA7Bs2TIAYmNjiYuLo379+i623vVJILeB5IXT+XtbnEs8R49fe9Dyu5acKnaKO2+6k409NvJWs7d8\nboEqpRTDhg1j9OjRnDp16rL35syZQ5UqVXjsscdwOBx06tSJqlWrMmvWrLSfjYyMpFq1ajgcDgID\nzVowr732GqGhoVSvXp2aNWvSqlUrypYtS/78+WnVqhXr15sRPQULFuTBBx8kd+7chIaGMnDgQJYs\nWZKlek+fPp3XX3+d8PBwbrrpJnr37n3NFElQUBC7d+/m1KlT5MmTJy1gW50KkkAuhE02n9hM3S/r\n8sXaLwh0BPJe8/dY1mUZVQtXdXtZWrvvlR01atSgbdu2vP/++5f1ZI8ePZrWy0518803c/Ro+jox\nGUeKpCpWLH174JCQkP8cnzt3DoALFy7QvXt3ypYtS1hYGE2bNiU2NjZLAfbo0aOXlX1lPTOaOHEi\nu3btolq1atSrV485c+Zken13kEBuA8kLp/PHttBaM3HdROpNqMfO0zupXqQ6a55dw51Jd5LL4fsD\nyYYOHcr48eMvy1GXKlWKAwcOXHbegQMHKFWqVNpxdlJMI0aMYNeuXaxevZrY2FiWLFmS5YemJUqU\n4ODBg2nHGb++UsWKFZkyZQr//vsv/fr1o0OHDsTHx1s+RFQCuRA5KCEpgciZkTwz+xkSkhLoUqcL\nq59ZTa1iteyuWo6pUKECHTt2vGwES6tWrdi1axdTp04lKSmJ77//nh07dtC2bdu0c7ISdDOek/Hr\nc+fOERISQlhYGNHR0QwdOvS6P5vRI488wnvvvUdMTAyHDx9m9OjR1yz/22+/5d9//wUgLCwMpRQO\nh4MiRYrgcDjYs+fam3ZkhwRyG/h7Xjgjf2qLY3HHaDq5Kd9s/IY8gXn4+oGvmdRuEnmD8gL+1RaD\nBg3iwoULaT3VQoUK8euvvzJixAgKFy7MRx99xK+//krBggXTfubKXu3VerlXPoRMPX7ppZeIj4+n\ncOHCNGzYkFatWmXpegCDBw/m5ptvply5crRs2ZInn3zymufOnz+fW265hXz58tGnTx+mTZtGcHAw\nefLk4fXXX6dRo0YUKFCA1atXZ6GVsk4WzRIiB6w5uoYHpj3Akbgj3Bx2M7M6z3JrL/xaCy4J75Gd\nRbOkR24Df8wLX4s/tMXMHTNp/FVjjsQd4a4yd7G629VTKf7QFsIaEsiFsND4teNp/0N7EpISeLrO\n0yx6chFF8xa1u1rCx0hqRQgLaK15e+nbDIoaBMDgpoMZ3HSwZaMXJLXi/SS1IoQH0Vrz0ryXGBQ1\nCIdyMK7NOIZEDPH7VQqnTJlC3bp1yZcvHyVLlqR169YsX76cIUOGEBgYSP78+cmfPz9VqlShV69e\nHD9+PO1no6KicDgc5MuXL+3Vrl07G/81nkUCuQ0kF5rO19rCqZ28MPcFRq0eRVBAENMfnk6Puj2y\n9LO+1hYZffzxx/Tp04c33niDkydPcujQIV544QVmzZqFUorOnTtz9uxZzpw5w4wZMzh+/Di33377\nZcG8VKlSxMXFpb1mzpxp47/Is0ggF8JNnNpJj197MG7NOIIDgpnZaSbtq7W3u1q2i42NZfDgwYwd\nO5YHHniAkJAQAgICaNOmDcOHD79sYk5AQADVq1fn+++/p0iRIowYMcLm2nsH359G5oH8abxwZnyl\nLZzaSbdZ3Zi0YRK5c+VmVqdZ3F3hbpeuYWVbqKHuS+vowa7l4leuXElCQgIPPvhgln/G4XDQrl07\n5s+f72r1/JL0yIXIptSc+KQNkwjJFcKcR+e4HMR92enTpylcuLDLGyqUKFGC6OjotOOjR49SoECB\ntNePP/7o7qp6LemR2yAqKspneqLZ5QttMWjxIEavHk1QQBCzOs+iWblmN3QdK9vC1V60OxUqVIhT\np07hdDpdCuZHjhyhUKFCacclS5bk0KFDVlTR60mPXIhs+GjFR7y97G0CVADTHppGi/It7K6Sx2nQ\noAHBwcHMmDHjqu9fbTSP0+lk9uzZNG7c2Orq+QTpkdvA23ug7uTNbfH1hq959fdXAZjUbhIPVst6\nDvhqvLktricsLIxhw4bxwgsvkCtXLu6++24CAwNZuHAhUVFRaXt4AiQlJbF7926GDBnCyZMn6du3\nr4019x7SIxfiBizcu5BnZj8DwKctP+XJ2k/aXCPP1rdvXz7++GPefvttihYtSpkyZRg7dmzaA9Dv\nv/+efPnyER4eTrt27ShSpAhr166lePHiadfw93H41yMzO23gC3lhd/HGtthycguNJjXi7MWzvNLg\nFT6850O3XDc7bSEzO72fzOwUIoccjTtK6+9ac/biWR6u/jDD7x5ud5WEkB65EFkVfymexl81Zu2x\ntTQs3ZCFTywkJDDE7moB0iP3BdIjF8JiWmuemf0Ma4+tpXyB8szsNNNjgrgQEsht4MtrarjKW9pi\nxMoRTNk8hbyBeZnZaSaF8xR2exne0hbC80ggFyIT8/6ZR7+F/QD4vwf/j1uK3mJzjYS4XKY5cqVU\nS2AkEABM0FoPv+L9wsC3QHHMuPSPtNaTr3IdyZELr7Mneg91x9clJiGGIU2HMDhisN1VuirJkXu/\n7OTIrxvIlVIBwE6gBXAE+BvorLXenuGcIUCw1npASlDfCRTTWiddcS0J5MKrJCQl0HBiQ9YfX88D\nVR/gp0d+wqE880OsjLH2DVY97KwH/KO13q+1vgRMA65czf0YkD/l6/zA6SuDuLic5ELTeXJb9P6t\nN+uPr6dCgQpMbjfZ8iCenbZIXQrWV16LFy+2vQ52vG5UZlP0SwEZV6k5DNS/4pzxwB9KqaNAPuCR\nG66NEB7i203f8uW6LwkOCGb6w9MJyx1md5WEuKbMAnlW/kQMBDZorSOUUhWA35VStbXWcVeeGBkZ\nSdmyZQEIDw+nTp06aTPZUnsj/nAcERHhUfWR48uPt/27jWdGPQNJMKrXKG4tcWuOlZ/Kk9rDjuPU\n73lKfXLyOCoqismTJwOkxcvMZJYjvxMYorVumXI8AHDqDA88lVJzgXe01stTjhcB/bTWa664luTI\nhceLvxRPvQn12HJyC4/XepxvHvhG8s/CVu7Ika8BKimlyiqlgoCOwKwrztmBeRiKUqoYUAXYe2NV\n9g9X9r78mae1Rb+F/dhycguVClZiXJtxORrEPa0t7CRt4Zrrpla01klKqZ7AfMzww4la6+1Kqe4p\n738BvAt8pZTaiPnD8JrWOvqaFxXCQ83ZNYfRq0cT6Ahk6kNTCQ0KtbtKQmSJrLUiBHD83HFqjavF\nvxf+5YMWH/Bqo1ftrpIQgKy1IkSWaK3pMrML/174l+blmvNyw5ftrpIQLpFAbgPJ/6XzhLb4cu2X\nzPtnHgVDCvL1A1/bNunHE9rCU0hbuEYCufBre8/s5eUFpgc+rs04SuUvZXONhHCd5MiF33JqJxGT\nI1h2cBkda3RkWodpdldJiP+QHLkQ1/HpX5+y7OAyiuUtxmetP7O7OkLcMAnkNpD8Xzq72mLHqR0M\nWDQAgPH3jadQnkK21CMjuS/SSVu4RgK58DtO7eSZWc9wMfkikXUiua/KfXZXSYhskRy58Dufrf6M\nnr/1pHhocbY9v40CIQXsrpIQ1yQ5ciGucDD2IP0X9Qfgs9afSRAXPkECuQ0k/5cuJ9tCa02PX3tw\nLvEcD1V7iPbV2udY2Vkh90U6aQvXSCAXfmPK5in89s9vhOcOZ0zrMXZXRwi3kRy58AvR8dFUGVOF\nUxdOMen+SXS5tYvdVRIiSyRHLkSKfr/349SFU0SUjSCyTqTd1RHCrSSQ20Dyf+lyoi3+PPgnE9ZP\nINARmONrjLtC7ot00haukUAufNql5Ev0+LUHAP3v6k/VwlVtrpEQ7ic5cuHThv85nP6L+lOhQAU2\nP7eZkMAQu6skhEskRy782oGYAwxdMhSAsW3GShAXPksCuQ0k/5fOyrbou6Av8UnxPFLjEe6pcI9l\n5biL3BfppC1cI4Fc+KQFexbw8/afyRuYlxH3jLC7OkJYSnLkwuckJidSa1wtdp7eyfvN36ffXf3s\nrpIQN0xy5MIvjfxrJDtP76RKoSr0adDH7uoIYTkJ5DaQ/F86d7fFkbNHGLZkGACjWo0iKCDIrde3\nktwX6aQtXCOBXPiU1xa+xvlL52lfrb1XPOAUwh0kRy58xopDK2g0qRHBAcHs6LmDsuFl7a6SENkm\nOXLhN5zaSe95vQF4teGrEsSFX5FAbgPJ/6VzV1t8s/Eb1hxdQ6l8peh/V3+3XDOnyX2RTtrCNRLI\nhdeLuxiXtpHy+y3eJ29QXptrJETOkhy58HoDFg7g/eXvc+dNd7L86eU4lPRPhO/ISo5cArnwavtj\n9lNlTBUSkxP5q+tf1L+pvt1VEsKt5GGnh5L8X7rstsWARQNITE7ksZqPeX0Ql/sinbSFaySQC6/1\n1+G/mLZlGrlz5ebd5u/aXR0hbCOpFeGVtNY0mtSIlYdXMvCugbzT/B27qySEJSS1InzW9G3TWXl4\nJUXzFvXa4YZCuIsEchtI/i/djbTFxaSL9F9ogvdb/3uLfMH53Fwre8h9kU7awjUSyIXXGbN6DPti\n9lGjSA2evvVpu6sjhO0yzZErpVoCI4EAYILWevhVzokAPgECgVNa64irnCM5cpFt0fHRVBhVgZiE\nGOY+OpdWlVrZXSUhLJWVHHmuTC4QAIwBWgBHgL+VUrO01tsznBMOfAbcq7U+rJQqnP2qC3F17y17\nj5iEGJqVa0bLii3tro4QHiGz1Eo94B+t9X6t9SVgGtDuinMeBX7SWh8G0Fqfcn81fYvk/9K50hb7\nY/YzavUoAD5o8QFKXbeT4nXkvkgnbeGazAJ5KeBQhuPDKd/LqBJQUCm1WCm1Rin1hDsrKESqN/54\nI23yz+0lb7e7OkJ4jOvmyJVSDwEttdbdUo4fB+prrXtlOGcMcBvQHMgDrATaaK13X3EtyZGLG7bu\n2Dpu//J2ggKC2NlzpyxTK/xGtnPkmLx46QzHpTG98owOYR5wxgPxSqmlQG1g9xXnERkZSdmyZQEI\nDw+nTp06REREAOkfpeRYjq923G10NzgKvR7tRdnwsrbXR47l2KrjqKgoJk+eDJAWLzOTWY88F7AT\n09s+CqwGOl/xsLMq5oHovUAwsAroqLXedsW1pEeeIioqKu1/oL/LSlv8vud37vn2HsKCw9jbey8F\nQwrmTOVymNwX6aQt0mW7R661TlJK9QTmY4YfTtRab1dKdU95/wut9Q6l1DxgE+AExl8ZxIW4UU7t\npP8iM/lnwF0DfDaIC5EdstaK8GjTtkyj80+dKZWvFLt77SYkMMTuKgmRo2StFeHVEpMTef2P1wEY\nGjFUgrgQ1yCB3AapDzb80fnzsHw5zJ4NP/8MgwZFMXMm7N4NycmXn/vl2i/Ze2YvVQtX5ak6T9lT\n4Rzkz/fFlaQtXJPZqBUhskVriIqC77+Hv/6CLVv+G7BThYRoqld30ratg4c6n+etpW8B8F7z98jl\nkFtViGuRHLmwxJkzMHkyfP457NqV/v2AAKhZE266CQIDweFI4siRS+zbF8yJExk+IDZ5C5oNolLu\nBmx7eTm5cvnWLE4hskr27BQ5LikJvvgC3ngDYmLM90qWhK5d4e674bbbIG+GTe4vXLjA1q1nCAsr\nRWwsbNwI0+ecZEG1ihAcB18toVbYXXz8sYPmze35NwlhJ3nY6aF8Nf+3bBncfjv07GmCeEQEzJgB\nBw7AsGHQuPHlQRxg6dKlaV+HhUGTJlCq43AIjqOyakWxhLvYtMlBixbQqpWTgwdz9t+Uk3z1vrgR\n0haukUAusi0pCfr1M0F40ya4+WbzIPOPP+CBByCXC+nto+cP8t3uzwD44N53WbDAQd++mrx5NfPm\nOahd28mMGRb9Q4TwUpJaEdly7Bh06gRLl5r898CBMGAAhGRxpGDG1ArAwFVd+WnvJNre3JkRDaek\nnXf6NAwY4GTJEtP3eP55zYgRity53f5PEsKjSI5cWGr5cnjoIThxAooXNyNTmjRx7RoZA/me2O20\n/e0WHDiY22Y7N+ereNm5WsM332g++ACSkhR33OFk7lwHhWUFfOHDJEfuoXwh/zd3LrRoYYJ4RASs\nX+96EIfLc+QjN72BUzvpUOGZ/wRxAKXgqacU33+vKFHCyd9/O2jY0Mn+/Tf+7/AkvnBfuIu0hWsk\nkAuXTZ0K7dpBQgI88wz8/rvpkWfHptN/s+Dwz+QOCOGFWwZd99xbboEffnBQpYqT3bsd3Hmnk40b\ns1e+EN5MUivCJZ9/Ds8/b9Icr70G779veso3KjW18tLaSFacWEi3av14pc77WfrZuDh47jnTM8+X\nT7N4seJ22W9C+BhJrQi3+uoreO45E8Tfew+GD89eEE/196nlrDixkHyBYXSr3i/LP5cvH0ya5KBF\nCydxcYp77nGydWv26yOEt5FAbgNvzP/98otJowCMHAn9+7vnukuWLGHszvcAeKbaa4QFFXDp54OC\n4JNPHDRp4iQ62gT1PXvcU7ec5o33hVWkLVwjgVxkKirKDDF0OmHQIOjd233XXnlkJVtjN1A4dzGe\nrHJjFw4KgtGjHdSr5+T4cQfNmjk5fOU+VkL4MMmRi+vasMGMRomLM7nxMWPck04BSHYmU3NsTbaf\n3s6bt4/m8co9s3W9c+egSxcnmzY5qFXLyfLlDkJD3VNXIewiOXKRLSdOwP33myDesSOMHu2+IA7w\n3ebv2H56OyVCSvNIhWezfb3QUBg/3kHp0iaYP/64E6fTDRUVwsNJILeBN+T/EhPNZJ9Dh6BBA/j6\na3C48W5JTE5kcNRg2AfPVnqZoIAgt1w3PBy+/NJBvnxOZs50MHCg93wK9Ib7IqdIW7hGArn4D61N\nGmX5crPc7M8/Q3Cwe8v4cu2X7I/ZT5mwMrQs9aBbr12+PIwa5SAgQDN8uOLrr916eSE8juTIxX+M\nHg0vvgi5c8Off+L2sdnnE89TYVQFTpw/wdR2U6lwqXHaWivuNHUqDBkCwcGalSsVt97q9iKEsJzk\nyIXLVq2Cvn3N11995f4gDjBq1ShOnD9BvVL1uK/Sfe4vIEXnzvDww04uXlR06OAkNtayooSwlQRy\nG3hq/i862jzUTEqCl14yQw7d7Uz8GT5Y8QEA7zZ7l2XLlrm/kAzefNNB1apO9u51EBnpxJM/FHrq\nfWEHaQvXSCAXgMmLd+liNoG44w4za9MKw5cPJyYhhmblmtG8vPVb/gQHmzHmefNqfvnFwSefeHAk\nF+IGSY5cAGa2Zp8+ZtTHunVQrpz7yzgWd4wKoyoQnxTPqmdWUa9Uvf+sR26VBQugVy/IlUuzfLmi\nXj1LixPCbSRHLrJk7VqzABaYvLgVQRzgraVvEZ8Uz4NVH6ReqZyNpPfcA088oUlKUjz6qJNz53K0\neCEsJYHcBp6U/4uPh8cfh0uXTI/1gQesKWdP9B7GrxuPQzl4u9nbad/PuB651V59VVGpkpM9exz0\n6eN5M4U86b6wm7SFaySQ+7l+/WDHDqhWzbq8OMCgqEEkOZN4otYTVC9S3bqCriM4GD7+2EFQkGbC\nBAe//GJLNYRwO8mR+7EFC+Dee83myKtWwW23WVPOphObqPN5HXI5crGr1y7KhpdNey+ncuQZff21\n5t13FQULOtmyxUGJEjlWtBAukxy5uKboaDNKBWDYMOuCOMDrf7yORtOjbo/LgrhdnnhC0bChWfa2\na1fPHpIoRFZIILeBJ+T/evWCo0ehUaP0B51W+PPgn/y661fyBubljSZv/Of9nMyRp3I44P33zXos\nv/3m4JtvcrwKV+UJ94WnkLZwjQRyPzRzJkyZAnnymMWwAgKsKUdrTf+FZgeKlxu8TNG8Ra0p6AYU\nKwZvvGFu/xdf1Bw5YnOFhMgGyZH7mehoqFEDjh+HTz81a6pYZfbO2dw/7X4K5ynMnhf3kD84/3/O\nsSNHnkpr6N7dyZIlDlq3dvLrrw63LtMrhDtIjlz8R58+JojfdRf0zN4+DteV7Exm4B8DAXi98etX\nDeJ2UwreesukWObO9ZwUixCukkBuA7vyf3PmwDffmFUNJ0507/riV5qyeQpbTm6hTFgZnqv73DXP\nsyNHntGVKZZjx+yri+SF00lbuEYCuZ84exa6dzdfv/UWVK5sXVkXky7y5uI3ARgWMYzgXG5ezNzN\n2rWDJk2cnD2r6NXL8yYKCZEZyZH7iZ494bPPzIJYK1da94ATYORfI+kzvw81itRgY4+NBDiuXZid\nOfKMjh6FNm00Fy4ofvnFBHchPIFbcuRKqZZKqR1Kqd1KqX7XOe8OpVSSUqr9jVRWWGflShg71kz8\nmTDB2iAemxDL20vNFPz3W7x/3SDuSUqWhL59ze/Kc8/J2uXCu1w3kCulAoAxQEugOtBZKVXtGucN\nB+YB8tw/EzmZ/0tMhG7dzAiNV16BWrWsLe/DFR9yOv40jcs0pk2lNpmeb3eOPKNHH4WaNZ0cO+ag\nf/+cT7FIXjidtIVrMuuR1wP+0Vrv11pfAqYBV/vQ2Qv4EfjXzfUT2fThh7B1K1SoAIMGWVvWsbhj\nfLzyYwCGtxiO8rKxfAEB8O67Zq/Pzz93sHy53TUSImsyC+SlgEMZjg+nfC+NUqoUJriPS/mWJMIz\nERERkSPl7NplHmwCfPEFhIRYW97QJUPTlqltULpBln6mSZMm1lbKRZUrw7PPmq+7dXOSmJhzZefU\nfeENpC1ck1kgz0pQHgn0T3mSqZDUikfQGp5/Hi5ehKeeguYWb8az89ROJqybgEM5eLf5u9YWZrHn\nnlOULu1k+3YHI0ZIv0R4vlyZvH8EKJ3huDSmV57R7cC0lI/RhYFWSqlLWutZV14sMjKSsmXLAhAe\nHk6dOnUt23EIAAAeH0lEQVTS/vKm5sT84Thj/s+q8t54I4pFi6BQoQg++sj6f1/XUV1JPphMt4e6\nUbVw1Sz/fEJCAoUK1WTVKnNcv755387j4GDo3HkpH3wAb73VlE6d4MAB97bX1Y43bNjASy+9ZNn1\nvel45MiRfh0fJk+eDJAWLzNz3eGHSqlcwE6gOXAUWA101lpvv8b5XwGztdY/X+U9GX6YIioqKu1/\noBXOnIGqVeHkSZg0KX2VQ6ssO7CMJpObkDcwL7t77aZEvqyvCztv3jwKFapp+/DDq+nTx8z4bNnS\n/NfqlL/V94U3kbZIl+3hh1rrJKAnMB/YBnyvtd6ulOqulOruvqr6F6tv0IEDTRBv3NikVayktebV\n318F4JWGr7gUxMHzcuQZDRzoIDRUM2+egx9/tL48CVzppC1cIxOCfMxff0HDhmYExoYNZoEsK03f\nOp1HfnyEYnmL8c+L/xAaFOrSz3vKhKBrmToVhgyB4sWd7NzpIL/nLRkjfJwsmuWhMubI3SkpCZ57\nLn3MuNVBPDE5kQGLBgAwNGKoy0EcPGsc+dV07GjGlh8/7mDwYGs7IlbdF95I2sI1Esh9yNixphd+\n883w5pvWlzfu73HsObOHaoWr0fW2rtYXaAOHA4YNc+BwaEaPho0b7a6REP8lqRUfcewYVKkCcXEw\naxbcd5+15UXHR1NxVEXOJJxhVqdZ3Fflxgr09NRKqrff1vzf/ynuvNPJ8uUOS1eOFCIjSa34kZdf\nNkH8/vutD+IAby15izMJZ2hWrhltK7e1vkCb9e6tKFzYyV9/OfjqK7trI8TlJJDbwN35v0WLzEO5\nkBCz64/Vdp3exZi/x6BQjLhnRLam4nt6jjxVvnwwYID5dXn1VSenTrm/DMkLp5O2cI0Eci938SK8\n8IL5+s03IYvzB7Kl38J+JDmT6FKnC3WK17G+QA/Rpg3Ur+/kzBl7FtUS4lokR+7l3n0XXn/d5Mc3\nbYKgIGvLW7xvMc2+aXZDk3+uxlty5Kn27oX77tMkJSlWrIAGWVtSRogbJjlyH7dvX/qiWGPHWh/E\nk53JvLzgZQD639U/20HcG5UvD11TBuj06OEkKcne+ggBEsht4a78X+/ekJBg1tFu1swtl7yurzZ8\nxfrj67kp/030bdDXLdf0lhx5Rs89pyhRwsmmTQ7GjnXfdSUvnE7awjUSyL3UrFkwezbkzw8ffWR9\nebEJsQxcNBCAD1p8QJ7APNYX6qFCQuDNN82vzuuv27thsxAgOXKvdOECVK8OBw6YUSovvmh9mS/P\nf5mP//qYRqUbsazLMrdtGuFtOfKMund3EhXloFMnJ1OnSp9IWENy5D7qnXdMEK9Tx6w5brWdp3Yy\navUoFIpRrUZ53c4/VnnzTQfBwZpp0xwsWmR3bYQ/k0Bug+zk/7ZvN9u3QfqGylbrM78PSc4kut7a\nldtK3ObWa3tjjjzVTTeZtW0Ann/eycWL2bue5IXTSVu4RgK5F0nd9efSJbOhck4MfZuzaw6//fMb\n+YPz807zd6wv0Mt07aooV87Jrl0OPvxQUofCHpIj9yLffgtPPAGFC8OOHVCokLXlJSQlUGNsDfae\n2cuIe0a4baRKRt6cI0+1ciVERkLu3Jpt2xTlytldI+FLJEfuQ86cMeupAHzwgfVBHODD5R+y98xe\nahSpQa96vawv0Es1aABt2jhJSFD07OlE+isip0kgt8GN5P9ef93s+nPXXdbv+gOw78w+3v3TbKL8\nWevPCAwItKQcb86RZzRggIPQULMl3IwZN3YNyQunk7ZwjQRyL7BqFXz+udn1Z9w4cmQJ1T7z+5CQ\nlMCjNR+ladmm1hfo5YoUgb59zf+YXr2cxMXZXCHhVyRH7uGSkqBuXbOhwWuvwfDh1pc5d/dc2kxp\nQ76gfOzouYOS+UpaVpYv5MhTJSfDI4842bLFwYsvaj79VIZpiuyTHLkP+PRTE8TLloVBg6wvL/5S\nPL1+M/nwoRFDLQ3iviYgAN5+2+wmNGYMrFtnd42Ev5BAboOs5v8OHEgP3p99BnnzWlenVG8vfZu9\nZ/ZSs2hNetbraXl5vpIjT1WtmhlZ5HQqunVzkpyc9Z+VvHA6aQvXSCD3UFpDr15mOv7DD0Pr1taX\nufXkVj5Y8QEKxRdtv7DsAaev691bUayYk3XrHIwZY3dthD+QHLmH+ukn6NDBLIq1fTuUtDjD4dRO\nmnzVhOWHlvNc3ecY28aNy/pdhy/lyDNauNBs+JE3rxlbXqaM3TUS3kpy5F4qJgZ6pmQ13n3X+iAO\nMHHdRJYfWk7x0OK82/xd6wv0cS1awN13Ozl/XvHcczK2XFhLArkNMsv/vfYaHD8ODRumr+VhpRPn\nTvDawtcAGHnvSMJzh1tfaApfy5FnNGiQg9BQzdy5Dn74IfPzJS+cTtrCNRLIPUxUFIwfb3b7GT8+\nZ8aMvzD3BWISYmhZsSWP1HjE+gL9RNGi8Npr5hNxr15OoqNtrpDwWZIj9yDx8VC7NuzeDUOH5sxw\nw5+2/USH6R0IDQpl6/NbKROWs8lcX82Rp3I64bHHzIPPyEgnX30lfSfhGsmRe5m33jJBvEYN6N/f\n+vKi46N5Ye4LAAxvMTzHg7g/cDjM2PKgIM3kyQ4WLLC7RsIXSSC3wdXyf2vWmMWwlEpPrVitz/w+\nnDh/giY3N6FH3R7WF3gVvpwjT1WhQvrD62eeufb0fckLp5O2cI0Ecg9w8SJ06WKmeL/0Us6sMz53\n91y+2fgNuXPlZsJ9E3AouRWs1LWronp1J4cOOejXz2l3dYSPkRy5Bxg0yKRVKlY00/HzWLyvcXR8\nNDXH1eRo3FE+vPtDXmn4irUFXoev58gz2rED2rfXJCcroqKgqaxFJrJAcuReYP16eO89k1KZNMn6\nIA7Qc25PjsYdpWHphvS5s4/1BQoAqlaFHikZrKefdnL+vL31Eb5DArkNUvN/iYkmpZKUZHKojRtb\nX/YPW39g6pap5A3MyzcPfEOAI8D6Qq/DH3LkGfXooahUycnevQ5effXyT6iSF04nbeEaCeQ2GjLE\npFLKlze9cqsdizvGc3PMDKOP7vmICgUrWF+ouExQEHz4oYNcuTTjxil+/93uGglfIDlym6xYYXrg\nWsPSpWbnHytprWk7tS1zd8/l3gr38ttjv6GU/etl+1OOPKNx4zQjRypKljTrlxcoYHeNhKeSHLmH\nOncudalTMx3f6iAOMGb1GObunkuB3AWYeP9Ejwji/qxbN0WtWk6OHnXQs6eMYhHZk6VArpRqqZTa\noZTarZTqd5X3H1NKbVRKbVJKLVdK1XJ/VX1H585R7N0LtWqZGZxW23h8I6/8bkamTLh/AqXye07v\n199y5Kly5YIPPnAQHKyZMsXBjz9KXjgjaQvXZBrIlVIBwBigJVAd6KyUqnbFaXuBJlrrWsBbwJfu\nrqiv+PVX8woKgm+/heBga8s7n3ieTj91IjE5kR6396B9tfbWFiiyrFw584kM4NlnzebaQtyIrPTI\n6wH/aK33a60vAdOAdhlP0Fqv1FrHphyuAm5ybzV9w9GjZpQKRPDOO1CzpvVlvjTvJXac2kH1ItUZ\nce8I6wt0UZMmTeyugq0efRQaNTrPmTMwdmyESzsK+bKIiAi7q+BVshLISwGHMhwfTvnetXQF5man\nUr4oORkefxxOnYJ77oG+fa0vc9qWaUxYP4HggGCmPTSNPIE5MEhduEQpGDToBMWLw5Il8P77dtdI\neKNcWTgny0NNlFL/A54GGl3t/cjISMqWLQtAeHg4derUSfvLm5oT89XjHj2iWLwYihaNoHv3KFJT\nw1aVN/mXyTz363NQGj659xNObz9N1PYoj2mP1OOEhAQKFarJqlXmuH59876/HNer15QCBZLp2zeK\n117bwODBL9G8OSQkmPft/v9j1/HIkSP9Kj5kPI6KimLy5MkAafEyM5kOP1RK3QkM0Vq3TDkeADi1\n1sOvOK8W8DPQUmv9z1Wu47fDD1euNEMNk5Pht98gd+70gGqFc4nnqDe+HttPbefRmo/y7YPfeuwo\nlXnz5lGoUE2/G36YSmvN2bP/cMcdlejYMYoffoigbFlYtw6/HpIYFWXt74g3ycrww6wE8lzATqA5\ncBRYDXTWWm/PcE4Z4A/gca31X9e4jl8G8tOn4bbb4OBBePll+Ogja8vTWvPoz48ybcs0qhepzupn\nVpM3KK+1hWaDv44jT5UxkF+8CI0awdq1cP/98MsvJvUi/JtbxpFrrZOAnsB8YBvwvdZ6u1Kqu1Kq\ne8ppg4ACwDil1Hql1Ops1t0nOJ0mL37wINSrZ/bftNro1aOZtmUaoUGh/PTITx4dxMXlgoNh+nQI\nD4dZs6z/oy98R5bGkWutf9NaV9FaV9Rav5fyvS+01l+kfP2M1rqQ1vrWlFc9KyvtLd5+G+bNg0KF\nzC9oUMoa41aNkV20dxF955unqBPvn0jVwlUtKced/HUc+dVERUVRrhx8/bU5HjAAli2zt052kXHk\nrpGZnRZZsMCspaIUfPcdlLF485090Xt4ePrDJOtk+jXqJ3tverH774dXXzXPVDp2NBtxC3E9staK\nBfbvh7p1TX58yBAYPNja8s5ePMudE+5k+6nttK3cll86/mL7qoZZJTny9Bx5RpcuQfPmpkfesCH8\n8Yf1k8eEZ5K1Vmxw7hy0a2eCeMuW8Oab1paX7Ezm0Z8eZfup7VQvUp3v2n/nNUFcXFtgoEnH3XST\nWWCtZ0+zwJoQVyOB3I2cTrMY1qZNULkyTJ1qNt+9krvyf1prXvztRebsnkPBkILM6jSL/MH53XLt\nnCI58nRX3hfFisGMGZA7N0yYAOPG2VMvO0iO3DUSyN1o8GAzZCwszIw6CA+3trwPV3zI2DVjCQoI\nYkbHGbK+uA+qWxcmTjRf9+4NixfbWx/hmSRH7iZTp5p1MxwOM+nnnnusLW/K5ik89vNjAPzQ4Qce\nrvGwtQVaRHLkV8+RX6lfP/jgA9M5WLECql25bJ3wWZIjzyFRURAZab4eMcL6IL5o7yIifzEFfnzP\nx14bxEXWvfsuPPAAxMRA69YykkVcTgJ5Nm3ebH7BEhOhVy/z8Tcz2cn/LT+4nPun3c8l5yV61+9N\nnwbevXmy5MjTXe++CAgww1jr1zejotq2xac3b5YcuWskkGfD4cOmdxQbCw89BJ98Yu2U6nXH1tF6\nSmsuXLrAk7Wf5ON7P7auMOFx8uQxz17KlzfT+Dt2NMMUhZAc+Q06fRoiImDLFrNV2++/m9EFVtl6\ncitNJzfldPxpOlTvwNSHppLLkZXFKz2b5MizliPPaNcuaNAAoqPNc5lvvjE9duGbJEdukZgYkwff\nsgWqV4eZM60N4ltObqH5N805HX+a1pVa813773wiiIsbU7myWfohNBSmTIHnn5cx5v5OArmLzp0z\n6ZR166BiRVi4EAoWdO0aruT/1h9bT8TkCE6cP0GL8i348eEfCQoIcq1ADyY58nSu3Bd33GG2DMyd\nG7780mwZ50vBXHLkrpFA7oILF+C++8z64mXKwKJFUKKEdeWtOryKZt80S+uJz+48m5DAEOsKFF6l\naVP46SezkfNHH8Hrr/tWMBdZJznyLDp71gTxpUtN8F661PTIrfLHvj94YNoDxCXG8WDVB5nWYZpP\n9cRTSY7c9Rz5lX78ETp1Mots9e1rgrqsY+47JEfuJtHR0KKFCd6lSpkFjKwM4lM2T6Hlty2JS4yj\n8y2d+b7D9z4ZxIV7dOhg1mUJDISPPzbrsjiddtdK5CQJ5Jk4cQL+9z/4+28oV86sRlc1m8t8Xyv/\np7Xmw+Uf8tjPj3HJeYk+d/bh2/bfEhgQmL0CPZjkyNNlJy/84INmXZbgYBg7Frp29e6hiZIjd40E\n8uvYtcssIbppE1SpYnrk5cpZU1ZiciLPz3me1xa+BsCIe0bw8b0f41Dyv0hkTZs2MHs2hITA5Mlm\nXfNz5+yulcgJkiO/hmXLzIzN6Giz5+Zvv0HRotaUdfL8SR6e/jBLDywlKCCIrx/4mk63dLKmMA8j\nOfLs58ivtGqVmfl56pS5d+fMgeLF3XZ5kcMkR36Dpk0zOfHoaPMLsWSJdUF83bF11P2yLksPLKVk\nvpIs67LMb4K4sEb9+mZhrQoVzDDZBg3Mp0rhuySQZ5CUZFaZ69zZrJ3ywgtmWdrQUPeWExUVhdaa\ncX+Po+HEhhw6e4g7b7qTNd3WUK+Uf213KjnydO7MC1eqZIJ5vXpmbZYGDeCHH9x2ectJjtw1EshT\nnDwJ995rlgoNCICRI2H0aGumPsddjKPD9A48P/d5LiZf5NnbniXqqShK5LNwULrwO0WLmpU5n3jC\nzIHo2BH69zfDFIVvkRw5sHy5GYd7+LDZleWHH6BJE2vKitofReQvkRyIPUD+4Px82fZLOt7S0ZrC\nvIDkyN2fI/9vGaZT0revCeJNm8L//R+ULm1ZkcKNJEeeiUuXzJ6aTZqYIN6ggVlVzoogfi7xHD3n\n9uR/X/+PA7EHuKPkHazvvt6vg7jIGUrBiy+a5SSKFzfPfGrXhp9/trtmwl38NpDv2gWNGsHbb5se\nS//+5mNoKQs6hgv3LqTWuFp89vdn5HLk4qmwp/jz6T8pX6C8+wvzMpIjT2d1XjgiAjZuNMMUz5wx\nSy8/84xZBM7TSI7cNX4XyBMT4Z13oFYtM8mnTBkTwN97D4LcPHnyUOwhHpn+CHf/393si9lH7WK1\n+bvb30TWiZSZmsIWRYuaseajRpnJQxMnmhU8Z8ywu2YiO/wqR758OXTvDlu3muOnnjIPNd29SXL8\npXhG/jWSt5e9zYVLFwjJFcIbTd7g1Yav+vQszRshOXLrc+TXsm2b6ZGvXGmOH3rITPEvUybHqyKu\nQ3LkKQ4cMAvw33WXCeIVK5qVCydPdm8QT3ImMX7teCqOrsjAPwZy4dIFOlTvwI6eOxjYeKAEceFR\nqlc3E99GjYK8ec1KilWrwuDBvr2NnC/y6UAeE2Ny31WqmF3ug4PNUp+bNkGzZu4rJ8mZxHebvqPG\n2Bo8++uzHI07yq3Fb2XhEwuZ/vB0yoRd3sWR/F86yZGns+O+CAgwe81u2waPPALx8TBsmNm8YuJE\n+9Zrkd8R1/hkII+OhkGD4OabYfhwuHjRTPLZscM83Axx05LeCUkJjPt7HJVHV+bxGY+z6/QuKhSo\nwNSHprLm2TU0L9/cPQUJYbEyZeD7700P/fbb4ehRk3apWtV8ck1KsruG4np8Kkd+8CCMGQPjxqUv\nFtSsmXmQWc+NEyYPxh7k8zWfM37deE5dOAVApYKV6NeoH0/WflJSKC6QHLl9OfJrcTrNMhVDh5rR\nXWAWi+vdG55+GvLls7d+/iYrOXKvD+Rap+f5ZsxIX4f53nvNGPFGjdxTTmJyIvP+mcek9ZOYvWs2\nTm0KurX4rQy4awDtq7UnwCE74LpKArnnBfJUSUkmJTlsGPzzj/le/vymp969u0m/COv5dCA/dMjM\nTps8GXbvNt/Llcvk+V56yexpmF1O7WTV4VVM3TKVqVumpvW+Ax2BPFzjYV644wUa3NQA5eJ2LFFR\nUURERGS/gj5g3rx5FCpUUwL5HZU89r5IToZZs+CTT0ynKVXjxqaH3qGDNesReWJb2CErgdyrtmI/\ncsTMRvvpJ7M2eOrfheLFoVs36NEDSpbMXhkXky6y7OAyftnxCzN2zOBo3NG092oUqcFTtZ/iidpP\nUDxU1gUV/iEgwGxc8eCDsGYNfPaZWcZi2TLzev55syH5I4+YyUZ589pdY//j0T3ypCRYvRrmz4d5\n88zXqYKDzXrhTz0Fd99teuM3QmvN9lPbWbxvMfP2zOOPfX9w4dKFtPfLhJWhfdX2PFH7CW4tfqvL\nvW9xfZJa8dzUyvXExZlg/tVXZn5Gqty5zVourVtDq1ZmqK/8ymSP16VWEhPNWiepf+n//PPy6cMh\nIdCypZm40LYthIW5Xo/ziedZd2wdq4+sZsXhFSw9sDQtZZKqZtGatKnUhoeqP8TtJW6X4G0hCeTe\nGcgzOnTIfEr+4Yf0yUWpSpc2axc1bWpSMZUrg8Mnx8pZxy2BXCnVEhgJBAATtNbDr3LOKKAVcAGI\n1Fqvv8o5lwXys2fN2NUtW0zwXrPGjO9OTLz85ypWNMH73nvN3plZ/djm1E4Oxh5k+7/b2XxyM5tO\nbGLjiY1s+3db2oPKVCVCS9Dk5ibcU+Ee7q1wL6XyWxtUJP+XTnLknp8jd8WJE+YT9G+/wYIFZihw\nRvnzQ9265lW7NtSsaeZ5XLk8hi+0hbtkO0eulAoAxgAtgCPA30qpWVrr7RnOaQ1U1FpXUkrVB8YB\nd17tei+8YB5M7txphgpeTbVq5i/3XXeZ/5Yte/XztNacvXiWw2cPc+jsIfad2ce+mH3sPbOX3dG7\n2X16N/FJ8f/5uQAVwK3Fb+WOkndQ/6b6NL25KeULlM/RXveGDRvkJk2xadMm/ve/mnZXwyP4wn1R\nrBg8+aR5OZ2mo7ZkiXmtWmVWGf3jD/NKFRhodjOqXNm8KlaEdes2UKJEBGXKuG/ehy/LLLNcD/hH\na70fQCk1DWgHbM9wzv3A1wBa61VKqXClVDGt9YkrLzZ2bPrXwcFmskGNGnDrrZpbbrtAhWpxOANj\nOZNwhpiEGFbGneHX1ac5deEUpy6c4sT5E5w4d4Lj545z7NwxziVef2fZYnmLUaVwFWoVrUWtYrWo\nWawmtYrVIk9gnsxbxkIxnrjcnE1iY2PtroLH8LX7wuEwi9PVqmVmjwIcO2YWq1u7FjZvNp/C9+41\nk/V27Mj40zF8+aX5qlAhM4ihVCkzsKFoUShSxLwKFkx/hYWZV548/peXzyyQlwIOZTg+DNTPwjk3\nAf8J5JXeak5A7nhUYDxJjvP8e+k8cxLPM+1CHM5lTlh25U9cX57APJTOX5pS+UtRLrwc5cLLUb5A\neSoUrEDlQpUJz+3m1bCEENlSogTcf795pTp/3oxT37XLfFrfuxcWLzbB+NAhOH3avDZvzloZAQFm\n0lJoaPorJMQE+JAQ8woONg9mg4PNKyjIvAIDzcCJwEDzCggwx7lymT9MAQGXvxyOy19K/fe/13vB\n5V+nyvheVmQWyLP6JPTK4q76c7uT/4BrLMYTkiuEfMH5yB+cnwK5C1AgpADhucMpHFKYwnkKUyhP\nIYrmLUrx0OIUy1uM4qHFCc8d7pUPIvfv3293FTzGgQMHcDoTiI09mvnJPkhrTWCguYf99b7Im9fk\ny2vXTv9eZOR+Jk826Zl//zVDj48cMVsypr7+/desqx4dbV6xsebZW3y8GSThYx9wruu6DzuVUncC\nQ7TWLVOOBwDOjA88lVKfA1Fa62kpxzuAplemVpRSnrnPmxBCeLjsTghaA1RSSpUFjgIdgc5XnDML\n6AlMSwn8MVfLj2dWESGEEDfmuoFca52klOoJzMcMP5yotd6ulOqe8v4XWuu5SqnWSql/MImTLpbX\nWgghRJocmxAkhBDCGjk6x0op9aFSartSaqNS6mel1A3MzfQNSqmHlVJblVLJSqnb7K5PTlNKtVRK\n7VBK7VZK9bO7PnZSSk1SSp1QSmVxXIZvUkqVVkotTvm92KKUetHuOtlFKZVbKbVKKbVBKbVNKfXe\n9c7P6cmyC4AaWuvawC5gQA6X70k2Aw8CfrdFToaJZi2B6kBnpVQ1e2tlq68wbeHvLgF9tNY1MJMK\nX/DX+0JrnQD8T2tdB6gF/E8pdde1zs/RQK61/l3rtPnxqzDjzf2S1nqH1nqX3fWwSdpEM631JSB1\noplf0lovA87YXQ+7aa2Pa603pHx9DjPxMJvrmXovrXXq6n1BmGeU0dc6187la54G5tpYvrDP1SaR\n+ediK+KqUkbK3Yrp8PklpZRDKbUBM7lysdZ627XOdft65Eqp34GrLdY9UGs9O+Wc14FErfUUd5fv\nSbLSFn5KnrCLa1JKhQI/Ar1TeuZ+KSV7USflWeJ8pVSE1jrqaue6PZBrre++3vtKqUigNeDzOxNn\n1hZ+7AhQOsNxaUyvXPg5pVQg8BPwrdb6F7vr4wm01rFKqTlAXSDqaufk9KiVlsCrQLuUZL4w/G2y\nVNpEM6VUEGai2Syb6yRspsx6GxOBbVrrkXbXx05KqcJKqfCUr0OAu4H/LA+eKqdz5KOBUOB3pdR6\npdTYzH7AVymlHlRKHcI8nZ+jlPrN7jrlFK11EmY28HxgG/B9xqWR/Y1SaiqwAqislDqklPLXSXWN\ngMcxIzTWp7z8dTRPCeCPlBz5KmC21nrRtU6WCUFCCOHlZNMlIYTwchLIhRDCy0kgF0IILyeBXAgh\nvJwEciGE8HISyIUQwstJIBdCCC8ngVwIIbzc/wNayl9cgSlPlgAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x = np.linspace(-2, 3, 100)\n", "x_min=0.0; x_max=1.0; mu0=0.0; sigma0=0.7\n", "plt.plot(x, norm.pdf(x, mu0, sigma0), label=\"Normal dist\")\n", "plt.fill_between(x, 0, norm.pdf(x, mu0, sigma0), where=(x_max>x)&(x>x_min), alpha=0.15)\n", "plt.plot(x, norm.cdf(x, mu0, sigma0), label=\"CDF\")\n", "plt.legend();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We want to look at the case where $f(x)$ is only defined on some interval $[x_\\text{min}, x_\\text{max}]$, that's hinted by the shaded area. That area is $A=\\text{cdf}(x_\\text{max})-\\text{cdf}(x_\\text{min})$. So the total integral is divided into $A$ + the area left of it ($\\text{cdf}(x_\\text{min})$) + the area right of it (not important). Since we know all areas, we'll solve $\\text{cdf}(x)$ for $x$:\n", "\n", "$$x=\\sigma\\sqrt 2\\cdot\\text{erf}^{-1}(2\\cdot\\text{cdf}(x)-1)+\\mu$$\n", "\n", "To generate our desired range, we'll evaluate that equation for different values of `cdf(x)`. We have to start at $\\text{cdf}(x_\\text{min})$ and step through $N-1$ parts of $A$. So we divide $A$ into $A_0=\\frac{A}{N-1}$ and can write down the result. The cases `i=0` and `i=N-1` are identical to `x_min` and `x_max`.\n", "\n", "Technical note: SciPy has the inverse error function and the pdf included. In C++ we need Boost." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX8AAAEACAYAAABbMHZzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XuQVeWZ7/HvA0gRacbWMpEgl47xMqhVtBMPMSraSSYe\nBBTNjApViYBzPNGoSdWEI3hBZEwKsMeEiplcJnps9GQ0OZrkgOIok6FLUULGCLngBQk0SFTUcBEI\nSZB+zh+9u23afVl7vfvy7t2/T9Uu9lrrXe/69Wr226ufXvvd5u6IiEj/MqDaAUREpPI0+IuI9EMa\n/EVE+iEN/iIi/ZAGfxGRfkiDv4hIPxQ8+JvZRDN7ycxeMbM5Odp8M7P9V2Z2RugxRUQkTNDgb2YD\ngW8BE4FTgelmNrZPm0nAie5+EvA/ge+EHFNERMKFXvmPBza5e4e7HwQeAqb2aXMxsBTA3dcCjWZ2\nXOBxRUQkQOjgfzzwaq/l7Zl1hdqMDDyuiIgECB38k84NYSn3ExGRMhgUuP/vgVG9lkfRdWWfr83I\nzLrDmJl+IIiIpODufS+wCwq98n8OOMnMmsxsMHAFsKxPm2XAlQBmdhaw2913ZOvM3VM/5s+fH7R/\nuR8x54s5m/IpX7UfsedLK+jK393fNbPrgSeAgcC97v6imX0hs/177r7CzCaZ2SZgPzAr5Ji5dHR0\nlKPbkok5X8zZQPlCKV+Y2POlFVr2wd0fBx7vs+57fZavDz2OiIiUTt28w3fmzJnVjpBXzPlizgbK\nF0r5wsSeLy0LqRmVkpl5LFlERGqFmeFV+INvNNrb26sdIa+Y88WcDZQvVAz5zEyPEjxKKbjmLyKS\nhH6zD1PqwV9lHxEpu0xpotoxalquc9jvyz4iIpJc3Qz+MdQ184k5X8zZQPlCxZ5PqqNuBn8REUlO\nNf8y6+yE1lYwg9mzYUAVftwWmyFJ+3xtcm0rZn3fdVB8m9DltDmKWZf2XIZ8f9PuE6Lea/47d+7k\nH/7hH1i5ciXHHnssCxcuZPr06Tnbt7S0sHbtWgYN6rrnZuTIkbz44ot5j1Hqmn/V56XoNT+F16NF\ni9yPOMJ98OCu57WQIUn7fG1ybStmfd91adqELqfNUcy6tOey2Dal2CdEvb6+u02bNs2nTZvm+/fv\n99WrV/tRRx3lGzZsyNm+paXF77333qKOkescZtYXPebWTdkn1rqmWdejs7OdEt+pVXSG7ud99T13\nhdoXapNrWzHre6/bvLm9YJvu56VcTnKM3vny7ZfmvBTalrRNttdGkn4robMTFi+GO+/sel6tPkLs\n37+fH//4x9xxxx0ceeSRnHPOOUydOpUHHngg735e5d+EdJ9/mc2eDe6wefN7ZYNqZehdughtn69N\nrm3FrO+97swzoaUlf5vudaVeTnKMK66AE04ovF+x56UU34tc0uxTDq2tMG9eVw53mJP1U8DL30cu\nU6ZM4Zlnnsm6bcKECSxbtoyNGzcyaNAgTjzxxJ5t48aNK3hBetNNNzF37lxOOeUUvva1r3H++eeX\nLngSaX5dKMeDOv+1UKQ/y/X6Xry4q/Q0eHDX8zRK0cfLL7/sl112mbe0tHhDQ4NPmTLFv/Od7yTa\n96mnnvLhw4cftu5f//VfvaWlJec+a9eu9X379vlf/vIXX7p0qQ8bNsx/97vf5T1OrnNIyrKPrvxF\npGpK8RtIaB87d+7kmmuuYcWKFQwZMoRLLrmEpUuXctRRRyXav6GhgXfeeeewdXv27GHYsGE59xk/\nfnzP8yuvvJIHH3yQFStWcP31lZsAWTX/Cok5X8zZQPlCxZxvwICuMs2NN6a/4yi0j3/5l3/huuuu\nY8iQIQD8+c9/ZujQoQBceOGFDBs2LOtj8uTJAJx88sm8++67bNq0qafPX/3qV5x++unpvqAK0ZW/\niPRre/fu5dRTTwVgw4YNnHbaaT23YD7++OP5dgVg6NChfPazn+W2227jnnvu4fnnn2f58uWsWbMm\na/s9e/bw85//nPPPP59Bgwbxwx/+kKeffpq77767dF9UAqnv8zezY4AfAmOADuByd9+dpV0H8A5w\nCDjo7uP7tsm087RZRCRuMd/nv2XLFpYtW8bIkSPZvn071113Xc/gn9SuXbu46qqreu7zX7RoEdOm\nTevZPmnSJM477zzmzp3L22+/zaRJk3jppZcYOHAgY8eO5Y477uDTn/503mOU+j7/kMH/TuBtd7/T\nzOYAR7v73CzttgAfc/edBfrT4C9Sp2Ie/GtFTBO7XQwszTxfClySp23Z7ySOua4JceeLORsoX6jY\n80l1hAz+x7n7jszzHcBxOdo58B9m9pyZXR1wPBERKZG8ZR8zWwkMz7LpFmCpux/dq+1Odz8mSx8f\ndvfXzeyDwErgBnd/Oks7nzFjBk1NTQA0NjbS3NxMS0sL8N7Vi5a1rOXaW1bZJ1z3OWxvb6etrQ2A\npqYmFixYUPGa/0tAi7u/YWYfBla5+18X2Gc+sM/d78qyTTV/kTpV74N/MRO7NTQ00PtTuQ4cOMAX\nv/hFvvnNb+Y9Rkw1/2XAjMzzGcBPs4Q60syGZZ4PBS4AfhNwzJxir2vGnC/mbKB8oWLPVw+63yfw\n5ptv8oMf/IBrr72WF154IWvbffv2sXfvXvbu3csbb7zBBz7wAS6//PIKJw4b/BcBnzGzjcCnMsuY\n2QgzeyzTZjjwtJmtB9YCj7r7kyGBRaR+dHoni1cv5s5n7qTT083KVoo+QqSd2A3g4Ycf5rjjjuPc\nc8+tQNLDaT5/ESm7XCWLxasXM2/VPMyMf2r5J+acW/ysbKXoI5ckE7utW7eOc889l/379/ds+/rX\nv057ezvLli3L2/+nPvUpWlpauO222wpmKXXZR+/wFZGqMbOe+nfvOnil+9i4cSO33norb731Fs89\n9xwtLS1MnjyZRx99tOC++/bt46/+6q8OWzds2DD27t2bd7+tW7fy1FNPcd9996XKHKpuBv/29vae\nOwtiFHO+mLOB8oWKOd/ss2fj7pgZs89ON7NbaB/VmNgN4IEHHmDChAmMGTOm6MylUDeDv4jUngE2\nILhME9pHoYndVq9enXW/8847j8cee+ywid265/RPMrHb/fffz80335w6dyjV/EWk7GK+1fPGG29k\n1qxZjB07lg0bNnDffffxz//8z0X1MX36dMysZ2K3KVOmsGbNGsaOHZu1/bPPPssFF1zAjh07en7Q\nFBLN3D6lpsFfpH7FPPhXemI3gGuuuYYDBw6wdOnSXF2+jwb/HGKua0Lc+WLOBsoXKoZ8MQ/+tSKm\nN3mJiEiNqpsrfxGJl678w+nKX0REgtXN4B/7/CUx54s5GyhfqNjz1YNvfetbnHnmmQwZMoRZs2YV\nbL9z504uvfRSGhoaaGpq4sEHH6xAysPpPn8RkUDHH3888+bN44knnuDAgQMF2/eeCG7dunVMnjyZ\ncePG9XyWcCWo5i8iZZez5t/ZCa2tYAazZ8OAFMWIUvRRIvPmzWP79u15p2zYv38/xxxzDBs2bOh5\nU9iMGTMYMWIECxcuzLmf5vYRkfrR2grz5nUN3O4wJ8U7dUvRRw5JJnbrLckF7MaNGxk0aFDPwA8w\nbty4ipfnVPOvkJjzxZwNlC9U1PnMuh7dz6vUx8aNG7n88sv55Cc/ybBhw7jooov47ne/y6OPPsqu\nXbuyPrLN2JlkYrm0E8GVmq78RaR6Zs/uulrvLtlUoY/Qid16S3Lln3YiuFKrm8G/2u9gLCTmfDFn\nA+ULFXW+AQPCyzSBfYRO7NZbkiv/tBPBlVrqwd/MLgNuB/4a+G/u/nyOdhOBJcBA4B53X5z2mCIi\npbZ3796eu2w2bNjAaaed1jO3z+OPP56oj0OHDnHw4EHeffddDh06xJ///GcGDRrEwIED39d26NCh\nfPazn+W2227rmQhu+fLlrFmzpnRfVBLunupB16B/MrAK+JscbQYCm4Am4AhgPTA2R1sPsWrVqqD9\nyy3mfDFnc1e+UDHkC319l9PmzZt9yZIl/vDDD/uSJUv84MGDRfcxf/58N7PDHgsWLOjZfuGFF/rC\nhQt7lnfu3OmXXHKJDx061MeMGeMPPvhgwWPkOoeZ9UWP4amv/N39JSj4a854YJO7d2TaPgRMBV5M\ne1wRkVL6yEc+wpe//OWgPm6//XZuv/32nNtXrFhx2PLRRx/NT37yk6Bjhgq+z9/MVgFf8SxlHzP7\ne+C/u/vVmeXPAR939xuytPXQLCISJ83tE66i9/mb2UpgeJZNN7v78gT9F/XdnjlzJk1NTQA0NjbS\n3Nzc88eq7tvVtKxlLdfespROe3s7bW1tAD3jZRrlvvI/C7jd3Sdmlm8COj3LH301n3/1xJwNlC9U\nDPl05R8u1lk9cx34OeAkM2sys8HAFcD73xkhIlLD0kzU9tBDDzF27FgaGho48cQTc95SWi6pr/zN\n7FLgm8CxwB5gnbtfaGYjgO+7++RMuwt571bPe9096+QVqvmL1K96v/KfPn06APfee2/PRG3PPvts\nzonaVq5cydVXX82PfvQjxo8fz+uvv467M2LEiJzH0Mc4ikjNyTVwdbrTum0bBswePZoBKaZnKEUf\nIdJM1Hb22Wdz9dVXJ5r+uVusZZ+q6/4DU6xizhdzNlC+UDHna922jXlbtjCvo4PWbduq1kcuU6ZM\n4eijj876uPjii4HcE7Vt2LAha5+HDh3il7/8JW+++SYnnXQSo0aN4oYbbuBPf/pTSbMXUjfTO4hI\n7THee69Q2uv1UvSxceNGbr31Vt566y2ee+45WlpamDx5Mo8++mjBfYudqG3Hjh0cPHiQRx55hNWr\nVzNo0CCmTp3KV7/6Vb761a+m/ApSSPPOsHI8iPgdgCISJtfr+1Bnpy/q6PDFHR1+qLMzVd+hffzh\nD3/wT37yk37gwAF3d586darv3r078f7PP/+8H3nkkYeta21t9Ysuuihr+507d7qZ+f3339+z7pFH\nHvEzzjgj73FynUMq/Q5fEZFQA8yYM2ZMVfsIndit2Inajj76aEaOHJk6b6mo5l8hMeeLORsoX6jY\n81VboYnd9u7dm/XRPaNn74na/vjHP7J69WqWL1/O5z//+ZzHnDVrFnfffTdvvfUWu3bt4hvf+AYX\nXXRR+b/YXnTlLyL92rXXXsuyZct44YUX2L59O4sWLSq6j29/+9tcddVVfOhDH+LYY4/lu9/9LmPH\nju3ZPmnSJM477zzmzp0LdH3c49tvv83JJ5/MkCFDuOKKK7jllltK9jUloVs9RaTs6v0+/0rQrZ4i\nIhKsbgb/WOuand7J4tWL+cLdX6DTO0vS153P3Jm4ryT7tLe3F2wXsj3XtqTru7+32doX03e52v7n\nqv8s6utLc05Ctmd7bZTr/5LUDtX8y6z1mVbmrZqHdzgnPHMCc85N/3Fz3X11//qXpK+k+xRqF7I9\n17ak6z/Ox3O2L6bvcrV96LcP0ba7LfHXl+aclOr7VGy70H0kXnVz5V/tWQtzMTPMjAEfGZDo8z2T\n9NX9vFT7tLS0FGwXsj3XtqTru7+32doX03e52p74NycW9fWlOSch27O9Nsr1f6m/KmZit4aGBoYN\nG9bzGDRoEF/60pcqmLaLrvzLbPbZs3F3zIzZZ8+ueF9J9ynULmR7rm2lWF9rbdN+7aXaXmy70H36\ni+73Cbz55ps9E7uNGzcu68Ru+/bt63m+f/9+hg8fzuWXX17JuF3SvDOsHA/0Gb5VE3M2d+ULFUO+\nXK/vzkOd3rGowzsWd3jnoXTv8C1FHyH27dvngwcP9ldeeaVn3ZVXXulz584tuG9bW5t/9KMfTXSc\nXOcQvcNXRGrNttZtbJm3pauM5DBmTvHv1C1FH7lMmTKFZ555Juu2CRMmsGzZspwTuyW5CWXp0qVc\neeWVpYpblLoZ/GOt+XeLOV/M2UD5QkWdz3r9/SBgZrfQPio5sVu3rVu38tRTT3HfffelCx2obgZ/\nEak9o2eP7vqkb8s8r0IfO3fu5JprrmHFihUMGTKESy65hKVLl3LUUUcl2r+hoYF33nnnsHV79uxh\n2LBhefd74IEHmDBhAmMC5zZKK/XdPmZ2mZltMLNDZvY3edp1mNmvzWydmf0i7fEKifU+/24x54s5\nGyhfqJjz2QBjzJwxjLlxDDYg3WV7aB+FJnbrfWdO78fkyZMBDpvYrVu+id263X///cyYMaPovKUS\ncuX/G+BS4HsF2jnQ4u47A44lIlIWhSZ2K6T3xG733HMPzz//PMuXL2fNmjU593n22Wd57bXXuOyy\ny0rzRaQQPLePma0CvuLuz+fYvgU4093/UKAfD80iInHqfmNYjLZs2cKyZcsYOXIk27dv57rrrusZ\n/JPatWsXV111FStXruTYY49l0aJFTJs2rWd734ndrrnmGg4cOMDSpUsTHyPXObSUc/tUYvDfTNcH\nvB8Cvufu38/RToO/SJ2KefCvFaUe/PP+eDOzlcDwLJtudvflCY9xjru/bmYfBFaa2Uvu/nS2hjNn\nzqSpqQmAxsZGmpube+5U6K5b5lpesmRJUe0rvRxzvt414RjyKF/95ZPSaW9vp62tDaBnvEyj7Ff+\nfdrOB/a5+11ZtgVd+be3t0f9Hy3mfDFnA+ULFUM+XfmHi7XsM9vdf5ll25HAQHffa2ZDgSeBBe7+\nZJa2KvuI1CkN/uFKPfiH3Op5qZm9CpwFPGZmj2fWjzCzxzLNhgNPm9l6YC3waLaBX0REKiv14O/u\nP3H3Ue7+AXcf7u4XZta/5u6TM883u3tz5nG6uy8sVfC+etc1YxRzvpizgfKFij2fVEfdTOksIiLJ\n6TN8RaTsNP9/aVTsVk8RkVLQhV186qbsE3tdM+Z8MWcD5QulfGFiz5dW3Qz+IiKSnGr+IiI1rOL3\n+YuISO2qm8E/9rpczPlizgbKF0r5wsSeL626GfxFRCQ51fxFRGqYav4iIpJY3Qz+sdflYs4XczZQ\nvlDKFyb2fGnVzeAvIiLJqeYvIlLDVPMXEZHEQj7MpdXMXjSzX5nZj83sqBztJprZS2b2ipnNSR81\nv9jrcjHnizkbKF8o5QsTe760Qq78nwROc/dxwEbgpr4NzGwg8C1gInAqMN3MxgYcU0RESqAkNX8z\nuxT4O3f/XJ/1nwDmu/vEzPJcAHdflKUP1fxFRIpU7Zr/VcCKLOuPB17ttbw9s05ERKoo7+BvZivN\n7DdZHhf1anML8Bd3/7csXVTsUj72ulzM+WLOBsoXSvnCxJ4vrbyf5OXun8m33cxmApOAT+do8ntg\nVK/lUXRd/Wc1c+ZMmpqaAGhsbKS5uZmWlhbgvW9AruX169fn3V7t5djzaVnLWq6N5fb2dtra2gB6\nxss0Utf8zWwicBdwvru/naPNIOBlun44vAb8Apju7i9maauav4hIkapR878baABWmtk6M/t2JsgI\nM3sMwN3fBa4HngBeAH6YbeAXEZHKSj34u/tJ7j7G3c/IPL6YWf+au0/u1e5xdz/F3U9094WlCJ1N\n969FsYo5X8zZQPlCKV+Y2POlpXf4ioj0Q5rbR0SkhlX7Pn8REakhdTP4x16XizlfzNlA+UIpX5jY\n86VVN4O/iIgkp5q/iEgNU81fREQSq5vBP/a6XMz5Ys4GyhdK+cLEni+tuhn8RUQkOdX8y62zE1pb\nwQxmz4YBOX7eJmlXqE3I9lJvq1R/seQO3Te0/zTt0rZPu085+pDUNX/cPYpHV5Q6tGiR+xFHuA8e\n3PU8pF2hNiHbS72tUv3Fkjt039D+07RL2z7tPuXoQzwzdhY95tbNj9po63JmYEZ7Z2fX8wLtep6n\naZNye3t7e/5902wrYX9580WQu+f/XtosSdsk6SNLu4KvjaT9hu6To4+Cr40qi3ZsCZR3Pn8pgdmz\nwR02b+56Xqhd96/AadqEbC/1tkr1F0vu0H1D+8/X7qmncrctpt/QfXL1Uei1IWWhmr+ISA3Tff4i\nIpJY3Qz+sdflYs4XczZQvlDKFyb2fGmlrvmbWSswBfgL8DtglrvvydKuA3gHOAQcdPfxaY8pIiKl\nEfIZvp8BfubunWa2CMDd52ZptwX4mLvvLNCfav4iIkWqeM3f3Ve6e2dmcS0wMk/zeO/jEhHph0pV\n878KWJFjmwP/YWbPmdnVJTre+8Rel4s5X8zZQPlCKV+Y2POllbfmb2YrgeFZNt3s7sszbW4B/uLu\n/5ajm3Pc/XUz+yCw0sxecveng1KLiEiQoPv8zWwmcDXwaXf/U4L284F97n5Xlm0+Y8YMmpqaAGhs\nbKS5uZmWlhbgvZ++WtaylrXcn5fb29tpa2sDoKmpiQULFqSq+Yf8wXcicBdwvru/naPNkcBAd99r\nZkOBJ4EF7v5klrb6g6+ISJGq8Savu4EGuko568zs25kgI8zssUyb4cDTZraerj8KP5pt4C+F7p+M\nsYo5X8zZQPlCKV+Y2POllfo+f3c/Kcf614DJmeebgea0xxARkfLQ3D4iIjVMc/uIiEhidTP4x16X\nizlfzNlA+UIpX5jY86VVN4O/iIgkp5q/iEgNU81fREQSq5vBP/a6XMz5Ys4GyhdK+cLEni+tuhn8\nRUQkOdX8RURqmGr+IiKSWN0M/rHX5WLOF3M2UL5Qyhcm9nxp1c3gLyIiyanmLyJSw1TzFxGRxOpm\n8I+9LhdzvpizgfKFUr4wsedLq24GfxERSS7kYxzvAC4GHPgDMNPdX83SbiKwBBgI3OPui3P0p5q/\niEiR0tb8Qwb/Ye6+N/P8BmCcu/+PPm0GAi8Dfwv8HvgvYLq7v5ilPw3+IiJFqvgffLsH/owGINuH\nuI8HNrl7h7sfBB4CpqY9Zj6x1+VizhdzNlC+UMoXJvZ8aaX+DF8AM/sa8Hngj8BZWZocD/QuBW0H\nPh5yTBERCZe37GNmK4HhWTbd7O7Le7WbC5zi7rP67P93wER3vzqz/Dng4+5+Q5Zj+YwZM2hqagKg\nsbGR5uZmWlpagPd++mpZy1rWcn9ebm9vp62tDYCmpiYWLFhQ2Zr/YZ2YjQZWuPvpfdafBdzu7hMz\nyzcBndn+6Kuav4hI8Spe8zezk3otTgXWZWn2HHCSmTWZ2WDgCmBZ2mPm0/2TMVYx54s5GyhfKOUL\nE3u+tEJq/gvN7BTgEPA74FoAMxsBfN/dJ7v7u2Z2PfAEXbd63pvtTh8REaksze0jIlLDNLePiIgk\nVjeDf+x1uZjzxZwNlC+U8oWJPV9adTP4i4hIcqr5i4jUMNX8RUQksboZ/GOvy8WcL+ZsoHyhlC9M\n7PnSqpvBX0REklPNX0SkhqWt+QfN6imFdbqzeOtWntq9G8w476ijMDNwpxNYvWcPEwqso9e/A4B/\nHDWKu159FU+4zvo8nz16NACt27blXB5gRqf7+9Z1f03Z1vf+mvNtz3Weit1Hapu+59VVN4N/e3t7\nzwx4MWndto1bOzroXL8empt5YteunlpbJ10fg/bvBdb1NtCMVbt387Ndu+hMuM76PO/+/Wreli2Y\nGZvWrOGjn/hEz7IDc8aMoXXbtvet6/6asq3v/TXn257rPOXaJ9bvbTflS6f7e+7r1+OXXpro/0k1\nxHr+QqnmX2aWeZS8zz5XSYXWve95geVs++TqN1++pF97mn2ktul7Xl2q+ZeZyj4q+0h2+p6XRsU/\nw7fU6nXwFxEpp37/Jq/Y78WNOV/M2UD5QilfmNjzpVU3g7+IiCSnso+ISA2r+H3+ZnYHcDFddyb+\nAZjp7q9madcBvEPXJ34ddPfxaY8pIiKlEVL2udPdx7l7M/BTYH6Odg60uPsZ5Rz4Y6/LxZwv5myg\nfKGUL0zs+dJKPfi7+95eiw3A23ma6x4uEZGIBNX8zexrwOeBPwJnufvuLG02A3voKvt8z92/n6Mv\n1fxFRIpUlpq/ma0EhmfZdLO7L3f3W4BbzGwu8A1gVpa257j762b2QWClmb3k7k9nO97MmTNpamoC\noLGxkebm5p63VXf/6qVlLWtZy/15ub29nba2NoCe8TIVdw9+AKOB3yZoNx/4So5tHmLVqlVB+5db\nzPlizuaufKGUL0zs+TJjZ9Hjduqav5md1GtxKrAuS5sjzWxY5vlQ4ALgN2mPKSIipZG65m9mDwOn\n0FXL/x1wrbu/aWYjgO+7+2QzOwH4cWaXQcAP3H1hjv48bRYRkf5Kc/uIiPRDmtsn8ntxY84XczZQ\nvlDKFyb2fGnVzeAvIiLJqewjIlLD+n3ZR0REkqubwT/2ulzM+WLOBsoXSvnCxJ4vrboZ/EVEJDnV\n/EVEaphq/iIikljdDP6x1+VizhdzNlC+UMoXJvZ8adXN4C8iIsmp5i8iUsNU8xcRkcTqZvCPvS4X\nc76Ys4HyhVK+MLHnS6tuBn8REUlONX8RkRqmmr+IiCQWPPib2VfMrNPMjsmxfaKZvWRmr5jZnNDj\n5RJ7XS7mfDFnA+ULpXxhYs+XVtDgb2ajgM8AW3NsHwh8C5gInApMN7OxIccUEZFwQTV/M/u/wB3A\n/wM+5u47+2z/BDDf3SdmlucCuPuiLH2p5i8iUqSK1/zNbCqw3d1/nafZ8cCrvZa3Z9aJiEgVDcq3\n0cxWAsOzbLoFuAm4oHfzLO2KupSfOXMmTU1NADQ2NtLc3ExLSwvwXt0t1/KSJUuKal/p5Zjz9a5p\nxpBH+ZRP+fLnaWtrA+gZL1Nx96IfwOnADmBL5nEQ6AA+1KfdWcC/91q+CZiTo08PsWrVqqD9yy3m\nfDFnc1e+UMoXJvZ8mbGz6HG8JPf5m9kWstf8BwEvA58GXgN+AUx39xez9OGlyCIi0p9U+z7/nlHb\nzEaY2WMA7v4ucD3wBPAC8MNsA7+IiFRWSQZ/dz+h+6rf3V9z98m9tj3u7qe4+4nuvrAUx8umd10u\nRjHnizkbKF8o5QsTe7609A5fEZF+SHP7iIjUsGrX/EVEpIbUzeAfe10u5nwxZwPlC6V8YWLPl1bd\nDP4iIpKcav4iIjVMNX8REUmsbgb/2OtyMeeLORsoXyjlCxN7vrTqZvCPlXc6Wxdv5Y0H38A7K1PW\n6j7m1ju3pjpmMfsnaVuoTb7t3um88W9v5NyWbb9i1oe2zZavXMdPsz5kW5LtxbbLtk8lXxvynryz\netaS7tmyKB4mAAAF+ElEQVTvYrOtdRtb5m3hw/Zhtp2wjTFzxlTsmGYGTsFj9j13xeyfpG2hNvm2\nb2vdxvC24XRYR9Zt2fYrZn1o22z5ynX8NOsBTvivE/Ke37TfmzTtsu1TyddGGrGOLaF05V9u1vUH\nme7nNXHMYvZP0rZQm3zb02wrZn1o23L0Wcr1IduSbC+2Xeg+UjJ1c+Xf3t4e5U/o0bNHg8OazWs4\nb/Z5FT0mlnleQN9zV8z+SdoWapNv++jZo1mzaQ2fOPETWbdl26+Y9aFts+Ur1/HTrAfYfOZmTvin\nE3Ke37TfmzTtsu1TyddGGrGOLcHSzANdjgeaz79qYs7mrnyhlC9M7Pmo5nz+paD7/EVEiqf7/EVE\nJLHgwd/MvmJmnWZ2TI7tHWb2azNbZ2a/CD1eLrHfixtzvpizgfKFUr4wsedLK2jwN7NRwGeArXma\nOdDi7me4+/iQ4+Wzfv36cnVdEjHnizkbKF8o5QsTe760Qq/8vw7cmKBd2W/k2r17d7kPESTmfDFn\nA+ULpXxhYs+XVurB38ymAtvd/dcFmjrwH2b2nJldnfZ4IiJSOnnv8zezlcDwLJtuAW4CLujdPEc3\n57j762b2QWClmb3k7k+nSptHR0dHqbssqZjzxZwNlC+U8oWJPV9aqW71NLPTgZ8Bf8ysGgn8Hhjv\n7m/m2W8+sM/d78qyTfd5ioikkOZWz1Tv8HX33wLHdS+b2RbgY+6+s3c7MzsSGOjue81sKF2/KSzI\n0afe4C0iUiGlus+/56rdzEaY2WOZxeHA02a2HlgLPOruT5bomCIiklI07/AVEZHKqal3+JrZRDN7\nycxeMbM5Odp8M7P9V2Z2Rkz5zKzFzPZk3vC2zsxurWC2/21mO8zsN3naVPPc5c1XzXOXOf4oM1tl\nZhvM7Ldm9qUc7apyDpPkq/L/vyFmttbM1pvZC2a2MEe7ip+/JNmq/f8vk2Fg5tjLc2wv7tylmRCo\nGg9gILAJaAKOANYDY/u0mQSsyDz/OPDzyPK1AMuqdP4mAGcAv8mxvWrnLmG+qp27zPGHA82Z5w3A\ny5H9/0uSr9rn8MjMv4OAnwPnRnT+CmWr6rnLZPhH4AfZcqQ5d7V05T8e2OTuHe5+EHgImNqnzcXA\nUgB3Xws0mtlxVEaSfFClmcu96/baXXmaVPPcJckHVZz13d3fcPf1mef7gBeBEX2aVe0cJswH1T2H\n3XcHDqbrYmlnnybVPH+FskEVz52ZjaRrgL8nR46iz10tDf7HA6/2Wt6eWVeozcgy58p37L75HDg7\n82vZCjM7tULZkqjmuUsimnNnZk10/Zayts+mKM5hnnxVPYdmNiBz88cOYJW7v9CnSdXOX4Js1f7/\n9w3gfwGdObYXfe5qafBP+pfpvj8VK/UX7STHeR4Y5e7jgLuBn5Y3UtGqde6SiOLcmVkD8DDw5cwV\n9vua9Fmu6DkskK+q59DdO929ma5B6Twza8nSrCrnL0G2qp07M5sCvOnu6yj4eWqHyXvuamnw/z0w\nqtfyKLp+uuVr0/3ms0oomM/d93b/eunujwNH5JoNtQqqee4KiuHcmdkRwCPA/3H3bC/+qp7DQvli\nOIeZY+8BHgPO7LOp6v8Hc2Wr8rk7G7g4836qB4FPmdn9fdoUfe5qafB/DjjJzJrMbDBwBbCsT5tl\nwJUAZnYWsNvdd8SSz8yOM+v60FIzG0/XrbbZaovVUM1zV1C1z13m2PcCL7j7khzNqnYOk+Sr5jk0\ns2PNrDHz/AN0zQa8rk+zqpy/JNmqee7c/WZ3H+XuHwGmAf/p7lf2aVb0uauZz/B193fN7HrgCbr+\nIHOvu79oZl/IbP+eu68ws0lmtgnYD8yKKR/w98C1ZvYuXVNjTKtUPjN7EDgfONbMXgXm03VXUtXP\nXZJ8VPHcZZwDfA74tZl1Dww3A6O7M1b5HBbMR3XP4YeBpWY2gK6Lzgfc/WeRvH4LZqP6//96c4DQ\nc6c3eYmI9EO1VPYREZES0eAvItIPafAXEemHNPiLiPRDGvxFRPohDf4iIv2QBn8RkX5Ig7+ISD/0\n/wHdvSs6ZZFWWAAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from scipy.special import erfinv\n", "\n", "def normalSpace(x_min, x_max, N, mu, sigma):\n", " area_min = norm.cdf(x_min, mu, sigma)\n", " area_max = norm.cdf(x_max, mu, sigma)\n", " A = area_max - area_min\n", " A0 = A/(N-1)\n", " seq = [x_min] + [np.sqrt(2)*sigma*erfinv((area_min+i*A0)*2-1)+mu for i in range(1,N-1)] + [x_max]\n", " return np.asarray(seq)\n", "\n", "fig = plt.subplot(xlim=(-0.2,4))\n", "parameters = [(1.75, 0.5), (1.75, 0.7), (1.5, 1.0), (0.0, 0.6), (2.0, 0.7)]\n", "for i, (mu, sigma) in enumerate(parameters):\n", " x_min=0.5; x_max=3.0; N=25\n", " points = normalSpace(x_min, x_max, N, mu, sigma)\n", " plt.plot(points, np.full(N, -i-0.3), \".\", markersize=5, label=\"$\\sigma$={:.1f}\".format(sigma))\n", " plt.legend()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Works like a charm! Example C++ implementation:\n", "\n", "```cpp\n", "#include \n", "#include \n", "\n", "std::vector normalSpace(double x_min, double x_max, unsigned int N, double mu, double sigma){\n", " boost::math::normal normal_dist(mu, sigma);\n", " double area_min = boost::math::cdf(normal_dist, x_min);\n", " double area_max = boost::math::cdf(normal_dist, x_max);\n", " double A = area_max - area_min;\n", " double A0 = A/(N-1.0);\n", " double area;\n", " std::vector x;\n", " for(int i=0; i