{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# iPython Cookbook - Curve Fitting" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\n", "import scipy as sp\n", "import matplotlib.pyplot as plt\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Generating some data to play with\n", "We first generate some data to play with. So in the first step we choose a functional relationship $func_0$ between our $x$ and $y$ that is parameterised by three parameters $a,b,c$..." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [], "source": [ "a,b,c=(1,2,1)\n", "\n", "def func0 (x,a,b,c):\n", " return a*exp(-b*x)+c" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "...we now choose the area over which we want to plot the function $x_{min}\\ldots x_{max}$ and the number of points we want to generate $N$ and generate the undisturbed curve $yvals_0$ over a random selection of $x$ coordinates $xvals$." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEACAYAAACnJV25AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X90XPV55/G3sYUU4yBBlZgUB02KIU5rsEJSx0lMdU0I\nMT9KbNYnJcQLSqCn6TlqmLNlqclJ0Ti7DS6H3Yyz2lOTuCAINeweEatYdZ2Goqsi1jYhQQKnxwQD\nI2NSAwqWiDE22Gj/eOabeyXrx5Xm98zndc4czcy9mlxu4Jmvnu/zfb4gIiIiIiIiIiIiIiIiIiIi\nIiIiIlLSPgx0A78A9gDfmOC87wHPA/3Ax/NzaSIiMhNnAY3p5/OA54CPjTnnCmB7+vmngF35uTQR\nEcmGTuBzY97bBPxJ6PVeYH7erkhEREY5ZRrnxrA0y+4x758NvBx6fQBYkNlliYjITEUN7POADuBm\n4PA4x2eNeT2SyUWJiMjMzYlwThXwMPAAlooZ6xVsktVZkH5vjHNH4IXpX6GISGV7AVg4nV+YasQ+\nC/h74N+B5ATnPAJcn36+DBgCXh3v2lpbR/jud0fo7h5hZKRyH62trQW/hmJ56F7oXuheTP4Azp1O\nUIepR+yfBdYCzwBPp9/7JnBO+vndWEXMFcA+4C3gqxN92Pr10NAAdXWQTILnTfdyRURkKlMF9l6i\n5eFbovyPtbaC79tzBXURkdyIkmPPmmQShofteSwGixbBunWVF+S9SvsHnoTuRUD3IqB7kZmx1Sy5\nNHLDDSP4PgwM2Og9lYLGRojH83gVIiIlZNasWTDNWD2dOvaMtbfbSL2pCRIJez40lM8rEBEpf3kN\n7ACrVllAD3N5dxERyVxeUzEjI5aK2bABjh6Fnh47MH++pWQqMd8uIjKZmaRi8h7YwxIJ++n7FtDd\naxERMTMJ7Hmtignz/SAF40buiYQFeI3aRURmLu859rCammDytKfHAv2ttyrnLiKSiYKmYhw3Sncj\ndgV2ERFT9OWOExlbJbN6dUEuQ0SkLBRFYJ83Dzo7bbTe0wPbtsHZZ9tKVRERmZ6CTZ6GtbUFzz0P\n+vqgqspKIEVEZHqKIscOln7p7obDh+HECZg1yx5XXw1bt+bxKkVEikjJ5tjBgvfatXD66RbQR0bs\nsXMntETqHSkiIlBEgR0sJfOjH8Gpp0JtrQX2RYtsFK9KGRGRaIoixx7W0QGnnGLBHCzfvns3/Oxn\n8Oyzhb02EZFSEGXEfg+21d1EYfUMYCvQD+wG/iCTC2prg+3b4aKLoLra+rdXV8NLL2nULiISRZTA\nfi+wcpLj3wR+DizB9j7dmOlFeR7ceSfMmWMtfoeH4fhxa/ur4C4iMrkogf1x4NAkxz8GdKefPwfE\ngA9kdlmWghkZgd5ee33sGDzwAFxzjSZTRUQmk43J037gmvTzpUADsCDTD43H4Z/+ySplamvtveXL\n4RvfgDVrMv10EZHylY3J0w1Y+uVpLA//NHBivBMTob68nudNua+hO7xjR9C/3fOsn3v4uIhIufB9\nHz/DnHPUovcYsA24IMK5L6XPOzzm/UkXKE0mmYQtW+CKK4JGYYsXj16xKiJSjgq1QKkWODX9/E+B\nHk4O6hmJx20yNayrS71kRETGE+Vb4EGgCajHyh5bgar0sbuBTwPtwAiwB7gRGB7nc2Y8YgcL4u3t\nUFdnKZnqauvnvnatRu4iUr5Kbmu86WppsZH6wYNWJVNdbSWRzc0K7iJSnkq6V0wUbW2WlqmpsaB+\n7JgF9t5epWVERJySCuxggX3tWgvuDQ22eKmuzrbY0+IlEZESDOxgdeznn3/y+319+b8WEZFiU5KB\n3bUcWLTIRu09PfZ+Z6dSMiIiJTV5Op7mZtszNZGwRypl1TMiIuWg7CdPx9PYaLn1RALWr7fnzc3K\nt4tI5Sq6fuzTFY/bz6Eh6wTZ0zN61K62AyJSaUp+xA7jb3qdSuX9MkREikJZBHbPs0csNnoyVf3b\nRaQSlXwqxnEpl74+C/Dh4B4+LiJS7kq+KmYs3z95lB7qFiwiUlJmUhVTNiP2sFTKHj09NqGaSEB/\nP2zdWuALExHJg7LIsYd5XlDb7qpkAH7yE+XbRaQylF0qJsylYBIJ6yfT2WmvlW8XkVJR9m17p2v1\naujutqA+MBDsnbpihdIyIlIalGMf4+abLZjHYrYq1dW7uwAvIlKOouTY78F2Tnp2guP1wA6gD9tB\nqTkrV5YFLt+eSsHs2ZZv7+mxXPvKlcq5i0h5ijK8vxjbw/R+xt/MOgFUA7dhQf45YD5wfMx5eU/F\nhC1aBGedZc9dtUwsZoFfOXcRKVa5SsU8DsQmOf4fwIXp56cDv+bkoF5wmzYFI/RwIO/rU2AXkfKS\njRz7D4DHgF8B7we+lIXPzDrPsyDe1ze6xv2ZZ+y5JlNFpFxkI7B/E8uve8C5wE+AJcBvxp6YCC0B\n9TwPL89D5XjcRu2uzUBPj+2d2t1tG3S4TpEiIoXi+z5+hhOAUfM2MWAb4+fYtwN/AzyRfv2vwF8B\nT405r6A59rBkMthpaWDAgntdHXz60xq5i0hxKdRGG3uBS9PP5wMfBV7Mwufm1OHDcOCAPT92DF57\nDbZtg5aWwl6XiEimogT2B4H/hwXsl4GvAX+WfgB8B/gk0A88CtwKvJH1K82ieBwGB+Guu2y0XlsL\nIyOwYIEFfJVBikgpK+uVp1Pxfbj8cnj3XThxwt6rrYX6ehu5K+cuIoWmlafT5Hlwxx2wbx/s2RO8\nB+PvyiQiUgoqesTu+L6NzoeGbDK1qcneX7wY2toKemkiUuEKNXla8lzrgUWLRrf67ehQvl1ESo8C\ne1o8DuvW2XMX3I8eHX9HJhGRYqZUzBgtLdDVZc/DaZlVqzSZKiL5p8nTLFizxkoeUykL7G4ydd++\nQl6ViEh0GrFPwPeDlr9gAT6RUMMwEckvTZ5mkefZZGqY72tlqogUP6ViJrFwYTBC7+mx0fvBg/Za\nZZAiUqyUipmCq4ppbw/2TW1s1GSqiOSHUjE54HnWV+bgQQvqw8PW0z2RUFpGRIqTRuwRtbRY24G+\nPgvuTU1QU2N7p2rkLiK5ohF7DrW1WfoFggVMe/dCZ6cWMIlIcdGIfZpaWqz7I1hA97xgglWlkCKS\nbTMZsSuwT1MyaaP3wUFLyUDQ6nfzZgV3EckuBfY88v2TFywlEsEoXkQkG3LVUuAe4ErgNcbf8/QW\n4Cuhz/sYUA8MTedCStHevUHrAbCgfuCARu4iUlhRvgUuBg4D9zN+YA+7CogT7IEaVlYjdrC0TGNj\nMHnqRvCqcReRbMllKiYGbGPqwL4F+Ffg78c5VnaB3XE9ZTwP1q+HhgZrR6BSSBHJVKHLHecCXwAe\nzuJnloTGRojFRr+3dy/s2FGQyxGRCpfNXjF/DPQySW49kUj89rnneXhlkoiOxy0t09lpFTKu9cCu\nXVYeqb4yIhKV7/v4GS6OyWYqZivwf4CHJjhetqkYJ5mEb38bjhyBY8fsvdpaOP98uO46pWVEZPoK\nmYqpBf4I+McsfV5Jisfh9tvhrLMszw6Wppk711oRJJOFvT4RqQxRvgUeBJqwEsZXgVagKn3s7vTP\nG7D8+nWTfE7Zj9idibbXe/55eOWVwl2XiJQeLVAqIi0tNnl6/HhQ5z57Nlx0Edx5p+rcRSQaBfYi\n4/uwYYNVyEAwqVpfb4FfOXcRmUqhyx1lDM+zXZjGGhzU5tgikjsK7DnW1gZXXQV1dZaKcY3Duro0\nmSoiuaHAngdtbbY6taoq2IUJrO599Wr1cxeR7FJgz5N4HG680UbuDQ2Wb/c8ePppuPVWBXcRyZ5s\nrjyVKbS1Wc69s9NG8K7NbyoVBHZVy4hIplQVUwCrV8POnbZnqiuFVLWMiIxHVTElYutWWLPG0jLV\n1fbe8LBVy7S1KS0jIpnRiL3AGhstFdPYaBtkNzVZp8jmZqVlREQj9pLU3Bxsjt3UZME9lYL2do3c\nRWRmNHlaYC6f7nZiCo/SNWIXkZlQKqZI+L6N0lOpICXjeUFwV5AXqUy52sxa8sAF7vZ2+9nTY++5\nckgFdhGJSiP2IuQ2mgr/TCSsBYFKIUUqi0bsZcLzbOTuAvv69fbTjeYV3EVkMlMF9nuAK4HXmHhb\nPA/4Lrb5xmD6tWRgvLy6W6U6NKR9VEVkclMN7y8GDgP3M35grwOewHZPOoDtsjQ4wWcpFTMDyaQ9\nYjHLuzc0wMGDcNNNCu4ilSAXqZjHsY2sJ3Id8DAW1GHioC4z5NIuQ0M2ieo269izx1oTNDUpNSMi\no2WaYz8PS8F0A+8HNgI/zPSiZLR43NIvBw8GbX/7+uDwYXjyyeAcERHIPLBXARcBnwPmAjuBXcDz\nGX6ujOHSLnv2WFAfHrYgX1UVHFNwFxHIPLC/jKVf3k4//g1YwgSBPeHKPADP8/BUnD0tbW2Wfnnz\nzWA3Jrdpxy23QEcH9PYW9hpFJDO+7+Nn2E8kSkI+Bmxj/MnTRUAbNnlaDewG/gT493HO1eRplrS0\n2J6pvm+tf4eHgy6RGzZYewJ9Z4qUh5lMnk518oNAE1bt8irQiqVfAO5O/7wF+CrwHvAD4HsTfJYC\nexYlk3DbbTAyAseO2XuzZ8PixXD0KOzdW9jrE5HsyEVgzyYF9ixbvRp+9jM4cMACPFje/Te/gauv\ntr7vIlLa1La3wmzdCvv3w2WXwfz5FtTjcXjvPTh0CFauVOtfkUqkwF4G1q2z3ZiOHg0CeSoFu3bZ\nhKqCu0hlUa+YMuB5sGkTbNwIS5ZYOeTAgK1S3bPH6t3deSJS/pRjLzPJJHzrW7YrkwvuAIsWwcKF\nakMgUmrU3VGIx60UsrfX2hAMDATHenstuGshk0h504i9TPm+BfC6Onvd02OTq42NVhK5Zo1SMyKl\nQOWOMkoyCVu2wC9/aa9dG4KjR+HCC4M+MyJSvBTY5SS+b5OpnZ1Bj5mGBmsDvGqVnaPUjEjxUmCX\nCbW0wPe/D3PnBv1lamvhlFPguuuUmhEpVpo8lQm5apj6ehvFuyDu+1YSWV9vrxXcRUqfAnsFaWuz\nvPszzwRpGbAGYnv32oImUHAXKXVKxVSgZNJKIjs6rCTSNRGrrbWR+wUXqM+MSLFQKkYicZOl9fU2\nqTo0FOyp2tgI3d2j0zUiUlrUK6aCeR40Nwe17k1NFtzBAvvy5TbpKiKlRakYIZm0gN7dba9dSeT+\n/baRx/btGr2LFIrKHSUjvm+P9nbrDulG8vG4BfaODvWaEcm3XAX2e4ArgdcYf3s8D/hH4MX064eB\n/z7OeQrsJWD5cnjiCZg1K9i8A2DOHHv88z9r9C6ST7naaONeYOUU5/QAH08/xgvqUiJ6e21F6vXX\n21Z7ztlnw4kTNpoXkeIW9VsgxsQbWnvAXwJ/PMVnaMReYubNg7ffhqqqoCSyutpSNB/4ADz7bGGv\nT6QSFGprvBHgM0A/sB34/Sx8phSBz38ezjvPerm7vu7Hjtnrl17SzkwixSobI/b3AyeAI8DlwEbg\n/HHO04i9BLng7ftWPROPQyJhde+pVHBMeXeR3MhlVUyMiQP7WC8BnwDeGPP+SGtr629feJ6Hp2hQ\nMnwfvvQleOcdS8WEd2f64AfVAlgkW3zfxw/9Obx+/XooQGCfj1XMjABLgf+bPn8sjdhLnPt3zfNG\nj9gTCXtAMKoXkezIVUuBB4EmoB54GWgFqtLH7gbWAH8OHMfSMddO5wKkdIT/wDp6NAjmNqAwnZ3W\nlkB/jIkUTpTA/uUpjv/v9EMqyJo1QWCH4Llb5ATWQVKjd5H808pTydjKlTaCB2tN4HLvc+bA5s0a\nvYtkQt0dpSDWrbOfvm9594EBayiWSgULmhTcRfJHI3bJGt+Ha66B994bvf1eTY0Ff6VlRKavUAuU\nRAAbld9+u7UkaGqylMzwsC1oGhrSgiaRfFFgl6yKx63Heyxmr8f2eHfBPZlUoBfJFeXYJetcPr2x\n0Ubqnje6ggasLLKvb/T5IpIdGrFLTnhe0Md9IqkUbNhgo3cRyR6N2CWnwoE9mbQqmaEhq5wBaw3s\n2hFoclUkO1QVI3nldmjq7BxdOROLWeDv6ICFCxXkRRxtjSclIZmE226Dd9+1zTvA+rwvWABvvAHv\nex/8wz8o9y4CKneUEhGPwx13wKWX2mi9ocH6vK9dazXwr7+uihmRTGjELgXj+xbkX3/dRuojI8FO\nTbNnw/z5sHq19aXR6F0qlVIxUnJcRcyWLbB/P7z2WrCJdm0tnHaa7eS0b5/txypSaRTYpWS5XZhi\nMTh40FoQhPu8b9hgjcbCPeFFKoECu5S8lha47z4491zbqcmtWp01C66/PjivuVnBXSqDujtKyWtr\ns59r1sC3vmW5drDqmYceguPHLT0Tiymwi0xkqm+Be4Arsa3vJtvv9A+BncCXgB9NcI5G7DJtiYSt\nUP3hD+Hii4MRfEODBfZ58zS5KuUtFyP2e4H/Bdw/yTmzgb8Fdkz3f1xkKi5gb9kS7LM6MGCPVAr2\n7oXDh4O+M1rYJDJ1YH+c8TemDvsLoAMbtYtklQvsS5faitVEIphoTSSCCdfOTmtV0Ng4+vdEKlGm\nOfazgS8Cl2CBXbkWyQlX6uh5sH07bNpkgX5gwNIz1dU2wep2bHLnilSiTAN7EliHBfRZTJGKSYR6\nt3qeh6f/8mSaPA/uvNOe9/VZa4LwwqaHHrIAf/CgNtOW0uT7Pn6GS6+j5MRjwDbGnzx9MfQZ9cAR\n4E+BR8Y5V5OnknUu9bJvH7z1VvD+/Pn2U1vySakrRK+Y3wM+kn50AH/O+EFdJCeam2316sKFVgbZ\n0GDvv/qqLWhqb9euTVJ5pkrFPAg0YaPxl4FWoCp97O4cXpdIJG403twMO3bY86EhW9w0MGDPEwnb\nh9VNsCpFI+VuqsD+5Wl81lczuRCRTMTjQVrGtR8AC+7NzRbQh4bsvc5O+9nTA1u35v1SRXJOLQWk\n7Lh0SyIBNTU2kRretam6GubOhXfegY98BG68USN4KV7qxy6CVc54nqVfli0LUi/h3u/vvWej+l/8\nwkbw2ndVyolG7FLW3GImsAnWV14J0jTOrFlw5pkwOGgBXqN3KSbq7igyiWTS8vBXXglHjpx8fO5c\n+NCHYPNmTbBK8VBgF4kgmbT0y+7dNnqvrbWNtZuaLKAvW2bvJxL2urFRq1ilcNS2VyQCNxI/4wzY\ntg3q64Njw8PWWGxoKKiuUYmklBqN2KWiXXCBVcW0t8NZZ1lO3rUnqK21/u/LlwcjeEcjeMkXpWJE\nZshNmjY2Whvg/fvh3XeD41VVcM45FuQPHrSJWPWBl3xQYBfJUDJpqRfftxr4XbuCVayuXcGbb8Lp\np8MHPwhPPlnQy5UKoMAukgW+bzn1oaFglerYJmNV6cYa//IvGrVLbimwi2SRC/CNjXDTTdZY7MgR\nW9zk1NRYnt61Et64EW6+WcFeskdVMSJZ5FawArS02Ah++3b4+c/tvRMnrD3wa69BR4dV13R326Sr\n+32RQlBLAZEI4vFgk4/Fi2HBAgvgqZQ92tqCc1Op0e2B1a5A8k2pGJFpSiYt597RYTs2ff3r8OKL\no6toTjkFzjvPtvCLxy2lE25vIBKVUjEieRBeqOR5tqDJ9+39PXssRfPee1YX395u5ZEwOmevAC+5\npMAuMkPh9IvrJplKWeCuqbFFTamUTbp6nj2fM8dq4cO/J5JtUQL7PcCVwGuMv+/pF4FvA++lH/8V\neCxbFyhSKjwP+vuDydNYzIJ5T4+9HhiwyVaXf4/FghG8WhZINkXJ21wMHAbuZ/zAfhrgKnwvALYC\nC8c5Tzl2qQgucLvR+Jln2qrVVMp60YAF//p6mDfPRvoPPAArV2o1q5wsVzn2x4HYJMdDyzaYBwxO\n5wJEys3YwLx0KaxbZyNyt9BpeNgep51mi6AOHoSuLmtnMN5niExHtnLsq4A7gA8Bl2XpM0XKwrp1\nFqiTSWskduCABfDhYQvy/f123pEjVmnjJls7Okbn8UWiijq8jwHbGD8VE3YxsBn46DjHRlpbW3/7\nwvM8PA1LpMIkk5ZPj8WsD81jjwVVNM7s2ZaLP3HCgrx2daosvu/jhxZCrF+/HnLUUiBGtMAO8AKw\nFPj1mPeVYxdh9OKlRMKC/JYtVvvuWgY3NNhWfV1ddk74d6SyFKqO/VzgRWAEuCj93tigLiJp4T9U\nV62yqpjeXkvPvPWWLXQaGAiOHzliLQ1cWkYLnWQqUb4FHgSagHrgVaAVSPe2427gVuB64F2seua/\nAD8d53M0YheZgFvNevhwkKbxfdumr6fHtu0DC/RDQ8GmHy7IK9iXr1yN2L88xfE70w8RmSGXQw+n\nXNyKVhe4nfBOTu54e7vl7js6bPQvlU0rT0WKSHjUPRgqHPb9ILjbXJpJpYKfqRQ89ZSN/tW2oLIp\nsIsUoXDL4FWrRr8OB/P77gtWt86fb3u0trcHvzN2sZRUBnV3FCkxiURQKXPTTSfvzwowaxZ85jNw\n6aXB77gmZCqdLC0zybGrH7tIiXGjb8+DtWvhnXfg3HNtFaszMmILoR54IBjh+76tclXpZPlTKkak\nxITTKu758uUWwHftgkWL7L26OkvRvPCCHRsasvfDFTRuExCN4suLUjEiZcClWTo6rDWw51n6pbnZ\nqmSOHw9q42trraSyuTnYrDuRUK/4YqVUjEiF8jwbdff2BhOnYAF87VqbUG1qshWtw8NBPXx4FO/S\nNG4Ury39SpdG7CJlLFwV4yZQa2psMdTgYNBGuLraJlybmmyDEN+30Xtzs31haAFU4cxkxK7ALlIh\nXLrGrVz1fbjmGlvtGq6qqaqCSy6xcxctsvMSCfsiUL/4/FNgF5EphUff7mdNjU28nnIKHDpkefjh\nYXv/U5+yn3v32gje/Y7aCueHAruITIvLow8NWcCPxeCRR+D000dPtoaDfCxmj/Z2q7ZZvdpSOJp4\nzQ0FdhGZEZemaWwM0jSxmKVoTpywwO9aCldXw+/8TtAvvq7Ofu/gQdi0ScE921QVIyIz4qpqPM8q\nZsDy6+edZwH72mttAdS8eRbgf/UrePVVez08DM8/b6N3N1nr/hLQYqjC0IhdRMYVbkHg+5Zfb262\nYw89BM89Zzn5sbs/XXqpfRn09dnov78fbr5Z7QxmSiN2EckaN4p3z6+6KkizXHut1cRffbW9npUO\nOydO2JfAc88Fo/Xu7qBOHux5S0se/gEqmFoKiEgkYytgBgftPZdj7+2FBQvsdX+/je4HBmwU7/uj\nF0N1ddnvqvtkbkQZsd+D7Zz07ATHvwL0A88ATwAXZufSRKQYuRbCLtCvWGF5+ZoaC+Z9fTa5Ojho\nI/kTJ6xnTX8/vO99NsH65pv2u+E+824Ur7x85qIE9nuBlZMcfxH4Iyyg/zfg+1m4LhEpEVu3Wsrm\n858PRt6LF9vCp9tvt3LJG26widbqapt0PXTIRvZ33QXbt9vvdHXZz3CgX7063/805SFKYH8cODTJ\n8Z1AemEyu4EFmV6UiJSerVtHbw4Co3vWfOITlo654QarsGlstM27f/pTO29wcPRovavL8vNOMqn+\nNVFlO8d+I7A9y58pIiUmPOm6YoX9rKuz92IxqK+3IL50Kfzyl7aq9a237AtheNhG8UND1pXScZOv\njY2wcaMqbSaTzcC+Avga8NmJTkiEduH1PA9PMyYiZW/rVvsZbmPg9nO94gp7JBIW8K+6ynLwTz1l\nm4WAlVTOmWOrYRcvti+E7m5YssSe79sHCxeWT4D3fR8/w4mGqLWRMWAbcMEExy8EfoTl4vdNcI7q\n2EVklHBVTCxmi5zc1n+xmOXjzz775E6U77wDl11mnShTKTvXtRx2rQ2SyfII9jOpY8/GiP0cLKiv\nZeKgLiJykvAf7VdddfJ71dU2Ivd9G9m//XbQ2uDHPw7OGxy0QN7ZaSkcz7PnjY12vNIalkX5FngQ\naALqsbLHVqAqfexuYDOwGtiffu9dYOk4n6MRu4hE1tICr7wyfirnO9+x+vhjx4KUTXW1rYK95BLY\nsSMoy4SgYRnYl8TGjcHnFjs1ARORshXuROnSLL4PzzxjqZn6+qAjZXW1Bf35862+/s034Y037Fgi\nYb8f/pxiTtsUKhUjIpJzLvD6vi12cpU24Ry7y9UDbNliOXpn3jx7nHNO8F5np32u+wn2l0Kpp20U\n2EWkpIRTLJ4XVMU4sZiNylOp4LxNm4IvARfsa2qsLfHYRVDhdgcdHcFnF+uIfjwK7CJS0sKj61Wr\ngglTt2E3BDtAAaxcCY8+an1tXnjBSieHhy2AL1gwut1BV1fwF4BL/fT1Ba+LlQK7iJSNcLB1gRiC\nunmAZctsG8B9+2xE7/tWUnn8uI3oDx2yYD40ZH1uXGCH0b1tHPdFUkyLpRTYRaRsjd3b1T3v7x99\n3nnnBUH7ssssXXP0qB3r6bGfp55qOfozz7SRvVsJG+5a6WroC921Uv3YRaSieF5Q6uh62oR72/zu\n71q9fGur9ZxfsgTOOMOC+qFDlr7p6YHHH4fdu23077iA7r4kli8f/V6++tArsItIxXKpk3AKxS2U\ncurq4MILrVzyC1+wfL3bOeroUVsotX69Bfu77rJRuwv2Tz01OrB3dQXB3nWvzEWbYqViRERC3GSs\n62kTrrhZtixI0TiuFcLevZbCGRqyVM+Pf2z96Ds7R6dkXLDv6rJza2vt/WyujtUCJRGRiMJVMS7H\n7oK2q6N/9FHYudNWxE4U8k47zYL+7NlBXn7s6liXq1+xQguURERyJlxD77iqmLo6G4H39lodPcCG\nDbBuXRDs3cbfb70V/H5Pjz3cFoIuoIcnYadLOXYRkRmIx+3hgr17Did/AfT2wl//tU3GVlcHm42E\nnThh+f2zzgp2lZopjdhFRLJkbHnlJz95ch/6zZuDzUbmzAlSMZs2wde/bumY++6z91yp5XQpsIuI\n5Ehvb/DcjeL7+oIg/8orweRpTU2QwnFtERIJq7iZLgV2EZE8csE+POkKo1fHZkpVMSIiRSbTqpgo\nk6f3YBtsPDvB8UXATuAo8JfT+R8XEZGTjdcKYTqiBPZ7sb1MJ/Jr4C+Au2Z2CZUn041qy4nuRUD3\nIqB7kZkSz3BpAAADQElEQVQogf1x4NAkx18HnsK2xJMI9C9tQPcioHsR0L3IjOrYRUTKjAK7iEiZ\niTrTGgO2ARdMck4rcBj4HxMc3wecG/nKREQE4AVg4ZRnhWSzjn2qL4lpXZiIiOTOg8CvgHeAl4Gv\nAX+WfgCclX5/GJtk3Q/My/9lioiIiIjIlFYCe4Hngb+a4JzvpY/3Ax/P03UVwlT34ivYPXgGeAK4\nMH+XlndR/r0A+EPgOHBNPi6qQKLcCw94GtgD+Hm5qsKY6l7UAzuAPuxeNOftyvJrqkWgUMC4ORub\nII0BVdj/GR8bc84VgGtI+SlgF+Upyr34NJBu/8NKKvteuPMeA7qA/5Svi8uzKPeiDvgFsCD9uj5f\nF5dnUe5FArgj/bweWwxZjv2tLsaC9USBfdpxM5vljkux/6NS2GKlh4AvjjnnauC+9PPd2L/E87N4\nDcUiyr3Yic1LgN2LBZSnKPcCbPVyB7bgrVxFuRfXAQ8DB9Kvs9gaqqhEuRf/AZyefn46FtiP5+n6\n8mmqRaDTjpvZDOxnY5OozoH0e1OdU44BLcq9CLuR4Bu53ET99+KLwN+lX5drt7go9+I84EygG1vR\n/Z/zc2l5F+Ve/AD4A6x4ox+4OT+XVnSmHTez+WdN1P8Yx5ZFluN/xNP5Z1qBVRp9NkfXUmhR7kUS\nWJc+dxb57TqaT1HuRRVwEfA5YC72l90uLL9aTqLci29iKRoPWwPzE2AJ8JvcXVbRmlbczGZgfwX4\ncOj1hwn+nJzonAXp98pNlHsBNmH6AyzHPtmfYqUsyr34BPanOFgu9XLsz/NHcn51+RXlXryMpV/e\nTj/+DQtm5RbYo9yLzwB/k37+AvAS8FHsL5lKUtC4OQe7+THgVKaePF1G+U4YRrkX52A5xmV5vbL8\ni3Ivwu6lfKtiotyLRcCj2OTiXGxC7ffzd4l5E+Ve/E9sRTtYTvkAlqYqRzGiTZ4WJG5eDjyHBazb\n0u+FFzMBtKWP92N/cparqe7FZmwy6On048l8X2AeRfn3winnwA7R7sUtWGXMs8A38np1+TXVvajH\nWpn0Y/fiunxfYJ5MtQgUKiduioiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiEih/X8NzGom9VMbfwAA\nAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "xmin,xmax = (0,1)\n", "N = 500\n", "xvals = np.random.uniform(xmin, xmax, N)\n", "yvals0 = func0(xvals,a,b,c)\n", "plt.plot(xvals, yvals0, '+')\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In a second step we generate the error term $err$ and add it to the previously computed values, storing the results in $yvals$" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEACAYAAACnJV25AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnWuMZVl1338FTGghzJRRS93i1aUM2IMtqLIt47YhqUuI\nNAM4phVLlnhEXX5GUTCeSJYz5oPr9gcnHclyKsYRM56EqcIjjT8MqGJHaEaOPLcNiCJGonsMziBm\n7NsabPeQJlSLRB4HlM6HfZfPurv2OWefe899/3/S1T11Hvvse2b6v9dZe621QQghhBBCCCGEEEII\nIYQQQgghhBBCCCGEEGKheS3wJPBl4EvAhxLnvB+4BjwFfBZ489R6J4QQojFnga3B9suBrwBvjM75\nUeDOwfa9wNF0uiaEEKINDoF3VBz/buBrU+qLEEKIMdkArhMs9zJ+GfidqfRGCCHEWLwc+AJwoeKc\ntwN/RrDahRBCzIiXZJxzB/AJ4BGCKybFm4GHCD72b6ZOuOuuu24/++yzo/RRCCFWmWeB1ze54EU1\nx9eA/0ywxPdKznkd8EngA8AzpT179llu376tz+3b7O7uzrwP8/LRs9Cz0LOo/gB3NRF1qLfY3zoQ\n7KeALw72fXgg5gAPAr9GcL98dLDv28BbmnZECCFEO9QJ+2eot+p/bvARQggxB9SJtpgAnU5n1l2Y\nG/QsCvQsCvQsxmNtive6PfAXCSGEyGRtbQ0aarUsdiGEWDIk7EIIsWRI2IUQYsmQsAshxJIhYRdC\niCVDwi6EEEuGhF0IIZYMCbsQQiwZEnYhhFgyJOxCCLFkSNiFEGLJkLALIcSSIWEXQoglQ8IuhBBL\nhoRdCCGWjJkLe6836x4IIcRyMXVhj4Vcwi6EEO0yc2EXQgjRLnWLWU+EXq8Q+EuXiv2dTvgIIYQY\nnakKe7d7UshtvxBCiHaYurD7b5BrRggh2mbmUTFyvQghRLtMXdhjIZewCyFEu6xN8V63b9++PcXb\nCSHE4rO2tgYNtXrmrhghhBDtImEXQoglo07YXws8CXwZ+BLwoZLzfgv4KnAN+IHWeieEEKIxdcL+\nbeBfAd8PnAf+JfDG6Jx3Aa8H3gD8AvDRcTulEEghhBidOmG/AVwdbP9v4H8Ar4rO+QngYLD9eWAd\nOJNqLFewJexCCDE6TXzsGwQ3y+ej/a8GnnN/fw14TaoBCbYQQkye3MzTlwOPAb9EsNxj4lCcZFxj\nrxeyTlM1Yax+TL8PBwfFftWPEUKIZuQI+x3AJ4BHgMPE8b8kTLIarxnsO8GVK10gCPiFCx22tjp/\nJ9om4N0u7O6qfowQYjXp9Xr0xnRv1AW9rxH8598gTKKmeBfwwcH3eWBv8B1ze3f39pBgd7vh0+ud\nLAgWC7s/RwghVoVREpTqLPa3Ah8AngK+ONj3YeB1g+0HgU8RRP0Z4P8AP92kAwD7+8NlfO+5Z9hl\nY24aCbsQQtRTJ+yfIW+C9YM5N/MiDUUJ334fdnaGhbtJBUiJvhBCFEy1bK8Ju7lfTLCvXCms9n4f\nNjbC/r09OD4O+69cGW7H2vLfQgghZrCCkom5nyyFIOjGpUvhbxPsK1dge3v4OmsrJegSeiHEKjOz\npfHiMMaNjULkbds+8XFrw6z/eFUmCbsQYpWZydJ429the339ZOx6vw9XB7muJti9Htx3X7HtXTid\nTmgvtuSFEGJVmbrFvr1dCPLxcSHI3iKPwx/7/ULYY0vfW/X2rQWyhRCrzExcMRcvVicg+eiZbjdY\n8zahmhJq21cWVSOEEKvE1F0x3pVipHziNrEaJy7F5/hvIYQQM1rzNFUnxu/zkTN1bTXZL4QQq8BM\nVlCqEl5v0ZcJfMrqz21fCCGWnZn42CGdgQrDCUr7++UWPkjAhRAixcyEPRZsmzD1oY9VNWJ8QpMQ\nQoiCmQm7x7JMIYRCnjsX/r5+vRD39fUQ397vF+UFzLqP68x4lKwkhFg15kLYvTD3++FjpQTMMt/a\nCq4ZKI+uSSFhF0KsGnMh7HFpARP6Xq/wt1sxMPtcuVIkL1VZ7EIIsWrMhbB7vO/dRD0+1usF692y\nTuOiYLmlgYUQYhmZS2E34f7MZ4bXPz08DC6ZnZ2Tog/DxcVSZQck6kKIVWDuhB2GhTleSg+K6Jm4\nsqMiZYQQYsZle0c9x7tkvJibdb+xMTwZmyrtK+tdCLGszLWwQ32WalynHYb/jgVeCCGWnbl0xRhW\n3bGKXMtbbhohxKow9eqOUO4WiaNZbEGOqnOgiG+H4cnWqggbIYRYVmYi7PG2YSK8vz+8IEevF7JO\n77uvvu56vGCHbVtNd4U9CiGWnblzxfjEJL/Y9Sj+cYU9CiFWkbkr22usr4/WTq5o55YjEEKIRWMm\nFntKfFP+9VOnimSkS5eKol+pJKS6tuOwR3PT1EXdyLoXQiwac+OKSQl1yv3StkvGlyGwc4UQYpGZ\nG2GHk0W+jDbE1q/K5C33w0PY2xt2u6QWAFFSkxBiUZgrYU8V+TKhPzqCJ55In5sjuCbo8QLZKYvd\nh1aOmtQkN44QYlbMlbAbPhrGhLbbhfPnh8UYmguor/64sxPeDDY2wqIe1p6VDh4HCbsQYlbkCPvH\ngHcDXwfelDh+GngEODto7zeA/XE7VpVQ1KQkQcqtsr4e2rd77OwMC3pc8x3qXTESciHEvJAj7A8D\nHwE+XnL8g8AXgV8liPxXCEL/naad8ZExlkFqkTAHB3DxYhDao6PyLFb7TkW9pNwqsSBXJUBViXfV\nQCL/vBBimuQI+6eBjYrjfw28ebD9CuAbjCDqkJdVCsHXfvZs+NvE3l/TRERPncpfMNv74VPX1PVf\nCCGmQRs+9oeAPwL+Cvgu4KdaaPMEZVUcY4s6jnqBaov57NniHGujiniytcxCF0KIWdGGsH8YuAp0\ngLuAPwQ2gW/FJ3adCnc6HToVCpiyhsvO8y4cqy+zvV0u6F6U+/3hwmHm+oHwXVaYzBbcloUuhGiT\nXq9Hb8y09zaE/ceAXx9sPwv8BfC9wBfiE7sNlC8l7CaysRj7KJpwn2qXzP5+Id5XrsDmZti2MgZ+\nkIDCJbO5Weyzb7+oR1X/hRAih9joveTdAZm0IexPA/8Y+CxwhiDqf95CuyeIy/BWjRNllnqnc/La\nw8Mg3o8/Ds8/X+w/c2a4GmS3CxcuFANHHHophBDzQI6wPwpsEyJengN2gTsGxx4E/g0hcuYaoajY\nrwD/q/WeZuIjYzw+SsYPgEdHYYFsS4iKiVdjit+Q/ALa/m8hhJgVOcL+3prjN4F/0kJfGlHlczdS\noYw+Ln1/H555JkTZ9PuFW8bi3OMkpbitFBJ2IcSsmcvM0xxyywjYt68/410y5mtP+eWrxNs+ioYR\nQswbCyvsufiSBPZtRb9M7ONIlypSiU9+YleJSUKIWbN0wh6HJ/oIl+3twr2ysVHsK4tsyUVhj0KI\neWLphL1MZONFsT110S2pBbRT96pCvnchxLRYOmGvoirevEp4m1jkmlQVQsyaF826A5OkbiJ0Ej5w\nibcQYtYstcWeE+FSViUyJ5xybw/uu6/8/m24cIQQoilrU7zX7du3b0/xds2IC4pV4ZOScks6NGlf\nCCGMtbU1aKjVS+2KmRQ5FSDHrOEjhBAjI2EfkOsa2dsLGaudTlFJstMJ+42UsDex7oUQYhyW2sfe\nhJzEJBPm69dDPHy/n7+wh0+SEkKISSJhb4DFu/ttHy65v1/UnLF9VnMmtwSChF8IMS4S9kzilZMg\nlPD1x0yUfaExO55TbmAUYddgIISIkbCPgAmphTqWiWvugtrjIGEXQsRI2CuIV2yyomH9/kkr3VO2\nfmq8OLa/B7QX6y6xF2K1kbBTLoSpFZv8ykkpF8v6Ohwfh+1USd9UjfimFn3dYCBhF2K1kbAzmhCW\nCXLclhdqE+R4MRBrrwnxGq9CCGFI2DOxha5zwiJj4fbWtdV+jzNXx60UqfIFQghjZYW9qRCae6XK\np57a710w29tFUlMqWanpoOGvGbUmvNw2QiwfKyvsbS+OUVZMzN+r2y0W9bDFs3Nj3ssW4x5XlCXs\nQiwfKyvsOeRa9VV/e+H02ac+5n1/v77cgB8c7Lqqc4UQq4uEnXYW2PDEGale2G3x7FxSg0tqEtbT\npDyC/PFCLB8SdtoVs5T/3AvpwUEQ96MjuHEj7Lt+vejDhQvDNd5Tg8u47pO23VBCiPlCwp5JUyvY\nL6Ld7w8vmB0LadPKj5O2quV3F2KxkbBnUid0ZZa1LxoG+QJeFQHTJmUTtRJ2IRYX1WOfICmXTEow\nrZhY1bVQLbaj1nqXgAuxfMhinwDmWrEQxoODsN9CFmOq1k216+KkpzYqQ8bXa0JViOVAwt4yXiAt\nGWl7ezgypqoaZEpc+/2TJYDbFltNqAqxPNQJ+8eAdwNfB95Uck4H+PfAHcDNwd8rSyqmPS4Wlirn\nm6oWmYpZt2qTtu0HgqMjOH++3MquGxDkWxdiOagT9oeBjwAfLzm+DvxH4B7ga8Dp9rq2HPikJCiq\nQ9a5UmzZPQii3e+Hz/ExXLtWxMN3OnD1Kuzulrt6/CIhOcJedo6EX4jFoE7YPw1sVBx/H/AJgqhD\nsNjFgNj9YsKb48P2pQUsXHJjI1x77tzwtXVrqZYlTNX1u6odIcT8Mq6P/Q0EF8yTwHcB/wH43XE7\ntSzELpnYci+r6R6LvcXA+0JiEEoRPPZYsNjNko8TneIBxVvlfsCxxUTK+iCEWBzGFfY7gB8E3gG8\nDPgccAR8NXVy1/kJOp0OnRVXjtwJSy/EVkjs5s0i6ubWrVBW+M47hwXd2vOTuLZgSNx2txvcOXEf\n9vbKFw5Z8f98QkyEXq9Hb9T45QHjCvtzBPfL3ww+fwxskiHsq4i3lnP93pCejIXgrul2YWur8N1b\n6KQXXr8kn5U0iM8p4/hYi3oIMU1io/eSt6gyGVfY/wvw28CLgZcCPwL85phtLj3eLRNPWOYIvU2k\nmkum34fLl9OJTlBY6T5U0u69s1P47i9ePFl9UgixeNQJ+6PANiHa5Tlgl+B+AXgQeBp4HHgK+H/A\nQ8CfTaSnC0pVbLrHC7v/22NuEUt62twM+7/ne+CFF8oTnVKTuF7c47BKvy9ez1UIMf/UCft7M9r4\njcFHJEi5RI6O4IknCqs7LhJWZrV7t0i/D48/Hs79kz8p7rWxEVwzcYVIv51b+92Q+0WIxUKZp1PE\nXCLdbkgk8jHtPjomXhe1rK2trdAOFBOkEKJkYqqqT0IYKEZdWFsIMV9I2KdImWB6Sz1eF9VCGWE4\nUenUqWD1m7BvbpZnqvqoGtvnq076Y3atd9+MiuLehZgNEvYpEces20Sl+a3jZKMy90ccoWKTopcv\np/3icZ2ZnFID8apPoyJhF2I2SNinRJOY9aYhrCkRjuvB7+3B4WE458qV4v7r6+GtwA8KNikbtzsJ\noZb4C9E+EvY5o0rk4mPr6+n6LlZnptcrfOiWzOQF/fg4fA4OwhtEvx9E/dq1EDq5vz983xwRblr+\nV8IuRPtI2GfAqLVa4v3Hx2lhtDozlqVqWHKSibwJNwz759fXCzdPHJZZRyzgiokXYvpI2GfAJIUu\n9uWfOxdE/MqVMDELwRq3Cdn4OgvB7HTC9vXrhbib9W+/Ied3lC0K0nRRD1n2QuQjYV8wer1gafvM\nU7PAzVL3cfPeejeBPj4eFtSbN4cTlq5eDduXLsGZM+lkphzK5gviPua0KWEXIh8J+4KRCk0sq8Hu\nI2JguBSwbysuCra1FbbPnQsWu28vhyYWeVNXjxCiHgn7kmJWvJES1TIr2AaAXi/42+3avb38omW5\nmau+jzFah1WI0ZCwLzBlUSbeVeMrOXrW10P4Y7zwhw9/NFfPzg7cfTecPRvcOLkDRBlesKuqTaq0\ngRCjIWFfYMomJU3UbbK0LOs0Xmw7bs/88imrumxVpljkU4IvwRZiskjYl4x4TdVYNFOiWrYG6+Fh\n2D44CIt43LpVHDs+Dq6ZOuu9qcukLAtWlSWFyEfCviRYxIpFs9hqSaOucdrpBNfL2bNBVP0bQLcb\n3D2pZf1GiXuPf0dK2OPwTCFEORL2JSGOcMmZYKyrB3P2bBBwGzCOjkLd914vWPG7u+G8e+4prmm6\nQlOuRV82YDRZiSq+pukxIRYFCfuCk4ocyRWnssnXeGLTRPXs2eGQSR9uOc7yefFviCtY+tIIVrc+\n9ulL2IUokLAvOCmf+TjiFIumX1WprN1xhbBqMtXPAfgSCEKIciTsS8g41mgcNbOzU0TGHBwU5/nJ\nzHEmS3P6YqGXcQGzfj/06fAwFC7zK1HlljFYXw8JWYqXF2I0bovJ8uST1cd3d5u1d/HieNdX8eST\nRX99v+PfYH9vbqbP2d0tPrn4c+Pr2vyNQrQBcLup2L5oAgIuZkTb1mVVVui4eAva15OJf4P97d8Q\nmtarr0IlDcQyIlfMkjNOWn6ZyObet+78OlGN3ULmijk6Gk7Gsnry8TJ/ZX2wssTQPIpHiEVgbYr3\nGrxViFkRR69M+16pypRnzsDzzwdhNv84nBTYnZ1CgC9dKkItq6pO7uzUT7iWVaAUYl5YW1uDhlot\ni11Mjbgypbeo40iYWNj96k+xkJcJc1Xse+pNQBa7WBYk7CvEpEUrdo/4+8a1anxsuu2v6p93n1SJ\ncY5o++3cevBCLBIS9hWi7VDEMh98rxdcJSnB9IJ6fBzE9+LF4DLZ3w+WeWouYGurmRjnToqOOkGs\nRCYxz0jYxUiUCVuOv9qXP7hwIe2LT+2z+5WJsV/Yu9sNtW3qBoBRxbltYddAIdpE4Y6iFUyMrSiY\nWdVVQj9OdmyKeKI0xxqfxuRpTvuawBVtIotdZFOVvWnVF60KZFkGqCflG7eImZ2d4ljKPZPCFgix\n/sUhkFVvGTmDzKiho3XtS9RF29QJ+8eAdwNfB95Ucd4PA58Dfgr4ZDtdE/NGkwUycickU7XhLfQx\nbqNsbVcTxmvXgmsHCr99fG58zyY0+f05lE0kKzpHjEudsD8MfAT4eMU5Lwb+HfA4042LFzOkanJy\nFFHyVm3Zqk3xeaPgq0FOsj5MTvt+PsCXXBZiXOqE/dPARs05vwg8RrDaxYozijCaCO7shDDJM2eC\nBZsq6pVaoamsGqVvG4oFSFJCGsfFj1tHvs66r+qXrHUxLuP62F8NvAf4RwRhV2rpilBmUTepie6/\nzRVx7hxcvz7sp89t1/rkhdmL/rlzxTqvVX7tqjVc2yIW8PV1uO++9u8jVpNxhX0PuJ8g6GvUuGK6\nzmzpdDp0ZJosFE3dF3ULWsQZpt1uUSbYW9Q7O3n3tQnXMqve3griDNiq9VTHnXC1+9ehpf+E0ev1\n6I05oz6usP8Q8HuD7dPAO4FvA7+fOrkrJ+JC03TysEm0SepeRjyRWnbf1L329oYjdux+flA5Pi4f\ntNqo/pjjtlFkjDBio/eS/x8yk3GF/e+77YeBP6BE1IWoEs+trSC2Fu5oa63GC3yUCa0X8LhtHx3T\n6RRZrhZzb4t028Dl7+GrP/qwzjYmXCcxgatEJwH1wv4osE2wxp8DdoE7BscenGC/xJxT5WKJxarf\nP+kCgWFXSyoO3q6JC4DZt2/r+Dht1cdL+0EQ683N0K/NzZN+99w3hHFfQCcVPilhF3XC/t4Gbf30\nOB0Ri0WZeKTEKlWt0dd7OTgoyvBevFhcV3XfKgHzbXtL+MKF4QlKE/NUWQO7xyxpKtK5Ga4S/uVH\nmadiJpg1bngLO3a3NHk7gODW8QIel/u1hKBeLx3d0+kUrp04O9ZH27Qpkm1kxOYkOo0yWGggWDwk\n7GIiWHmAbjftP66y+GPL01vpZT5pyCv45ZOTzLe/tXXy3MPDk/2MhT0enOJzmjCOT92HjMbzBeMi\nYV9MJOxiIjTxH+dmgla1WeZfTglvrxeE+9q1IgQy7ke/X16fxu6Tmsgd1yJuOqGa2l+XDFXXplh8\nJOxi5qSSkJpOJHrXSCpByZ9nomzx6+ae2dsLgg8hSercufD32bNhn1WttPZzV1+qi+eP5w3GfQ6p\nfU3DVDUQLDYSdjFxJpW5Gf9dFjHj6fXg8uWwIPatW2Hfxga88ALcf/9woTCfJOX3X7067Ju3kEoL\nn/R9auKLH9ftMa7w+oGlzWgdMX0k7GLiNBWbnPOr/O5lbhnbf/48PPFECHU8Pi5EOU7p99fbItxQ\nVJG88848V0hMbjLUJJ9biv19WeTLgoRdzB1NxCWeEIUi4iV2bcQibG6VlKV94cLwuRZXD0XC061b\n6ftAvTvDzi1LhhrVtVNH1XWpOQMJ/WIiYRcLTew+iBOTqtjZGQ6N9NfYYOHXYTWhvn69SG6KBThn\noteTmwzlaTNSJR4Q4zmDWQq7InJG50Wz7oAQ4+BDI1OhkmWYYFRVVOx0CuHd3S1E0BYCia16366n\n7dr1KWadTDUJlvE3TQtZ7GJpMIs919fsz626JhYY7y5pk7o+pFw7Vu43FRGUe08/X6GJ0uVAwi4W\njrpEpXEjQ/y3uWA2Nops1H6/KEyW8q/v7QUXj+2Ls1tTIY51/S5z7aRi1kf9/TmLf08ahVq2g4Rd\nLBwpv3bZJGhZZAycjEixdlJWfyygqTVZe70QCnl4GL7LBKot37H9Fp/dO04hML+AeO79J/HWolDL\n8ZGwi4WniRjEESnxpKsR119JLdMX4weNsklRX2oBRrNK19eL669cCfc6d64oeVwWqVNHU5HW5Ob8\nImEXC804wnL16nDN9zIRN7GOo2Y8e3vwwAPBVfO3f1vUlL/77pC5Gg8UvpqlZb6m3DMp4glfq6Bp\nv2VaVq5NCu/tTWZZPw0aoyNhFwtNTlRKGVtbhRDat308qegML8KXL8PTT8Pzz4d9L31pManpJzb9\nhG3sI/fJQT7LNZf19frM2zYsbP9WYnH3DzwwPKdQde0k3yBEgYRdLBU54mJZpH6pvLhEr+0314lN\nlvpFQ+x+tt3tBv/6+nrYlyodbCs2mcDHiUopca4SRNtvy/tV+fDbEPZUFE1uxuokXDdyB6WRsIuV\nIhZi+04JhD/XW/RVHB8XiU+pdiAI/9ZWcAWZINsqTv1+SIDy4p8j7P7vuvjvuL0m4miD3dFR4d+/\nfr34vTs7k3HLVPVHwn4SCbsQNAuVLBOTTicIXp2wmZibL9/Y2Egv11eW4JR6+3jggeASsrZT67TG\nkT9NxNHPO9y4EZK0Ll0K34eHJ2vbTyp8UYJejYRdrCxNrFS/aEjsSvHi9cQT1aV8vUVtbwG+dMH1\n6+Hb1mS1id1UDZnU2wcML+AdC2DZUoVN6XSKJQVt4fG6tx7fz3GySv3ANOlY90UdQCTsYmVpYqVC\n+Ee+vX1y8euccEsT76tXQ3VIw6zcu+8eFmyfBerDI3PwFr73u/uY935/eODwv7XJc+n14CUvCe2W\n1adPiWOOYFa9GaUmoHOvb4KEXYglpiyqpen1UB7fXhYNU1Vrxq61cw4OgrUfW+W+7/F9q+Lz6/z7\nvR687W31/YuFPofUXEDs1kktYVh2/SohYReiAXWTk6MKSTxJ6tupS/VPvVFY+GO8uHXMuCGSqSiZ\nKvxi4qk2qki9GVmSWRl1kUJl1yx6WQMJuxANqPuHnStQRkpEer3hSc+cOu1lbxSpSJ5RImnqMMvf\n/10mjvZd5rKqu86TqsPjrzdBrxsAPDmutXlHwi5EQ0a12srS/GNfdKr9XHHJFWmzdlMZtzdvwunT\n4byqSpKeXq88xNPX4en1ymP54+usn37bRwr5Gjn+2rhvOWGqKcYd7KranbTlL2EXYkSa/gNtOmE4\nqgDEQpk67kXXXCixlR+XUvB++VS79oYR42vypPrQlDqL2gaQw8MwUR0PXLn3npQAS9iFmGPaTtE3\nv3OcmNR08EitBFXnh05hIZhV1/nMXG+F270htDHKs+p0ilLJkO/v9pa7PY+dnSIEdBrMeuJWwi7E\nBKmbiEtNfHp3g507bh/K/PoWonjqVBA/WwbQzvvxH69226yvh+2470dHxT0ODoat9iaTmcfH1UsH\n1k3sWh/qIouMMh+9uaFyyYnosT5OYgCQsAvRgKb/QHMm4vzEp7eQmyxuPU6ffL98W9evF752XzAt\n9ll7F0xq0tZEzl+TisZJiWEduZPZuYuItOmjr2p30hOyOcL+MeDdwNeBNyWOvx/4FWAN+BbwL4Cn\n2uqgEPPEJP+BetdBk/ZTfTKBtslKw0+A2seHHq6vp+9x9Wp11Ms994QKl97qtySop58urPw4Kcpn\nkcbJVLl9qwrVtE88SPrvqut8obgcl9a8hEnmCPvDwEeAj5cc/3PgHwK3gHuB3wHOt9I7IZaIpqGS\nue6DsrZSA1CcGWpi5AcQK8Hr93e7hbDGLpVeD154IdSoMfG8556wvbERxHF7u2in0xm2/qvcRHZd\n3RtL07elqrIKdt3eXrEvdwETfyw1KO/tTUfkc4T908BGxfHPue3PA68Zp0NCLAqjTAa23WbZ9bk+\nbO+Hjq87dWq43PDmZlHe2Aucd3dUxc/7+9mgYX5sKE+mmlUc+X33Dfv46/phv6PquR8eTqf6Zds+\n9p8FPtVym0LMJZO2vHL9wp6UsPul9Orivv119947LEKxBW/7bGL13Ll0bLrdywtfqrxBPBhUDU5N\nE5n8gJJ7Tdm+uB+5100qLj5Fm8L+duBngLe22KYQK8Uk/LSppfQ8Ze1aeWHvb/aFyWJrPVWszIQ/\n9vXnvk1UnVc2F1HmZhnFNVPnj9/fTy/c4gdlW+Dcau1bWxcu5K08NQptCfubgYcIPvZvlp3UdU+y\n0+nQmfaMghBzzjiTs6MOCmXXeR+/+cq928TE3M41YfNJQHW/xx9LWdll5K7aZIwbV172VuPdTzZx\nHA8Y3W4YXP2ktj8e96vX69Eb07xvQ9hfB3wS+ADwTNWJ3Vk5y4SYEyaZuJIzKNRZsl7MfSTJzk5h\neVvbsdskzlTN7XNV32Ls+fl+1rlZUu1akbRRBsLLl9PXxZFDEOL5vRvM3FVVcfHB4O38XR8u+Ztk\nkiPsjwLbwGngOWAXuGNw7EHg14DvBj462Pdt4C2NeyLECtBE2CcxANS1mbI4/bGqdsuMzFEyaONr\nbdssXO9q1QtjAAAItklEQVQWMpH0fagK+ex0hte4rYpiSfHCC+lnFNfEgbDwyvlBjODFi8H1srUV\nrjk+Lh9QfOjnKOQI+3trjv/c4COEaJG2XAdtXbezM3xeHHfuF/22c+IImqbkDIRXrxbbvgxC1RtM\nKkGqrh/eDWVtfcqFiqQSyuIJZ7PSbQI57lfqvqOgzFMhJswsElfGEfYyIY6P+QiYlM+5KXWJRlCU\nLThzZrgfRr8/PACVtePLJqQyXlPRRSle+cqTBdPKJpxv3jwZneTvX9a/UZCwCzFhxpkQnQVVFvYo\nFm7uYGYTonVhid0uPPLIsJvFC38qkcgPQL1eiFLx5/t7pIQ9dpOkwj7j32hvMak3hbhPRjzB6mP7\nmyBhF2KJmOTkLAwLVVVM/CiDmZ8QLfPzm+C97W3Dvm1zgZQJYTxJ62P7c/pWlbhlpRIuXDh5TWoh\nkFSb/u/4TWSEuVMJuxDTZNIRvqMKe52VbMcODmB3N2xfvNis2mLdfc1t4SdJ44W5bWDx/uwbN4qS\nvN6ariJn5SXrkwm39en++0/OLZiLKE7QKnse/vqySKVxkLALMUUmLeyjUmdhp47lupTa+s3exeHv\n733qTerHWxx+7LZJhW5abL71I+XqKYvE2d8Pg88LLxRvFPv7IVnpVa8qLPt44tm+ZbELsYLMYnLW\n2q+j7g3Ci7VZ3d5PXtamx4t52cSpv1e3m37bSN0n5evudEJfL1+Gs2dPLjLi200NmOY6ssnVn//5\nsK+NNyBDwi7EgtP25GydEOecZ1RFusRtmSiWJUn531lV/yaHuHJmLOKWkHXqVIhFtyqT6+vDE6MW\nXlk2iJZlyHa7xSDk22sLCbsQYohcYR+HlOCb0KWSpEx49/eH/fw+Uqbufv6tZnPzZLhm2eAYT2ha\nYpFlkVpiUkynczK7tNMJln63G+Lvr10r7nXqVPVvaIKEXYglYh58+KO6hqomHL2w7u6myxnktN3r\nFcv4XbhQ9LPuudn9T50K15w7F3zkvV4Q+r29QsT97/fC7d1ONiDcuBGOWe35suqYTZGwC7FEzIOw\nl1m/PmKlSvCbunvi2PO6fnlfeOo6H7ECw0lGV66EgaXXK+rN5GBuG3vjODwMon/qVGjzzJmisFru\nYFOFhF0IMRVGmQuIRfbSpWKSMVUjJkXKn19VUiD15hDXjM8NqzS8H/7GjVA/Zn19eJWoe+8tJpDH\nRcIuhJgY475BpLJgfZhjzhtASthz7llFp3MyISl1/8PDkwOAFRHzJXvrBqemSNiFEBOjKkEnhyqR\nrapEWddmHGtufUq5dcq2c5a429o6uc9H+VSVbxgHCbsQYupUWdRVlIU5VoUv1i0HWDcgNPX5eyzS\nx/fHh3Ba+OQobVchYRdCzJRUzZSUKKcWp0gV4YLm4j2pJK+UVd7UPz8KEnYhxFxRJspVYjiq+Hpr\nuc0kL8+kC7OlkLALIaZOGxZy7iToJBKumop1WfjkJPzrIGEXQsyAXAu5qnRAG8Je1rc66oS9rlpm\nWQx9W0jYhRBzS5lPPaZNkWyjnaqBqyqGvi0k7EKI1mkiWm2Im91vkpUu22p7Gv52CbsQonUmIey5\n5QUmNQk6atupiJ9UDH187jhI2IUQC0EserOqQ9+UOA4+HhxsEPSZqOMiYRdCtMK0hbbOgm4atTIr\nd8okfO0SdiFEK0zSDTIK8y7sfk6g2213MJSwCyEWnnlyveQQx7V3OqHKY1tvNxJ2IUTrTFtoR7nf\nLH30qbcb+diFEHPNIljQi+w6quNF7TUlhBCiKZOo7pgj7B8Dngf+tOKc3wK+ClwDfqCFfgkhxNSY\n5RvGJO6dI+wPA/dWHH8X8HrgDcAvAB9toV9LTa/NpVIWHD2LAj2Lgmk/i0VwHTUhR9g/DXyz4vhP\nAAeD7c8D68CZMfu11OgfcIGeRYGeRYGexXi04WN/NfCc+/trwGtaaFcIIcQItDV5uhb9fbuldoUQ\nQjQkFuQyNoA/AN6UOPYA0AN+b/D308A2YcLV8wxwV+MeCiHEavMsYR6zdTYoj4p5F/CpwfZ54GgS\nHRBCCNEejwJ/Bfxfgi/9Z4B/PvgYv02wyK8BPzjtDgohhBBCCCEyuZfgY/8q8K9LzlmVZKa6Z/F+\nwjN4Cvgs8ObpdW3q5Px/AfDDwHeAfzqNTs2InGfRAb4IfIkwf7Ws1D2L08DjwFXCs9iZWs+my1wn\ngb6Y4I7ZAO4g/Md4Y3SO98f/CMvrj895Fj8K3DnYvpfVfhZ23h8B/xX4yWl1bsrkPIt14MsUIcOn\np9W5KZPzLLrAvx1snwa+wXLWt/oHBLHOmcfM0s02a8W8hfAfqg98mxAl857onFVJZsp5Fp8Dbg22\nP8/yxv7nPAuAXwQeA/7n1Ho2fXKexfuATxDyQQBuTqtzUybnWfw18IrB9isIwv6dKfVvmrSeBNqm\nsKcSlV6dcc4yClrOs/D8LMWIvGzk/n/xHopyFMuaB5HzLN4AvBJ4EvgC8M+m07Wpk/MsHgK+nxC8\ncQ34pel0be5orJttvtbk/mNchWSmJr/p7YRIo7dOqC+zJudZ7AH3D85dIz+/YtHIeRZ3ECLL3gG8\njPBmd0Twry4TOc/iwwQXTYeQA/OHwCbwrcl1a25ppJttCvtfAq91f7+W4nWy7JzXDPYtGznPAsKE\n6UMEH3vVq9gik/Msfogiwe008E7C6/nvT7x30yXnWTxHcL/8zeDzxwQxWzZhz3kWPwb8+mD7WeAv\ngO8lvMmsEjPVzZcQHv4G8Peonzxd5mSmnGfxOoKP8fxUezZ9cp6F52GWNyom51ncDfw3wuTiywgT\nat83vS5OjZxn8ZvA7mD7DEH4Xzml/k2bDeY4CfSdwFcIgvWrg32rmsxU9yz+E2Ey6IuDz3+fdgen\nSM7/F8YyCzvkPYtfJkTG/Cnwoan2brrUPYvThFIm1wjP4n3T7uCUUBKoEEIIIYQQQgghhBBCCCGE\nEEIIIYQQQgghhBBCCCGEEGJ++P/y3s+HPFBaYgAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "sig = 0.05\n", "err = sig * np.random.standard_normal(N)\n", "yvals = yvals0 + err\n", "plt.plot(xvals, yvals, '+')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##Fitting the curve\n", "We now import the relvant module, and we now choose the functional relationship $func$ between our $x$ and $y$ that is parameterised by three parameters $a,b,c$ of the curve we want to fit. We here use the same relationship as above, but obviously we can fit any curve we want (eg, linear, quadratic etc). We also have to provide an initial guess for the paramters, here $a_i,b_i,c_i$." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# import the curve fitting module and standard imports\n", "import matplotlib.pyplot as plt\n", "from scipy.optimize import curve_fit \n", "\n", "# choose the function to be fitted... \n", "def func (x,a,b,c):\n", " return a*exp(-b*x)+c\n", "\n", "# ...and provide initial estimates for the parameters\n", "a0,b0,c0 = (0.5,0.5,0.5)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We now execute the curve fit and plot the results. The coefficients will be in the tuple $coeffs=(a,b,c)$ and $fiterr$ will contain an error estimate. If $fiterr$ is NaN (or the algorithm terminates with too many iterations) this means that the algorithm did not converge, which means we need to adapt the initial values (see below)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "a=1.00755718879, b=1.94018195497, c=0.985806031536\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEACAYAAACnJV25AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnX18XGWZ97+TJm1oCz3pC7BIO9MtUIqQFCgQ2tQZRR5L\n1bXqAw9Y16YKXWEViwsL+rhmsq5r2dWH4BsIStNdt7Ct2ChCLRY6aYMUN5WkokARmVgqoS9kUlra\nvM7zxz1nzj1nzsycmUwmyeT6fj7zmZlz7nOfe06b37nOdV/XdYMgCIIgCIIgCIIgCIIgCIIgCIIg\nCIIgCIIgCIIgCMKYZjawA/g98AJwq0OblUA7sBd4Bqgs2OgEQRCErDkTWBj7PBV4GVhga3MlMC32\neRmwuzBDEwRBEPJBE3BVmv0VwOsFGosgCIIwRHxAB8pyT8XtwAMFGY0gCIIwJKYCrcCKNG3eC/wB\nZbULgiAII0SpizZlwKPAj1GuGCcqgQdRPvYupwbz5s2Lvvrqq7mMURAEYTzzKnBONgeUZNjvAX6E\nssQbUrSZA/wU+CTwx5Qje/VVotGovKJR6urqRnwMo+Ul10KuhVyL9C9gXjaiDpkt9iUxwd4LPB/b\n9uWYmAP8APgqyv1yX2xbH3B5tgMRBEEQ8kMmYW8hs1V/Y+wlCIIgjAIyibYwDAQCgZEewqhBroWF\nXAsLuRZDw1PAc0Vj/iJBEATBJR6PB7LUarHYBUEQigwRdkEQhCJDhF0QBKHIEGEXBEEoMkTYBUEQ\nigwRdkEQhCJDhF0QBKHIEGEXBEEoMkTYBUEQigwRdkEQhCJDhF0QBKHIEGEXBEEoMkTYBUEQigwR\ndkEQhCJDhF0QBKHIGHFhD4VGegSCIAjFRcGF3S7kIuyCIAj5ZcSFXRAEQcgvmRazHhZCIUvg6+ut\n7YGAegmCIAi5U1BhDwaThdzcLgiCIOSHggu7/g7imhEEQcg3Ix4VI64XQRCE/FJwYbcLuQi7IAhC\nfvEU8FzRaDRawNMJgiCMfTweD2Sp1SPuihEEQRDyiwi7IAhCkZFJ2GcDO4DfAy8At6Zo923gFaAd\nuDhvoxMEQRCyJpOw9wG3Ae8GqoG/BxbY2iwHzgHOBdYA9w11UBICKQiCkDuZhL0TaIt9Pga8CJxl\na/M3wIbY5+cAAzjDqTO3gi3CLgiCkDvZ+Nh9KDfLc7bt7wL2a99fB8526kAEWxAEYfhxm3k6FfgJ\n8AWU5W7HHorjGNcYCqmsU6eaMGb9mHAYNmywtkv9GEEQhOxwI+xlwKPAj4Emh/0HUJOsJmfHtiVx\nWfNknj25hlDIYMWKAAsXBuKibQp4MAh1dVI/RhCE8UkoFCI0RPdGpqB3D8p/fgQ1ierEcuBzsfdq\noCH2bkelJ519NuxXnptgUL1CoeSCYHZh19sIgiCMF3JJUMpksS8BPgnsBZ6PbfsyMCf2+QfAEyhR\n/yNwHFidqrNBSiiZNw8iETCM+PbGxsQyvh/4QKLLxnTTiLALgiBkJpOwt+BugvVzbk5WwiA0NzNY\nUUGH109DRxNgEA5DbW2icGdTAVJEXxAEwWJEFtooAeZ2NPPYjNU8xRaamy2rPRwGn0+1a2hQxn0o\nBM3N1vG6Ja+/C4IgCIUW9okTobc3/nXGzGjcMg+HrWb19eq7KdjNzeD3q316lEwqQRehFwRhPFNY\nYd+3j7558ykb6IHp05n35m6YPp2/rbiUzddu5q6g8rv7fNbEqin85jYT0+/utCqTCLsgCOOZwq6g\ntN5Lw0AnW2atYV75X5iz/xkA5nVtZ+7dawiWbyIchrZYrqsp2KEQrF1rfTZ97s3NSsD9/mRLXhAE\nYbxScB/7Qr/B+5o3se+c5QnbP85PKQ0FoKmJUJuREP4YDlvCbk9Y0q16810WyBYEYTwzIpOnq1bB\nuQ0bVSjM449Dfz+lDCgTfM0aAps2JbhaNmywJlSdhNrcliqqRhAEYTxR8MWs464Uw4CmJpg1Cw4f\nVg0qKuAvf4Hlywls3EgwaCQlLunoWauCIAiCouAWe5IIt7bSfeFipr1nIbS2wjPK787q1QS2bMmu\nrwzbBUEQxgMjsoJSgvB6vdzzDweUS6avz9re1MSe77QktydxAjVj/4IgCOOMEfGxQ6I4m5Odf2tc\nyryu7fE2l9y6FPztBAKVSceCCLggCIITIybs9snOQAA2n9zMnXdXxKvdeIBoVRWeXbugpibheD2h\nSRAEQbAYMWHXMbNMKTdYwi6eYWmiuC9dyv3/2sVzL6u6MmZ5AbP8gL3OjI4kKwmCMN7IqhTkEIlG\no47rbySIb20teF7Yy0N7qhIGN4iHkvY2qKxMjK7JgB7nLgiCMNbIpWzviEye2tEtap8PVn2zEs+u\nXQxqbUqIEq2q4oefbokXBQsG1Y1AMk0FQRAsRoUrRsfyvdfwuRu6+O7D0ymJrbTnAW5cv5TzH+vi\nhz4jbo3bi4I5Tcw6lQYWBEEoRkalsIdCqoxvy28MqmhjL4lumfO+tQaff1PSsaawpyo7IKIuCMJ4\nYNQJOyQKczBYyY/+vIvPrI9NqE6axPTfbOOLz87inkgr9fd648dJpIwgCMIICbubSBW9zetzarj7\nzi7u+tMa2LaN0qNHOQ1Ye6+Pwx9uB1Sc+4YNqr3PZx0bDjuX9hXrXRCEYmVUCzvo7QwIbFK1ZWJ4\ngH95rArPUb+qO4Oq565HwdgFXhAEodgZFVExqXD0mbe2JrTxgAqRqa1N25e4aQRBGC8UvLojpHaL\n2KNZ/H5r0jPe5jUvf/hsOzffnzih2r3tWQ5OjnDfw4Zj32bZX0EQhGJnRBKU0iUNmRExZoZpXZ3a\nbhjWYhsAdHRYy+yZzJrFK9MWce7/bATDSLpR1NVJ2KMgCGOLXBKURl1UjBnuaE6A2ldHiuP18q3b\nO7lrwwLo7ISpU+HQIc49tBXOOAP27SMQ8ErYoyAI446RL9ubAsPI3KZ6mQEvvgjXXgvV1daO3t6k\nomF23GSrSkarIAhjkVFRKwac/evl5XDmmcp6r69XS+qZlnzSzSESUZZ6b6/6Pm0aXHYZbN5MqM2g\nsdHqp67OWnYv3U1G6swIgjDSjGlXjJNYO4lqSqE1DNi3D847T4l7dzds3x5fQzWVS0YvQ2COQxAE\nYSwzaoQdLKvdLPJl4lpsvV646irYulV9Ly2FJ5+Eq6+m5Qub6Z9qxC11MzKnqQkaGhLdLvY6M+YY\nRPQFQRgLjBpXjI4ZGdPYaAn97t2wbZsVJWMKbVKyUySiwl5aWuDIkfjmFy64lgt/vyne3m61O03S\n6pZ8tkgdeEEQ8sGYLdtrJxCwEopMEa6uVqKeFNcesh1sGMoMv/zyhM3v/sNmKCkhUNoSv1nU1qon\nA9P3Xl+v+s1HKWCZeBUEYaRwI+wPAW8Cv0uxfybwS6ANeAGozcfA0iUUuRHNlls2cqhmBX+Y/xGi\nxG530SjRpUt58dkIPp91jtpaNVm7apW6cfh8ljvItOwznVOEXBCE0YIbH/t64DvAf6TY/zngeeBL\nKJF/Gfgx0J/tYPTIGLOgl7n83YYNlvDu3p06i9V8395qULNrC7MASkogatV0v/nbC7j5xRfjSUz2\ntVftE636+NItwae/i39eEISRwo2w7wJ8afa/gVleEU4DjpCDqEN6UQXLwt62zQqDNMVePyZJRHfu\nJLrUWkeVzk44/XTw+5kxYTOhkOEovPb6MnrkjJPIZxq/IAhCIchHVMyDwNPAX4BTgevy0GcSPl+y\nUNrjzE2LP7lMbw2lj3VRc1MsSxWgrw+2b+dH5Qu4/7IXCQSMeB/psIdHprLQBUEQRop8CPuXUf71\nADAP+BVQBbxtbxjUVDgQCBBIo4BO1nCqdroLp7lZbfP77RZ0LEt1gSbugHGyk7V3n8E9J/fRbXjj\nrh+wfO1OhcnCYbHQBUHIP6FQiNAQJ+3chtD4gMeAixz2PQF8HXgm9v0p4E6g1dbOdbhjKkyRDYeV\nz90e+giWFZ8yazQSYU9lLee99SynHj9oDQ748R3tPHWwMu6CMYuQmd/b2qC9XYl7c3NiJqweMikI\ngpAvRirz9CXg/ShhPwOYD/wpD/0mYS/Dm05InSx85UYxeOzTTQTXRqCiIr7PA6z89yq+Qpg/o5bb\nO+OMxGqQwSCsWJF445B4dUEQRhtuhP1hwI+KeNkP1AFlsX0/AP4VFTnTjgqf/EfgrbyP1CV6ZIyO\nXhtG+cMN7iivYMrJrnibEiCMj/AcP/99QxMny42MCUv6YiD6d0EQhJHCjbDfkGH/YeDDeRhLVqTz\nuZs4hTKaQh0IwFdfeJ47f3o5M6MH4wH9HmDun5t53yNreKdxU8q+04U9irALgjCSjKpaMdngRjz1\nyBW9/owpvqde6OX7F77JtEgHt317bjzWfdBTwsVHfkXZ16+GhZuTagjrNW0kGkYQhNHGmBV2t5hu\nEtMNEwxaRb9Msff7vXDra6x6YDHTJ71DSSRCybGIqg4ZW7QDr/K7O5X71SNzJDFJEISRpuiE3R6e\nqIdB+v3Wmtc+n7XN54OLV3jZu+IAgX9bblWHBFUCePFiOHAg5Tkl7FEQhNFE0Ql7KpG1Fw/TSYhu\nWbhRqf/Pfhbff3jOQr4b6ydXi1x874IgFIqiE/Z0OAlrUjSLWR0yEFAm/UUXMfPwPoINBpSVMe0L\nrdwW9GZ1joT+BUEQhplRWbY3X2SKYklrcTc1qbVUd+6EQ4fUikyHD/OFe+fC1Veruu8ZzikIgjAS\njMqFNgqF3R/vlMkKwKxZcPhw4sEzZ3LfF1/h5i+lXnXbdf+CIAgpyCXzdFwLu07ahas7OmD+fOjp\nSdh8pOxMZhx8MSkcMuv+BUEQUlA0KyiNOrxeVTisrCxh84y+ThUO2dGRsF231AVBEAqNCHuMjK4R\nw4BXXuHt087isGemtb23l56587nvG5bP3UnYh7J+qiAIQjaIsMfIJOyhEATXe/nWbQc4J/oKfSUT\n4/smRXu4ue70jJOqIuyCIBSCcRXuOFSstU8Nylbvo2fufCZFY3732MIdO85bQ/0hVWMmFFLJT2Z1\nyExISKQgCPlAhN0l9pWT8Hp56Oud3PxtbeGOSZN4b8823p48iwdvauXiFVYZAjflBnIRdrkZCIJg\nR4Q9B0whvflLBtz8In+oruWCBR54+mk4epSpwNp7fYSb/Mxta4ovuwf5j4wRYRcEwY4IexrsKzaZ\nRcPCYb0Ou8Gm65uUYM+aFT/WA8ztaFaLZr/yCni9SYtj6+eA/BUQE7EXhPGNCDuphdBpxSZ95SS7\ni2Xayla+cO9cStDi9fv6GPTNJXLJVUyev5lQyEjKiM22gFimm4EIuyCMb0TYyU0InQXZy7PVr3Fl\nbWIyUwlRpv92O3e/vIDvzXkxwTVjCnS259dvAJL4JAiCjgi7S8zk0kwCvO0lL1d2dsLKlfRte4qy\nAUvgTz3eyZ13V/DDg7s451M1CZOxQ60UORwuHUEQxibjVtizFUIzPD1dhcg4hgGPP85zv4hQ86PV\n9Oz4NZO6DwLK9/6Z9UtpOC0MeNP34zBmext9jddcasKL20YQio9xK+z5XhzDPD7xJmHAli1MikSg\noiK+3YOKmmlr9NPQ3UQ3RsaY9+TFuK3fMRRE2AWh+Bi3wu4Gt1Z9uu9KOA3YtYvo0qXxSj4e4OLu\nZlouWMOHjm/KmJWqL/EH6W9EItSCML4RYSe1EObDvRH/XFODJxxm0KdFzUycyPkvbmFftAwCS1QN\n+BQLZ4N1czG3pRt3pvGJP14QihcRdvIrZk41YSwh9bKe19g7dTEvli/k3UdaODV6VDVqbqZr1rk8\n8s+JNd6dbi5DdZ/IGq2CUNyIsLskWytYX0Q7HFb+c9WHl2nBA1RDbAEPq4+K/sOqRMHN6Wu8D7dV\nLX53QRjbiLC7JJPQpbKsTYE3reIEa761FebOBX0Bks5OOP103qr0M3375gSBHw6xTTVRK8IuCGMX\nKds7jDi5ZBIE0+uF117j7dPOgplajfe+Pqbv2Q4LFiSUAU4ntrmWBBYBF4TiQyz2YcD0s5shjBs2\nqO1myGICXi+ndh9QAr5AqxQJ6vOaNYRu2ZQgwKkSlIYi0jKhKgjFgwh7ntEFsrlZiaLfn5hI5CjC\nhkHLgy8y/fZa5rz+LFOPH+TAWYs49PIpnHl9AC6ZDBs3EmozhsVVIhOqglA8ZBL2h4APAgeBi1K0\nCQD3AGWoqcBAnsY2JnGKabcXC9OTjUxCIQh8yIAPNUEkwgtL1nDhMw/wrhUrYG8zbAVmzWLRRIPH\nrmkFvElW9u7dUF2d2srOdEMQ37ogFAeZfOzrgWVp9hvA94APAxcC/ztP4yoazKSiYBDq6ixBd3Kl\nxDEMvnnZJoINBq8cmAzAACXQ38/Udw7zL4/OZ91dkbgQG4bq++TJ5BuG3ncmP3ymgmSytJ8gjA0y\nCfsuoCvN/k8AjwKvx74fTtN23GF3v5iWen19YvlfJ2pr1f5z/2cjv/Fey8lTrJIEp9DDP949neWz\n9xIIqLymdOjncCPOIuyCMLYZqo/9XJQLZgdwKnAv8J9DHVSxYHfJ6NZ0qpruZtv4sYbBE7WbGDi/\ng0WfnB+vFllClMturKL/pjJ6p7TSFKmkvd06bsUKWLs2+YZiLRCSeMMxFxNxHIMgCGOKoQp7GXAJ\ncBUwGXgW2A284tQ4qM3IBQIBAuNcOdxOWAYCcGXAC8s6GayYHi9J4AFKo308c6yKoy9UUHPq86xY\n4Y0LutmfPolrLhhiH4PuKtJpaLAiLiVaRhCGn1AoRGiIj8dDFfb9KPfLidhrJ1CFC2Efj+jWsr44\ntpvkJwAMg9YftrHoxospYTC+3wNMG+hi79s+PO9rByoThFdfkm/DBiXuZr+Zzh2JyKIeglBI7EZv\nvW5RuWSoCUo/A2qACSiL/QrgD0Pss+jR68nYJywz3ajfmVfJ5244wnPTrtakXeEBqKqCvXsTtuvL\n+tXVJVaJNH359fXqBpDO7y8Iwtggk8X+MOAHZqKs8zqU+wXgB8BLwC+BvcAg8CAi7AmkSvyxL2xt\nF3YnS9p0i5x+nkH1w0+ybEEHP3/5PMoGexMbVlVBWZkqWVBZ6TiJqz852EsB69v0MacpXyMIwigi\nk7Df4KKPb8ZeggNOLpHdu2HbNvU9uUhYaveM7hYJh+GXv/Ty0Svf5HPPXMcH+BUeIErMcu/rUwI/\nYwaBPXswV2tyqj6Zbswg7hdBGGtI5mkB0V0i1dXWZ3sSk9+fOt5d72vhQlhUbbCu9EnuaN7LntIr\nKOnvo5QBq+GRIzB/Pi0/6WR7qzK57dUnQd0ocl1YWxCE0YUIewHJtDBGKKRE3YxiCYWgrS0xKiUc\nVq/ycmX1V1erfROqKpnYdgL27mWwKnFylZ4eaj52OjV+P2zeTChkJFSd1G8gpnWuu29yRTJZBWFk\nEGEvEPaY9VWrlIiafms97FCf3LRjj1AxJzvXrTP3VdLAEXafs5J5rz0Vj3unrw+2b4fp05l7SRtQ\nmTC2VEXGRNgFYewhwl4gsolZzzYqJVmEDc4PPk7LLyJUXb+AU49rFSOjUR7aU8XxyTN4YM0eDJ+X\ntrbEydKqKud+h0OoRfwFIf+IsI8y0omcfZ9hJEa4mFhhiwZtx1/k1/NrOf/lnyckNk05cYTPfXc+\n37q9kw0bDFatUsdVVUF7u8pcbWxMPK8bEc62/K8IuyDkHxH2EcB1QlKG7ZGIszDW1lruHJ/PYJOv\niY5f7OXBPRdTqvneywZ6uOvuCj5/ygym1O8BrzfuHjLdPPawzEw4VbcU4RaEwiLCPgIMp9DZffle\nr4qgad5TyaElR/jWmyv5a933jrLeoz4fYa+fg4ubaG42CASUqHd0WOJuRtCYv8HN70i1KEi2i3qI\nZS8I7hFhH2OEQspFYsbAm9t8PstS1+Pmdes9HDZ4+P2PU34ywqfuXsBZWL53DzC3o5nvvn46fzvV\nT3P1Zr7UbHDGGc7JTG5INV9gH6ObPkXYBcE9IuxjDKfQRCdhNCs26mJoijwYrONFbgytZOZvtiYs\npl0y0MeVx7Zzxd3T2fFXbTz5RmVCf27IxiLP1tUjCEJmRNiLFNOKN7GL6slyg+8ue5zgf3fQfeFi\nppX3wmGrnH4JUX75hgqPiX6tlJLn99DwdKXromVuM1f1MTr9BlmHVRCyR4R9DJMqykR31eiVHHUM\nQy3QEcRL/bEDfOPvI/xNkxU9Ey9NAHgG+xmouph9NxwhEjEcfebZCK0u2OmqTUppA0HIDRH2MUyq\nSUlT1P1+q50esui02Hb1MoML1jWpypBXXIHn5MmEvicwyPceriCKB96/k1B/jWMYpF3knQRfBFsQ\nhpehlu0VRhlm5qpentcutGnXYK2shBMn2PPtXY5lgUuIMrh0Kf/0+YjjxKh9W7YuE/vx5nepLCkI\n7hGLvUgwI1YCAauQmM+X5UIeGpd+voaPN3Xx5VdqmfvGs0zvPxjfVwLsfKGCgavK+P6aVurvt8oT\n5BL3bv8dTuUNzHo5giBkRoS9SNBrzejfMx2jv9s5dbbBpU81se6uCNc+8lHmdoQowSoNXDrYx40/\nvIKff+BE/JhsV2hya9GnumFksxKV/Zhs9wnCWEGEfYzjFDniVpxSTb7aJzbDnQZfC+xgcm+E4Btr\nOP2ZJlVUzONh4t8s5/49AXwlk2HjRsDI2mdu/w32CpahkBWzb9att/v0RdgFYWSICsNLXZ1637Ej\nP/35/dbnHTus/qPt7dGe0vJotL1dNVKR8NHo1KnRE5OmRaMzZ0aj4XBO54yfw/ZdH4vTfvtx2ZzD\n7T5BGAmAaEZ1tSEWexEyFGvUHjVTW2ut8LRhg9mqEuPfT7C2Epg8WW1atAj27aO8pxt6gLlz4aqr\nYPPmnGc+zbEEg2os5pOAYSifezisxtTUpAqX6StRuS1jYBhqwRKJlxeE3BjpG1/Rk8lSz9YaXbXK\nxfFdXdHotdeq95kzLetdf/n9ar9trOZ49XHbf4P5varKuU1dnfVyi9421ROCIIwWyMFil3DHIiLf\n1mW6rNA4hgGbNqn31laYNCm5TXMzLFiQENqiW9B6iKP9N5jfdaM/23r16ZCSBkIxIq6YImcoafmp\nRDYlXi90dkJtLQOPb2VCf6+1r7MT1qxRN4EYmUTV7hYyXTG7dycmY23YYK1I5TSxascsSwzZR/EI\nwljAk7lJ3og9VQgjhZmYVAjuWdvBbZsXq+iZQ4dg0SJa6n7FD39iUL/pfKad6KSXMq6Z0cpvj3hZ\ntcryj4NzjXlTgOvrVXKV2S5V1cnaWmuxkFSkqkApCKMFj8cDWWq1WOzCsNBteOHAAeV+WbMGHniA\nGsOg5kNAUyec6Aag9YiP8Bw/cxua4v4We7YsKFE3xdsu5KmEOV3su9OTgFjsQrEgwj6OGG7RsrtH\nFAaBWzYRMKw2i/rKmBrb6wHm/rmZaEUFHr9fhbiQHEWju0/SibEb0c6lHrwgjCVE2McR+RR2J/+1\n7ts269A4jeHZH7VyxQ1zKdEm+z0Azc0MTJ/FLVGDeyKtdBveuAgvXJidGLudFHU1QeyAJDIJoxmJ\nihFyIpX7w42/+srrvZSEX2PQwW04IdrP6Rzmtnt9BH+xiMDCSFK/qcTYFFu11quql2POK6QS4VzF\nOd9+efHzC/lEhF3IC6bf21xr1bSqUwqW18tzD79G95Qz4fTTYflyKLU9QO7ZA+eeC8uXU37SCpVM\nJcb2iVI31nghJk/d9C/CLuQTccUIrkmXvWmGqJt14FNlgOpceb2X0JlvxMXV+/sAczusFbN7POVM\nOnwYtm7lzq0VRO+G//zEE8y56RrH/traLBdNfX1yCGQq94lbt0quoaOZ+hdRF/JNJmF/CPggcBC4\nKE27y4BngeuAn+ZnaMJoI5sFMtxOSOr9zW1rIlz5Yebsb6Fj9hLmzj8Ftm8HrFivT21cDlsr4NJL\nYfNmQm1GXBjb22HFCvV51apkC95sl6v7Jd8LhJg3Cr3ImdN5BCFbMgn7euA7wH+kaTMBuBv4JYWN\nixdGkHSTk7mIkrJqDRo/vcvyla+NqIzVzs7Exl1dSvDnz6f0Ry/jFEWT/jzDv56qm/71Ust6yWVB\nKAQ+4Hdp9q8FbkHdBD6ept1Il1wQ8oi9joxJLpUld+xQ5WRWrVKlZc44Q72vWhWNfuPOruhbFy6J\nDkI0eumlyXVoPJ5otLQ0GvX7ox9c0pXQp/lu1pIxy9bU1SXXm8n2N2TzO53qz7gZlyBEoyNT3fFd\nwEeA96HcMZJaOk5INTGZTU10/d10RXi90NFh+emrAwbty1pUv5EIb589n1OPx1ZzKimBwUHo74fm\nZh6jAqZX8Oz3nydwvTdhPKGQ6ttc5zWdXzvdGq75wv5kYBiwdm3+zyOMT4Yq7A3AXViL6qR1xQS1\nZ81AIEBAHIljimzdF5kWtLCHIQaDVplg3TVRW2ue16Dh+Mv8en4tHjzMe2MnE4++FT/eA9DVRfUN\nPvjcDBVV47Vi4UMh1Zc+3mAwfVXhoU64grt2svSfYBIKhQgNcUZ9qMJ+KfBI7PNM4BqgD/i5U+Og\nOBHHNNlOHmYTbeJ0LhO9nAAYXBBsUh87OmDxYlWLpq8v3t4DcOQI+Hy8M3kGP7hpD92GN/5UYL+p\nRCKpb1r5qP7optiaRMYIJnajt17/D+mSoQr7X2uf1wOPkULUBSGdeC5cqMQ2HFbuksZGK9HIKk9g\nE1pvrB5NRwd9f30eZYNaNckYk985wm33z49NwBoEAqrvxkYr5t7vV23NG5d+Dr36ox7WmY8J1+GY\nwJWMWAEyC/vDgB9lje8H6oCy2L4fDOO4hFFOOheLXazC4WQXCCS6Wpzi4M1j7AXAzPd4X14v37rj\nTe56eTUcP07f0zspG+ixTtTTQ/9Zs/m/J3rw/HOUaHQJX72wiXDYoKoq2e+e+ISQ+slkqA+gwxU+\nKcIuZBL2G7Loa/VQBiKMLdKl6NvFyimlX6/3smGDVYZ31SrruHTntQvYyXIDtmwB4AufiHDr/6xk\n3mtPUTbQw4GzFjHj8D7K6YMovJdmdr5QAcCAp5QHPruHW4KVSecaafdItiLtNsNVhL/4kcxTYUQw\nrXETU8hdjT8iAAAZoUlEQVRNd4xONk8HAOddbnD+xschEuGFJWu48JkHVGmCw1Ybc5a/NNrPzfdX\nwY+nwgsvxCdbGxrUU4S9cqQeg55PkcxHRqybRKdcbhZyIxh7iLALw4JZvjcYdPYfp7P47ZanbqWn\n8kmDg5VvGPzk2k1caKCW7Vu8GA4eVOGRGh6AY8fUAtyvvQZeL01NziWBdWG335zsbbJhKD51PWTU\nPl8wVETYxyYi7MKwkI3/2G0maLo+U/mXdT88Bw7w7CMdvPvGaiYdP8QkBhIbR6NEfT6iePgZp/Jc\n8+WsO7mZk+WGo8A7RcwM1SLOdkLVabv9Wg93lq0w+hBhF0YcXbhNsp1I1F0jdsta58rrvXD9G9zy\niQiffWwZFx17Lu6WsZIxokzjKP+L7bz/7goGmcD+xhoaOtSEq9m/29WXMsXz2+cNhnodnLZlG6Yq\nN4KxjQi7MOwMV+am/btjxIyNUAjWrYPduw3uO7YbgMXv6mDTGzVMn3ScU050JbQvAUoYYG5HM4cm\nnEHrsn20tXkTMm/D4cTwSX1M2fjih+r2GKrw6jeWfEbrCIVHhF0YdrIVGzft0/ndU7llzO3V1bBt\nG1RVqQnSrzd6eTC0n+DqWMJTT49KcLJRNtBL9Q0+ngdoh1sxePn8Nq4MehPa5WoRu500TkU2182J\nxkaxyIsFEXZh1JGNuOhWsT6JqIt7Kh99OGy1UwLrhZsO8P5FEc780mrO+avj8LvfKZGPZbbqNTOm\nE6H6Bh8/fHIX53yqJitfudk2VTJUrq6dTKQ7zmnOQIR+bCLCLoxp7O4DPWwyE7W1KuPVLL6lly3g\nQ1sIhZQVW2V08Nl753MKPUl9eIAbN/ih9DMQ3MfnT0yGtRsJBAxX7gy3yVA6+YxUsd8Q7XMGIyns\nEpGTOyUjPQBBGAp6aKRTqGQqTMFIV1ExEFDCe1uDl3vv7ISPfAQCAQbtDZubYd8+aG5mxm+2QkUF\nTJumatmQ/9r1Tox0MtVwUIy/qVCIxS4UDabF7tbXrLdNd0woFMtsbVLFx/7trgjX/ffH8IVDlOza\nCTU1MHly4kFHj8L8+TBnDt8PH4THy1QsvdebfAJtTOnG4OTaMcv9OkUEuUG3ys3iaMLYR4RdGHNk\nSlQaamSI/l5frz77fFY2arjTYF74aVUGYTsE+iGwcaOy1HV6euCPf2RyNKqyXs87D666ipZbNlLz\nISN+HjfjThWp4hSznuvvd7P493AjoZb5QYRdGHPY/8j1ydNMYpAuIsXsx8nqtwuo3TcOBnu+vYuL\nb/XTTwkT6efAWYs48402JhDLdO3tha1bWbK1AiZMgJoayqubyGZpP6ffomf3DqUQWG1t9ufPt9hK\nqGV+EGEXxjzZZrnqESn2SVcTe/2VcFiJeTrL8bEjNfyDf4Cp/RHWdcVq1KxYoTqYOlWVLSAWWTMw\nAM3N3NlcQf+/lbJ/9hIOfK8pbsmnwzCs8TY3q3F5vVbJY6eIIDdkK9IyuSmArHkqDAP2NUKd1hdN\nRVVV4lqrfr+11qq5/qi+Nmm6895zj1qrddIkazlWrzca/dj7uqLPea+N/vrhcLS3ZGLymq0Or96p\nRjQaDrv6Dea47O+FwFz39p57hqd/Wf9VwQiseSoII0rK2jAuWLjQcqmY7+ZLxyk6Q59wXLcOXnoJ\n3nxTbZs0yZrUXLvWIBTapNpfuY+jF1YTHRxk2pSoWvnJgbJjEaI+n7Lsy2KTrpXJZYV1DCNz5m0+\nLGzdlWXG3d9/v7qWmfrOZWJXyA0RdqGocCMujY2W28LcZi/Ra243q1Saqzjpi4aY5zM/B4MqcMYw\n1DY9lFLdLLzUH3sDvx+WVUe4sXklp+77LZPe6kwaZzwRqq+PaFUVnpkzHaNqzHOby/uli4zJh7A7\nRdG4zVgdDteNuIOcEWEXxhV2ITbfU1WGNLfpFn06IhEr8cmpH1DCv3ChwY/PfJxIBM7+cwu16/30\nAx4mMAlr/VazMBmHD6uwymuuUTHzkyfDRpUIZR9zpvhvp4qSbsXRvNnt3m359zs6rN9bW5s+NyDf\niLA7I8IuCGQXKplKTAIBJXiZhM20rs0J2dfn1FDGAKtWwW+fjtDQvZrBoxHeM2E3/aXlTO6JKCFv\naVFLTMUeNU5ccAkdg3OI9E6moWsj3Rjcf79yCZl9O63Tao/8yUYc9RIOnZ1qbri+Xr03NSXe0My+\nhyN8UQQ9PSLswrglGytVXzTE77eE0V6nZtu29KV8dYvafAow3UM+H2zYb/DP/i1EIkokq4wOPn5v\nDY/e1EL3ei+fPzGZGQCLFnHKpEmc/4wS+S5UDP3xboO6j7fxrUZvfOz2iCG3SVzpCASUmJuuGDdP\nPeb5zXHlin5jGu5Y97F6AxFhF8Yt2VipoP7I/f7kxa/dhFua4t3WBu3t1nbTyj3//MSnBisL1Ett\nZD+NDbED1m6ENWvggQfgE59Q20pK8AyqQgdTT0b45qM+mFACS5eyu7qJQMBIinkPhy2r3pw/cPot\nma5LKASlparfVPXpU1XazHWy1R6imso9lq/JYhF2QShSdDFxint3ezw4F/3SrX67kCTUmjEM2LRJ\nfd6oRL5329NMPGqVGfYADA5CczP/0HImPDCZwKWX8i+3b47HyevnTRefn07YTGGvqXHeb+9DF3o3\nOM0F2N06TksYpjp+PCHCLghZkGlyMlchsUez6P2kTPWPifzEjg56Lqlm4ludCWWFjzOZKQPvQFcP\nbN/O4qdmQHSQf/KUQmgJNDURCiWXNtDJJI5OUTLpMF02mRbdznQusy8zySwVudTQKYayBiLsgpAF\nmf6w3QqUiZOIhEKJk54Z67R7vUw68gZs3crg8uWcKK9gyvuuZEpvL2zfrtpMnUpJLPN1QrRfKevs\n2Xz12HE8RDlebsDqtrRFytJhWv7pfpf+21M99WRbB8h+07Mfbwp6phuAjhvX2mhHhF0QsiRXqy1V\nmr/dF+3UvytxueYadu6IWpZzJMIfqmu5YIFHlTMwRR544/SLmPZWB5NjSY1TT0biC3mfKDuVztmX\ns/nazdTfbYVT2itJ2n9bqhBPvQ6PWc/HPgHtdJz9d5vt7ZPZ9mPtY3MTpurEUCZ4M/U73Ja/CLsg\n5Ei2f6DZThjmKgDxYwyDg/c3cUEAZf6vXAm//S1ccgl/9V//Beeeq6pOxjAX8p7Sd5R5f9rOnXdX\ncCfgaZgG7e0E1ztH2piYTxh29Jo8+hiH4rZKZ1GbN5CmJjVRnU29H3s/wyHAIuyCMIrJd4q+6XfW\nLdNsBdDu57b6MeDxxxMbt7Zy+NxqKgaPMGGgDztxf313N/h81AHUw+rZS2DtLwi1GQmZuboVbp4b\nlJjmcq0CAatUMrj3d+uWu3k9amtVVFKhGOmJWxF2QRhGMk3EpQulNKNh8h2yZ43JS33fG3xwSYSv\nvLSS+cd/Sym9nHryLSvjFRhELbVmfvfuf4bBigr8QE1JGTdd3IrfX5kw9t27E2vK6FZ7NpOZkUj6\npQMzTeyaY0i3ipVOKh+96YZyi5uIHnOMw3EDEGEXhCzI9g/UzUScHkqpTwZms7j1UMYEBtUtyppv\n+UWEvk+v4dUP3crK9e+nd/I03pq1gLkdzQl9lMTeSwf7eGhPVXy7px544gmCXJMQUmmfLDXHm871\n5MbH7XYy2+0iIvn00afrd7gnZN0I+0PAB4GDwEUO+1cC/4i6ob8N3AzszdcABWE0MZx/oLrrIJv+\nncZkir05WWmiT4CaLz300DAMmi7YRGAOTOYkdXdA+ckIH/yP6zj/4E7Kpk2Bt96Kt9ct+/i25cu5\njdM4um4iNy1s5aWwN+7jtidF6Vmkemii89icf3+6UE3zZb9J6u/pjtMLxWV6yhhNYZJuhH098B3g\nP1Ls/xPwHqAbWAY8AFTnZXSCUERkGyrp1n2Qqi+nG5A9M1TPSDUxS/Ba2w2C5U+qm0JthEMfXs2s\n0ojyt5w8mXRuDzCNo9ADjzzni28f3OBh4LTpTNywB5/P61jewD7n4Pc7/x6dTBmq9mtgfs7ko29o\nsLa5XcBE3+d0U25oKIzIuxH2XYAvzf5ntc/PAWcPZUCCMFbIZTIw332mOt6tD1v3Q9uPKy+3xKm+\nHqqqIBIxCHxtS/xYT0sLg0v9nJw4lVN6jyrrfcIEtUJUiVXqAGACUSYcPcJgbCK2/2tlPDrYSiik\nas3bk5ZMRiqOfO3aRB9/pnGkyhzWaWoqTPXLfPvYPwM8kec+BWFUMtyWVy6LSzsJu76UXqa4b/24\nZcsSRUi37C3/eQ3hVQNs2AAXzY6w8dQ1vPPFr1D5pQ9SfsG8RLUGKCmhJCb2pYN97KUKz94KOHqU\nQaKU/OuE+MIi6W5O2SYy6U8nbo9Jtc0+DrfHDVdcvBP5FPb3Ap8GluSxT0EYVwyHn9ZuIbqNLDHL\nC+v+ZvNYe2y6zwfhsMGFjbE6Np/ZD5EIb733o0xr28nr77qMCQf20zf73ISJWA9AVxcQm5DtG1S+\noAkTqOvvV5OxALt2JRWlSTUXkcrNkotrJpM/vrHReeEW/abc1qYs9XBY1a43+1qxwt3KU7mQL2Gv\nBB5E+di7UjUKalcyEAgQGMlAT0EYhQxlcjbXm0Kq43Qfv89n+ZnBWmTDbKMLm5UEZDD9+R0AeGO/\nJbg2Atddp4T60kvhmWeSBxSNQn9/4qTs0qXg8ah9JSWcMWspvNSUekbV4TcORW5SPdWYpZdBxcrr\n38GKqFm7NnFSW99vH1coFCI0RPM+H8I+B/gp8Engj+kaBkfKWSYIo4ThTFxxG1qZ7jhdzPVIktpa\ntc/UG3sIoC5Yaf/MDQOefFJ9jkRg9Wol1rffDldfnVrsp0yB48fV58FBFrzZDBWqBn0dMPi1Uu5f\ns4eDZ1S6uqmZa8TmciNct875uLa2ZPfZ7t2JbjCvV31PFxevDN5AfAz1+klc4kbYHwb8wExgP+o6\nlsX2/QD4KlAB3Bfb1gdcnvVIBGEckI2wD8cNIFOfThanvi9dv6mMzJSRJIYBW7ZY30+cUGK/Zg3c\neitcdRWDvX2U7NoJH/2oJew2PIBnsJ9b7q9icEIpX/FAaWgJLbc3sT1kJIV8BgKJa9ymi2Jx4uRJ\n52tkr4kDauGV6liM4KpVyvWycKE6JhJJfUPRQz9zwY2w35Bh/42xlyAIeSRfroN8HVdbm9jOHneu\nL/ptttHbuRqTXm++p4d/DkKwBmhtpeeSao51DzJ94GBS7LxJyUC/8tU3N1Ozczo1ANEodYDn66Ww\nZ098Yhbc32jtYZimmD+hhYo4JZTZJ5xNKz0QcBdtM5zCLgjCEBiJxJWhCHsqIbbv06szOvmcsyV9\nopEXPv8G9fXw4Gq1+PexBZdhHP4T0UOHkoXe9MebXwH6+xmougQPg9QRJVoPdzCF3oZT+MD9raiZ\ngPTZsHamT7d855CcHKbfBA4fTo5O0mPjnf6fiLALwihlKBOiI0E6CzsXC9ftzayx0V1Y4o3BGs7Z\nOcAnr1Pf19d38D8TawgueJgvt/8fjp2zkNkHnmPKCWtVKZMJly5UVnuMKRyH7uNU3+Cj97MzmThw\nkvccPwGlJdCqYuztNy/nsM/k32g+xTiVUtDb2f9v2J8OckGEXRCKiOGcnIVEoUoXE5/LzUyfuE3l\n5zcFr6ZG9217+b5vP6cDs9sPULcSpkU6uPnHiym/5N3Q3w+//jX85jdw112O5/YAE7tVDWMVdjkA\nVVXxipZRwFNRAc8/T3m5N+EamKUSVqxI7NMe9qiP3/zsFG2jPxGZ58gWEXZBKCDDHeGbq7BnspLN\nfRs2QF2d+rxqVXbVFjOd13Rb6PVjdLHXLWDdn93ZaZXktaxpLzQcSD7hxo3w0Y8yENrJhE2P0P2Z\nLzKt9Dh0ddFbNpmS/l5Ko/3xGjimiyceb+/zsRXwbFNi/1Xg7R9XUDnwPJGAN2nxkFTXQ38yySUx\nKhMi7IJQQEZr6kYmC9tpn1uXUr5+s+7i0M+vT+qmi9xRQmoQ9u5gA1D3ewh/7Fr+blkHV95Rw8SW\nFuju5uTCKyhf7BB2GfPb62LvAaYNdBHGh6ceBktK6Whcwrrrm+KrTzU2qpvPyZOWa6WxUSUrnXWW\nZdnbJ57Nd7HYBWEcMlJVBd362tO108XatLr1glup+tTRxVwX+VTnCgbtTxteuH5/vN2y95wg1BTh\nyNXXMbV1J+9MmMo7A2X0zpmfVL7YxBT7ksF+5nY0c1fDmdxJD9TH9s2YAV1dRBmkn1IevGYPm15U\nNexvukk9ieTjCchEhF0Qxjj5npzNJMRu2pmki3Sx92WKbqokKf13pqt/4wZ75Uy7S6h2rUF4ypOU\nf0DFovv9sKw6wnU/X8lfH2nl7WNR+nqiTJ9wFHp7EzuLJVMlROocURO5HqAMFXN/s7mvHk5MOo1/\n+tBezMicoSLCLghCAm6FfSg4Cb5pbTslSZnC29iY6Od3m1xkf6qpqkoO10x1czT3VQcMPr37cQLX\nWVmkn70+wtUbV3Pp+cfhd7+DSy5RQq8tHA7AKaeoBCwNXfgn9xzlm4/6YGIZfPjD6kbw6qvqJDkg\nwi4IRcRo8OHn6hpKN+GoW/R1dc7lDNz0HQpZSxCuWOGu1K5+/vJydYzXq3zkv9xt8Mj0LdQuh7Wx\nSgktv4gwY/9qphzr5OwDuwmfvYSfffy/+FTLTczY86uU5/AA9PUR/ekWPMRi8F9/Pf3A0vVVIKJR\nLWFAEITiR68pYxd83epOJax2wXaauDX7deuC0vtwuiHoWanmWM0FP5qb1bhDIVUSoK3Necz232m6\nfi54fSu3P7WcP5RfzLyTL1BGH/2UUk4/AIN42D97Md79z8C0adDdbYp0VlotFrsgCAUhl7kAp4Ql\nc5LRLKSlR8o44eTPT1dSwOnJQb8Z6CGZbrEqXl7DsmVRqqutzNQPLonw1d9fR9XRXUx6/jm8c+ao\nejn//u9wxx2webP7E8UQYRcEYdgYqmvIKQtWD3NMlyxl78PtmNzkAgQCyQlJTudvakq+AZhFxKyS\nvQbv/MuTfCMEwcpYI7NezqZNKswyS0TYBUEYNtIl6Lghncimq0SZqU/z5XRDSHcj0D+7WeJu4cLk\nbXqUT1YF0rJAhF0QhIKTzqJOR6owx3Thi5mWA8x0Q8g2xFPHjPTRx6OHcBpG5izVXBBhFwRhRHGq\nmeIkyk6LUzgV4YLsxXu4krycrPJs/fO5IMIuCMKoIpUopxPDXMV3KEXL3DLchdmcEGEXBKHg5MNC\ndjsJOhwJV9mKdaqFOoar7IMIuyAIBcethZyudEA+hD3V2DKRSdgzVctMFUOfL0TYBUEYtaTyqdvJ\np0jmo590N650MfT5QoRdEIS8k41o5UPc0i0vZ54jHzH1+ei7EP52EXZBEPLOcAi7m3bDOQmaa99O\nET9OMfT2tkNBhF0QhDGBU02XkahDny32OHj7zcG8CVqZqENHhF0QhLxQaKHNZEFnG7UyUu6U4fC1\ni7ALgpAXhtMNkgujXdj1OYFgML83QxF2QRDGPKPJ9eIGe1x7IKBKA+fr6UaEXRCEvFNooc3lfCPp\no3d6uhEfuyAIo5qxYEGPZddRJkry15UgCIKQLcNR3dGNsD8EvAn8Lk2bbwOvAO3AxXkYlyAIQsEY\nySeM4Ti3G2FfDyxLs385cA5wLrAGuC8P4ypqQunW8RpnyLWwkGthUehrMRZcR9ngRth3AV1p9v8N\nsCH2+TnAAM4Y4riKGvkDtpBrYSHXwkKuxdDIh4/9XcB+7fvrwNl56FcQBEHIgXxNntpXW43mqV9B\nEAQhS9wuf+0DHgMucth3PxACHol9fwnwoyZcdf4IzMt6hIIgCOObV1HzmHnHR+qomOXAE7HP1cDu\n4RiAIAiCkD8eBv4C9KJ86Z8G/i72MvkuyiJvBy4p9AAFQRAEQRAEQXDJMpSP/RXgzhRtxksyU6Zr\nsRJ1DfYCzwCVhRtawXHz/wLgMqAf+FghBjVCuLkWAeB54AXU/FWxkulazAR+CbShrkVtwUZWWEZ1\nEugElDvGB5Sh/jEW2Nro/vgrKF5/vJtrcSUwLfZ5GeP7WpjtngZ+AXy8UIMrMG6uhQH8HitkeGah\nBldg3FyLIPCN2OeZwBGKs77VUpRYu5nHdKWb+awVcznqHyoM9KGiZD5iazNekpncXItnge7Y5+co\n3th/N9cC4PPAT4BDBRtZ4XFzLT4BPIrKBwE4XKjBFRg31+IN4LTY59NQwt5foPEVkrwngeZT2J0S\nld7lok0xCpqba6HzGaw7crHh9v/FR7DKURRrHoSba3EuMB3YAbQCf1uYoRUcN9fiQeDdqOCNduAL\nhRnaqCNr3cznY43bP8bxkMyUzW96LyrSaMkwjWWkcXMtGoC7Ym09uM+vGGu4uRZlqMiyq4DJqCe7\n3Sj/ajHh5lp8GeWiCaByYH4FVAFvD9+wRi1Z6WY+hf0AMFv7PhvrcTJVm7Nj24oNN9cC1ITpgygf\ne7pHsbGMm2txKVaC20zgGtTj+c+HfXSFxc212I9yv5yIvXaixKzYhN3NtVgMfD32+VXgNWA+6klm\nPDGiulmKuvg+YCKZJ0+LOZnJzbWYg/IxVhd0ZIXHzbXQWU/xRsW4uRbnA9tRk4uTURNqFxRuiAXD\nzbX4f0Bd7PMZKOGfXqDxFRofozgJ9BrgZZRgfSm2bbwmM2W6Fj9ETQY9H3v9ptADLCBu/l+YFLOw\ng7trcTsqMuZ3wK0FHV1hyXQtZqJKmbSjrsUnCj3AAiFJoIIgCIIgCIIgCIIgCIIgCIIgCIIgCIIg\nCIIgCIIgCIIgCIIgCIIgjB7+P3eMUMf7yeSmAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# exectute the curve fit...\n", "coeffs, fiterr = curve_fit(func, xvals, yvals, p0=(a0,b0,c0))\n", "\n", "# ...and plot the results\n", "print (\"a=%s, b=%s, c=%s\" % (coeffs[0], coeffs[1], coeffs[2]))\n", "plt.plot(xvals,yvals, '+')\n", "plt.plot(xvals,func(xvals,*coeffs),'r.')\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to find the initial values we might want to plot the function with some initial parameters $a_t,b_t,c_t$ to allow for a rough manual fit" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "a=0.5, b=0.5, c=0.5\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEACAYAAACnJV25AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnX2YXVV97z/zFkII5KQGEl4zlir4gknBYASuOV61gFSN\nVWkRagar1NtrNfepTy9679NMfJ5e6lOf3ilyaxGFiWKx+ML4UggVy5kAZQIpzARqoxKcMUQSDM1E\nQpgwycz9Y80ve501e++zz8ucmTnz/ficZ87Ze+211zmG71r7u37rt0AIIYQQQgghhBBCCCGEEEII\nIYQQQgghhBBi1nMmcD/w78CTwCdiylwNDADbgYeAN9StdUIIIcpmGbBy4v1C4CfAa4IybwYWTby/\nDOirT9OEEELUgh7gbSnnFwPP1KktQgghqqQdGMKN3JP4FPClurRGCCFEVSwEtgFrU8q8FfgxbtQu\nhBBiGmjNWK4N+DZwO86KieMNwC04j31/ePLss88e37lzZyVtFEKIucxO4LfKuaA5Q5km4Cu4kXhX\nQpmzgO8A1wBPxbZs507Gx8f1Gh9nw4YN096GmfLSb6HfQr9F+gs4uxxRh2wj9osnBHs78PjEsc9M\niDnAzcBf4OyXL04cGwUuLLcxQgghqieLsD9I6ZH9RyZeQgghppksVoyoMfl8frqbMGPQbxGh3yJC\nv0V1NNXxXuMTfpEQQoiMNDU1QZlarRG7EEI0GBJ2IYRoMCTsQgjRYEjYhRCiwZCwCyFEgyFhF0KI\nBkPCLoQQDYaEXQghGgwJuxBCNBgSdiGEaDAk7EII0WBI2IUQosGou7AXCumfhRBCVIeEXQghGgxZ\nMUII0WBk3cy6JnR2wsaNMDjoXu3tsGlTdD6fdy8hhBCVU3dh9/+CE3f/sxBCiOqQFSOEEA1G3YU9\ntFpkvQghRG3RnqdCCDGD0Z6nQgghSgr7mcD9wL8DTwKfSCh3I/AzYAD47Zq1TgghRNmUEvZR4H8A\nrwNWA/8deE1Q5p3AbwGvAq4DvliLhmnhkhBCVEYpYd8D9E+8Pwj8B3BaUObdgEWjbwVywNK4ysoR\nawm7EEJURjkeezvOZtkaHD8d2OV9fgY4I64CibUQQkw9WRcoLQS+BXwSN3IPCWdsY8Nfurs7KRTc\noqSOjjyQLwp3LBSgu9ud37gxOq4VqUKIuUKhUKBQ5Sg4SwhNG/AD4B6gK+b83wMF4BsTn3cAa4C9\nQblxGGfDBvchn3dC3tnp/ppwd3YWvwy/jBBCzBWmItyxCfgK8GPiRR3ge8CHJt6vBoaZLOoAbNgQ\nCXY4UveJ66y6u0u0VAghBFDairkYuAbYDjw+cewzwFkT728G7sZFxjwFvAhcm1ZhoRAJt9ktPT1R\nUrCNGzlm1/ij9MHB+Lps5K/RvBBCOEoJ+4Nkm2D9eJabhV65CfzAAKxc6d6vWxclBisUoKPDfe7t\njawZq0fCLoQQk6lrdkcT3+5u98rnI7G2EfngYJTKd3AQ+vudsJeDhF4IMZepq7AbJuK++PqjdJtU\nhciP9ydaLXpm0yb3vrc3snI0ghdCzHWmJR+72SqDg9GkqL/hRl9fVNY89+HhSLBNtK0ziIugEUKI\nuUpdhT0UXBth5/PFG27EjbhLxbLHTcpmuU4IIRqNaRF231v3jxtmp9irtzcq5wu1/zcUcO3KJISY\nq0yLxx5OhiZ54nGeu08o8EIIIaYpH3tHx+RjJs42ejdBjxNt33aJQ0IvhJjLTMuI3bdaoNgTHxws\nDouM88gLheJycfULIcRcZVqEHeI9cbNcOjud2C9fPjn9gBG3EjUOhT4KIeYa0ybsPoODxaGQ4ER9\naCgKifzWt+DgQffZylg0TUdHsnhL2IUQc40ZIey+MNtCo0IhEnaIBNzsm97e6Jr+fom3EEIYM0LY\nk1aghjaM2Tf+qlT/r587JvTvBwfTR/ZCCNEozAhh9zFh7u6G7dsj2wVcFsiVK90rLn+MH+ce+vdJ\nXr0QQjQa0yLsab63f9y3XqA4lr2rK8r8aKPyQsEdX7++Zk0VQohZx4wT9qxl1q8vThAG0eh+eLjY\ns7coGz9RmEbvQohGZcZZMT5ZOoDQd4fi9+H1SjUghGh0piW7Y1KSLpv0tJzsa9Ykl/Hr6etztowf\nChleZyj8UQjR6JS1QWqVjI+PjwOT0+z6+Hlhens5tvl1LjfZO49L1xvmlYnrLOImWIUQYiZSyWbW\nM86KMcG1lAFhSGOpa+MyRfr527PWJYQQs5VpyxUTElosK1Y4cc8aQZOlbj8FcJa4dtk2QojZyLRY\nMaWwEXUu5yJcwIn9hg3VLTTyrRp7+YuaID7+XQghpotKrJhpSdubhL+v6caNkajn807UOzujcMVK\n6o67zrbmi0sFnDXRmBBCzCRmlMdufritPPVHyybAUCzSWe0SE20/rt2O+4Luf960KfLlK5lslZUj\nhJgOsgj7rcAVwHPAeTHnlwC3A8sm6vs80F1No+K2zuvvjwTZTx1Qjnj611ikzNCQs3aGhqJ7+QnH\nqrFiJOxCiOkgi7DfBnwB+GrC+Y8DjwOfxon8T3BCf6Sahvm5YPJ5WLs2fhFSaJ+EnnncZh65XPE9\nLArHF3SbaLX7lNqGTyIuhJgpZBH2B4D2lPPPAm+YeH8S8DwVinpog0AkuJs2FY+ye3qcQIcLknxh\nDxN/+R3DypXuXE+Pqz/Ocunri+9M4p4U0jJLhvUKIcRUUguP/RbgX4BfAicCV1ZaUVxWRsNP52vC\n6Ue3GGl7oYb3ivvsXz8yMvk6X8CTwi2TvoMQQtSDWgj7Z4B+IA+cDfwQWAG8EBbs9FQun8+Tr2AI\nG7cIKVxdajZKUtKvNJvGroFoww+/3Jo1Tqwtxl6jcyFELSkUChSyjlATqIWwXwT85cT7ncDPgXOA\nbWHBzjKGr3Ej6lDAIQpJtKptZJ8Ug24RNybeAwNuGz6IvHc/W2R7e3SPFSsi68f++uGX4f0k7kKI\ncgkHvRv9UWNGaiHsO4C3Aw8BS3Gi/nS1lSZZJX56gEpsDot68a/t6XHivXkz7N0bHV+6tHhBVGdn\nNImbtMDJt2gk7EKI6SCLsN8BrMFFvOwCNgBtE+duBv4PLnJmALfg6c+B/6x5SzMSimo4qenbNP39\nboPsnTudYJ97Llx2WXFETpgCOCkKp9y4eiGEmCqyCPtVJc7vA95Vg7aURan8MeFkqNkvvb3OJzfx\n3rMHzjjDvffzyFjoY1zdceItMRdCzBRm1MrTcihHSP3wR8OslDDjY7h5dlxdtiWfpT4YHCwOy/TL\nSvCFEPVm1gp7GnHRKv7E6549TnBthG6ZJCGbEIcefVxYphBCTBcNKexh6gBwIm7vr79+cq73UtFF\nSaGNWROFyXsXQtSLhhR2I9wPNW0kXSovTNLCo7iNPZLql7ALIepBQws7pItptUKrSVQhxExkzgh7\n2kjaJkH9a7J0CHGbc8TVDVqZKoSoHw0v7MZU5XWJG7X7gq+8MUKIejOjdlCaDfhb61lu+HDyNW43\nJjsuhBBTjYR9gnKsEYt9X7MmCplMEnO/fgm7EKIezBkrphSlhL2ry+3BaikJNmxwx9eti6617fxM\n7Ht6XN53yzWTJaRS3rsQolok7BkZHo4sGAuh9PPP+Jt6xOWXyTJBW4mwqzMQQoRI2DPgi7efk91f\nCJWEL+C2kXat2yZhF0L4SNhJFsdw5Wo+70R97VpnsSTtp2qrUcO8MXbcv99UhERK7IWY20jYybbN\nnY22Ozth/frovOGvRLXzYd4Ym3QN7+ffJ8uI3vx+iO8MJOxCzG0k7BnJmiQsbpNrf4Nuf7s938bJ\nIsR2jfn9huLjhRA+c1bYy7VAbLSdNLKPw+oKt/WzhGT22cS+3E4jPFeJpaPRvRCNx5wV9iyrQuNG\n23HX+tebqNqWev4OS3bNunVRDncLjezuTk9/EE7gGrZPa6WrXCXsQjQec1bYs5BVLOPK+WGRvmCH\ne676+6aW6lz8CVzlmxFCJCFhZ+oFMrRJ1q1zIp7LRXZMkn2S1mkkkcXSUXIyIRoXCTvZxCzr5KZZ\nMJs2RatU+/pg9eriyBcTV38fVrtPqXvV4nwWG0oiL8TsRMKekazib+Xa24ujXsKNOcKyfhkrl5Tr\nvR6CK2EXYvaiJGBTSLgiNU0sk1L/ppUpVT4LEm8hGg+N2KeIXM4lAWtvj1ag+il+kxYoxRGuVLVr\nksqUQ1LMvbx3IWYvWYT9VuAK4DngvIQyeeD/Am3AvonPc5K0KBaYHCljx+yapA2zfWH3P4f3rjZm\nXQufhJj9ZBH224AvAF9NOJ8D/h9wKfAMsKQ2TZudhKPbuOyO/spTI8wQCcmbcPs5ZwoFNzl7773u\nCWHt2slt8Fe5ViP88t2FmB1kEfYHgPaU8x8Evo0TdXAjdjGBH5/uj+R9cY8TSz8T5MaNTrSHh2Fo\nKKq3vT1aBLVhQ9Q5hMR1IqU8/UqeCIQQM4NaeOyvwlkw9wMnAn8LfK0G9c56wigWX+TDfO6hp22b\nc/jHrKxZPIOD0N/v/tpxqz9MZeCvWPWfDvw22CpYv81CiNlHLYS9DTgfeBuwAHgY6AN+Fhbs9IaT\n+XyefIOrR/j1bPk/FI/Ww5DIpLr8jiGXg82bnbAPDcHIiDtuMfP9/aW9/rDeuNWvhYJLd9DerglV\nIepBoVCgUOU+mrUQ9l04++WlidcWYAUlhH0uYql2Q7/b//8wLX7d3g8OutTAPT2u/MqVLpVvZ6cT\n+zhCr998eohEO9xAxL8uye8XQtSWcNC70R9RZaQWwv5d4CagBTgOeBPwNzWot2HxR+n+MUgXdv8a\n21u1s9OJ9F/9VbTQKemecZOz/v00CheiMcgi7HcAa3DRLruADTj7BeBmYAewGdgOjAG3AD+ueUtn\nKWFseKHgRu4DA8VWiZU1vzxt8jIU8CVLYM8euP764nJG2orVUNiT2g3pbRNCzByyCPtVGcp8fuIl\nApKSeK1dG9kaHR3R+U2bojj1uBG078t//vPu786dUd1h3vjQ8vHDI/1jYfks+WSEEDMTrTydAYRb\n6EFxiKQ/crY87u3t8OKL7vjy5W4C1YTY98/jrJaOjijpmIVe+rHxVrZaNLoXYnqQsNcJE08LO7TU\nvWFYYlzYYdpq0O5ud6ynJzpmm4L4Am9tMCsnzAkfhj+mef5ZyXqtOgAhaouEvU5ksTbCuPMsmCj6\nC5pM0G2HJhP5Bx+Ea64p7ljMN/dj3VescNfncsX7q06VAEvYhagtEvYZRlzEjH8uJK4TsNF4GJse\nvgc34vefEpYvdxO7a9e6p4D584vvVUqA46wjS2Es8RaiPkjYp4FSApd0PjzuR8eEScL80be/oYeF\nOA4OFu+zClFn4Ns69947+UmgVNvDdpba8i/rwieN7IXIhoR9GqhU2CFaCepvhG3CbWkI/EnUlSuj\nhVGG3yGYuFpn0NfnRuy5nKv/+OOj+ktt6F0OYTuzWE8SdiGyIWGfZfiCmLYa1B+VGzYSt5Wm7e2T\nOwWAyy5zaQkg8ujb250vn0XQbTQ+OBjlog/b7pPlSUAIkR0JewNiwhouZFq5cnKETei5gxN1O7dj\nR+TZWybJLE8cViYulDMkbcWsNv8Qonwk7LOYUpOpNtq2srYQqlCI8szE2Tn79rnyw8Owd68b9ff0\nwLx5UX1JIZHlhDiaYKdZPFooJUT5SNhnMaGAdnVFfnpvb7QIKU4su7qKo2R8wfTDLvfsiXLAQ3Lq\n3yRhT5sIjhPsMCJIvroQ5dM83Q0QtcFG4ZAu6kZWO2XZMldXe7sTeKvXQhhLUa4oxwm7/ZXAC5EN\njdgbBBPicORdTQROd7ezX8B1FmefHaU0sJw2di5MaJb1/qXa4LfFnwg2KrWBFFIpGhkJ+ywnbnKx\nUIj2Pi1FWhmLhglXxFpagjAlQfi+HPxtAP18OBZVE6Y6MCTsQkxGwj7LifOqS6X9TSPsKMx68VeQ\n1po0v92eCoaGisXdLy8xFqIYCXsDUs1o1B/9rlkTWSuWQ96iV/xt/rJMlqa1p1xh/ta3ijsf27ov\nzQZKCpvM5VwYqEIqhaiMcTG13H9/+vkNG8qrb8OG4mvKvT6N+++P6vPrDb+DfV6+PL4dfhuzts8v\nt25d8jkhZgLAeLliqxF7A5Fkv1Q6Gk1LSFYtWRKd+amMh4Yii6avLz7NcSUph7XqVTQiEvYGp5oF\nPqEoVhK6GDfZaaGZAwPR39DDt3bb576+qI57741SHtgkcdzcQpjkzMjlihdx+atvZb2IRkDCLlKp\nxD9P8s79p4eBAZciOJdzOWj8UMZw022I0hxYPXGbj4TXJI3GfU89RMIuGgEJ+xyiHqKVFJZo988S\nzRInyPv2ZRtl+51HUjn/fTmbmggxW5CwzyFqKexJomwphS0fvH9vG8Xb5h5DQ058V6xwVsuyZcnp\ngZcsKc4Ln5TRslzSEpCloRBLMZORsIuKSLJZLE98Pu/CJZOSevkj61CkS42ik8TYX0AF0X6waWTN\nWBl3r3KtqVrWKUQaEnZRNaH9kSWfTBhxE3YMlonS6gljz+MILZy4cnFJyuI8/VpSSrArjecXIoks\nwn4rcAXwHHBeSrlVwMPAlcB3qm+amGmkhU6GI2OLUkkjFDFLR3DppfF53EttsWfpByzCZmBgssee\nFCmThanKDS9BF7Umi7DfBnwB+GpKmRbgc8BmoKkG7RIzEF/A0raz8ycpS9Vn+OI2MpJ8TdyI258I\n9aNrws7BT1FQiUiXEzqatX4/B49WvYpakUXYHwDaS5T5U+BbuFG7mAMkhRJWIkgmbh0dbtS9eLGz\nY2zkncslLzpKi4IJy9hiprhImbADqFZUS3UCYbuT5iSEqIRaeOynA+8B/itO2Mte/ipmH0kTmOVM\nKPp/bfXo8uVR3nc/s2SWUbQ9RZQ7sRl2UllWr1bricdNKvt/haiGWgh7F3A9TtCbSLFiOr1hSz6f\nJ69/xbOK0M9O2s7OL5/0f3EYwQLRvqqDg8W+vR/dkmZXhNv2hdE3vb2u4xgcjDbvLhU1kybsWSdd\nFREjyqFQKFCoMpdHLYT9AuAbE++XAJcDo8D3woKdpf4rEjOacjxmKE+s+vpcff39btJz/nx3rcW1\nb9jgyq1blzyJmiTA/l+zd3za25M98Vrkksky11ArYVcHMfsJB70b/X+QGamFsP+m9/424PvEiLqY\nu/jhfEkTisuWRYK9bBls3hydDzf18Os0LHVvoTC5br9sR4cru3698+7jNvP2xTx8MrGycd+h0rmF\npLrKFem0Vb9ibpFF2O8A1uBG47uADUDbxLmbp6hdYoaTZrGEYuXbLnEjfvO5LTJk+fL48v59QwEL\nI2DCyVC7R6HgxNo8fHuF5f1wzfAJIUtkUBayTLCWK+xZykj4G58swn5VGfVdW2lDxOyinJDAJC86\nKTJk/vzJe7fGvU+zT9LO2QYiw8PJ36GrK/68L4z9/fHX15sw6seI6xgr6SzUEcw+tPJUTAml4rP9\nFxTbH/feGx+SmGZdQPEoO2003NnpRunmt69cWdx2f2I0Thj9e4ZUKpxZvmOWztFW/dYqbFLCPjuR\nsIspIYwaSbIzDLND/A4hrs4ksfbDJ9OEyLdkhobcMetU1q8vtm18K8nubfZLPh+/srUaYQ9FPvyO\nccQJeFrMvBZBzQ0k7GLKKEc0/OiVSiK9Qosmyb4xbHWrhUB2d8NTT8GTT7rjQ0PueKEAr399lJFy\n0yZYurS8tmUV+1otjIo7lqWzyDLJrY5gdiBhF1NOubHe5Zb3PeaNG50d4Vspoc3R3g5798KiRS4y\nZmjIHbPUwP51JoKWlMz8+Y99zHnsYd1ZoluyhlaWK6LVCq/l0Sn3qUHMPCTsYsqpRKCylglXiZrw\n+iN//2nAr7uvz43c29snXxO2w6wiPxzywIF4T78c+yQttLJcH7/cDjFE+782DhJ2MavxfWnf188y\nyrz++uI6/GvWrnV/OzqitMGFgouUGRpyIrhiRbZJzVKrZZPCNJOIE/asNk7cdf7kazhnMJ3WiyZu\nK0fCLhqKXC45jNHHH8XDZF9//froeDgSN7GPyzcf93SQNBls7a0G3xefChGUsM9OJOxi1pE2Gh4e\nzj5y9f10i0lPGg37hCtRw3q7uqIIm6TRvB33O6FS0TxJvnzaxHEpsnZAYnYhYRezjlJhj5WM8kxg\n/bBJiCZj29ujRUuDg9GmHmG4IrhJyJUro44jtDgq8cjjvrMJsb9WwF/lWy7l7v86FSNqReTUBgm7\nmPWUIwZhWVuk5E8cxk3GtrdHTwP5fPwOT1avXy4ux42fiqAS8fKvD3O5++0vF387wixMhbArIqc2\nSNjFrKbc8Dy/bE9PJOibNkX7rZrYh0KZVm9Xlxupj4y4UMr2djfJapE3Ntq3Zf9htsquruyWiqUc\n9r9TWrhlVsoVafvtPv5xuOmm0uXlmdcPCbuY1VQjFLYzk6UYMFsjaSQed99CwQl6f3+0ktXqzuWi\nyBu/c4gbyff0FPvyaXuzxk3Yrl0bpUZIi9SpVlj9Jx4Lz7z9dnj/+0vXXYn/LypDwi4aiizi0t0d\njc57e6MwxiTbxvfUBwedgFrUTBi/bukK+vsjEfcnJTdudGGSVk9S2KJvDaUJoj8JWyqtQrnHk+4X\nTrZaJ1TrEXmWuvQUEI+EXTQUpf4jD4XY/sYJhF/WH9GnhVMODsK55xa3JRxBm9WzebMbqYPrYFau\ndOf6+oonXMPIl7CNPlnCPeNWwmYVR+vs7r4bHn3U3e/AAZdD/+BBZxP5tsxUTYZOdZjnbEfCLgTV\nTzzatbbtXlpdFn2yenXUqXR0OGEfGIjy2MTt9uQTPn2AS5Nw4EB07b59LlUC1CZyxp8Utk1LbPI5\nn5/stSfNf1SSD8inXoI+WzsOCbuYs5QzSu3ri0IX7ZjtuGSfIfKdw0yNJsIPPgg7d0Z19/Q4Qc/l\n3HmIwiVN1JNSDcQ9ffj3suPhytbQ3690JG3fc3jY1RW3ctXak2Q5pZFmHVnnlNb2Ws0pSNiFmEVU\n4iuHohhXV6mUw3EpBLq63F/z1i0OHkovHArTAvg7UIVlQjHM5VyETiUTrnb81a9Obptfhy/0WUhL\nnmZhnqXSOvjtnEtI2IWYAdhkrI3My4kpt6gYWz1rE8IWWmmsXRvdp9QeskYWYX/kkahsqbhzs55K\n7fSUdn1YV9q1lXjxjbBISsIuRBlkmZwtt44w+sY8a/OvkywOH1+ww3DKLNkm/VF0JfaDWVP+5yRx\ntL9x7Uq7Lq5d4RxEeL2lcE6bgA6x39jqmY2LpCTsQpRBKQ83y+gxKfqmszMSopCs4pIlKiaXi8Ru\n06ZoBG0Trn198eKayyXnwCkUircY9L+niaq/iCrMmR93HSTvkmWrb5PmH/xVwzapXW6HVe0Eb1q9\nUz3yl7ALUSGVCMVUlofJ8fhxI10b3UMkkL6NYqN835e340ntCiN4/DI2b+CLdNYoJH9Tcb+T8Vff\nxs135PPFdpZvT5VjqVSafTPt/zsJuxANTlzIYj6fbcONJEJhLxXamJS/xiJe7H3cdX400IoV0T19\n8Zo/f7KYxbUrLGNPH+Xmqzc6OqKOKYv/b20IrZzh4ehpJSvTHU0jYReiDMqdWCtVPi3ixvfVp6pN\n1plYbhyzLYyBAXdsxw6XAye0bUx8ly93ZdeujcJDrf333uti8+OsF18AQzHMmlc/SUQrEda43yzL\nPEUp6j0hm0XYbwWuAJ4Dzos5fzXw50AT8ALw34DttWqgEDOJNA+4FuUryfwYdw8TEvO0DX/kGfrQ\nhYKLrDHR7u2NskbaxiJ+Ph1frIaG3CbfNol6771uNergoDtuGTItvt7aFWeN+G33o2fibJG0J5Kw\nnrCjTLNKwvv7Nk4SSbaR34F3dbnw0qmekM0i7LcBXwC+mnD+aeAtwAHgMuBLwOqatE6IOUQoOL7I\nV1JXXIcSRteEUTQQ5bAZHCzOfnn22VG7wpH20qVuRP+xj7ljl14adQbWQdj38dtl9ZvQ2z396+K+\nTzW/gb0v1VHa2gLI5tH76ZrD+xk9PfWxaLII+wNAe8r5h733W4EzqmmQELOFarzvLOcr2Vw6ydZI\nK58W6ugLYE9P8UjdH4kWCm6EHhdm6U+s+tExNiL2UxjbiB6ia8oNi/T/ZiXu91q/vtjjT2pH2jxB\nWLaStlVCrT32PwLurnGdQsxIai3sIeXuaOTfwxcYC2+EbBuRhJOchlkhoTjD5BWvfjv8UX9477io\nm7jvE0fSpGrcaDyXKz9GvtT9wXVEcXHzvm308Y/Dk0+632FoKGqrLRibionWWgr7W4EPAxcnFej0\nfv18Pk9+OqeNhZiBhJEm1UTHGGE0R1o6BJ/2dmdH9PS4kevAQNQeyz1vom1ePTix8xOhmdCm7ToV\n146075yWrz6O0CaB+Bj5kLgnAF+IbdK5UHBRQTZfYWmU83mXhM2PzgnvNXnUX6BQZRB9rYT9DcAt\nOI99f1KhzrguWYg5RKnRWbmTrWHdlURehNf5o3A/j4wJk7UrbpRsx7OKbtLTQdxna2sYueM/Qdh3\n8M/F1WPzBFmv8T93d0ergv3rYPLEaFdXJOobN0ZPNWnhk27Amz/2vTaGN8lALYT9LOA7wDXAUzWo\nT4iGZSrjm7N0Ckl2Q9JEa6nkY/6INmmQmTbhWOq3CCdprdPwo1xMJMOok6TIoHAUnjbZGYfv/fsT\n3OHmKP48Qj7vnmpe/3qOpTxOs8fse1c6cM8i7HcAa4AlwC5gA9A2ce5m4C+AxcAXJ46NAhdW1hwh\nhDEVHUCWOpM6H4tRD8XHMJGLE9Gs987aFh9LfgbFUSdZrJe4+pOOxXnpDz4Y1btpU7SHbfj04oeJ\nwmRbqtZmRhZhv6rE+Y9MvIQQMVRqkVQj7JVe6wt2XJSIPzka+uOhiCeNlpPa1tUVb0+Ek7T5vAuv\n9C2i8DvEEf7/YKGVaZEt4W8Rxv4vXDi587D7dHc7sbdOYN++aO7B32rR//7Dw5Nz+FSCVp4KMcVU\n45tXQjV2jy/YSd50lgibSr6zv6F3nIftWya3317ccViYZFzUid8ev8OxNA5pE7i+JWLX+zloBgai\nlBDhRHc4QKHIAAASYUlEQVRn5+Rc99Y2s4s6OqL6/NBKfzFYBRa7hF2IRqMSYS8VF+5H6vhx5+WG\nKGbBhDT0sP12XHJJsbdtwp7kSyd1Ov6kcNwThj/Jatf4OWgMi45J+j5JlJNOuBwk7ELUkamaOK2W\nUiPspMnVpLrSPvtYOCVE9ga4xGEweULS9pT1Q0F9j70U9jQAxTaJb+lYe6yMvzrWfqf+/skWEUR1\n29PGpZe69u7Z4/LlWLmlS919LBePXwcUrxnQiF2IGc5UCXulPn41VDMRa6xfH9kllsoA3Hfo7CyO\nBw8tjrholLVrS7fZ6vMTn/nYgq5QuO1zLueEffNmV4d1EFCcTtgIJ0n7+txG5oOD0dqA9nbXmWWd\nLyiFhF2IBqCWPn7WCd1qhD3ueJygDgwUWyF2jW/RhAu5shDOI4Qi3tHh7mELs/ycNYXC5F2uIPk3\nD7+rZbrs7HRx7ean+xuVVIuEXQhRxFRF6vjECbuNtsNol7hIk3XrittjZbJGxKxZ4wQ1fBow7N5+\n7Lw/8eqnUDB7KLyXRbf09LjvZnVYSuNCwU322j0OHoxveyVI2IVoMKbbx6/UFvJDHeN8+jgbxu5n\nf0tNYIbhiv6ka9p1lmumr89du2iRE2Xz4/v7iydk7b3feZjY33uv89wHBqCtzdU3PAw7d7oOIski\nKgcJuxANxnQLe6nRL2Rfwp+2NV3aitakctY52N+49Af+yByKJ3Wvv97541Zv2mpbf3Rv5Wxk393t\n7t3T475jPl+ci6daJOxCiLpQyTyATZ76TwAWZulHq6QthApH8kkhkVnbWs5yf8sr79PfH1kxvn+/\nYkXt1jhI2IUQU0YtImeSRNb3vsNzpUhKf+BHzZRqT1oETtJEL7hJUn/ewMpWkn8/CQm7EGLKyLJo\nJ25Jf5xHn8ulr3z1KeXzl9MZJNlEWTe3Np/dX+Rlx/0FSv4EbLVI2IUQdafSEXoYJ562SKoc8Z7K\ndQAWnx8X356WzqAaJOxCiLozFYnR4kbtpcQyLpIFsnUEWYXYNrhOu3+tkbALIepOViEtZ1u/kCxR\nM5UKazX5eCy+3RKH1eKpIETCLoSoOdVkmPQpta1flvuV245atLvUE0mp1arVImEXQtSccoS9VkIa\nxp9X6pVnXb2aVnel0Tq1QsIuhJhWatkBTKWg1qLu0I6B4pQEtbJkJOxCiJow1Rkmyw1pLJda2Uch\n4eIoP91AuH2ehF0IMaOot/1Q6n6VTG5OhX00HSkeJOxCiIZkKgW12rrDTUMsp0ytImUk7EKImlPv\nUWo1YYv13qAkrN+3ZJQrRggxY5ktwj7d0StTRfN0N0AIIeYy/urXWpFF2G8F9gJPpJS5EfgZMAD8\ndg3aJYQQdWU6Jjn9+9Zb2G8DLks5/07gt4BXAdcBX6xBu4QQoq5Ml7BPBVmE/QFgf8r5dwMTofZs\nBXLA0irbJYQQokJqMXl6OrDL+/wMcAbOvimieWMz44wD0NbcxrbrtvGGpW+oQROEEEIYtYqKaQo+\nj8cVGr8/OjzaPsqKv19BE02MM04LLcxrmcfI0ZFj4g9wwakXcN+H7iM3P2XzQyGEaBAKhQKFLPvu\npRAKchLtwPeB82LO/T1QAL4x8XkHsIbJI/ZxOstu3yTmNc/j5bGXj31uaWrhkjMvoeeqHom/EKLh\naGpqguxaDdRmxP494OM4YV8NDBNjwwC845XvYPToKH27+7jgtAt4aNdDReebaWaMsdSb+aIOcHT8\nKL2/6OUVn3sFbc1tHB47TBNNnDTvJEbHRmltbmVe6zy2fXQby3PLK/+WQggxS8jSC9yBG4EvwQn2\nBqBt4tzNE39vwkXOvAhcCzwWU8/4+HhksQyPDNPR08Gh0UM8sfcJzj/1fG54+w1c+rVLOXTkEL8+\n/OvKvlFG2prbGB0bPfa5pamFx/74MXn+QogZRSUj9rIKV0mRsJdiaHiI1V9ezdHxo4weHWX48DB3\nvv9Ouge6eeSZRzgwcoDc8TmGR4YZHRuluamZsfHi0X5LUwtHx4+W3dBWWjnCEcB1AKctPI3hw8O0\ntbRp5C+EqCsNJexZGRoe4pJbL+Gfrv4nPn3fp3lsz2Ocd/J5nDDvBDa+dSNXfP0K9hzcw5HxIzW9\nr0362vsT553IhWdcyDc/8E15/UKImjEnhT0L2/du58JbLmTV6at46vmnOO+U82hrbePXI7/mwV0P\npl7bRBMnHXcSBw4fyHQvX/CN1qbWoo5lYetCVp+1Wp2AEKIkEvYKMdtnjDFWnLKCtpY2Do0eou+Z\nPrZ+dCuLjlvEOTedw+Gjh1PryTL56zO/ZT4vH3256JrWplZOmHcCq05fJeEXQkjYp5LhkWGu/s7V\nbNu9DZpgxSkroAme2PsE553irJ8tv9jCf770nzW757ITlvE7Z/8OQweGeHT3o4yOjRZN+IJ7Qthy\n7RYuOeuSmt1XCDFzkLBPM0PDQ1z0lYs4Z8k5PLH3CV54+QXefMabWdC2gAOHD/DwrocZn/gfwMJ5\nCzn48sHUOk9ecDK/OvSritt0zP8//UK+eaWeAISYbUjYZwEm/itPXcnfvfPv+Oj3P0rvYC+L5i9i\n5dKVLvLm2W089+JzvPG0N5I7Lsd9P7+v4gifrLQ0tbBo/iIeu+4xRf0IMYOQsDcIwyPDXPf96/jS\nu74EwHXfv47//Zb/zeW3X85vLv5Ntv1yG1+84otc+71rp6V9bzrtTWz+w80a/QtRByTsc4yh4SEu\nvOVC9o/sZ2xsjNaWVr723q/xobs+xAWnXcDO/Ts57+TzWDBvAd/7yfcmRetUS1tzG+865108f+h5\nFrQt4GfP/4yhA0OMj4/T1tJGW3ObVv0KUSUSdpHIg794kDW3reEb7/8Gn9j8CRiHW99zK+/9x/ey\n6vRVPPLMI5PSNWTBD+9sbW7lyNjk9QJNNPG2V77tmMd/7k3n8tPnf3rsumaaWXz8Yv7tun9TByBE\ngIRdVEyY4uHcJecysHeA/SP7WTx/Ma9+xat5dPejvOWst/DQMw8xOjZKE01cdOZFPLTrId542ht5\nev/TqVFBH3jtB7jzA3eS+6tcpnUBrzj+FbQ0tfDcoeeOHZMNJOYaEnZRF7bv3c6bvvwmtn5kK2ct\nOuvYfMCBkQNuYnjZSg4cPlCU5G3lspXcv+5+cvNznPzXJ7Pv0L6q2tBCC/Na5zF6dJSj40cn2UyL\n5y/m8T9+XE8AYtYjYRczBnsCGB0bpa25je613cdG2ccWhI2PMTo2yv6R+A26FrQu4NCRQ1W1o4km\nzj/1fPa+uJd9h/YxcmQktpzy/ouZioRdzDps4Vf/s/18+d1f5j3feA+jY6NcfObFfP33vs76zev5\n1Yu/4uFnHmbR/EVFncBUhICG6R/8eykcVEwHEnbR0Fgn8NgvH+P8087nhrfdwOW3X87rTnkdo0dH\nefiZh1nYtpC2ljb2HdpXJNCLjluUOd9Pudz9wbu5a8dd3Pr4rUUdzQWnXsApJ5zCP7zvH/QkICpG\nwi7EBDYP8MM//CE3br2Rv37HX/Mn//QnbHt2G2NjYxwZO8Lw4eGa3W/N8jX0DvUmnm+iibbmtqLI\nIz+iSJFBIgkJuxBlYn4/TfDND3yT3//m7/PqV7yaHft2cO6Sc9mxbwevOfk1/GTfTzh05BDDI5M7\ng7s/eDdfeOQL3PPUPVPSxoVtC5nfNl/rAeYoEnYhphBbEDY8Mszxbcdz4PAB7v7g3Vz+qssZHhnm\nyjuvpPcXvTAOP1r3I27ceiM9O3omJW6rlDNOPIPLX3U5P33+pyxoW8DJC04+liDupSMvTYoMOuWE\nU2iiiZePvswFp12gbKGzFAm7EDMM3xK64YEb2PbsNhiHFctWwDg88dwTHD56ODEyyFjQtoAf/8mP\nWdez7pjlU26CONtlrLWplYvPupjdv97Nzv07izoEO9fzB9ocfqYgYRdiFjI8MsyV37ySws8LNDc3\nc9fv30XHdzs48NIB7vqDu7ju+9fx4IcfZHluOe/8+ju556l7MiWIW9i6kINHXPbQhW0LOThanEk0\naaVwHP7Wky20cBR3v7bmNrZdt017BU8hEnYhGpykBHGX3n4pI6MjHDh8gNbmVvLL89zy7lv45OZP\n0kQTB18+yH0/v+9YPeedch67X9hd0/0DjLamNsYZp7mpedJk8byWeaxctpLfOP43FC2UEQm7ECKW\n4ZFhrv721Tz27GOcf+r5fP19X+fAyAFWf2U1Y2NjvPbk17Jj345JYaJTTVtzGyfMO4EXDr/A+Pg4\nLc0tvOucd/HdHd8tegp54NoH5uxmMhJ2IURV2Obwd7z/Dt535/sYfmmYH637EZ8tfJYtv9jCm09/\nMwuOW8Chw4fo293nUjoQv0gsbv/fLNtHxl2Xdjw839rUysVnXkzPVY0xTzBVwn4Z0AW0AF8GPhec\nXwLcDiwDWoHPA90x9UjYhWgwtu/dzqqbVzHW5NYG3Pn+O+ke6Kb/2X7uueYePnnPJ9kytIVVp69i\n14Fd3HPNPaz60qrETKJ+YjmfE1pP4MUjL5bVtqYJeTuWfbSplXef++6ipwGzh0467iSOjB3hglMv\nmHE7jU2FsLcAPwHeDuwGHgWuAv7DK9MJHAd8GifyPwGWAuHznIR9gkKhQD6fn+5mzAj0W0TMld/C\ndhF73cmvo63V5e3/1EWf4h1fe8exxHJv/+zb2T5/O0fGj7Dl2i289x/fW1biuEpH/X65JpqKni4W\nz1/M4vmLeXr46aKyU51wrhJhby1x/kLgKWBw4vM3gPdQLOzPAjYlfhLwPJNFXXjMlf+As6DfImKu\n/BbLc8vZ/We7Jx1/6X+9dOz97877Xbb9xbZjn7d9dNux+YBXLn4lj+5+lFWnreLp4addVtEzLqKt\npY0nfvUE5y87n627t/L8S89Pukfc00CcPeTvTWzsH9kfG5a6f2Q/7X/bDringiUnLGFsbIx9h/YV\n1Rt2KlOZgrqUsJ8O7PI+PwO8KShzC/AvwC+BE4Era9Y6IYTAdQbP/tmzmcv7TwVHxo/wr7v+lUc+\n+ghnLTqLK++8ki2/2MIbT30jPx/+Ofdccw8b7t/AOOPc/bO7ExeUzW+Zz/Ftx6euOTgyfoQ9B/fE\nngs7iq2/3Mrizy0uOtZEE1e8+gpeOPwCO/fvZPmiyp4CSgl7Fu/kM0A/kAfOBn4IrABeqKhFQghR\nJUlPBQD//KF/nnTsrj+4CyheUPbZ3s+yZWjLscRyfR/pAzi2+vjE407khcMvTJoviFszkJVxxvnB\nT39w7PMzv36monpKsRrY7H3+NPA/gzJ3Axd7n38EvDGmrqdwHYVeeumll15ZXycRv4lAFbQCO4F2\nYB5uZP6aoMzfABsm3i/F2TW/UeuGCCGEqB2X4yJdnsKN2AH+eOIFLhLm+8AA8ATwwXo3UAghhBBC\nCFEGlwE7gJ8x2Y83bpw4PwD8dp3aNR2U+i2uxv0G24GHiMJGG40s/yYAVuFCZX+vHo2aJrL8Fnng\nceBJoFCXVk0PpX6LJbg5vn7cb9FRt5bVn1uBvTjXI4lp080WnGXTDrQR78m/EzfhCi50sq9ejasz\nWX6LNwOLJt5fRmP+Fll+Byv3L8APgPfVq3F1JstvkQP+HThj4vOSejWuzmT5LTqBGybeL8GtkSkV\nyTdb+S84sU4S9rJ0s7l27QKKFzSNEi1o8nk3sGni/VbcP+SlNW7HTCDLb/EwYBtxbiX6j7mRyPI7\nAPwp8C0ge4Lx2UeW3+KDwLdxQQgA2Zdbzi6y/BbP4hY9QuMvfnwASEvKX5Zu1lrY4xY0nZ6hTCMK\nWpbfwuePiHrkRiLrv4n3AF+c+Dxeh3ZNB1l+i1fhosruB7YBf1ifptWdLL/FLcDrcIsfB4BP1qdp\nM5KydLPWjzVZ/4MM8x404n/I5XyntwIfpng9QKOQ5XfoAq6fKNtEfbOO1pMsv0UbcD7wNmAB7qmu\nD+etNhJZfgstfiwms27WWth3A2d6n88keqRMKnPGxLFGI8tvAW7C9Bacx56+P9rsJMvvcAHuURyc\nl3o57vH8e1PeuvqS5bfYhbNfXpp4bcGJWaMJe5bf4iLgLyfe7wR+DpyDe5KZa0yrbmZZ0ORPAqym\nMScMIdtvcRbOZ1xd15bVlyy/g89tNG5UTJbf4lzgPtzk4gLcZNpr69fEuqHFj5NpJ9vk6bToZqkF\nTQA3TZwfwD12Niqlfosv4yaEHp94PVLvBtaJLP8mjEYWdsj2W3wKFxnzBPCJurauvmjxY8QduLmE\nl3FPbR9m7uqmEEIIIYQQQgghhBBCCCGEEEIIIYQQQgghhBBCCCGEEGI6+f8HcIfdVfbC9gAAAABJ\nRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# manually fit the curve to obtain a viable set of starting parameters\n", "at,bt,ct = (0.5,0.5,0.5)\n", "plt.plot(xvals,yvals, '+')\n", "plot(xvals,func(xvals,at,bt,ct), 'g.')\n", "print (\"a=%s, b=%s, c=%s\" % (at,bt,ct))\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Licence and version\n", "*(c) Stefan Loesch / oditorium 2014; all rights reserved \n", "([license](https://github.com/oditorium/blog/blob/master/LICENSE))*" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3.4.0 (default, Apr 11 2014, 13:05:11) \n", "[GCC 4.8.2]\n" ] } ], "source": [ "import sys\n", "print(sys.version)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.4.3" } }, "nbformat": 4, "nbformat_minor": 0 }