{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": true }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "#Lmfit provides a high-level interface to non-linear optimization and curve fitting problems for Python.\n", "#It builds on and extends many of the optimization methods of scipy.optimize.\n", "#The advantage is the ease on parameter settings, change of fitting algorithms,\n", "#estimation of confidence intervals, and multiple model combination and a set of ready-to-use fitting models.\n", "#I'll show some examples like line deblending and continuum estimation." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": true }, "outputs": [], "source": [ "#let's start with a simple one\n", "def parabola(x, a0, a1, a2):\n", " return a2*x**2 + a1*x + a0" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[-8, 11, -10, 60]" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGUJJREFUeJzt3X+MJ3V9x/Hn+764SBCFW/A8gRVMqcQfUWFDXbXt1kWl\n1HjYpBfatJyVeDEpxkvayG0M1uYS92xTcza1LYeCd61ViJZyof6CrRvT8AXZo4DAQQ8t4JE9Dr9S\n0TRhvb13//jO9zp8b+b7ne93Zr4z853XI9ns9zvf+X7ns7Pzfc9n3p8fY+6OiIiMv3VFF0BEREZD\nAV9EpCYU8EVEakIBX0SkJhTwRURqQgFfRKQmMgn4Zna6mX3NzB41swNmNmNm683sDjM7GPw+I4tt\niYjIcLKq4X8O+Ja7Xwi8GTgAbAcW3f0CYDF4LiIiBbG0A6/M7BXA/cBrPfRhZvYYMOvuK2a2EVhy\n99el2piIiAztpAw+43zgWeAmM3szsB/4GLDB3VeCdQ4DG6LebGZbga0Ap5566sUXXnhhBkUSEamP\n/fv3/8Tdz+q3XhY1/GngbuAd7n6PmX0OeB74qLufHlrvOXfvmcefnp725eXlVOUREakbM9vv7tP9\n1ssih38IOOTu9wTPvwZcBDwTpHIIfh/JYFsiIjKk1AHf3Q8DPzazTn5+DngE2AdsCZZtAW5Luy0R\nERleFjl8gI8CXzazCeBHwB/TPpncYmZXA08CmzPaloiIDCGTgO/u9wNR+aO5LD5fRETS00hbEZGa\nUMAXEakJBXwRkZpQwBcRKUiz2WRhYYFmszmS7WXVS0dERAbQbDaZm5tjdXWViYkJFhcXmZmZyXWb\nquGLiBRgaWmJ1dVV1tbWWF1dZWlpKfdtKuCLiBRgdnaWiYkJGo0GExMTzM7O5r5NpXRERAowMzPD\n4uIiS0tLzM7O5p7OAQV8EZHCzMzMjCTQdyilIyJSEwr4IiI1oYAvIlITCvgiIjWhgC8iUhMK+CIi\nNaGALyJSEwr4IiI1oYAvIlITCvgiIjWhgC8iUhMK+CIiNaGALyJSE5nMlmlmTwA/B9aAo+4+bWbr\ngZuB84AngM3u/lwW2xMRkcFlWcP/LXd/i7tPB8+3A4vufgGwGDwXEZGC5JnS2QTsCR7vAa7IcVsi\nItJHVgHfgTvNbL+ZbQ2WbXD3leDxYWBD1BvNbKuZLZvZ8rPPPptRcUREyqHZbLKwsECz2Sy6KJnd\n8eqd7v60mb0SuMPMHg2/6O5uZh71RnffDewGmJ6ejlxHRKSKms0mc3NzrK6uMjExweLi4kjvcNUt\nkxq+uz8d/D4C3ApcAjxjZhsBgt9HstiWiEhVLC0tsbq6ytraGqurqywtLRVantQB38xONbPTOo+B\n9wAPAfuALcFqW4Db0m5LRKRKZmdnmZiYoNFoMDExwezsbKHlySKlswG41cw6n/fP7v4tM7sXuMXM\nrgaeBDZnsC0RkVJqNpssLS0xOzt7PG0zMzPD4uLiCcuLYu7lSZtPT0/78vJy0cUQERlI0bl6M9sf\n6hIfSyNtRURSKluuPo4CvohISmXL1cfJqlumiEjthPP2ZcrVx1HAFxEZQlTefn5+vuhi9aSUjojI\nEKqStw9TwBcRGUJV8vZhSumIiAyhbH3sk1DAFxEZ0szMTCUCfYdSOiIiNaGALyJSEwr4IiJ9lGlO\n+zSUwxcR6aHoeXKypBq+iEgPVexvH0cBX0Skhyz72xedGlJKR0Skh6z625chNaSALyISofuGJoME\n56iboUSlhhTwRUQKlqY2HvfeTmqos7yIqRiUwxcR6ZKmoTbuvZ3U0I4dOwrr6aMavohIoJOKmZyc\nHLo23qsmX/RUDAr4IiKcmIrZtWsXrVZr4IbaMk+qpoAvIsKJqZhWqxV5Q5OoBtluRdfk4yjgi4jQ\nOxXTUYaulWlk1mhrZg0z+08zuz14vt7M7jCzg8HvM7LalohIGlEDoJI0qlZ91G2WNfyPAQeAlwfP\ntwOL7r7TzLYHz6/NcHsiIgPrVUvvl4opQ9fKNDKp4ZvZOcDvAF8ILd4E7Ake7wGuyGJbIiJppKml\nl6FrZRpZ1fB3AR8HTgst2+DuK8Hjw8CGqDea2VZgK8DU1FRGxRERiZa2ll7WBtkkUgd8M3sfcMTd\n95vZbNQ67u5m5jGv7QZ2A0xPT0euIyKSlTJ3m8xbFjX8dwDvN7PLgZcCLzezfwKeMbON7r5iZhuB\nIxlsS0QktSrX0tNIncN393l3P8fdzwOuBP7d3f8Q2AdsCVbbAtyWdlsiIjK8POfS2Qm828wOApcG\nz0VEpCCZDrxy9yVgKXjcAuay/HwRERmeZssUEakJBXwRGStF30awzDSXjoiMjarPdZM31fBFZGxU\nfa6bvCngi8jY6IyibTQaiUfR1ikFpJSOiIyNQUfR1i0FpIAvImNlkFG0USmgcQ74SumISG0NkwKq\nMtXwRaS26jaRmgK+iNRanSZSU0pHRKQmFPBFRGpCAV9EpCYU8EWkFuo0wCqOGm1FZOzVbYBVHNXw\nRaSSBqmxa46dNtXwRaRyBq2xdwZYddYf9wFWcRTwRaRyBp0SoW4DrOIo4ItI5QxTY6/TAKs4Cvgi\nUjmqsQ9HAV9EKkk19sGpl46ISE2kDvhm9lIz+76ZPWBmD5vZXwTL15vZHWZ2MPh9RvriiojIsLKo\n4b8AvMvd3wy8BbjMzN4GbAcW3f0CYDF4LiKSOY2iTSZ1Dt/dHfhF8PQlwY8Dm4DZYPkeYAm4Nu32\nBtFsNtWoIzLmNIo2uUwabc2sAewHfgX4vLvfY2Yb3H0lWOUwsCGLbSWlg0CkHup2m8I0Mmm0dfc1\nd38LcA5wiZm9set1p13rP4GZbTWzZTNbfvbZZ7MoDqCh1CJ1UbfbFKaRabdMd/8fM/sucBnwjJlt\ndPcVM9sIHIl5z25gN8D09HTkSWEYGkotUg/qk5+ctSvfKT7A7Czgl0GwPwX4DvAZ4DeBlrvvNLPt\nwHp3/3ivz5qenvbl5eVU5QlTDl+k+vQ97s/M9rv7dL/1sqjhbwT2BHn8dcAt7n67mTWBW8zsauBJ\nYHMG2+qr++CIOkB0AIlUg9rispVFL50HgbdGLG8Bc2k/fxBJDg4dQCLVkbZBVpW7FxurkbZJGmrV\nmCtSHWkaZDuVu+uuu465uTn10WfM5tLp1VDbOdNPTk6qMVekItI0yKq75onGKuDHHRzdaZxdu3bR\narV0mSdSAcNOkqaeeicaq4AP0QdH95m+1WoxPz9fUAlFpJ+43PsgOXl11zzR2AX8KDrTi1RHXMeK\nYTpcaArlFxurRts4nTP9jh071CtHpOTiOlaow0V6tajhg870IlURd0WuK/X0ahPwRaQa4nLvysmn\nl3pqhSxlPbWCiFSHBkkNb5RTK4iInGCQAK4R8KOhgC8imRs0gGuQ1GjUopeOiIzWoD1qNKf9aKiG\nLyKZG7RHjRpkR0ONtiKSi3AOH1Awz5EabUWkUJ2xL2qQLQ/l8EUkVxohWx4K+CKSKzXIlodSOiKS\nq14NshpsNVoK+CKSu6i5rJTbHz2ldESkEMrtj54CvogUQrn90VNKR0QKocFWo5c64JvZucBeYAPg\nwG53/5yZrQduBs4DngA2u/tzabc3KmpMEsmf7lMxWlnU8I8Cf+ru95nZacB+M7sD+CCw6O47zWw7\nsB24NoPt5U6NSSIyjlLn8N19xd3vCx7/HDgAnA1sAvYEq+0Brki7rVFRY5KIjKNMG23N7DzgrcA9\nwAZ3XwleOkw75VMJakwSkXGUWaOtmb0M+Dqwzd2fN7Pjr7m7m1nkLG1mthXYCjA1NZVVcVJRY5KI\njKNMZss0s5cAtwPfdvfPBsseA2bdfcXMNgJL7v66Xp+j2TJFyiurjgzqEJG9kc2Wae2q/BeBA51g\nH9gHbAF2Br9vS7stEYmXZyDNqiODOkQUK4sc/juAPwLeZWb3Bz+X0w707zazg8ClwfPSaTabLCws\n0Gw2R/pekSx1Aul1113H3NxcqmMy6rju1ZFhkO+BOkQUK3UN393/A7CYl+fSfn6e0tQ2VFORMsnq\nnrBxx3XcHay619+1axetViv2KmPQO2FJtmo90rb7S7J3797El8S66bKUSVaBNO64juvIEF7/hRde\n4JprruHYsWOxlSB1iChWrQN++EvSaDS46aabOHr0aKIau2oqUiZZBdJex3XUqNjw+mbG2toax44d\n61kJ0uja4tT+nradhq6nnnqKG264gbW1NRqNBjt27GB+fj7Re5N+wdQ7Qapg2ON6cnKSbdu2Kc1Z\ngKS9dGof8DuS5uSHDdrK+UsdqFJTDN3EfEBJLonTBG3l/KUOlK4pNwX8kH4Ha5qgnTTnrxqSiORF\nAX8AaRpq876CEBklVUyqSQE/gfDBnaYnRJ5XECKjoopJdSng9xF1cEf13gmfFIChTgrq6ilVoIpJ\ndSng95Hk4A6fFBqNBmZ2vD9/v5GHYRqUIlWgikl1KeD3keTgDp8Ujh07BoC7Jxp52E29HGRUhs3D\nq2JSXQr4fSQ5uLtH7HZq+ElHHoqMWto8vCom1aSAn0C/g7v7pABEjjzUpa+UhfLw9aSAn5Huk0Ln\n8Zve9KZSXPqqG52EaVxIPWlqhRpQN7p6iwva/YK5jpvq0NQKcpwu3+urV9DWuJD6UcCviGEurcOz\nGKobXb2EZ4HNezoQqQ4F/BIY9NI6Sd/+Yd4j46F7XMhJJ7W/5nlMByLVooBfsCR50vClddK+/d2X\n461Wq+/8/jIewv97gA9/+MNMTU0dD/YLCwuJA7i6X46XLG5iLilE5Um7dS6tG40G69atY21tre9N\noMPv0eV4vXT/76+66qrjJ/usbnQu1aQafsGS5EnDl9ZJ+/brcry+ktx/Vo2w9aRumSWgWyXKKKib\n5fjSLQ5F5ASqLIynkfbDN7MbgfcBR9z9jcGy9cDNwHnAE8Bmd38ui+2NA33xpAhxjbA6Hushqxz+\nl4C/BfaGlm0HFt19p5ltD55fm9H2Ki3vG6aLDEKpnvrIpJeOu38P+GnX4k3AnuDxHuCKLLY1DpL0\nzOl8CdWjQvKW5HiU8ZBnL50N7r4SPD4MbIhaycy2AlsBpqamcixOeQw6x36vHhVZ3Gkriq4u6kMj\nautjJN0y3d3NLLJ12N13A7uh3Wg7ivIUbdA59uO+hL3utJXmslyX+PWiLrz1kWfAf8bMNrr7iplt\nBI7kuK3KGXSO/X6jacN32krbx7pK/bV1JZINjaithzwD/j5gC7Az+H1bjtsaS/2+hHF32kp7WV6V\nS3xdiYgMJqtumV8BZoEzzewQ8Oe0A/0tZnY18CSwOYttyf+Lu9PWMDNqht9T9kv8LGaCLJKuSqQo\nGnhVY0XOqDls0Mur3SIPUX+juuRKHnQDFIkVVUNOOgtnVtsfNhXTaybIMgXGuL+xV/tI+P4F4fmS\neu0fnRhkEAr4NRM3V7qZHW/8zTs9kqZRuLt94aqrriploIv7G+PaR8L/FzPj2LFjff8XasOQQSng\n10xcDTnpLJxZSNMoXJb2hX4167i/MardZWFh4UVXW+vWrTueruq1f6rUm0pKwt1L83PxxRe75Ouu\nu+7yU045xRuNhp9yyil+1113vei1T3/60y9almc5RrWtrPXah93r9fobw58zMTHhJ5988vHPvP76\n6/vun6TlkPEHLHuCGKsafs30qiGn6YudJJfcvU4VaqNRf1fSmvUgNwmH5O0R4TKV4WpHKiTJWWFU\nP6rhl19UrTVJTbOKtdG4Mmf1twzzOVXcj5I/VMOXrA3T86SjivnmuDJn1Y4wzOdUcT9KeSjgS2KD\n9jwJq8ro3bBeZc4qJTXo51RxP0p5aOCVJNarG+AwOfy8y5rFtsrYz72MZZJi6RaHkosqBBv1T5e6\n0UhbyUUVeteUMc9dhROljD8FfBk7Zctz64pDykIBX8ZOWUbjdpTxikPqSQFfSiWr1EeZUk9lu+KQ\n+lLAl9IY19RH2a44pL4U8KU0yp76SHP1UaYrDqkvBXwpjTKnPtKOQRApAwV8KY0ypz7irj7GNQ0l\n40kBX0qlrKmPuKuP7hPB3r17S3nCEgEFfElgXFMWg/xdcVcf4RNBo9Hgpptuyvweu+O6/2X0FPCl\npyJz13l+/jCpmKirj/CJ4KmnnuKGG27ItNFZKSPJkgK+9FRU7jrvz88yFdM5ETSbTfbs2ZNpo3PZ\ney5JteQe8M3sMuBzQAP4grvvzHubkp2kueusA1Hen59HKiaPRucy91yS6sk14JtZA/g88G7gEHCv\nme1z90fy3K5kJ0nuOo9AlNfnR90eMMtUTNaNzmXuuSTVk+v0yGY2A3zK3d8bPJ8HcPeFqPU1PXK1\nxOXYyzoXfVyaqHv5rl27aLVaCrBSGWWZHvls4Meh54eAXwuvYGZbga0AU1NTORdHshRVm80y9551\nbTkuTRSuRU9OTrJt2zY1kspYWld0Adx9t7tPu/v0WWedVXRxJKWoxtCFhQWazWZhZWo2mywsLDA5\nOcnExASNRiPyloXz8/O0Wq0TTgoi4yLvGv7TwLmh5+cEy2RMjaJf+iAGTdeokVTGWd4B/17gAjM7\nn3agvxL4g5y3KQXKu1/6oLqvOFqtFvPz87Hrq5FUxlmuAd/dj5rZNcC3aXfLvNHdH85zm1K8PPul\n9xLVyDtMjb2s0zuIpKWbmEuuRjUtgGazlDorSy8dqblR1ZZ7DdRSjV2krfBeOiJJdHraxPX26aRu\nonrgiEibavhSekn69quxVaQ/BXwpXL8ce9J5dZS6EelNAV8KlaT2rr7xItlQwJdCJam9K10jkg0F\nfClU0tq70jUi6SngS6FUexcZHQV8KZxq7yKjoX74IiI1oYAvIlITCvgiIjWhgC8iUhMK+CIiNaGA\nLyJSEwr4IiI1oYAvIlITCvgiIjWhgC8iUhMK+CIiNaGALyJSEwr4IiI1kSrgm9nvmdnDZnbMzKa7\nXps3s8fN7DEze2+6YoqISFppp0d+CPhd4PrwQjN7PXAl8Abg1cCdZvar7r6WcnsiIjKkVDV8dz/g\n7o9FvLQJ+Kq7v+Du/w08DlySZlsiIpJOXjdAORu4O/T8ULDsBGa2FdgaPP2FmUWdQMLOBH6SuoSj\no/LmS+XNl8qbvyzK/JokK/UN+GZ2J/CqiJc+4e63DVqqbu6+G9iddH0zW3b36f5rloPKmy+VN18q\nb/5GWea+Ad/dLx3ic58Gzg09PydYJiIiBcmrW+Y+4EozO9nMzgcuAL6f07ZERCSBtN0yP2Bmh4AZ\n4N/M7NsA7v4wcAvwCPAt4E8y7KGTOP1TEipvvlTefKm8+RtZmc3dR7UtEREpkEbaiojUhAK+iEhN\nlD7gm9nNZnZ/8POEmd0fs94TZvaDYL3lUZczVI5PmdnToTJfHrPeZcG0E4+b2fZRlzNUjr8ys0fN\n7EEzu9XMTo9Zr9D9229/WdvfBK8/aGYXjbqMobKca2bfNbNHgqlHPhaxzqyZ/Sx0nHyyiLKGytPz\n/1uy/fu60H6738yeN7NtXesUun/N7EYzO2JmD4WWrTezO8zsYPD7jJj35hcb3L0yP8BfA5+Mee0J\n4MwSlPFTwJ/1WacB/BB4LTABPAC8vqDyvgc4KXj8GeAzZdu/SfYXcDnwTcCAtwH3FHgMbAQuCh6f\nBvxXRHlngduLKuOg/98y7d+IY+Mw8Joy7V/gN4CLgIdCy/4S2B483h71Xcs7NpS+ht9hZgZsBr5S\ndFkycAnwuLv/yN1Xga/Sno5i5Nz9O+5+NHh6N+0xE2WTZH9tAvZ6293A6Wa2cdQFBXD3FXe/L3j8\nc+AAMSPNK6Q0+7fLHPBDd3+y6IKEufv3gJ92Ld4E7Ake7wGuiHhrrrGhMgEf+HXgGXc/GPO6056k\nbX8wXUORPhpc9t4Yc9l2NvDj0PPYqSdG7EO0a3FRity/SfZXKfepmZ0HvBW4J+LltwfHyTfN7A0j\nLdiJ+v1/S7l/aU/SGFcJLNP+Bdjg7ivB48PAhoh1ct3Pec2lM5CE0zf8Pr1r9+9096fN7JXAHWb2\naHCWzVyv8gJ/D+yg/QXaQTsN9aE8ypFUkv1rZp8AjgJfjvmYke3fcWFmLwO+Dmxz9+e7Xr4PmHL3\nXwTtPP9Ke4BiUSr3/zWzCeD9wHzEy2Xbvy/i7m5mI+8TX4qA732mbzCzk2hPw3xxj894Ovh9xMxu\npX1plMsB26+8HWZ2A3B7xEsjnXoiwf79IPA+YM6DRGLEZ4xs/0ZIsr9KNZ2Hmb2EdrD/srv/S/fr\n4ROAu3/DzP7OzM5090Im/krw/y3V/g38NnCfuz/T/ULZ9m/gGTPb6O4rQTrsSMQ6ue7nqqR0LgUe\ndfdDUS+a2almdlrnMe2GyIei1s1bV17zAzHluBe4wMzOD2opV9KejmLkzOwy4OPA+939f2PWKXr/\nJtlf+4Crgt4kbwN+Frp8HqmgvemLwAF3/2zMOq8K1sPMLqH9XWyNrpQvKkuS/29p9m9I7FV/mfZv\nyD5gS/B4CxA1+WS+saGoVuxBfoAvAR/pWvZq4BvB49fSbs1+AHiYdqqiqLL+I/AD4MHgH7Wxu7zB\n88tp9974YcHlfZx2zvD+4Ocfyrh/o/YX8JHOcUG798jng9d/AEwXuE/fSTul92Bov17eVd5rgn35\nAO3G8rcXWN7I/29Z929QnlNpB/BXhJaVZv/SPhGtAL+knYe/GpgEFoGDwJ3A+mDdkcUGTa0gIlIT\nVUnpiIhISgr4IiI1oYAvIlITCvgiIjWhgC8iUhMK+CIiNaGALyJSE/8HEZ53kghQj+UAAAAASUVO\nRK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#and generate some fake data\n", "Xf = np.linspace(-7, 10, 100)\n", "Yf = parabola(Xf, 3, 0.1, 0.5)+np.random.randn(100)*3\n", "plt.plot(Xf, Yf, 'k.')\n", "plt.axis([-8, 11, -10, 60])" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "['a0', 'a1', 'a2']\n", "['x']\n" ] } ], "source": [ "#Lmfit works with Model objects that is the law of the fitting function\n", "from lmfit import Model\n", "\n", "#first, let's create a model from scratch using our parabola function\n", "pmodel = Model(parabola)\n", "\n", "#LmFit automatically recognizes what is a parameter and an independent variable\n", "#note: the Model object params are not the ones that will be minimized, they are model properties\n", "print(pmodel.param_names)\n", "print(pmodel.independent_vars)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n" ] } ], "source": [ "#the central part of LmFit is the Parameters object. it behaves like a dict and contains the parameters\n", "#and values that are going to be adjusted. we create them like this:\n", "params = pmodel.make_params()\n", "print(params['a0'])" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Parameters([('a0', ), ('a1', ), ('a2', )])\n" ] } ], "source": [ "#we can set each parameter value separately, and their range\n", "params['a0'].set(value=1, min=-5, max=5)\n", "print(params)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[Model]]\n", " Model(parabola)\n", "[[Fit Statistics]]\n", " # function evals = 15\n", " # data points = 100\n", " # variables = 3\n", " chi-square = 1047.238\n", " reduced chi-square = 10.796\n", " Akaike info crit = 240.874\n", " Bayesian info crit = 248.690\n", "[[Variables]]\n", " a0: 2.54746403 +/- 0.478801 (18.80%) (init= 0)\n", " a1: -0.01636483 +/- 0.080041 (489.11%) (init= 0)\n", " a2: 0.51167689 +/- 0.014953 (2.92%) (init= 0)\n", "[[Correlations]] (unreported correlations are < 0.100)\n", " C(a0, a2) = -0.697 \n", " C(a1, a2) = -0.560 \n", " C(a0, a1) = 0.219 \n", "\n" ] } ], "source": [ "#so, let's perform the fit giving initial values for every parameter (the fit will fail otherwise)\n", "result = pmodel.fit(Yf, x=Xf, a0=0, a1=0, a2=0)\n", "#the Results object has several attributes, including a pretty report\n", "print(result.fit_report())" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[-8, 11, -10, 60]" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VOX1+PHPyZAAIosCRUBxqdS9IKaWUdFoFEWtSm1R\nWxAVxI1W6lbiUvUHX6KibWiLLaAC7qBWQcU1EjcCCogiIkWLIAiCAQS3JCTn98czkSGZyUwyy53l\nvF+veWVm7p17n9xMzn3uuc8iqooxxpjMl+N1AYwxxiSHBXxjjMkSFvCNMSZLWMA3xpgsYQHfGGOy\nhAV8Y4zJEnEJ+CLSQUSeFJGPRWS5iPhFZE8ReUVEVgZ+7hGPfRljjGmeeNXwJwAvqurBQC9gOTAa\nKFXVnkBp4LUxxhiPSKwdr0SkPbAEOECDNiYiK4ACVV0vIl2BMlU9KKadGWOMabYWcdjG/sAmYKqI\n9AIWAVcDXVR1fWCdDUCXUB8WkRHACIA2bdocdfDBB8ehSMYYkz0WLVr0lap2jrRePGr4+cB84FhV\nXSAiE4BtwB9UtUPQeltUtdE8fn5+vi5cuDCm8hhjTLYRkUWqmh9pvXjk8NcCa1V1QeD1k0Af4MtA\nKofAz41x2Jcxxphmijngq+oG4HMRqcvPFwIfAbOBoYH3hgKzYt2XMcaY5otHDh/gD8AjIpIH/A+4\nGHcymSkiw4DVwKA47csYY0wzxCXgq+oSIFT+qDAe2zfGGBM762lrjDFZwgK+McZkCQv4xhiTJSzg\nG2OMR8rLyykuLqa8vDwp+4tXKx1jjDFNUF5eTmFhIVVVVeTl5VFaWorf70/oPq2Gb4wxHigrK6Oq\nqoqamhqqqqooKytL+D4t4BtjjAcKCgrIy8vD5/ORl5dHQUFBwvdpKR1jjPGA3++ntLSUsrIyCgoK\nEp7OAQv4xhjjGb/fn5RAX8dSOsYYkyUs4BtjTJawgG+MMVnCAr4xxmQJC/jGGOOB6mqYOBG+/z55\n+7SAb4wxHpg0CUaOhLlzk7dPC/jGGJNkW7fCbbfBiSfCgAHJ268FfGOMSbJx42DzZrjnHhBJ3n4t\n4BtjTBKtWgUTJsDQoXDkkcndtwV8Y4xJotGjoUULGDs2+fu2gG+MMUkybx7MnAnXXQfduyd//xbw\njTEmCWpr4U9/gq5d4frrvSmDDZ5mjDFJ8Oij8M47MH067L67N2WIS8AXkc+A7UANsENV80VkT2AG\nsB/wGTBIVbfEY3/GGJNOvv3W5e7z82HwYO/KEc+Uzomq2ltV8wOvRwOlqtoTKA28NsaYrDN+PKxb\nByUlkONhIj2Ruz4bmB54Ph04J4H7MsaYlPT553DXXXDeeXDssd6WJV4BX4FXRWSRiIwIvNdFVdcH\nnm8AuoT6oIiMEJGFIrJw06ZNcSqOMcakhmHDvqK6uppBgxZ7XZS4BfzjVLU3MAC4SkSOD16oqoo7\nKTSgqpNVNV9V8zt37hyn4hhjjPf+/e8PeeWVTtTU3MngwcdRXl7uaXniEvBVdV3g50bgaeBo4EsR\n6QoQ+LkxHvsyxph0UFMDY8Z0Aj5HtZiqqirKyso8LVPMAV9E2ohI27rnQH/gQ2A2MDSw2lBgVqz7\nMsaYdDF1KnzxxV7k5t6Mz1dJXl4eBQUFnpYpHs0yuwBPixsBqAXwqKq+KCLvAjNFZBiwGhgUh30Z\nY0xKKi8vp6ysjIKCAg491M+NN7qbtHfddTmvv34wBQUFSZ2wPJSYA76q/g/oFeL9CqAw1u0bY0yq\nKy8vp7CwkKqqKvLy8jjrrE/46qtuvPACHHWUn2OO8TbQ17GhFYwxJkZlZWVUVVVRU1NDZeWBPPHE\nXgwfDkcd5XXJdmUB3xhjYlRQUEBeXh45OT5gAm3a1DJunNelasgCvjHGNFN5eTnFxcUAlJaWct55\nM6mtPZHi4hZ06uRx4UKwwdOMMaYZ6uftn39+Lm+//Wt69YLLLvO6dKFZwDfGmGYIzttXVVUxZswO\n1qyBhx92E5ykIkvpGGNMM9Tl7X0+Hy1aHMJbb/kZPBj69fO6ZOGl6HnIGGNSm9/vp7S0lLlzy5g9\n+wqWL89h/HivS9U4C/jGGNNMfr+fDRv83HSTG/p4r728LlHjLKVjjDHN9O23MGoUHHEEXHWV16WJ\nzGr4xhjTTGPHwpo18MYbqXujNpjV8I0xJoK69vbBwxt/9BHcfTdcdFFq36gNlgbnJGOM8U799val\npaX07evnyiuhbVs3m1W6sIBvjDGNqN/evqysjJUr/bz+OkyaBOk0b5MFfGOMaURde/u6Gn6fPicz\nZAj07QvDhzdtW8FDKHsxVLIFfGOMaURde/u6QD1t2i+oqICXX4acJtwFDZUaSnbQt4BvjDEh1K+N\n+/1+3n4bJk+Ga6+F3r2j/yyETg1ZwDfGGI+Fqo0fdZSfyy6DHj3gttua9lm/398gNeTFdIcW8I0x\npp5QtfGyMj/LlsGzz8Luuzfts3VXCMGpIcvhG2OMh+pSMR07dtylNn7ggady4YVw7rlw5pmNb6Ox\nmnxd4PeKBXxjjKFhKqakpISKigpOOKGAv/ylD7m5MGFC5O2kQk0+HAv4xhhDw1RMRUUFRUVFTJ8O\npaVw773QvXt0TSu9rsmHYwHfGGMInYrZuBGuuQaOPdbNYpUKTStjEbexdETEJyLvichzgdd7isgr\nIrIy8HOPeO3LGGNiEWpsnLpUzJgxY34M5H/6E2zf7ppi5uSEviGbTuJZw78aWA60C7weDZSq6h0i\nMjrw+s9x3J8xxjRZY7X04FTMCy/Ao4/CrbfCoYe6z6ZC08pYxKWGLyJ7A2cA9wW9fTYwPfB8OnBO\nPPZljDGxiKaWvn27S+EccggUFe18P9RVQDqJVw2/BLgBaBv0XhdVXR94vgHoEuqDIjICGAHQo0eP\nOBXHGGNCi6aWPno0rF0Lb78NLVvuuixVb8hGI+aALyJnAhtVdZGIFIRaR1VVRDTMssnAZID8/PyQ\n6xhjTLxEajb55puuRc6oUZCmcT0sUY0txopIMTAE2AG0wuXw/wP8AihQ1fUi0hUoU9WDGttWfn6+\nLly4MKbyGGNMc33/PfTqBTt2wNKl0KaN1yWKjogsUtX8SOvFnMNX1SJV3VtV9wPOB15T1cHAbGBo\nYLWhwKxY92WMMYl0662wciVMmZI+wb4pEjnF4R3AKSKyEjg58NoYY1LSggVwzz1ujPvCQq9Lkxgx\np3TiyVI6xhgv/PAD9OkD33wDH34I7dpF/kwqiTalYz1tjTFZ7/bbYflyePHF9Av2TZHIlI4xxiRd\nqF60jXn3XTcR+bBhcOqpCS6cx6yGb4zJGE0d6+aHH2DoUOja1eXvM53V8I0xGaOpY93ccotL5dx/\nP7Rvn5wyeskCvjEmY9T1ovX5fBHHunn7bVerP+ecDSxeHH0KKJ1ZKx1jTEaJZrz6b791Hay+++4H\ntmzpQXX15rQc7rhO0jpeGWNMKvH7/RQVFTUauP/8Z/j0UzjjjKeort6ctsMdN5UFfGNMVnnpJZg4\n0Y2Vc8klB0SdAsoE1krHGJM1Nm+Giy9249uPGwetW6fu/LOJYAHfGJM1rrwSNm2C55+H1q3de+k8\n3HFTWcA3xmSFxx6DGTNg7Fg48kivS+MNy+EbYzLe6tVwxRVufPs/Z/FEqxbwjTEZraYGLrwQamvh\n4YehRRbnNSzgG2My2vjx8MYbcPXVnzBjRnZ0sAoni891xphMt2iRGz7hpJO+4u67f051dXRj7GQq\nq+EbY9JSpFExv/kGfvc76NIFjjnmIaqrox9jJ1NZDd8Yk3aiGRXzj3900xW+9hq0bNmXe+7J+3H9\nTO9gFY4FfGNM2gk1KmZwwJ8xA6ZOhZtuAhfbs6uDVTgW8I0xaaduVMxQNfZVq2DECOjb101KXieb\nOliFYwHfGJN2/P7QNfbqape3B3j0UcjN9bCQKcgCvjEmLYWqsd98M8yf71I6++/vUcFSmLXSMcZk\nhBdfdHPTjhgBgwZ5XZrUFHPAF5FWIvKOiLwvIstE5PbA+3uKyCsisjLwc4/Yi2uMMQ198YXrTXv4\n4VBS4nVpUlc8aviVwEmq2gvoDZwmIn2B0UCpqvYESgOvjTEmrmpq4Fe/+pqtW6spKlry4yiYpqGY\nA7463wRe5gYeCpwNTA+8Px04J9Z9NVWkjhnGmPQ3fPhaFi9uz44dlzF8+DH2/96IuNy0FREfsAg4\nEJioqgtEpIuqrg+ssgHoEo99RSuajhnGmPT24oswbdreiDyA6lSqqnwN2uSbneJy01ZVa1S1N7A3\ncLSIHF5vueJq/Q2IyAgRWSgiCzdt2hSP4gChO2YYYzLH55/D4MHw059+S8uW12fNNIWxiGuzTFXd\nKiJzgdOAL0Wkq6quF5GuwMYwn5kMTAbIz88PeVJojsY6Zhhj0ltVFZx3HlRWwvPPt2Hz5ueyvhdt\nNGIO+CLSGagOBPvWwCnAncBsYChwR+DnrFj31RThOmYYY9JLeXl5g//j66+H8nJ4/HE46CAA60Ub\njXjU8LsC0wN5/Bxgpqo+JyLlwEwRGQasBpLSMrb+lyPUlyDUF8gYk3pC3YtbvdrP3/8Oo0a5Wr6J\nXswBX1U/ABrMEKmqFUBhrNtvimhu1NrNXGPSR/17cTNmLOW++/wce6zrZBWJVe52lVE9bcvKyqis\n/Ak1Nd3C3qi1m7nGpI+6e3E+n4/c3E48/fQQ2rSBmTMjj5NTV7m75ZZbKCwstOaaZNhYOscdV4Dq\nb4Bt5OaessuN2rozfceOHe1mrjFpou5e3GuvlfHii5cxf35rSkuhW7fIn400hHI2yqiA36+fn7vu\nWs4NN/ThhBNW0LdvZ6BhGqekpISKigq7zDMmDfj9fl591c9bb8GECXD88dF9zlrqNZRRAR/guusO\n4fvv4S9/6cyECe7GTv0zfUVFBUVFRV4X1RgTRnDuffNmP7fe6trc5+eXU1wcXU7eWuo1lHEBH9ws\nN4sXw3XXQa9edqY3Jp0EX5G3aHE4LVospHfvFlxyyQJOPrlpDS5s0pNdZdRN2zo5OfDgg6597m9/\nC3vt5c70Y8aMsVY5xqS4nVfku1NZOQOo5JlnYP7816zBRYwysoYP0LYtzJoFRx8NZ50F8+bZmd6Y\ndFBQUEBubitqah4DDmD8+BX06HG4XanHQUbW8OsceKBrvrV8uRsru7bW6xIZYyLx+/2ce+4KYAA3\n3LCGK644/Mf37Uo9NuLGNUsN+fn5unDhwrhvt+7m7c03w5gxcd+8MSYO6m7Ufvfdbxk79kCuuALu\nvdfrUqUHEVmkqvmR1svYlE6wP/4Rli6FsWPh4IPh97/3ukTGZL6m9HKtu1FbWfkLamuvJT9/KxMm\ndEhSSbNHVgR8EVdT+OQTGDYMDjgA7GrQmMRp6hAmrpd8V2prnwRWMWDAHHJz/5S8AmeJjM7hB8vL\ng6eegr33hnPOgdWrvS6RMZmrqUOY5OcXojobaEHLlr9lwIC+SSlntsmagA/QsSM895wbQ/uMM+Dr\nr70ukTGZKXgMnEgtaqqrYfz4o8nJOYRhw15g7txJdkM2QbLipm19paVw2mlQUABz5kQehMkY03TB\nOXwgZD5fFS6/HCZPhgcegIsv9qiwaS7am7aoaso8jjrqKE2WqVNVQXXYMNXa2qTt1pisM2/ePG3d\nurX6fD5t3bq1zps378dld93l/g9vvNHDAmYAYKFGEWOzKqUT7KKLXDPN+++H4mKvS2NM5gqXz58x\nA264wU1iYs2lkyMrWumE8//+H6xa5cbe2WcfGDLE6xIZk3lC9ZB94w3XGbJfP5g2zQ2HYhIvK3P4\nwSor4fTT4Y03XD7/lFOSuntjskJwPr99ezdj1V57QUnJuyxe/KqNZhmjaHP4WR/wwbXWOf54V9t/\n4w3o3TvpRTAmK6xdC8ceC1VVMHHiYgYPPs6mG42DaAO+XUgB7du72n2HDjBgAPzvf16XyJjMs2WL\nax23ZQs8/zysWPGSjX6ZZBbwA7p3h5decjWP/v3hyy+9LpExmeP77+FXv4KVK+GZZ6BPn6a11Tfx\nYQE/yCGHuJrH+vWuJmIds4yJXXU1nH8+zJsHDz0EJ53k3rfRL5Mv5oAvIvuIyFwR+UhElonI1YH3\n9xSRV0RkZeDnHrEXN/H69oUnn4SlS2s58sjVlJUt8LpIxqSt2lo3ftXs2fCPf8CgQbsu9/v9FBUV\nWbBPknjU8HcA16rqoUBf4CoRORQYDZSqak+gNPA6LXToUE5OziWsWrUPhYVf8eab870ukjFpR9UN\nS/7QQ66d/VVXeV0iE3PAV9X1qro48Hw7sBzoDpwNTA+sNh04J9Z9JUtZWRm1tQ8DV1JbewYjR7al\npsbrUhmTXm6/3dXq//Qn19fFeC+uOXwR2Q84ElgAdFHV9YFFG4Au8dxXIu28mXQfLVrczAcfHMZV\nV7kaizEmsvHjXcC/+GK45x43RLnxXtx62orI7sBTwChV3SZBf2FVVREJGS5FZAQwAqBHjx7xKk5M\n6m4m1XUUefZZN/xC69bw17/al9eYxkycuHPIhClT7P8llcQl4ItILi7YP6Kq/wm8/aWIdFXV9SLS\nFdgY6rOqOhmYDK7jVTzKEw9+/85Jz/v2he++g5ISaNUKxo2zL7HJPtHMYDV1KowcCWef7XL3Pl/z\ntmMSI+aAL64qfz+wXFX/GrRoNjAUuCPwc1as+/KKCPztb24YhjvucEH/1lu9LpUxu0pkII1mBquH\nHnItcvr3dwOjhRp2vKkzYZn4ikcN/1hgCLBURJYE3rsRF+hnisgwYDUwKMznPRXtP4mIu1T94Qe4\n7TY32NPJJ1tNxaSGeAbSUP8ToUa8rFtWXl5OSclGnnzyLE48UXj6aWjZMvS2G9uOSYJoxlBO1iOZ\n4+GrNj5Odzg7dqheeKEbw7tFi9ub9FljEmXcuHHq8/kUUJ/Pp+PGjWvWdsL9TzT2fm7uEIUdmpNT\npv/4x306bty4sP8PzfmfM5ER5Xj4WT08cv3axoMPPhixxu7zuZl5li5dynvv/QWopqqq2GoqxlOh\nhiBujnA18PoNGeq+6yUlm6iungq8heqZXHNNJbW1tWGvMsJtxyRJNGeFZD28rOHn5eVpy5Yto655\nvPnmPPX5Hg7U9Mfp229bTcV4a968eY3WrqPdRrQ18GnTVEVqNSenTHNy2mqLFi00Jycn5qsM03RY\nDT+y4NrGmjVrmDJlStS5xeOO81NWBqNGvc+iRUXMmgV+f+Otd6x1gkmk4JZlsWwjmhr4/ffDpZdC\nYaFw442tmD+/iI4dOzJq1KiYrzJMAkVzVkjWI9k1/GDR1mzq16JqalSvuMLl9K++Ovz8uJa7NJli\nwgT3fT/tNNXvv991WTyuMkzTYTX8pommZhOuJcTEia6p5t/+Bt98A5MmNWx/bK0TTCYYN84NkzBw\nIDz2WMPWOPG4yjCJYwE/SKQva7igLeK6j7dt6+bJ3bYNHn4Y8vJ2fjbam2qW9jGpSBVuvNH1Qxk8\n2HWwamHRI+3Yn6wJGgvaIm7skPbt4dprYft2N8xymzZueSxXEMZ4qaYGrrwSJk+Gyy6De++FBQus\nYpKOLOBHIbjWHSloX3MNtGvn/jFOPhmeew46dnTLmnsFYYxXKivh97+Hp55yNfyxY2H+fKuYpCsL\n+BGEqnUXFRWFXK/uRHDYYXDBBV/wxBMD6dcvh5degn32ibyveLWlNiYevv4afv1reO01d39q1Cj3\nvlVM0pcF/Aii+XIHnxR8Ph8iwo4dO/D5TmLVquc4/PAfmDhxFYMH92p0X9YpxaSKdevg9NPho4/g\nwQdhyJCdy6xikr4s4EcQzZc7+KRQW1sL1DV3LSUn5xh++GEWQ4bsx5Yty/jDHw5rdH/WysEkS7gG\nAsuWwYABsGWLm+O5f/9dP2cVk/QlrglnasjPz9eFCxd6XYwGIrWcCVfDFxFqa2upre0OzMHnO4Rp\n03wMHpz838GYYOEaCLz2Gpx7rmtmPGcOHHmk1yU10RCRRaqaH2k9q+FHIVKtu36NB1ytf2fPwy/I\nzT2FQw/9mCFD2vPJJ254ZRtT33glVKry44/9jBgBP/uZq9nvt5/XpTTxZgE/TuqfFOqeH3HEET+e\nCI46qj2XXeaab/73v24QtlatklM+a99vggWnKnNzW/Lxx0O48UbXsuyJJ6BDB7eefW8yTDTdcZP1\n8HJohWSprVUtLnZd0/1+1fXrE79PG9Yhu4Ub7mDevHl6663j9YQTvlJQvfRS1aqqXZfb9yY9EOXQ\nCnGdxNxEJgKjR7tOWe+/D7/4BSxalNh9hrp8N9mhLld/yy23UFhYSHl5+Y/LunXz8/TT1/Hmmx35\n29/ckCDBs1TZ9ybzWErHI+eeCz/9qZv787jjXHrnggvCr9+cS+u6z3Ts2NGa0WWZur/9mjVrQjYr\nLiuDQYOgqsrdnD311IbbsOaXmcda6Xhs40bo338b77/fjvPP/4KHHurWYIyS+i0qSkpKqKioaDT4\nN+czJjOEazWWl5fHq6+W8s47fq67Dnr2hGeegYMOanxblsNPfdZKJ018+mk5K1achsj/8fjjI1mx\n4mteeKE9XbrsXCf40rqyspKRI0c2OqtQ/c9UVVVRUVERsoewyTzBf3uASy+9lB49evDLX57E7bcf\nwMsvw/HHb+bZZ/ekXbvGt2X9QjKL5fA9VlZWRnX1t6j+AZGhLF26G336wFtv7Vyn7tLa5/ORk5ND\nTU1NxLxq8Gfscjy71P/bX3jhhQwcWMTw4Ufw8ssdEbmFd97Zh2XLyiNvzGQUC/geC/7nbNXqCe67\n7yN22w0KCuDOO6G2dmc7/zFjxjBx4kRatmwZMZAHf8YGt8ou9f/2n33mJz8fNm5UcnJOR3Us1dWV\ndhM2C1kOPwXUz5Nu2wbDh7v20GecAdOmQadO4dc3JpTvvnMDnk2ZAsceCzfcsIjzz+9no1xmoGhz\n+BbwU5SqG3f8mmvc8MoPPwwnneR1qUy6WLoUzj8fli+HP//ZTcyTm2uVhUwVbcCPS0pHRB4QkY0i\n8mHQe3uKyCsisjLwc4947CtTlJeXU1xcvEu76GAicNVVsGCBm1Tl5JOhqMg1ozMmHFX4xz/g6KOh\nogJeegmKi3e2r/f7/RQVFYUc8bWx76PJENH0zor0AI4H+gAfBr13FzA68Hw0cGek7WRDT1vVpk+Y\nXlo6X4cPd71z+/RRXbYsyQU2aWHdOtX+/d33ZMAA1Q0bovuc9ahNfySzp62qvgFsrvf22cD0wPPp\nwDnx2FcmiKYHY3APyTPPPJFLLinnP/+BNWugTx83IUVgJGaT5VTh8cfhiCPgzTddKvD559mlaW9j\nrEdt9khkK50uqro+8HwDEPLrJyIjRGShiCzctGlTAouTOqJpMhnqn3DgQPjwQzjlFJfbLyiAmTPf\n+/FSPJ6X5XaJnx42boTf/tb10j7wQHjvPbjiiqaNxGpNeLNINJcB0TyA/dg1pbO13vItkbaRLSkd\n1fADWgUvD3eZXVur+sADqrvvXq3wrYpcq7m5rbRly5ZxuSy3S/zUV1ur+sgjqp06qeblqd5xh2p1\ndfO3F+n7aFIbUaZ0EtnT9ksR6aqq60WkK7AxgftKO00dYz94XRG4+GJYseJf3Hnn/qjeTXX1IOAy\nYEnM84ym05yl2djqZM0aV4ufM8cNvjd1KhzW+ERqEVmP2iwRzVkhmgcNa/jj2fWm7V2RtpFNNfx4\nmDdvnrZq1VpFfqewQaFaRe7WVq06ZkUNP13KGS9VVar33KPapo17lJSo7tjhdalMKiCZNXwReQwo\nADqJyFrgVuAOYKaIDANWA4PisS+zk5uSzl0F9OmzhokT4dlnr6V9+5GsW9cS1ci53FA15FSfszTS\nSJCprjlXJW+95Wr1H34IZ54J//wn7LtvggtqMk80Z4VkPayGH7u33lLt1cs1zTvllMabcNavIU+a\nNClpedzm5oyDy5yXlxe3+xaJEOp3bGqT3GeeWai//737e/booTprVrJKb9IJKZDDNx449lhYuBD+\n/W+4+Wb4+c9dzfC221yPXQhdQ452FM54CDeBdjTCjQSZalci4X7Hxu6PBM9fcPXVRVRWjkT1YPLy\narnpphyKiqBNm4b7SdUrMZOCojkrJOthNfz42rRJ9corVX0+1Q4dVO+6S3Xu3Pkha8gtWrTQnJwc\nBdTn8+m4ceMSVq5x48apz+dr1r7SJW8f7ncMV/6693NycjUnZ7jC5wqqIjP1+usnhtxHuhwLk3hY\nDd906gQTJ8KVV8L118MNN0D79ofwww+/R3UqUPNjDbljx46MGjUqKbMbxTKTUqrcX4hUsw73O9Yv\nP0BxcTGrV6+hsrI/tbVjgCMQWYDIYFq2fIeBA0tDliGdWlOZFBHNWSFZD6vhJ9bcuaqHHrpNXd/M\nFZqbe5G++eau+eVUz+Gngqbm4Rtb3qpVa83JGaAi7/z4d8nLu0D//e/I91Oshm/qEGUN30bLzDKq\ncOedK5gwoQMbNnTh4IPhppvcyIr1p1ZsimhyyemYbw5V5uLiYm655RZqamrw+XyMGTOmybOJqcJF\nFz3Bgw/uA/QFVlNQ8DqFhV9QWHhCo8cnuExA2h1TE3/Rjpbpea0++GE1/OSpqVF98knVI45QBdX9\n91e9917V775r/HPNbXmSjrXRSPn25vwu1dWqjz2m2ru3BnL0n6nIFdqqVfuotpOOx9EkHskcPM2k\nn5wcOPdcWLIEZs+Gn/zE5fp79IBbb4Uvv2z4meAB3QoLC38cZyeawbfScYCucGVuzmxi27bBhAlu\nvJsLLnCTk0ydCq+/vp7/+799eO21F6LaTjoeR5NCojkrJOthNXzv1NaqlpWpnnWWqogbn2XIENXy\ncrdMtektT4KlY800HmVevlz1qqtUd9/d1ej79XNt6WtqvCuTyTxYDt8018qV8Pe/w/TpsH079O4N\nI0bAT3/6DuecUxCy/Xyq5fDjta/mbOf77+Gpp+C+++D11yEvz90jGTnSjX0Tq3S8F2ISy6Y4NDH7\n5ht45BE3vvoHH0CrVnD88Zvo1q2USy7Zl379UjPYxNKxq7lUYd48eOghmDEDtm6FAw5wcxMPG+ZS\nZsYkSrRB7byaAAAMLklEQVQB39rhm7B23x0uu8zV7t97Dx54AB55pDNbt57PnDkwaJB7HHMM+Hxe\nl3anZLVPV3XH5YknXJBftQpat4aBA12gP+EEd68ErFZuUoMFfBORiJtlq08fuOceeOEFePRRl7L4\n5z9d7fXss+Gss+DEExt2/0+2WDp2RVJd7QYye/ZZmDUL/vc/d7I76SQ3fMXAgdC27a6f8eKKw5hQ\nLOCbJmnZEs45xz22b3fB/z//gccegylT3PKCAujfHwoL3bR7OUluCxbP3riq8Omn8OqrOx9ff+3y\n8oWFrg/D2WfvHKcoFOsRa1KF5fBNXFRWwhtvuEk55syB//7Xvd+pE/Tr59I+xxwDRx7p0h7heJ36\nqK52QxCXl7v5Yd98E9atc8v22cedyM48E04+2aW8omE1fJNodtPWeOrzz2HuXCgthbffdrVkcOmP\nww5z6aGf/9w9P/RQ6N4d5s9PXmBUhU2b4OOPYelS91iyxD0qK9063bq5k9Xxx7sA37Nn0+aKDeb1\nicxkNgv4JqVs2OBqzYsWucfixW4C7jpt2kDbthvZsKEc+AyRDfzmN34uvfQs9tzTpUzat4fddnPp\nlHCBd8cOl2ravh22bHFBfdMmWL8e1q51J6JVq1zT023bdn6uQwd3AvrFL3Y+9t9/1/1Y0DapygK+\nSXmbNsGyZe7xySewcOFm3n57Pap7A+3Dfs7nc/cKcnJ2tg6qqnLpmB07wu9vt91cWmbffV1tvWdP\nOOggd5+hW7fGa++NpWXsRGC8Zs0yTcrr3Nnd4N3ZiGZPystXUFZ2L0cffRJ77/1LNm2CzZuhosLV\nyL/91j0qK6GmBmpr3Sfz8iA3190faNvW5dc7dHD76NwZ9toL9tij+SmZcDdeLT9v0okFfJNS/H7/\nLgHzoIM8LEyQcE09658IHnzwQavtm5RlAd9ElKkpi6b8XuGaegafCHw+H1OnTmXHjh1xre1n6vE3\nyWcB3zTKy9x1IrffnFRM/auPuvfqTgRr1qxhypQpcW1vbykjE08W8E2jvMpdJ3r78UzF1J0IysvL\nmT59elx7+FqnLRNPCQ/4InIaMAHwAfep6h2J3qeJn2hz1/EORInefiJSMYmYbzeRw0SY7JPQgC8i\nPmAicAqwFnhXRGar6keJ3K+Jn2hy14kIRInafnCaKBGpmFBpn1ikyqTtJjMktB2+iPiB21T11MDr\nIgBVLQ61ftza4YcKDoMGuSmdvvsOTj+94fKLLnKPr76C3/ym4fIrroDzznM9d4YMabj82mvhV7+C\nFSvcEJP13Xyz6665ZAmMGtVw+bhxbuyBefPgxhsbLi8pcQPTv/oqjB3bcPmkSa5Jy7PPuhHO6nvo\nIdcIfcYM+Ne/Gi5/8kk3DsK0ae5R35w5riH7vffCzJkAfL1tG1u3bqVDhw60X7zYrXf33Xz96KM7\n32/XzrWVfOEFt3zMGNf9NljHjm4AeYCiItdDK2j7u/XsSeeXXnLLR41yxzDYz34Gkye75yNG7BzX\noU7v3lBSQnl5OZ/160e3mhokJ4devXrRvl071vXoQc8nn6SqqoonVTn6wAOprq7eWf7CQrjlFret\nAQPcgPfBzjwTrrvOPbfvXsPlCfju7aJu1q+774bnntt1WTO/ez/ae294+GH3PIbvHgCDB7vef8H8\nfigOGQ6bJFXa4XcHPg96vRb4ZfAKIjICGAHQo0ePBBfHxFP7du1cQAzy2erVrHn/fbS2dmdQbWzw\nnGi237lzPIpLWVkZPWprUYDaWrZu3Ur7du3o3r37j7XoXk88wSf1yx+XvRuTAqKZFqu5D+A3uLx9\n3eshwD/DrW9THKa/+tMgXn755Q0mPU+2uonXJ02aFHF6wHDTOBqTyohyisNE1/DXAfsEvd478J7J\nUMlol94U9Vv7lJSUUFFRETYfbjdJTSZLdMB/F+gpIvvjAv35wO8SvE/joUS3S2+q+q19KioqKCoq\nCru+3SQ1mSyhAV9Vd4jISOAlXLPMB1R1WSL3abyXyHbpjQnVUas5NfZ4t7QxJlXYaJkmoZI1LICN\nZmmyWaq00jFZLlm15cY6almN3RgnybONGtM85eXlFBcXU16/nXRAXerG5/PZzVZjwrAavkl50Yyr\nYzdbjYnMAr7xXKQce7Tj6ljqxpjGWcA3noqm9m5t442JDwv4xlPR1N4tXWNMfFjAN56KtvZu6Rpj\nYmcB33jKau/GJI8FfOM5q70bkxzWDt8YY7KEBXxjjMkSFvCNMSZLWMA3xpgsYQHfGGOyhAV8Y4zJ\nEhbwjTEmS1jAN8aYLGEB3xhjsoQFfGOMyRIW8I0xJktYwDfGmCxhAd8YY7JETAFfRH4rIstEpFZE\n8ustKxKRT0RkhYicGlsxjTHGxCrW4ZE/BH4NTAp+U0QOBc4HDgO6Aa+KyM9UtSbG/RljjGmmmGr4\nqrpcVVeEWHQ28LiqVqrqKuAT4OhY9mWMMSY2iZoApTswP+j12sB7DYjICGBE4OU3IhLqBBKsE/BV\nzCVMHitvYll5E8vKm3jxKPO+0awUMeCLyKvAXiEW3aSqs5paqvpUdTIwOdr1RWShquZHXjM1WHkT\ny8qbWFbexEtmmSMGfFU9uRnbXQfsE/R678B7xhhjPJKoZpmzgfNFpKWI7A/0BN5J0L6MMcZEIdZm\nmQNFZC3gB54XkZcAVHUZMBP4CHgRuCqOLXSiTv+kCCtvYll5E8vKm3hJK7OoarL2ZYwxxkPW09YY\nY7KEBXxjjMkSKR/wRWSGiCwJPD4TkSVh1vtMRJYG1luY7HIGleM2EVkXVObTw6x3WmDYiU9EZHSy\nyxlUjvEi8rGIfCAiT4tIhzDreXp8Ix0vcf4eWP6BiPRJdhmDyrKPiMwVkY8CQ49cHWKdAhH5Ouh7\n8hcvyhpUnkb/vil2fA8KOm5LRGSbiIyqt46nx1dEHhCRjSLyYdB7e4rIKyKyMvBzjzCfTVxsUNW0\neQD3AH8Js+wzoFMKlPE24LoI6/iAT4EDgDzgfeBQj8rbH2gReH4ncGeqHd9ojhdwOvACIEBfYIGH\n34GuQJ/A87bAf0OUtwB4zqsyNvXvm0rHN8R3YwOwbyodX+B4oA/wYdB7dwGjA89Hh/pfS3RsSPka\nfh0REWAQ8JjXZYmDo4FPVPV/qloFPI4bjiLpVPVlVd0ReDkf12ci1URzvM4GHlRnPtBBRLomu6AA\nqrpeVRcHnm8HlhOmp3kaSZnjW08h8Kmqrva6IMFU9Q1gc723zwamB55PB84J8dGExoa0CfhAP+BL\nVV0ZZrniBmlbFBiuwUt/CFz2PhDmsq078HnQ67BDTyTZJbhaXCheHt9ojldKHlMR2Q84ElgQYvEx\nge/JCyJyWFIL1lCkv29KHl/cII3hKoGpdHwBuqjq+sDzDUCXEOsk9DgnaiydJoly+IYLaLx2f5yq\nrhORnwCviMjHgbNs3DVWXuBfwBjcP9AYXBrqkkSUI1rRHF8RuQnYATwSZjNJO76ZQkR2B54CRqnq\ntnqLFwM9VPWbwH2eZ3AdFL2Sdn9fEckDzgKKQixOteO7C1VVEUl6m/iUCPgaYfgGEWmBG4b5qEa2\nsS7wc6OIPI27NErIFzZSeeuIyBTguRCLkjr0RBTH9yLgTKBQA4nEENtI2vENIZrjlVLDeYhILi7Y\nP6Kq/6m/PPgEoKpzROReEemkqp4M/BXF3zeljm/AAGCxqn5Zf0GqHd+AL0Wkq6quD6TDNoZYJ6HH\nOV1SOicDH6vq2lALRaSNiLSte467EflhqHUTrV5ec2CYcrwL9BSR/QO1lPNxw1EknYicBtwAnKWq\n34VZx+vjG83xmg1cGGhN0hf4OujyOakC95vuB5ar6l/DrLNXYD1E5Gjc/2JF8kq5S1mi+fumzPEN\nEvaqP5WOb5DZwNDA86FAqMEnExsbvLqL3ZQHMA24vN573YA5gecH4O5mvw8sw6UqvCrrQ8BS4IPA\nH6pr/fIGXp+Oa73xqcfl/QSXM1wSePw7FY9vqOMFXF73vcC1HpkYWL4UyPfwmB6HS+l9EHRcT69X\n3pGBY/k+7mb5MR6WN+TfN1WPb6A8bXABvH3QeylzfHEnovVANS4PPwzoCJQCK4FXgT0D6yYtNtjQ\nCsYYkyXSJaVjjDEmRhbwjTEmS1jAN8aYLGEB3xhjsoQFfGOMyRIW8I0xJktYwDfGmCzx/wFpmojJ\nYr1VxAAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#let's see how it is\n", "plt.plot(Xf, Yf, 'k.')\n", "plt.plot(Xf, result.init_fit, 'r--')\n", "plt.plot(Xf, result.best_fit, 'b-')\n", "plt.axis([-8, 11, -10, 60])" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XlcVXX+P/DXm0URUVDcNUUU2QQDUYRUlHHcIhvNptXM\nLNNxGpussaap5msz32Zaxqa+U2lW6ky2ambm8nPfURaRVRERRHFXXFBE4PP7Q7gDsl0u595zl9fz\n8fCR95zP+Zz38dqL41k+H1FKgYiI7IuT3gUQEZH2GO5ERHaI4U5EZIcY7kREdojhTkRkhxjuRER2\niOFORGSHGO5ERHaI4U5EZIdc9Npxhw4dlI+Pj167JyKySUlJSeeVUh0ba6dbuPv4+CAxMVGv3RMR\n2SQRyTemHS/LEBHZIYY7EZEdYrgTEdkhhjsRkR1iuBMR2SGGOxGRHWK4ExHZIbsN96tXr4JTCBKR\no7LLcC8uLkZoaCjmzJmjdylERLqwy3DPyMjAE088gdzcXL1LISLShV2Ge1paGgYMGKB3GUREurHL\ncE9NTUVISAgA8Lo7ETkkuwz3o0ePwtfXF15eXrh8+bLe5RCRDbt58yaefvppLFmyRO9SmsRmw728\nvBwlJSU1lt24cQMXL16Ei4sLnJ2d0aNHD5w4cUKnConIHqxcuRL+/v5YtmwZbt26pXc5RrO5cFdK\n4aOPPsKoUaMwfvx43LhxAwBw7NgxhIWF4d5778Wzzz4LAA2Ge2FhIUaOHIkXX3wRFRUVFqufiGxL\nfHw8JkyYgAkTJmDdunV6l2M0mwt3AOjatSs2b96M3/3ud3jnnXcAAPPmzcPatWsxa9YsjB07FkDD\n4b5gwQK8+eabCAgIwPz58y1WOxHZluzsbPTr1w9xcXFYv359jXV13dP79ttvMWrUKKxZs8ZSJdbJ\n5sJdRDBx4kQ4OTlhwoQJ2L59Oy5evIiysjL4+vriiSeegIgAqD/cS0tLcfDgQQwdOhTTp09HUlIS\nTp8+belDMSgqKsLXX3+t2/6JqG4lJSVo2bIlRAR9+/bFkSNHDIE+bdo0REdHG64eALffsfn000+x\natUq/Otf/9I1V2wu3KtzcnJCREQEHn30Ufz617+utb6+cN+9ezdGjhwJ4PYPi/nz52PmzJkoKysz\ne811+eGHH/DnP/8ZW7ZsafK28+fPR3FxsRmqInJMly5dMvz+wIEDCA8PN3wOCAhAdnY2Dh8+DHd3\nd7z88sv429/+Zli/dOlSTJs2DR4eHpg3bx4++OCDWv3n5eWhtLTUvAcBGw93AHj11Vfxv//7v3jo\noYdqrevUqRPOnDlTa/natWsxbtw4w+ewsDA88cQTiI2Nxb59+8xab13WrFmDdevWYeHChU3a7tat\nW1i2bBn++Mc/mqkyIsczdOhQQw7Ex8djyJAhhnVjx47F+vXr8eOPP+L+++/HhAkTsGfPHsPDHatX\nr8bkyZMBADExMUhKSsLVq1cN25eXl+Ppp5+ucbZvLjYf7m3btkV4eLjhUkx1Tk5Odd4sTUlJqfWS\n06RJk7Bq1Sq8/PLLFr0jfuXKFVRUVKB3794oLi7GzZs3jd42MzMTkydPRnZ2thkrJHIc586dQ48e\nPfCXv/wFALBv3z4MHjzYsH7EiBFYt24dfv75Z4wYMQIigoceegjLly9HdnY2evXqhRYtWgC4fVXg\n6aefxqJFiwzbf/zxx/jVr34FT09Psx9Lo+EuIp+LyFkRSa9nvYjIByKSIyKpIhJeVztrkZ+fj169\netX5w6B9+/YYMmQI0tPrPFSz+PnnnxEXFwfg9k/6bdu2Gb1tYmIiIiIi0K5duxr/lCQi0+zatQsT\nJ06Et7c3CgoKcPnyZXh5eRnWt27dGu3atcNrr71mCPEpU6bgyy+/xNy5cw1P6lWZNGkStmzZgp9+\n+glvvfUWEhMTMWvWLIscizFn7ksAjG1g/TgAfpW/ZgD4uPllaadNmza4cuWK4fO6deswfvz4ettH\nRkZa9NJM1T/vAGDkyJHYuXOn0dtWhXtwcDAyMzPNVSKRw9ixYweGDx+OCRMm4KWXXkJgYGCtNl99\n9RVGjRpl+NyyZUssWLAAL7zwQo3r8wDg7OyMpUuXIi0tDf7+/vj888/h7Oxs9uMAjAh3pdQOABcb\naHI/gGXqtngAXiLSVasCm6tHjx44efKk4fPGjRtrfDF3Gjx4MPbv32/2upRSmDx5Mvz8/NChQwcA\nQGhoKFJTU+tsX1ZWVuPaHfDff4UEBQUx3Ik0kJWVhcDAQIwePRq5ubn4n//5H6O2Cw0NNTykcacO\nHTrgj3/8IyZNmgQnJ8tdCddiT90BFFT7fKJyWS0iMkNEEkUk8dy5cxrsunHVn5gpKSlBaWkp2rZt\nW2/7bt26obCw0Ox1xcfHo2/fvnjzzTcNy1q0aIFbt27VuE+wdetWfPHFFxg9ejQmT55sWFf1OJaI\nIDg4GBkZGWavmcieXblyBR4eHhAReHh4YP/+/WjTpo3eZZnMojdUlVKLlFIRSqmIjh07WmSfd911\nF44dOwbg9j+5YmJiGt3G3d291tAGWvvss8/wzDPP1Fru7++Pw4cPAwDOnDljCP81a9Zg/PjxWL58\nOQCgoKAAPXv2BAD4+voiJyfHrPUS2bu9e/ciOjpa7zI0o0W4nwRwV7XPPSqXWYWRI0ca3ipbuXJl\ng9fbq/j5+eHIkSNmrevUqVPo06dPreURERFITk4GACxfvhxz5szBtGnT4O7ujmnTpmHlypUAbo9Z\nHxwcDABwcXFBRUUFh1EgaoYtW7Zg+PDhepehGS3CfTWAJyqfmhkC4LJS6pQG/Wqiffv2aNmyJZYu\nXYri4mIEBQU1uk1gYCCysrLMVlN5eXm9195CQkKQlpYG4Pb9gerP47dt2xbl5eW4du0a0tPT0b9/\nf8O6gIAAwxk/ETVNWVkZEhISMHDgQL1L0Ywxj0J+BWAvAH8ROSEi00VkpojMrGyyFkAugBwAnwL4\njdmqNdFrr72GwsJCLFiwwKj2Wob7xYsX8cMPP9RYlpeXh969e9fZPiAgAIcOHUJpaSlcXFwMj1tV\nGTt2LH788Ufs3r0bd999t2G5pW4EE9mj9evXY+zYsXU+Im2rXBproJR6pJH1CsBszSoyg6CgIKPO\n2KsEBATg/fffb/Z+y8rK8Pjjj6NFixZwdXU1PM+elZWFgICAOrdp2bKlYeyb0NDQWuunTp2KmJgY\nhIWFoV27doblUVFRmD17NnJycpCfn4+FCxeiVatWOHHiBHr06NHsYyGyZ8uWLatzqABbZvNvqJqD\np6cnLl5s6OlP4yxZsgRxcXFYsmQJ/vOf/xiWHzp0qM7nZ6t4eXnh559/RmRkZK117u7ueOONNzBv\n3rway3v16oWnnnoKLVu2xMSJE/Hqq68afhCcP3++2cdCZK+qBvfq0qWLzpVoi+Fej+pP2Zjqq6++\nwjPPPAMvLy8UFRUZljd05g4As2fPxqpVq+oMdwCIi4ur82bspEmT8Kc//QkTJ07EhQsXEBMTg5df\nfhlbt25ttNbi4mLD8MlEjuDy5cv47LPPMHPmTMydO1fvcjTHcK/HxIkTa10rb4rjx4/D19cXrq6u\nAFBjiIBTp041eJZwzz33ICUlBZ06dTJ5/0uXLsVTTz2FX/ziF9i8eXOj7b///nv84x//qDE+9e7d\nu2sMa5CRkYENGzZwFEqyCXl5eTVOWCoqKvDcc8+hsLAQSinMmDEDJSUlmDdvXr0nUraM4V6PX/zi\nF1i1ahUuXLhg0vY7d+7EsGHDDJ8HDhyI5OTkGi8fWYKvry+OHj0KpRTOnz9f79AK3333HUaNGoW8\nvDwAt/918dxzz+Hvf/87gNuXkp5//nmkp6fjvvvuq3GpZ+XKlXj44Yc5GTlZlX//+9/YvHmzYa6E\nBQsWwM3NDbNmzcLKlSsRGBiI2bNnIyoqSudKzYPhXg83Nzf885//xAsvvGDS9neG+6BBg7B//36c\nO3euWWfkpggNDUVCQgKmTJmCV155pdb6kpISuLi4YNSoUYbw//zzz/HFF18gLS0NZ8+exV//+ld8\n+OGHmDt3Ll5//XW89957hu0//fRTdOzY0aIDrhE1ZseOHfj222+xfPlyVFRUYO3atXj77bfx5JNP\nYt68eXjxxRf1LtGsGO4NCAsLw5kzZ0w6I83Pz4ePj4/hc1W4N3Yz1RweeOABTJ48GVOnTkVoaCgO\nHDhQY31GRgb69++PyMhIxMfHAwDS09MREhKC119/HRMmTIC3t7fhPsGwYcNw8OBBALefCHJxccG4\nceOMuvxDZAkZGRkICAhA27Zt4e7ujuXLl2Po0KGGmdyys7Ph4eGhd5lmxXBvhK+vb5NvrFZUVMDZ\n2bnGpZeqIQ2q/tJZ0pAhQzB37lw89NBDmD59Oj777LMa6w8cOICwsDD069cPhw8fxvXr1+Hm5gYn\nJydERkZi1qxZeOuttwztnZ2d4ebmhuLiYsP8ksOGDWvSiJZE5vDGG28gISEB33zzjWECnzlz5mDb\ntm01hvuw5ABeerH/I2wmU0KrarTGOwUHB2Px4sUYOnSoVuUZxcnJCXPmzIGIICQkBEeOHEFJSQmu\nXr0KpRRSUlJw9913w8nJCS1atMCOHTtqTFAwdepUtGrVqkafUVFR2Lt3L1JTUzFgwAC0adMGN27c\nqHOqwldeeQXh4eEcc57Mqry8HNu3b8c//vEPrFq1yjBOTFRUFBYvXuxw73sw3BsxaNAgJCUlNWmb\nzMzMOi+9TJ8+Hf/+97/h7e2tVXkmmTJlCsaMGYNf/vKX+P3vf4+kpCTDG7NhYWF46aWX6py2sLqY\nmBhs374dBw4cMLxsFRERgcTExBrtTp06hSNHjuCFF15o1tNHRI1JT09HZGQkli9fjnXr1jnE2XlD\nHPvojeDj44P8/PwmbZOVlVXnG7GBgYFNelPWXB5//HFs2bIFe/fuRXBwMJYsWWL4HyEqKgr+/v7w\n9fVtsI/w8HAkJycjOTnZEO6xsbG1Jvn+4IMP8Nvf/hYTJkzA6tWrzXNAZPfquu9VXFyM48ePGz7v\n2bMHUVFREBF0717nqOMOheHeCBcXlzovNTSkvjN3a1J1T+CZZ56Bv7+/Yfno0aOxbNmyRrd3cXGB\niMDPz6/GD4aNGzca/rxOnz6NtLQ0xMTEoG3btlBKmX0oZbI/SimMGzcOTz75JPbs2WMY9+kvf/kL\nJk2aZHhBMD4+3m4fazRFo2PL0O0gu3XrluGFpMacPn3aZl9lFhG4u7sb1faBBx6o8UOsZcuWmDJl\nCn7961+jVatWOH36NN577z3DjeXIyEjs37/froZVJfP74YcfMHz4cMTFxWHFihVIS0tD9+7dcfLk\nSbz++utYvHgxXnzxRZw9exadO3fWu1yrwXA3QtWlmb59+zba1tIvKelp2rRpdS4bNmwYOnbsCE9P\nzxp/DsOGDcP27dvh7++PxMRE3HvvvZYsl2zUjz/+iLfffhudO3dGaGgoysvLER8fj0GDBqG0tBRL\nly5FaWlprRFUHR3D3Qh9+/ZFTk6OUeF+6tQpdOvWzQJVWaeqSzV1GTRoEP7whz8gOTkZJ0+eRHR0\ndI2RLYnqcucZubOzM+655x4At6emvH79OjIyMqz+Uqil8Zq7Eaqe/zZG1QS7VJubmxs+/fRTPPDA\nA3jttdewePFivUsiK3ft2rVGXzbq378/Fi5ciAEDBlioKtvAcDdCaGio4Y3MxmRmZlrFEzHWqn//\n/njssccQExPDyUWoUcnJyQgPD2+wzfPPP48LFy5gyJAhFqrKNjDcjdClSxecOmXczIEpKSl1TrJB\nNbVp0wbXrl3TuwyycikpKQgLC2uwTffu3fHdd9/VO7uZo2K4G0FEDE/MNEQphRMnTvAZWyN16tQJ\nZ86c0bsMsmK8zGk6hruRjJmA+tixY42+/EP/NWjQoFpvtBJVV1BQgLvuukvvMmwSw91IEyZMwPz5\n81FaWlpvmy1btiA2NtaCVdm2uoYrIKpOKeXwwwiYin9qRho2bBimTp2KuLi4emci2rp1K0aOHGnh\nymzXgAEDkJKSoncZZKWKiorg6empdxk2i+HeBPfeey/GjBmDvXv31lqnlMKFCxfQoUMHHSqzTa1a\ntUJJSQlncKI6NTbXMDWM4d5E0dHR2L17d63lGRkZCA4O1qEi23bXXXfhxIkTepdBVqhqKGoyDcO9\nicLDw+scAnjv3r2Gt+bIeGFhYbVmhtLC1q1bsWbNGs37JctJTExERESE3mXYLIZ7E7Vs2RJlZWUo\nLy+vsbxqWjpqmuDgYGRkZGjaZ2pqKt5//318+OGHmvZLllVYWOjQQ3k0F8PdBMHBwcjMzKyxLDc3\nl49BmsAc4f7OO+/g3XffhZeXF2d/slHXr1+vNfsXNQ3D3QTR0dHYs2dPjWXl5eVwdnbWqSLb5e3t\nrWkAHz16FEop+Pn5YejQodi1a5dmfZPlJCQkYODAgXqXYdMY7iaIioqqcVP1/Pnzuk+dZ8ucnJxq\nXeYyRXl5Od544w384Q9/AHD78VWGu23aunUr3xlpJoa7Cbp06YL8/Hz89NNPUEohISGBd/WbwdfX\nF7m5uc3qo6KiAg8//DBiY2MNY/sEBgYaPZonWZeEhATeTG0mhruJVq5ciaSkJMyYMQPffPMNJk+e\nrHdJNsvPzw85OTnN6mPx4sWIiYnBU089ZVjWsmXLBt8oJut0/fp1uLq6Gj3zGdWN4W4ib29v/PnP\nf0bXrl1x/Phx9OrVS++SbFbVZCjNsWHDBkyfPr3W8tatW3P0SRuTmZnJJ880wHBvpj/96U/47LPP\n9C7DpvXt2xdHjhxpVh8lJSV1Pl0RGBhomFCZbENmZiZfCNQAw72ZWrRowXGkm6lXr17Iy8szefuL\nFy/WO11fcHAw0tLSTO6bLC8jI4MT3mjAqHAXkbEiclhEckTk5TrWe4rITyJyUEQyRKT2zMlE9XB1\ndUVZWZnJ2x88eLDeG9ojRozApk2bTO6bLC87Oxv+/v56l2HzGg13EXEG8C8A4wAEAXhERO78sTob\nQKZSagCAEQDeExFORU5Gc3Z2NjngU1JS6p0/s3PnzigqKkJJSUlzyiMLunnzJlq2bKl3GTbPmDP3\nwQBylFK5SqlSAF8DuP+ONgpAGxERAB4ALgIw/VSMHE6fPn1w9OhRk7ZtKNwBYPTo0fjxxx9NLY0s\n6ObNm3Bzc9O7DLtgTLh3B1BQ7fOJymXV/R+AQACFANIAzFFKVWhSITmE5gxDcO7cOXTq1Kne9c8+\n+yw+/vhjXL161dTyyERKKezZswc3b940qn1+fj6fPNOIVjdUxwBIAdANwN0A/k9E2t7ZSERmiEii\niCSeO3dOo12TPTA13G/evIkWLRq+AtiqVSvMnj2bTzVZWHZ2NsaPH48lS5Zg9uzZRm3DMZq0Y0y4\nnwRQfRLDHpXLqpsGYKW6LQfAMQC1RtlXSi1SSkUopSI6duxoas1kh4KCgmoNxmYMYx+bmzhxIlav\nXo2KCv6D0hIKCwvx7LPPYtGiRVi0aBFKS0tRUFDQ6HYMd+0YE+4JAPxEpHflTdKHAay+o81xAL8A\nABHpDMAfQPPeJyeH4uXlhcuXLzd5u82bNyM6OrrRdi4uLrjnnnsQHx9vSnnURCtWrMALL7xgmNx6\nzJgx2LJlS6Pb5ebm8tFijTQa7kqpMgC/BbABQBaAb5VSGSIyU0RmVjZ7E0C0iKQB2AxgnlLqvLmK\nJvvUpUsXo1442r17NzZt2oSSkhKsXbsWY8eONar/8ePHY926dc0tk4wQHx+PqKgow+fY2Fijwv3Y\nsWPw8fExY2WOw6hr7kqptUqpfkqpPkqpv1Yu+0Qp9Unl7wuVUqOVUiFKqf5Kqf+Ys2iyT6+88gpm\nzJiBd999t942K1euxEcffYTExEQ8+OCDmD59utFDLQ8ePBj79+/XqlxqwJ3zCXft2hVnzpxpdLsb\nN27A3d3dnKU5DBe9CyCq4ufnh8WLF+O1115DXl5erTM4pRQ+/PBDbNiwodGbqHVxdnaGp6cnLl++\nDE9PT42qpjudOnUKXbp0qbW8Xbt2KCoqgpeXV53bVVRU4PbT1KQFDj9AVsXf3x/PP/88Fi5cWGvd\n0aNHERQUZFKwVwkPDzfLnK30X6mpqQgLC6u1PCgoqMHLbgUFBejZs6c5S3MoDHeyOlFRUUhJSam1\nfMuWLc2ewGHgwIF1TnBO2jl06BACAwNrLW/scdfMzMw6tyPTMNzJ6ogIWrVqhevXr9dYvm3bNsTE\nxDSr7/DwcCQnJzerD2rYoUOHEBBQ60noRsM9KyuL4a4hhjtZpYiICCQmJtZYdunSpRo36Uzh7e2N\nixcvNqsPalhBQQF69OhRa3mfPn0aHLc/KyuLo0FqiOFOVunOSci1vAnq7e0NviFtPkopODnVjhYX\nFxeUlZVBKVXndoWFhejWrZu5y3MYDHeySpGRkdi7d6/hc0pKimbz1A4bNgw7d+7UpC+qqaioCG3b\n1hp5xKC+KRWVUlBK8WkZDTHcySq1atUKFRUVhgGnDhw4UOcTGKaIiYnBjh07NOmLajp8+HCd19ur\nDBkypM63hM+ePdvg4G/UdAx3slrVgyAhIQEDBw7UpF9/f3+kpKSYNNwBNay+m6lVhgwZUuNfZFV4\nvV17DHeyWvfddx++//57lJWV1XrjsTlEBG+99Rbi4uKwfv16Tfqk2xoL9969e9c5Xy6flNEew52s\nVmhoKLKzs40eHKwpoqKisG7dOixYsEDTfh1ddnY2/Pz86l0vIvDz80N2dnaN5Qx37THcyao999xz\nePvtt/Hggw9q3reHhwd8fHzqPJMk0xgzNkxcXBzWrFlTY1ldw01Q8zDcyarFxcVh8+bNZjurmzRp\nEqfg08iNGzeMmvt0+PDhNR5zBYDy8nK4uHCoKy0x3MmhRUZGcjgCjezbtw+RkZGNtvPw8EBxcbHh\n861btxjsZsBwJ4dm6iQhVNuOHTswfPhwo9p6eXnh0qVLAIDjx49zwDAzYLiTw6trHBtqusTERERE\nRBjV9u677zYMDnf06FH07dvXnKU5JIY7ObyQkBCkp6frXYZNu3LlCpydnY0ejjksLMww9HJOTg76\n9OljzvIcEsOdHN6wYcPw008/Yc6cOUhLS9O7HJu0evVq3HfffUa3rx7uR44c4Zm7GTDcyeHFxsYi\nJSUFmZmZWLZsmd7l2KTVq1fjV7/6ldHtO3bsiPPnb0+z3Niz8WQa3qImhyciWLZsGVq1aoV7770X\nFRUVdY5qSPW7du0a2rdv36Rt3NzcUFxcDKWU0fPgkvH4N5gIt+f3dHNzQ0BAAPLz8/Uux6acPn26\nzjlTGxMSEoIVK1bA39/fDFURw52oGj8/P76x2kQHDhwwaTjmsWPHYv78+ZoNCEc1MdyJqunXrx/D\nvYlMHY45OjoaW7ZswcSJE81QFTHciaqpa1AralhqaipCQ0NN2rZnz55o3bq1xhURwHAnqsHHxwd5\neXl6l2FTrl69qtkUiKQdhjtRNa6urigrKzN8dqQ3V3NycvDll1/WO8cp2RaGO9EdvLy8cPr0abz/\n/vsYM2aM3uVYxI0bNzBlyhTk5eVh5syZRgf85cuXG5wzlfTDcCe6w+TJk/Hqq69i3759CAwMxKFD\nh/QuyexSUlIwbtw4vPrqq+jevTu+//57o7Y7duwYevfubebqyBQMd6I7jBs3Dhs3bsR7772HqVOn\n4quvvtK7JLPbv38/Bg8eDACYPn260dMP5ubmwtfX15ylkYkY7kR3cHNzQ05ODrp164aIiAgcPHhQ\n75LMLiEhAYMGDQIAdO/eHYWFhUZtxzN368VwJ6pD1eiGLVu2xM2bN3WuxvwuXrwIb29vw2dPT0/D\neOsN4YiO1ovhTtQILy8vFBUV6V2G2Vy8eBFeXl41lg0ZMgTx8fGNbpufn8+JNqwUw52oEaGhoXY9\nFHBiYqLhkkyVESNGYOvWrUZtz0HWrBO/FaJGhIaGIjU1Ve8yzKb6zdQqxhzzuXPn0LFjR3OWRs1g\nVLiLyFgROSwiOSLycj1tRohIiohkiMh2bcsk0o+9h3tSUlKtsWGcnJzg5eWFgoKCerfLyMhAcHCw\nucsjEzUa7iLiDOBfAMYBCALwiIgE3dHGC8BHACYopYIBPGiGWol00aNHjwZDzpZdvXoVZWVlcHd3\nr7XulVdewe9+9zvMnDkTt27dqrU+LS0NQUFBtZaTdTDmzH0wgBylVK5SqhTA1wDuv6PNowBWKqWO\nA4BS6qy2ZRLpR0QgIqioqNC7FM2tWLECkydPrnPdgAED8MMPPyA8PBwLFy6stX7Tpk2IiYkxd4lk\nImPCvTuA6qctJyqXVdcPQDsR2SYiSSLyhFYFElmD3r172+WAYitWrMCkSZMabDN9+nSsXLmyxpg7\nhYWF8PDwQJs2bcxdIplIqxuqLgAGArgXwBgAr4lIvzsbicgMEUkUkcRz585ptGsi8xswYACSk5P1\nLkNTx48fR7t27RoNaGdnZ4waNQrbt//3VtrSpUvx2GOPmbtEagZjwv0kgLuqfe5Ruay6EwA2KKWK\nlVLnAewAMODOjpRSi5RSEUqpCN5lJ1syYsQIbNu2Te8yNLV8+XKjA3ry5MmG8WZu3bqFTZs2YezY\nseYsj5rJmHBPAOAnIr1FpAWAhwGsvqPNjwCGioiLiLgDiASQpW2pRPrp27ev3c3QtHPnTsTGxhrV\ntl+/fsjJyUF5eTnWrl2LuLg4Pt9u5Rr9dpRSZQB+C2ADbgf2t0qpDBGZKSIzK9tkAVgPIBXAfgCL\nlVLp5iubyLJExK6emrl69Src3Nzg6upq9DbDhg3Drl278N133+Ghhx4yY3WkBRdjGiml1gJYe8ey\nT+74/A6Ad7Qrjci6jBo1Cps3b8aTTz6pdynNtm3bNowYMaJJ2zz44IP4/e9/D1dXV3Tr1s08hZFm\n+O8qIiPFxsZi8+bNepehiY0bN2LUqFFN2iYwMBBvvvkm3n//fTNVRVoy6sydiIDOnTvj7NmzUEpB\nRPQup1kOHTqEgICAJm935xg0ZL145k7UBMHBwcjMzNS7jGY5ceIEunfvbvM/oKhhDHeiJhg1ahQ2\nbdqkdxl9aVA6AAAOaElEQVTNsnHjRvzyl7/UuwwyM4Y7URMMHz68xss8tmjTpk1Nvt5OtofhTtQE\nHh4euHnzZp0DadmCiooKnD9/Hp06ddK7FDIzhjtRE0VGRmL//v16l2GS9PR0hISE6F0GWQDDnaiJ\nbPm6+759+xAdHa13GWQBDHeiJho0aBD27dundxkmSU5ORnh4uN5lkAUw3ImayNXVFU5OTjY5vnte\nXh569eqldxlkAQx3IhP06NEDJ0/eOTiqdSsrK4OLiwufb3cQDHciE/Tr1w/Z2dl6l9Ek6enpnBbP\ngTDciUxgi+G+fft2TovnQBjuRCawxXDfuXMnhg4dqncZZCEMdyIT9O7dG7m5uXqXYbSKigpcu3YN\nbdu21bsUshCGO5EJXF1dbeot1V27dmHIkCF6l0EWxHAnMlHr1q1x7do1vcswyjfffMPZkxwMw53I\nRMHBwcjIyNC7jEYppZCdnY3AwEC9SyELYrgTmah///5IT7f+qYJPnTqFnj176l0GWRjDnchEISEh\nSEtL07uMRqWmpnKwMAfEcCcyUZ8+fZCamgqllN6lNCgtLQ2hoaF6l0EWxnAnMpGLiwuGDx+OjRs3\n6l1Kg9LS0njm7oAY7kTN8Jvf/AYfffSR3mU06Pz58+jYsaPeZZCFMdyJmqFTp07o3LkzUlNT9S6l\nXhwozDEx3ImaadasWVi6dKneZdTpypUraNOmjd5lkA4Y7kTNFBoaarXPuxcUFPAxSAfFcCdqJicn\nJ7Rt2xaXL1/Wu5Rajh8/znB3UAx3Ig0MHToUu3bt0ruMWhjujovhTqSByMhIJCQk6F1GLfn5+Qx3\nB8VwJ9JASEiIVQ5FcPz4cc6Z6qAY7kQacHd3x/Xr1/Uuo5aLFy/Cy8tL7zJIBwx3Io14enqiqKhI\n7zIMqoZF4HPujonhTqSRAQMG4ODBg3qXYXDmzBl06dJF7zJIJwx3Io3079/fqp53z8zMRFBQkN5l\nkE6MCncRGSsih0UkR0RebqDdIBEpE5HJ2pVIZBuCgoKQmZmpdxkGWVlZnKDDgTUa7iLiDOBfAMYB\nCALwiIjUOh2obPd3AP9P6yKJbEGvXr2Ql5endxkGGRkZPHN3YMacuQ8GkKOUylVKlQL4GsD9dbR7\nDsAKAGc1rI/IZjg7O6OiokLvMgyOHj0KHx8fvcsgnRgT7t0BFFT7fKJymYGIdAcwEcDHDXUkIjNE\nJFFEEs+dO9fUWomsXvv27XHhwgW9y8CZM2fQqVMnPinjwLS6ofo+gHlKqQZPW5RSi5RSEUqpCI4v\nTfYoKCgIWVlZepeBrVu3YuTIkXqXQToyJtxPArir2ucelcuqiwDwtYjkAZgM4CMR+ZUmFRLZEGu5\nqbpz507ExMToXQbpyJhwTwDgJyK9RaQFgIcBrK7eQCnVWynlo5TyAfA9gN8opVZpXi2RlbOWcM/N\nzYWvr6/eZZCOXBproJQqE5HfAtgAwBnA50qpDBGZWbn+EzPXSGQzfH19kZubq2sNt27dgouLC6+3\nO7hGwx0AlFJrAay9Y1mdoa6UerL5ZRHZJhcXF5SXl+taQ0ZGBvr3769rDaQ/vqFKpDG9J+5ITExE\nRESEbvsn68BwJ9JYv379cOTIEd32n5SUhIEDB+q2f7IODHcijfn5+SE7O1u3/efl5XEMd2K4E2nN\nz89PtzP30tJSuLq68mYqMdyJtKbnmXt6ejpCQkJ02TdZF4Y7kcbat2+PS5cu6bLv/fv3Y9CgQbrs\nm6wLw53ITKpmQrKkPXv2ICoqyuL7JevDcCcyg169eiE/P9/i+z179iw6d+5s8f2S9WG4E5lBWFgY\nkpOTLbrPs2fPggPyURWGO5EZhIeH48CBAxbd54EDB/h8Oxkw3InMoH///khLS7PoPtPT0znsABkw\n3InMwM3NDSKCY8eOWWyfHFOGqmO4E5nJO++8g+nTp+P48eMW2d/p06d5M5UMGO5EZtK3b1988skn\nePrpp83+WGRFRQVEhG+mkgHDnciM+vXrh8jISGzdutWs+8nPz+d4MlQDw53IzGbNmoUvvvjCrPvg\n9Xa6E8OdyMy6deuGs2fPoqKiwfnjmyU9PR3BwcFm659sD8OdyALM/dx7RkYGw51qYLgTWcCoUaOw\nadMms/V/4cIFdOjQwWz9k+1huBNZQEREBJKSkszSd2lpKZydnc3SN9kuhjuRBXh6euLKlStm6XvP\nnj0YMmSIWfom28VwJ7IQb29vXLhwQfN+N2zYgDFjxmjeL9k2hjuRhZjr0kxycjLCw8M175dsG8Od\nyEJCQkKQnp6uaZ9FRUXw9PSEkxP/V6aa+DeCyEKCgoKQlZWlaZ+7d+/GPffco2mfZB8Y7kQW0rVr\nVxQWFmra586dOzF8+HBN+yT7wHAnshARgVJK00HE0tLSEBoaqll/ZD8Y7kQW1KVLF5w5c0aTvqp+\nUPAZd6oLw53IggIDAzW77n7u3Dl06tRJk77I/jDciSwoKCgImZmZmvR16NAhBAQEaNIX2R+GO5EF\naXnmznCnhjDciSyoV69eyMvL06SvQ4cOITAwUJO+yP4w3IksyNnZWbNx3XNycuDr66tJX2R/jAp3\nERkrIodFJEdEXq5j/WMikioiaSKyR0QGaF8qkX1o164dLl261Kw+lFIoKyuDq6urRlWRvWk03EXE\nGcC/AIwDEATgEREJuqPZMQAxSqkQAG8CWKR1oUT2Qovr7jk5OfDz89OoIrJHxpy5DwaQo5TKVUqV\nAvgawP3VGyil9iilqk5F4gH00LZMIvsRGBjY7Cdm9u3bh8GDB2tUEdkjY8K9O4CCap9PVC6rz3QA\n65pTFJE902KMmR07dnAMd2qQi5adichI3A73ofWsnwFgBgD07NlTy10T2Yy+ffviyJEjJm9/9epV\nnDhxAn369NGwKrI3xpy5nwRwV7XPPSqX1SAioQAWA7hfKVXnjARKqUVKqQilVETHjh1NqZfI5rm6\nuqKsrAxXr17FsGHDMHfu3CZt/+WXX+Lxxx83U3VkL4wJ9wQAfiLSW0RaAHgYwOrqDUSkJ4CVAKYo\npbK1L5PIvvj4+ODJJ5/E888/j6ysLJSVlRm97erVqzFx4kQzVkf2oNFwV0qVAfgtgA0AsgB8q5TK\nEJGZIjKzstnrALwBfCQiKSKSaLaKiezAggULEB4ejokTJyIsLAwpKSlGbZefn4/OnTujVatWZq6Q\nbJ1oOfxoU0RERKjERP4MINq4cSMOHjyIF198sdG2f//73xEWFobRo0dboDKyRiKSpJSKaKwd31Al\n0ll0dDR2795tVNvNmzcjNjbWzBWRPWC4E+msdevWuHnzZqPX3TMyMuDv7w8XF00fciM7xXAnsgLh\n4eE4cOBAg22++uorPPLIIxaqiGwdw53ICsTGxmL16tX1rldKIT4+HlFRURasimwZw53ICowcORIZ\nGRn1PjWzZcsWREdHQ0QsXBnZKoY7kRUQESxcuBAvvfQSqj9FVlxcjLS0NLz77ruYM2eOjhWSreGd\nGSIr0bFjR3z33Xd49NFHERISgqKiIhQUFCA0NBSzZs2Ct7e33iWSDWG4E1kRLy8vrFmzBgcPHoSH\nhweH9SWTMdyJrIyTkxPCwsL0LoNsHK+5ExHZIYY7EZEdYrgTEdkhhjsRkR1iuBMR2SGGOxGRHWK4\nExHZIYY7EZEd0m0mJhE5ByBfl53XrwOA83oXoQMet+Nx1GO3h+PupZTq2Fgj3cLdGolIojHTV9kb\nHrfjcdRjd6Tj5mUZIiI7xHAnIrJDDPeaFuldgE543I7HUY/dYY6b19yJiOwQz9yJiOyQ3Ye7iOSJ\nSJqIpIhIYuWyB0UkQ0QqRCSiWltXEVla2T5LRF6ptu6vIlIgItf0OI6m0uK4RcRdRH4WkUOV2/1N\nr+Mxlobf93oROVi53Sci4qzH8TSFVsderc1qEUm35DGYQsPvfJuIHK7sJ0VEOulxPJpRStn1LwB5\nADrcsSwQgD+AbQAiqi1/FMDXlb93r9zWp/LzEABdAVzT+5gsddyVvx9ZubwFgJ0Axul9bBb6vttW\n/lcArADwsN7HZqljr1w2CcByAOl6H5cFv/MabW39l0POxKSUygJQ10zyCkBrEXEB0ApAKYArldvE\n17ONzWjqcSulrgPYWrltqYgkA+hhuYq1YeL3faWyjQtu/2CzyZtTphy7iHgAeAHADADfWqxYDZly\n3PbG7i/L4PaXuUlEkkRkRiNtvwdQDOAUgOMA3lVKXTR3gWai6XGLiBeA+wBsNkexGtLsuEVkA4Cz\nAK5WtrV2Wh37mwDeA3DdbJVqS8u/60srL8m8JrZ8JgfHmEN1qFLqZOX1s40ickgptaOetoMBlAPo\nBqAdgJ0iskkplWupYjWk2XFXnuV8BeADG/iz0Oy4lVJjRMQNwJcAYgFstED9zdHsYwfQFkAfpdTv\nRcTHEkVrQKvv/LHKftrg9qW4KQCWWeIAzMHuz9yVUicr/3sWwA+4/eXW51EA65VStyrb7wZgk68q\na3zciwAcUUq9b656taL1962UKgHwI4D7zVOxdjQ69igAESKSB2AXgH4iss2cdTeXVt95tX6u4vb9\nhob6sXp2He4i0rrypzBEpDWA0QAauvt/HLfP0KraDwFwyNx1ak3L4xaRvwDwBPC8OWvWglbHLSIe\nItK1crkLgHth5X8PtDp2pdTHSqluSikfAEMBZCulRpiz9ubQ8Dt3EZEOlctdAcQ10o/10/uOrjl/\nAfAFcLDyVwaAVyuXTwRwAsBNAGcAbKhc7gHgu8q2mQBeqtbX25XbVFT+9896H5+5jxu3b54qAFkA\nUip/Pa338VnguDsDSACQitv/g38IwEXv47PU3/VqffrAyp+W0fA7bw0gqfI7zwDwTwDOeh9fc37x\nDVUiIjtk15dliIgcFcOdiMgOMdyJiOwQw52IyA4x3ImI7BDDnYjIDjHciYjsEMOdiMgO/X+EOgB5\nB3uArAAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#ok, now let's turn adventure mode on\n", "#we're going to perform a high-resolution spectrum deblending. let's have a look at it:\n", "spec = np.genfromtxt('mg.dat', unpack=True)\n", "plt.plot(spec[0], spec[1], 'k-', lw=0.7)" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": true }, "outputs": [], "source": [ "#for that, we'll use a constant model for the continuum, a voigt for the Mg line and gaussian for additional lines\n", "from lmfit.models import VoigtModel, GaussianModel, ConstantModel\n", "\n", "#on LmFit we can combine several models to fit data, so to make things easier, let's separate them\n", "#calling a the model to make the continuum:\n", "cont = ConstantModel(prefix='cont_')\n", "#on built-in models we can guess the parameters instead of setting each one\n", "pars = cont.guess(spec[1], x=spec[0])\n", "#this will create the Parameters object that will store every parameter of each model" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": true }, "outputs": [], "source": [ "#now the Voigt model for the main Mg line\n", "mg = VoigtModel(prefix='Mg_')\n", "#we can update the parameters using the guess\n", "pars.update(mg.guess(spec[1], x=spec[0]))\n", "#but if we take it like that it's going to be a mess, so we set some constraints\n", "pars['Mg_center'].set(5183.5)\n", "pars['Mg_sigma'].set(0.2)\n", "pars['Mg_gamma'].set(0.2, vary=True)\n", "#now, we have 4 lines to fit a gaussian, and will be easier to create a function" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def call_gauss(cen, no, pars):\n", " \n", " #since we can have a lot of parameters, it is good to identify each one\n", " label='g'+str(no)+'_'\n", " gauss = GaussianModel(prefix=label)\n", " #so they would be called g1_sigma for example.\n", " \n", " #now we're going to include this model parameters inside the existing Parameters object\n", " #I'll do it separately, because using guess would use the entire data\n", " pars.update(gauss.make_params())\n", " #and give them initial values\n", " pars[label+'center'].set(cen, min=cen-0.05, max=cen+0.05)\n", " pars[label+'amplitude'].set(-1, max=-0.001)\n", " pars[label+'sigma'].set(0.2, max=0.3)\n", " #now we return this Model object, since we already modified the Parameters object\n", " return gauss" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n" ] } ], "source": [ "#here we have the approximate lambdas\n", "lbd = [5181.2, 5182.7, 5184.2, 5184.6]\n", "#so let's sum the models\n", "finmod = cont + mg\n", "for n, l in enumerate(lbd):\n", " finmod = finmod + call_gauss(l, n, pars)\n", "print(finmod)" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[Model]]\n", " (((((Model(constant, prefix='cont_') + Model(voigt, prefix='Mg_')) + Model(gaussian, prefix='g0_')) + Model(gaussian, prefix='g1_')) + Model(gaussian, prefix='g2_')) + Model(gaussian, prefix='g3_'))\n", "[[Fit Statistics]]\n", " # function evals = 297\n", " # data points = 341\n", " # variables = 17\n", " chi-square = 0.169\n", " reduced chi-square = 0.001\n", " Akaike info crit = -2560.961\n", " Bayesian info crit = -2495.819\n", "[[Variables]]\n", " cont_c: 1.03932300 +/- 0.003784 (0.36%) (init= 0.7928076)\n", " Mg_amplitude: -1.20464721 +/- 0.027784 (2.31%) (init= 10.18333)\n", " Mg_sigma: 0.00066706 +/- 3.622285 (543020.52%) (init= 0.2)\n", " Mg_center: 5183.61619 +/- 0.002489 (0.00%) (init= 5183.5)\n", " Mg_gamma: 0.42473265 +/- 0.017813 (4.19%) (init= 0.2)\n", " Mg_fwhm: 0.00240229 +/- 13.04497 (543020.52%) == '3.6013100*Mg_sigma'\n", " Mg_height: -720.449419 +/- 3.91e+06 (543022.53%) == '0.3989423*Mg_amplitude/max(1.e-15, Mg_sigma)'\n", " g0_sigma: 0.08432131 +/- 0.005053 (5.99%) (init= 0.2)\n", " g0_center: 5181.21352 +/- 0.004771 (0.00%) (init= 5181.2)\n", " g0_amplitude: -0.03739167 +/- 0.002145 (5.74%) (init=-1)\n", " g0_fwhm: 0.19856152 +/- 0.011900 (5.99%) == '2.3548200*g0_sigma'\n", " g0_height: -0.17690805 +/- 0.008844 (5.00%) == '0.3989423*g0_amplitude/max(1.e-15, g0_sigma)'\n", " g1_sigma: 0.08534003 +/- 0.005706 (6.69%) (init= 0.2)\n", " g1_center: 5182.65022 +/- 0.005175 (0.00%) (init= 5182.7)\n", " g1_amplitude: -0.03475291 +/- 0.002326 (6.70%) (init=-1)\n", " g1_fwhm: 0.20096042 +/- 0.013438 (6.69%) == '2.3548200*g1_sigma'\n", " g1_height: -0.16246078 +/- 0.008907 (5.48%) == '0.3989423*g1_amplitude/max(1.e-15, g1_sigma)'\n", " g2_sigma: 0.07452672 +/- 0.004988 (6.69%) (init= 0.2)\n", " g2_center: 5184.24999 +/- 0.004427 (0.00%) (init= 5184.2)\n", " g2_amplitude: -0.03397825 +/- 0.002432 (7.16%) (init=-1)\n", " g2_fwhm: 0.17549701 +/- 0.011747 (6.69%) == '2.3548200*g2_sigma'\n", " g2_height: -0.18188597 +/- 0.010123 (5.57%) == '0.3989423*g2_amplitude/max(1.e-15, g2_sigma)'\n", " g3_sigma: 0.06498473 +/- 0.004030 (6.20%) (init= 0.2)\n", " g3_center: 5184.55809 +/- 0.003775 (0.00%) (init= 5184.6)\n", " g3_amplitude: -0.03240719 +/- 0.001961 (6.05%) (init=-1)\n", " g3_fwhm: 0.15302735 +/- 0.009491 (6.20%) == '2.3548200*g3_sigma'\n", " g3_height: -0.19894829 +/- 0.010362 (5.21%) == '0.3989423*g3_amplitude/max(1.e-15, g3_sigma)'\n", "[[Correlations]] (unreported correlations are < 0.100)\n", " C(Mg_amplitude, Mg_gamma) = -0.961 \n", " C(Mg_sigma, Mg_gamma) = -0.956 \n", " C(cont_c, Mg_amplitude) = -0.905 \n", " C(Mg_amplitude, Mg_sigma) = 0.872 \n", " C(cont_c, Mg_gamma) = 0.834 \n", " C(cont_c, Mg_sigma) = -0.747 \n", " C(g2_sigma, g2_amplitude) = -0.679 \n", " C(g1_sigma, g1_amplitude) = -0.664 \n", " C(g3_sigma, g3_amplitude) = -0.639 \n", " C(g0_sigma, g0_amplitude) = -0.638 \n", " C(Mg_amplitude, g2_amplitude) = -0.557 \n", " C(Mg_gamma, g2_amplitude) = 0.519 \n", " C(Mg_gamma, g1_amplitude) = 0.514 \n", " C(cont_c, g0_amplitude) = -0.493 \n", " C(Mg_amplitude, g1_amplitude) = -0.481 \n", " C(Mg_sigma, g1_amplitude) = -0.475 \n", " C(Mg_gamma, g3_amplitude) = 0.469 \n", " C(Mg_amplitude, g3_amplitude) = -0.440 \n", " C(Mg_sigma, g3_amplitude) = -0.433 \n", " C(Mg_sigma, g2_amplitude) = -0.410 \n", " C(cont_c, g2_amplitude) = 0.401 \n", " C(Mg_center, g2_amplitude) = 0.395 \n", " C(Mg_amplitude, g0_amplitude) = 0.379 \n", " C(Mg_gamma, g1_sigma) = -0.336 \n", " C(Mg_gamma, g0_amplitude) = -0.319 \n", " C(g1_amplitude, g2_amplitude) = 0.318 \n", " C(Mg_amplitude, g1_sigma) = 0.316 \n", " C(cont_c, g0_sigma) = 0.313 \n", " C(g1_amplitude, g3_amplitude) = 0.311 \n", " C(Mg_sigma, g1_sigma) = 0.307 \n", " C(Mg_amplitude, g2_sigma) = 0.306 \n", " C(cont_c, g1_amplitude) = 0.292 \n", " C(Mg_sigma, g0_amplitude) = 0.272 \n", " C(cont_c, g3_amplitude) = 0.272 \n", " C(Mg_gamma, g2_sigma) = -0.270 \n", " C(Mg_center, g2_sigma) = -0.263 \n", " C(g2_amplitude, g3_amplitude) = 0.255 \n", " C(Mg_gamma, g3_sigma) = -0.249 \n", " C(Mg_amplitude, g0_sigma) = -0.240 \n", " C(Mg_sigma, g3_sigma) = 0.234 \n", " C(Mg_amplitude, g3_sigma) = 0.228 \n", " C(cont_c, g2_sigma) = -0.218 \n", " C(g1_sigma, g2_amplitude) = -0.209 \n", " C(Mg_amplitude, Mg_center) = -0.204 \n", " C(g1_sigma, g3_amplitude) = -0.203 \n", " C(Mg_center, g3_amplitude) = 0.202 \n", " C(Mg_gamma, g0_sigma) = 0.202 \n", " C(Mg_sigma, g2_sigma) = 0.197 \n", " C(cont_c, g1_sigma) = -0.191 \n", " C(Mg_center, Mg_gamma) = 0.186 \n", " C(Mg_sigma, g0_sigma) = -0.172 \n", " C(g1_amplitude, g3_sigma) = -0.168 \n", " C(g1_amplitude, g2_sigma) = -0.168 \n", " C(cont_c, Mg_center) = 0.146 \n", " C(Mg_sigma, Mg_center) = -0.142 \n", " C(cont_c, g3_sigma) = -0.136 \n", " C(g0_amplitude, g2_amplitude) = -0.135 \n", " C(g2_sigma, g3_center) = 0.124 \n", " C(g2_center, g3_sigma) = -0.122 \n", " C(Mg_amplitude, g1_center) = 0.121 \n", " C(cont_c, g1_center) = -0.112 \n", " C(g1_sigma, g2_sigma) = 0.111 \n", " C(g1_sigma, g3_sigma) = 0.110 \n", " C(Mg_gamma, g1_center) = -0.110 \n", " C(g2_sigma, g3_sigma) = -0.101 \n", "\n" ] } ], "source": [ "#now we perform the fit\n", "result = finmod.fit(spec[1], pars, x=spec[0])\n", "print(result.fit_report())" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Mg 5183.62\n", "g0 5181.21\n", "g1 5182.65\n", "g2 5184.25\n", "g3 5184.56\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XdYFWf6//H3Q28iCBYUBQsKgmUVsYtRY4vRVONq4qYa\nN8VsevLbkmySzabuN5uiiSbZJLrRTTUmMfbeBbsidoMVFRWVzrl/fxwkKCoHOIcDh/t1XVx4Zp6Z\n+YwHbuY8M/OMERGUUkq5FjdnB1BKKWV/WtyVUsoFaXFXSikXpMVdKaVckBZ3pZRyQVrclVLKBWlx\nV0opF6TFXSmlXJAWd6WUckEeztpwaGioREZGOmvzSilVIyUnJ58UkfpltXNacY+MjCQpKclZm1dK\nqRrJGHPQlnbaLaOUUi5Ii7tSSrkgLe5KKeWCtLgrpZQL0uKulFIuSIu7Ukq5IC3uSinlglyyuGdk\nZPD222+zePFiZ0dRSimncNpNTI70xRdf8NRTTwGQnZ2Nj4+PkxMppVTVcskj9y1bthT/+/Dhw05M\nopRSzuGyxd3NzbpraWlpTk6jlFJVz+WKe0FBAdu3b2fo0KEAHDp0yMmJlFI13d69e2vcgaLLFfc9\ne/aQk5PDkCFDAC3uSqnK+fHHH2nVqhW9evWioKDA2XFsViOL+7lz55g8eTKTJk1CRADrEfs777zD\nG2+8AUD//v0JCgq6anHPzs5m3rx5Ne6vsVKqar3++usA/Prrr/zvf/9zcppyEBGnfHXu3FkqYsmS\nJRIeHi6AAPLxxx+LiMi4ceOKpw0aNEhERNq1ayfDhw+/4nqeeeYZAcTb21vmzp1boSxKKdeWm5sr\n3t7e8vjjj0vr1q2lf//+Ni2zfv16ycnJcUgmIElsqLE17sg9KCiI+vXrs2zZMvr27cszzzxDamoq\nkydP5o477qB37968+uqrAISHh1/xyP306dNMnDiRTp06ER0dzR133MHp06erelcuIUWfQJRS1cfm\nzZvJzc2lR48eDB8+nOXLl3P+/HkANmzYwEcffYTFYiluf+7cOWJiYujSpQsJCQmcOHHCWdFr3pG7\niIjFYhERkZ9++kkASUhIEEAOHjx4SbsHHnhAGjRoUGr5jz76SABJTk6WLVu2iDFGnnjiiQrnqaxJ\nkyZJs2bNSuW3RVpamgMSKaVERP79738LIGlpabJw4UIB5IcffpDly5cX9xRMmjSpuP3//d//CSAv\nvviiuLm5ybPPPmv3TNh45F4ji/tFeXl5EhwcLID06dOn1PyXXnpJgFIfj4YPHy6RkZHFfyQeeOAB\nAeS1116TzMzMSucqrwYNGgggvXv3LtdyW7duFUD+9re/OSiZUrVLbm6uPPbYY7J3714REfn9738v\n4eHhxfMCAgLkj3/8o4wfP178/f2lU6dOUrduXTl79qwUFhZKZGRk8e/xyJEjJTAwUM6cOVO8fovF\nIjfddJNMmzatwhldtrjPnSvSrp1I9+4iAweKdO16SOLikmTs2LPy7LMiL78s8s47Ih9/LDJ+/GKB\nITJ9+mFJThbZtUvk+PEc8fX1k4ceeqh4nVlZWdK3b18BpF+/flJYWFihbBWRmpoqgDRp0kQAOXDg\ngM3LTpo0qfjo4fDhww5MqVTtMG/ePAGkcePGIiISGRkpt912W/H8ESNGSEREhISFhcltt90ma9as\nEUDef/99Wbp0qQDFhTspKUkAeeONN4qX//rrrwWQKVOmVDijyxb35ctFbrpJpH9/ka5dRWJjRSIi\nROrVE/Hysu5R2V8XpEmTczJ6tMhHH4lkZFj/ol4sllOnTq1Qtop49dVXBZAFCxYIIO+9957Ny953\n333Fxf3bb791YEqlaoe//vWvl/xOAfL2228Xzy95QPXf//5XLBaLxMfHS1RUlNx5553i7+8v58+f\nL27fv39/CQsLk9OnT8v+/fulSZMmEhcXJwUFBRXOaLfiDnwKpAPbrjLfAO8Ce4AtQCdbNmyPbpkr\nyc0VOXVK5OBBkVmz9ggkyPPPL5CZM0W++EKkd+/vxd39HRk+PF8aN7b+D9SpI/L22yIFBRapX7++\n3H333Q7JdiWdOnWSrl27iohI69atZfDgwTYv26FDB+nZs6cA8vLLLzsqolK1Rt++faVNmzbi5eUl\nzZs3F0BWrlxZPP/gwYMCSNeuXSU/P19ERH755Zfigj9+/PhL1rd69Wrx8PCQmJgYqVu3rgQFBcnG\njRsrldGexb0P0OkaxX0o8EtRke8GrLVlw44q7iVlZmYKIK+//nrxtKioKBkyZIiIiFgsIsnJIjfc\nYP2fGDNGZOjQGyUmJsbh2UREpk79ViBannrqU9m7V+S++x6TwMBAm7qFsrKyxN3dXf785z9LZGSk\njBo1qgoSK+W6cnJyxMfHRx5//HG5/vrrBZB69epJdnb2Je32799fXNgv+vOf/yz33HPPFS9/nDFj\nhnTt2lVGjhwpKSkplc5p124ZIPIaxf0j4PclXqcCYWWtsyqKu4hI3bp15ZFHHhERkV27dl2x68Ni\nsfbVg0ivXivFGHPJSRB7On9e5KWXRBo3Pn+VLqPz0rRpriQmitx9t8iLL4p8/HGePPHE9/LMM+9K\nQYH1JPDF/ryvv/5abrjhBmnfvr1D8ipVW6xYsUIA+f77753SRWsrW4u7PYb8bQKUvM3zUNG0o5c3\nNMaMA8YBNGvWzA6bLlvJa91nz54NUDzuzG+54C9/gbQ0mDKlOyI9Wb9+PQMGDLBrlt27YcgQ2LsX\n6tdPxd9/Hh988CS+vp5kZcG2bcd5++3/0rjxTRQUtOCXXwo4ftwD8ARuAuBf/yokMtIdP7+GwNvs\n2NGLOnVOkJLyD/LzC/D0dMlRnJVyuGXLlgHQq1cvgoOD6dGjB+3bt6/UOgsLYd8+2LEDdu6EXbus\ndWDUKHjoIXukvgZb/gJw7SP3n4BeJV4vBOLLWmdVHbkPGTKkuJtl4MCBEh0dfdW258+LhIcXCKyV\nt9/+l11zpKeLNG0qEhoqMnPmaXFzc5Pnn3/+kjaFhYVSt25defDBB0VE5KmnnhLwll697pVXXlkv\nTZq8LMHBH8ntt1ukQYNDAhcuOeoPDMyXbt2sR/xvvSWyYIHIyZN23Q2lXNbgwYOlbdu2FVq2sFAk\nNVXkm2+sn8xHjRJp317E2/vST+YNG4r07i3y6acVz0kVHrkfBpqWeB1eNK1aGDp0KI8++ijfffcd\nS5Ys4dFHH71qW39/ePlld+65J4E5c9bzxBP2y3HffXDiBKxcCefObcZisZCYmHhJGzc3Nzp16sSm\nTZsoKChg2rRpjBgxmJkzPwGgXr31PPTQg7z4Yi+efPJJjh49zqxZG5g/P43773+TTp3+iDExzJkD\nn33223qbNoWOHX/76tQJIiKsn1iUUnD+/HmWL1/O2LFjy2xbUAApKbBhw29fmzZB0Y2rADRvDm3b\nwqBB1u9t20J0NAQGOnAnLmfLXwCufeR+A5eeUF1nyzqr6sj95MmT4unpKYAEBgaWeUdnfr6It/cx\nqVt3vd0y/PRToYDIm29aX1/sz7vSHamPPfaY+Pn5yapVqwSQGTNmFM9LS0sTQF599VUJDw+XO++8\nU0R+O+IfN25ccdsTJ0Tmz7duc8wY6yWjbm6/HUGEhYnccot1/sqVIpedM1KqVvn4449LXRlzUUaG\nyM8/i/y//yeSmCji6/vb75Gfn0iPHiKPPGI9Gk9KsvYAOBJ2vFpmOtb+83ys/en3AeOB8UXzDfAB\nsBfYig1dMlKFxV3Eei15v3795Mcff7SpfZcu3wuIpKRYKrVdi8UiTzzxpHh4JEvTpjmSm2udPmHC\nBPH39y++Q7akTz75RAB59NFHBZBDhw5dMj8hIUFCQ0MFkIkTJxZPv/7666VDhw6yevVqmT9//hXz\nZGWJrFsn8sEH1oLfosVvP6ReXtYbw/78Z5FFi0QcNOaRUtVSt27dpG3btmKxWOTkSZH//U/kwQet\nB0UXf0fc3UXi40Uee0xk2jSRHTtEKnG5eoXZrbg76qsqi3t5vfjiRwKF8tRTlfsT/Nlnnwn0LOpr\ne6F4+sCBA+Vq+79+/XqhaLTKJk2alJr//fffCyBBQUFy7ty54ukvv/xy8XKAvPvuu/Ldd9/JzTff\nLKmpqVfNePSoyPffizz9tEi3btYfYLAenVx/vcjrr4ts2WK9okgpV5ScvE2gr/Tvv1Y6dxYx5uI5\nLJHBg61X0i1a5Pgjcltpca+E+fPnCyyU8PBzFS5qFy5ckCZNmkhIyGzx9c0W8C0eIqBp06YyZsyY\nqy5njBFAbr755lLzLRaLjB49Wv79739fMj07O1tuu+02iYuLkyFDhoiHh0fxjRWXn7i9kvz8fNm4\ncaOcOSMya5bIhAmXHrW0aCHypz+JLF5s7bpSqibLzLQenY8aZRFPzywBEQ8Pi/TubT0hunp19f05\n1+JeCdnZ2eLt/aiAyNatFVuHdQwJf/HxKZDhw48KILNmzZJz584JIK+88spVl7041vwHH3xQ7u1a\nLBY5evSoBAcHS9++fSUhIUESEhLKXO7iePgXB0y6uK7Dh61DNAwd+tuZ/3r1rFfkzJ/vnI+lStni\nwoULsn///uLXZ86I/O1v+2XIkIISV7EcF5gsw4ZNFieMGVghWtwr6cYbxwuIvPpqxarXo48+Kl5e\ndwuIzJ2bJW5ubvK3v/2t+OajssaCSU9Pr9T4EydPnpT8/Hx54YUXxM3NTTIyMkRErtjPX1hYWHyU\nf3HQo2PHjklMTEzxiJO7d++Wnj0HScOGf5TY2CQJDLQUn5h96KFsmThxteTna6VX1cfo0aMF3OX1\n1zfLqFEinp75AiL+/ukyYsRegd7Sp8918sorr1Tqd62qaXGvpB9//FFggzRqdPX+6mvp2LGjNGw4\nX+rXt14DGxcXJ0OHDpWpU6cKINu3b7dz4iu7eNfd9OnT5cEHH5R27dpJXl7eJW1Kjk198W7e4cOH\nF/fhHz58WG6++WapU6eODB06tGgMjT/J11+LDBuWJ5ArIBIZeVI+/dR64lYpZ9q1K1vgHwJHBESC\ngwvF0/MjqVNnQPHPeps2bST34lUONYgWdzvo0mWuQL5s316+B2KcOXNGwF18fS9I0dWKcs8990hI\nSIg8//zz4u7uXmU/VAUFBdKwYcPiH2hAZs6ceUmbN998UwCJi4uT+Ph4yc7OFi8vL7n55pvFw8ND\nEhMTBawPIBARufvuu8XX11dycnKKhjkNFnhI6tZNK+62efZZkePHq2QXlRIR60n/xYutl/i6uVkE\nCqRJk2SBm+T++x8WQJYuXSojRowoHmagJtLibgcTJ+4SEHn66UXlWm7ZsmUCCQIiX35pnXbxEseY\nmBhp06aNA9Je3R/+8AcBpH///tKoUSO58cYbL5k/ZswYCQ8Pl+eff148PDyKx7SeOXOm/OlPfxJA\nIiMj5ezZsyIi8t133wkgy5cvl3fffVfA+rCURo3CZNEii9x2m/Waej8/66WVeiSvHKmgQGTGDOsd\nodaDC4tER/8gwcG/u+Q+l5YtW4rFYpGCggK7DODlLFrc7eDs2QKBPOnYcU65lvvwww8FXhA3N0vx\n7f8XH8oByIgRIxyQ9uq2bNki/fv3l8OHD8tzzz0n7u7ucuTIkeL5bdu2lWHDhskPP/wggNxwww0C\nyIkTJ+T06dNy1113yfr1v93UdfLkyeKTwvfdd5+EhITIlClTLulu2rlT5Pe/l6J++QsyZsyHl2xT\nqcoqKBD5739FYmKsP2fNm2cL3CtRUe3Fz8+vePjdZ555Rtq3b1+ppx9VJ1rc7SQoaLv4+pZv/OUJ\nEyaIm9s66dbtt5OXFouluLg78+PgxZEx7733Xhk5cqR88cUX4ubmJn/961/l2LFjxRnLGva4Xbt2\nMmDAAImPj5d+/frJvn37iq+vL+nTT/cL7BLIlyFDvnLgnqnaZNmy347UY2OtlzW+995EAYofW7lk\nyRJnx3QILe520qPHUoE8yciwvY+8T59bBArl73+/dPrkyZPlo48+snPC8rvY51jya+7cuSIixa9f\ne+21a67j0UcfFT8/P/Hx8Sl+uHjz5s1LfSp54IEHxNu7gQQGrhYQ+fxzvRtKVVx+vsjjj1srV9Om\n1u6Yi48/GDNmjDRu3FhOnTolP/300xWvDHMFWtzt5Mkn5wuIfPqp7SdV69W7X8A6Zkt1lJ+fL9u2\nbZMDBw7IyJEj5dNPPy3+RUhISBCg+NLJq/nmm2+K/xAsXLhQRETuv/9+qVu3bvGDDA4cOCBeXl4y\nfvx4+fjjqQILxMOjUJKSHLt/yrVYLBZ54YUXZOrUmdK/v/US3EceEXnppbdkyJAhxQ/TaN68+SXP\nO3VVWtztZN681QL5cscdu2xqb71S5v/E0zNPauBVVnL8+HHZtm1bme3S09MFkJCQkOJiPmvWLAHk\nhRdekNmzZ8uAAQPE29tb0tLSitoHS2DgWYmO1rFrlO2mTZsm4CYwWyBfIiJekGnTphWfKJ0wYULx\nz+Nbb73l7LgOp8XdTqw/NGulefNDZTcW6zMTIVnat093cDLnGzhwoDz33HPFry0WS/F18BevkX/n\nnXeK57dt21bi4/8iIPLSS9VkoA5V7bVq1UoaN35DQKRjx0nStGlTAaROnToycOBACQ0Nlblz5wpY\nHzTv6rS424nFYhEvr/fE3T3PpqPN9977QqBQJkw45fhw1VBGRoZ88803smLFilJdO+PHj5c6depI\nWNhGgdOye/cxJ6VUNYX1JH+oeHnlytCh1mvZDx8+LO+//74cOHBAPv30UwHknnvuKb7Cy9XZWtzd\nUNdkjCE8/ACFhZ4kJ5fdftGiXMCNYcPqOjxbdRQcHMytt95Kz549CQ4OvmTeoEGDOHfuHEePPggE\ncffdy5wTUtUYSUlJwATy8z156y3rA2YaN27Mww8/TEREBH369AHgP//5D02aNCE0NNS5gasRLe42\n6NgxC7A+RaksmzcHAvn07Onu2FA10IgRI/jiiy94+OEuhIamsHp1F6zD2ih1ZWvWbAAeYtiwQmJi\nSs9v0aIFTZo0ASDmSg1qMS3uNujevQWwm8WL88pse/hwC0JD9+Pn5/hcNY0xhrvuuov333+fESOO\nYrFEMmNGtXkio6qG5swRIITx46/8RFBjDFOnTiU2NpaRI0dWbbhqTou7DTp27AisYNUq6yChV7Nn\nz1Fyc9vTrt3pKstWU913XwiQyeTJZf/BVLXXjh0d8PY+x/XXX73Nddddx7Zt23jggQeqLlgNoMXd\nBh06dABWcvasF7t2Xb3dlCk7AC9uvTWkqqLVWPHxbXFz+4G1axuRn+/sNKo6Onv2PFlZicTG7sXT\n09lpah4t7jaoX78+DRrsBsA6gu6V/fJLDpDPXXe1qKJkNZenpydRUdvIzfVjyRJnp1HV0cyZaUAQ\nffrkOjtKjaTF3UbjxiUCJ/nnP1dw4cKFUvNFhNTUMEJC9hIYqP+ttujfvxC4wHffWZwdRVVDP/9s\nLeq33BLk5CQ1k1YhG7300t+JiUln795mvPHGm6Xmb9q0l7y89nTpUrrwqyvr3r0jsJhfftF+GVVa\ncrI/sIWuXZs7O0qNpMXdRsYYHn64LRDBnDmppeZ/8skuwINRoxpWebaaqkuXLsBCDh705tAhZ6dR\n1YkIpKWFUbduCl5eXs6OUyNpcS+Hfv2s3zdurEf+ZWcB588vAHK5/fYmVR+shoqKiiIgYD0Aixfb\nd91nz55l79699l2pqjJ79kB+fgCtW591dpQaS4t7OURHQ0hIFvn517Np06ZL5h040JrQ0D34+Rkn\npat53Nzc6NjRHQ+PMyxaZL/1Ll26lKCgINq2bcuZM2fst2JVZRYssBb1Pn18nJyk5tLiXg7GwLBh\nFmAgS5asK56+ZUs2eXnRxMenOS9cDRUX1xZjlrBokVzzHgJbiQjPP/88AHl5eSxfvrzyK1VVbs6c\nU0AWN9yg/e0VpcW9nO68MwDwY+bMrOJpU6acAmDECL3qo7zi4uLIz5/Lr78a9u2r/PqWLFnC6tWr\nefvtt/H29mbp0qWVX6mqclu3CrCNhIROzo5SY2lxL6fERPD1zSA5OR4RobAQvvoqEFjFdde1dHa8\nGic2Nhaw9snYo2vm1VdfpWHDhjz00EN069aNxfbuzFdV4tChYOrVO4K/v7+zo9RYWtzLydMT+vff\nRW7udXz11S6+/RbS0wPx8Hifli21uJeXtbjvIiDgApXpQSksLOT+++9nwYIFPPnkk/j4+NCjRw+2\nbNlCXp4OcVCTHDyYRX5+Pdq31/NXlaHFvQL+8pf6wDFGjfJh7NhsvL13kZBwGA+PKw9upK6ufv36\nhIaGEhqaatOom1fz2muv8cknn/D000/zpz/9CbD+4SgoKGDXtcaMUNXO11+nADBgQJiTk9RsWtwr\noGvXlnz3nRAcbMjN3URu7gDGjLnD2bFqrKioKNzd17JvHxw7Vv7l8/Pz+cc//sEtt9zC66+/jmfR\nQCTWTwWwfft2e8ZVDrZw4XEARo5s6+QkNZsW9wq6+eYwTp1qSosWdwJpOtxoJbRq1Ypz534BbBsz\n/3KpqalkZ2dzyy23YMxvH+Wjo6Nxc3Njx44d9oqqqsDu3W64uWXSqlWAs6PUaFrcK8EYQ3JyMsnJ\nyfoEmEpo1aoV6elz8PERVqwo//IbN24ELg7N/BsfHx9atmypR+41zLFjdQgMPI7RLvdK0eJeSUFB\nQXTqpJdrVUZUVBSQT1xcVoWO3Ddt2oS3tw9JSdEkJkKvXvD++1BQYB2ued26ddYHBqtqLy8vjwsX\nGtOkiY7RVFk2FXdjzGBjTKoxZo8x5rkrzK9rjPnRGLPZGLPdGHOP/aMqV9WqVSsAIiIOsXEjXGHQ\nzWvatGkT9epN5O673Tl1CrKz4dFHYdAg6N9/BGlpaWzYsMEByZW9bd26G4igTRt9TGVllVncjTHu\nwAfAEKAt8HtjzOVnOh4GdohIB6Av8LYxRkf7UTaxHrlDYOBWCgpg3boyFihBRFi3LpSjR+/h/vth\nyxZIToZPP4WlS+HLL+/Azc2Lb775xkHp1ZWIwJtvwqhR5Xs/ly5NA9yIjw90WLbawpYj9wRgj4js\nE5E8YAYw4rI2AtQx1rNZAUAGUGDXpMplBQUFERYWRlbWQqB8J1X37k3j/Pk3CAs7xfvvg1vRT/Q9\n98DHH8Py5Z60bPk5EydO5MiRIw5Ir67kww/hmWdg1iwYMaKQ06cLbVpu48bzAPTqpaOrVpYtxb0J\nUHLQlENF00p6H4gBjgBbgcdERO/FVzaLjY1lz571xMaWr7i//fZpIIJnnjmOt/el8+6+2/q1Z88d\nZGd35y9/+YsdE6uryc+Hl16CHj1y6dnzrxw7Zujd+3ubznvs2mVt066dDhhWWfY6oToI2AQ0BjoC\n7xtjSn2uMsaMM8YkGWOSTpw4YadNK1cQFxdHSkoKPXsKq1ZBoQ0HeiLw9ddNgDXcd1+zK7Z5/31o\n3drg4zONadO+Iz093b7BVSk//WS9XyElZRwrVrxFSMhWtm9vz/TpM8pcNi3NDw+P0wTpw5cqzZbi\nfhhoWuJ1eNG0ku4BvhOrPcB+IPryFYnIZBGJF5H4+vXrVzSzckGxsbFkZWXRuvUJMjPBlqsX16yB\nU6dCadjwB+rUufI10f7+8NFHcO5cKPn5T/Cf//zHzsnV5aZPB2/vDCyWn1m/fj2vvtoeaM3UqWW/\nqRkZIQQH6x9ge7CluK8HoowxzYtOko4CZl3W5legP4AxpiHQBrDDGH+qtoiLiwPA3X01cO2umaNH\nj3Lw4EEmTboAXOCWW67dA5iYCGPGgDHP8dVXG+0VWV1BYSEsXCiIzOH2228mLi6OW24xgIUVKwKv\n2TWTl5dHbm4zGjfWyyDtocziLiIFwCPAXCAF+EpEthtjxhtjxhc1exnoYYzZCiwEnhWRk44KrVxP\nfHw8UVFRPP74TQQGnmf58isXgcmTJ9O4cWMiI6OZOjUf+JYnnxxX5vrffBM8PYUNG0Zy+vRpO6dX\nF23cCBkZhry8n+nRowcAoaEQHn6a8+c7kZKSctVlU1LSgMa0aqX3JNiDTX3uIjJbRFqLSEsR+UfR\ntA9F5MOifx8RkYEi0k5E4kRkmiNDK9fj4eHBK6+8AkBm5i8sWJBdqs2BAweYMGEC/fv3Z+zYb4Eg\n7r3Xw6bROMPCYOzYdOAWJk3Sa94dZeHCi/9aQPfu3Yun9+1rgB6sXp181WWXL7cOLNS+va/jAtYi\neoeqqjZGjhxJSkoKgYFbOXHCj8OXndn5+eefyc3NZdKkSaSnD6VpU5gyZbTN63/jjTAgnUmTIuzy\n1CdV2tq1EBiYTlBQHq1bty6efuONdQE/Fi3KvOqy69dbH4nYp08jR8esFbS4q2olOjqa0aMjAPj5\n50uff7po0SIiIiLw82vFvHkwduxv17XbIjjYg4iIqRw61KrEEaayp3XrwNt7M7GxsbiVeHN697be\ncbpx49Xvbdy+3frQ+YSEeo4NWUtocVfVzqhRMcAFZs3KKJ5WWFjI4sWL6devH//9r8FisRb38ho8\n+FeMOcKrr+qhu70dOQKHD0N29lJiYmIumRcWBr6+pzh4sMFVl//1Vx+8vNLx83N00tpBi7uqdrp2\n7YQx61i//re7klJSUjh9+jSJiX35/HPo0QNKfOq3WUJCB0TeZvFiw5o1dgytWL/e+v38+YVER5e6\nEprIyJNkZbXl/PnzpeaJCBkZIYSGZpSapypGi7uqdnx8fGjSZBfp6U2K+90vDuvr7d2THTvgD3+o\n2Lq7dOkCfIS/fy7//Kd98iqrTZvAGAG2XLG4d+5sAaJYsaL09e7p6ekUFrYgIiLf8UFrCS3uqloa\nMMDa3z5jhvX5pxs3bsTX15c5c5rj7w93VPDBV7GxsQQHexEVNYdZs2DbNnslVlu3QoMGmUBWqW4Z\ngBtusHbJzJyZVmremjW7gRDatfMuNU9VjBZ3VS2NHNke2MKUKdmIwIYNG4iJ6c3//ufGXXdB3boV\nW6+bmxu9e/fmzJmX8feH11+3a+xabcsWCAz8FW9vbyIiIkrNHzQoBICVK0s/sHz5cuuj9bp3D3Fs\nyFpEi7vnvhC8AAAcKklEQVSqlnr37o2b22RSU+uyaJGFTZs2kZ//HHl5MGFC5dadmJjIgQPJ3H77\naaZPh0OH7JO5NsvKgj17ALbSunVr3N1Lj8ceHAwBAUfZuze41LwNG3IB6N1br5SxFy3uqloKCAig\nR49UPDzSeeCBHM6evZUdOxJ54AG4wif+crntttuoW7cuq1aNxmIRJk+2T+babMcO60BumZkrr9gl\nc1Hr1plkZ8exf//+S6bv2eOFm1sOkZH6bD170eKuqq0xY26loGAUBw96AJ/QoYOF116r/HqbNWvG\n9OnTOXNmAyKz+fDDQvJK9xSocrg4qkB6+uIrnky9aNCgekBTpk9ffMn09PQQgoKOc4UDflVBWtxV\ntTVmzBgCAtZjscTQrt2TLF/uYbehYIcMGcLChQuB9zlxwp3vvrPPemurnTvBw0MQ2X3N4j50qHU0\n2K+/PlA87ezZs+TmtiQiovQlkqritLiraqtOnTpMmTKFO+/swTvv3GD3m1tiY2Np0+Ygvr6H+eAD\n+667ttm5Exo0OAcUXLNb5ne/A2MsbNniTV7Rx6VNm/YBzWh7+cM7VaVocVfV2qhRo5g6dSr9+vWz\n+7qNMQwdOpj8/HdZsQI2b7b7JmqN1FTw9T2Ij4/PNYu7vz80bZqJxdKJ7UWD9s+fbx1vpk8fvTXV\nnrS4q1rtd7/7HQUFk/HxsTBxorPT1EwFBbB7N5w/n0y3bt3wvvx5h5dJSDBAF5KTrTemrV0LYGHY\nsKsPTaDKT4u7qtU6duwInKFLl4PMmAHZpUcaVmU4cADy8iA9fTl9+vQps/2wYXWAUObPtz5qc+fO\nIDw8dtO4sb9jg9YyWtxVrRYdHY2XlxeNG88nMxNmznR2oppn507rd5EdJCYmltl+yBA3wMLq1cEU\nFMCxYxGEhOxxbMhaSIu7qtU8PT2Ji4tj//7/4Oubzj//efnjgVVZUlOt34OC0undu3eZ7Rs0gIYN\nf+Xw4Xh++qmQgoIgYmP3Ojhl7aPFXdV6Y8eOZd26NWRnf8TWrY1KPSREXdv27QUYc4Jbb70OT09P\nm5YZNOgoFksnbr7ZHTjO9dfrgGH2psVd1XoPP/wwvXr1Ar4A3HnrraPOjlSjJCdfQCSF4cOH27zM\nhAl+wEaMKQRepFu3zg7LV1uZaz2N3JHi4+MlKSnJKdtW6nIiwsmTJ2nQYBcNGrTh2LFQjN4Jb5M6\ndbI4f34ahw8Po3HjxjYtk5+fj79/GPn5AmRw9uxZAgMDHRvURRhjkkUkvqx2euSuFNZr3uvXr0/d\nujNJTw9Fjztsk5EB58/7ERBwmLCwMJuX8/T0pEOH5kAGjRo10sLuAFrclSohLi4FY/KZMcPZSWqG\ni2PKtG5twZTzo85bb71FUFAQd1R0cH51TVrclSohJiYMT8/FfPONdZRDdW1btxYAkJBQp9zLJiYm\nkpGRwb/+9S97x1JocVfqEq1btyYv77/8+iusW+fsNNXf+vWZQDZdu9reJVOSMQY3Ny1DjqD/q0qV\nEBUVBfyAh4eFr74Szpw54+xIVWbVqlW89dZb5OTk2LyM9cg9lVatmjsumKoQLe5KlWAdrvYs0dFp\nTJqUTsOGjUhPT3d2LIdbtGgRvXr14umnnyYxMZELFy7YtNy+fd7ATlq0aOHYgKrctLgrVUKbNm1o\n2bIl27a9SHZ2Q/Ly2jF16lRnx3K4n376CW9vb6ZOnUpSUhITbHiWYXY2nDpVBw+P3TRq1KgKUqry\n0OKuVAnGGG677TbgB4wpICzsMT777DNnx3K4tWvX0rlzZ+68807uvvtuvv32WywWyzWX2bULwI2G\nDU9rv3k1pO+IUpf54x//SJ8+7ejVK5ecnGFs27aNs2fPOjuWw+Tn57NhwwYSEhIA6NOnD2fPniXl\n4nWOV3FxdvPmuY6OqCpAi7tSl4mIiGDp0qWMHevP6dNBQDu2bdvm7FgOs23bNnJycoqLe48ePQDr\nCdZrSUkRoJC4OC9HR1QVoMVdqau44YaL/xrGli1bnBnFodYVXfPZtWtXAFq1akVoaGiZxX3jxmxg\nP3FxrRwdUVWAFnelriIsDOLjBXf3m1y6uK9du5bQ0FD27o0kNhaaNjU0bfoCS5YsueZy27YVAjuJ\njY2tkpyqfLS4K3UNw4YZCgvjSU7+1dlRHGbdunVER49kxAhDYSFERMCmTQ9z4EBb9u3bd8VlcnLg\n4EE/YLMW92rKpuJujBlsjEk1xuwxxjx3lTZ9jTGbjDHbjTFL7RtTKecYNgzAja1bm5Z59UhNdO7c\nObZv30Fa2tP4+cGSJbBgAbRunQtM5Ntvf7rictu2gcXiTp06e6hfv36VZla2KbO4G2PcgQ+AIUBb\n4PfGmLaXtQkCJgLDRSQWuN0BWZWqcr/7HdSte4GcnP4cOHDA2XHsbs6cOcD1HDwYyd//Do0aga8v\nvPeeNxDBM8/sok+fPuzfv/+S5TZssH6PibH9blZVtWw5ck8A9ojIPhHJA2YAIy5rMxr4TkR+BRAR\n17+lT9UKbm6QmHgBGERy8lZnx7G7L774Am/vZ2jUSLj//t+mDxhg6NatgJCQV9i4cQvjx4+n5LMf\nVq7MBs7Qt29E1YdWNrGluDcB0kq8PlQ0raTWQLAxZokxJtkYM9ZeAZVyttGjA4FAfv4509lR7Co9\nPZ3Zs1PIzb2OBx80eJW4otEYeOIJD06dCuKOOz5n3rx5LF36W2/rokXngWTuvHNM1QdXNrHXCVUP\noDNwAzAI+KsxpvXljYwx44wxScaYpBMnTthp00o51rBhPhiTy4oVwc6OYldffvklFssdgBv33lt6\n/k03QXg47N8/DH9/f6ZPnw7AyZPCoUMhNGiQQrt27ao2tLKZLcX9MNC0xOvwomklHQLmisgFETkJ\nLAM6XL4iEZksIvEiEq8nYVRN4e8PYWGpHDjQ2qVOqn7++Rf4+NxH797QrFnp+Z6e8NBDsGiRO717\nP8j3339PQUEBkyenAm7ceeflH+BVdWJLcV8PRBljmhtjvIBRwKzL2vwA9DLGeBhj/ICuwLXvXVaq\nBklMzKewsDXz5u1wdhS7OHToEJs25ZCT04JRo67e7oEHwNsbLJaHOHHiBMuWLeOzzw4BF/jLXwZW\nWV5VfmUWdxEpAB4B5mIt2F+JyHZjzHhjzPiiNinAHGALsA74WERc935tVevcd5/10PbTT13jevf5\n8+cDNwIwfPjV24WGwujRsGJFC3x9G/PJJ1+ye3cnWrRIITjYv2rCqgrxsKWRiMwGZl827cPLXr8J\nvGm/aEpVH9ddVx9391MsX+7j7Ch2MW/ePDw9HycuTggPv/azTx97DD77zNC8+Wd8+eX3QD0efdQ1\n/si5Mr1DVSkbuLlBq1YHOHYsluzsmn1tt8ViYd68DRQUxHPjjWU/1LpDB3jqKdi373pgIgEBG3jk\nET2RWt1pcVfKRoMHewANmTZts7OjVMrmzZvJyEhAxI0bb7RtmX/8Az74AB5+GI4c6YSHh7tjQ6pK\n0+KulI3Gj28JwPTpJ52cpHLmzZsH3EjDhoV06mTbMhevnHn/fahTx6HxlJ1ocVfKRtHRAfj6HiA5\nuWZf7z5v3kLc3AYzbJg7+gAl16VvrVLl0K5dOpmZHTly5JSzo1RIVlYWy5dnYbEEcv31zk6jHEmL\nu1LlcMstdQA/Jk+umVf6Ll++nPz83gD06+fkMMqhtLgrVQ73398KKOTHH887O0qFzJ8/H2MGEhdn\nQW8Sd21a3JUqh5AQT+rW3c+OHaHOjlIhK1duAHowcKD+6rs6fYeVKqdOnc6Rk9OBlJT9ZTeuRgoL\nC9m0yQ8Rb/r3d3Ya5Wha3JUqpxEjggAfpk/f6+wo5bJr1y5ycnri7m6hTx9np1GOpsVdqXK6/fZG\ngIUVK2rWjTwbNmwAEmnXLpuAAGenUY6mxV2pcmrc2BdPz1RSUmpWv/uqVZuALgwc6OvsKKoKaHFX\nqgLCwlJJT29Ffr6zk9hu7tzTgCeJifprXxvou6xUBbRtexKLxZf166XsxtVARkYGe/eGY4zQo4ez\n06iqoMVdqQrobb0PiNmza8b17tbnn/aiVasLBAU5O42qClrclaqA+PhmwA7mz89zdhSbzJz5M9Cd\n66/X/vbaQou7UhUQFxcHLGXz5joUFDg7zbXl5ubyzTf7AH/69q1ZV/ioitPirlQFhIWF4e+/kdxc\nL7ZudXaaa1u4cCFZWdaxfS92JynXp8VdqQowxhAXlwnAqlVODlOG5ORkoA8tW1po1MjZaVRV0eKu\nVAV17hyKMUdZvbp6XzGzefNW3Nz66CWQtYy+20pVUHx8Z0RWsnhxrrOjXFNycg4WS5B2ydQyWtyV\nqqDbb78dX99NHDniw9Gjzk5zZVlZWRw8GA5Ar15ODqOqlBZ3pSooICCAkSOthfO776pndU9JSUGk\nC4GBObRs6ew0qippcVeqEl58cQSQy5Qp250d5Yr2798PdKNjx1yMcXYaVZW0uCtVCZGRYdSvf5Bt\n2wIoLCx0dpxSdu48DsTQu7eXs6OoKqbFXalKSkgopLDwd6xbt9nZUUpJSrL+ivft6+PkJKqqaXFX\nqpJuvrkh4M2XX+50dpRSUlODAAsJCdonU9tocVeqkm64oR4ACxdmOzlJaUeONCMg4FcCA52dRFU1\nLe5KVVKjRhAQkM6+fdXr9k8ROHeuLeHhh5wdRTmBFnel7KBNm1Pk5nbi6NFjzo5SbPPmLESCiYnJ\ndHYU5QRa3JWyg1693IEw5s1LdXaUYrNnZwDQrZuTgyin0OKulB3ccou1S2bOnLNOTvKbJUtygEwG\nDGji7CjKCbS4K2UHPXoEYsx5kpK8nR2l2Nat/sA6YmJaOzuKcgIt7krZgYcH1Ku3l0OHmjo7CgBZ\nWXD8eAOCg3fh66tPX6qNbCruxpjBxphUY8weY8xz12jXxRhTYIy5zX4RlaoZoqJOkpPTmsxM5z+a\nKTkZRNyJijrl7CjKScos7sYYd+ADYAjQFvi9MabtVdq9Dsyzd0ilaoKuXQE8+PHHI86OwooV1j8w\n1kyqNrLlyD0B2CMi+0QkD5gBjLhCu0eBb4F0O+ZTqsYYOtR6M9O8ec6/9PCXX04DexkwoIOzoygn\nsaW4NwHSSrw+VDStmDGmCXAzMOlaKzLGjDPGJBljkk6cOFHerEpVaz16RAF7WL/ew6k5RGDjRm9g\nLYmJiU7NopzHXidU3wGeFRHLtRqJyGQRiReR+Pr169tp00pVDwEBAfj7b2X//oaIE5+8d+gQnD8f\nSETEUerWreu8IMqpbCnuh4GSlwCEF00rKR6YYYw5ANwGTDTG3GSXhErVIBERR8nJCSYtrey2jrJ8\neR4AiYnV57JMVfVsKe7rgShjTHNjjBcwCphVsoGINBeRSBGJBL4BHhKRmXZPq1Q116mT9Xmqq1Zd\n80OsQ1nvTM1hyJDGTsugnK/M4i4iBcAjwFwgBfhKRLYbY8YbY8Y7OqBSNUmfPsFANnPnOu9O1TVr\nBEime/fOTsugnM+mMz8iMhuYfdm0D6/S9u7Kx1KqZmrfPgZIYvXqWKdsPy8PDhwIwcfne5o16+GU\nDKp60DtUlbKjtm3bAqvZsyeQ3Nyq3/6WLVBY6EWbNmcw+tDUWk2Lu1J2VKdOHYKCUiks9GDjxqrf\n/rJl1pOpffroM1NrOy3uStlZdPQZANasqfptz5t3FjhCv35RVb9xVa1ocVfKztq1C8HN7ZBTivuG\nDZ7AGrp0ia/6jatqRYu7UnYWFRWFxbKSVasKq3S7J0/CiRNBBARsp3FjvQyyttPirpSdtW7dGlhD\nWpo7R49W3XbXrrV+b98+S0+mKi3uStnbxeIOVdvvvmhRFlDAoEGhVbdRVW1pcVfKzqKiovDx2Ymb\nW0EVF/cLwFYSE7W/XWlxV8ruPDw86NgxGn//1Cor7hYLpKTUwZi1dOnSpWo2qqo1Le5KOUCnTp3I\nyVnK+vVCfr7jt7dzJ+Tm+tCkySH8/Pwcv0FV7WlxV8oBOnfuTH7+MrKzDVu3On57Fz8hdOyY4/iN\nqRpBi7tSDnDdddfh6ZkMwOrVjh/cfenSPCCDbt1CHL4tVTNocVfKAZo3b85LL90PHOWf/1xCYaFj\nr3lfvrwAWEP79nEO3Y6qObS4K+Ugzz77DDExmRw+HM7MmY57vMHZs3DggC+wmthY54xGqaofLe5K\nOYgxhrvuigKiePXVyQ7bztq1IGLw9t5IZGSkw7ajahYt7ko5UM+e1l+xDRs8SHPQs/dWrwawEBt7\nATc3/ZVWVvqToJQDde4M7u4CdGfu3LkO2caqVeDhsZMOHZo7ZP2qZtLirpQD+ftD+/bg7Z3IL7/8\nYvf1WyywZo2FgoLlxMXpyVT1Gy3uSjlYt24GiyWeZctWImLfyyJTUiAz0w1YpcVdXUKLu1IO1r07\n5Of7cvJkqN373a397QCradeunV3XrWo2Le5KOVi3bhf/1Z3169fbdd2rV4On51liYjwJCwuz67pV\nzabFXSkHa9UK6tUTjOlBUlKSXde9cqWFgoKVDBky2K7rVTWfFnelHMwYa7+7t3cfkpOT7bbe06ch\nNdUNkZUMGjTIbutVrkGLu1JVoHt3yMlpybZth+y2zlWrrN/d3NbSs2dPu61XuQYt7kpVgYv97keP\nhnP27Fm7rHPZMjAmn86d8/H397fLOpXr0OKuVBVISABjBOhGSkqKXda5ZIkFWEdiYoJd1qdcixZ3\npapAYCBEReUBPexS3C9cgORkEFlKjx49Kh9QuRwt7kpVkf79PYFebN26s9LrWrMGCgvdgGV685K6\nIi3uSlWRfv3cgADWrKn82O7W/nYLnp5JNG+uY8qo0jycHUCp2qJvX+v3lJRGlV7X0qUQGLiP8PBG\neHjor7EqTY/claoioaHQsOFxzpzpyIULFyq8nuxsa7eMMcuJjo62Y0LlSrS4K1WF4uMzgZ5s3bqr\nwutYsQJycyEz81tiYmLsF065FC3uSlWhwYO9AV9+/PFEhdcxdy54elqwWBbrYGHqqmwq7saYwcaY\nVGPMHmPMc1eYP8YYs8UYs9UYs8oY08H+UZWq+UaObAQUsnRpxfvJ582D5s2PAFl07drVbtmUaymz\nuBtj3IEPgCFAW+D3xpi2lzXbDySKSDvgZcBxD4xUqgZr0MALX9/NbN/erELLHz0KW7dCQMAqGjZs\nSLNmFVuPcn22HLknAHtEZJ+I5AEzgBElG4jIKhE5XfRyDRBu35hKuY4WLVI4c6YVR46Uf9n5863f\nMzKmk5CQgDHGvuGUy7CluDcBSj5h4FDRtKu5D7ji88SMMeOMMUnGmKQTJyre56hUTdazp/U46Icf\n8su97Ny5EBJSwIEDP5CYmGjvaMqF2PWEqjHmOqzF/dkrzReRySISLyLx9evXt+emlaoxEhNDgF+Z\nNi2Du+66i5dffhmLxVLmcnl5MHs2hIYm4+XlydixYx0fVtVYtpzVOQw0LfE6vGjaJYwx7YGPgSEi\ncso+8ZRyPbGxbYGfWLVqLKtWfQ3kMmTIEOLj46+53KJFcOYM5OX9H7fffjt6gKSuxZYj9/VAlDGm\nuTHGCxgFzCrZwBjTDPgOuEtEKn4Br1K1QNu2bQkP3wIE0LXr/wNg8eLFZS739dfg65tPVtb3etSu\nylRmcReRAuARYC6QAnwlItuNMeONMeOLmv0NCAEmGmM2GWPs+ywxpVyIp6cnK1f+HW/vLMLCniY6\nOrrM4p6fDzNnQmjoaho2DKZfv35VlFbVVDZdbCsis4HZl037sMS/7wfut280pVxXs2YNufde+M9/\nYPToIXz11RQKCgquOk7MkiWQkQGZme/y0EN36Hgyqkx6h6pSTvKHP0BODnh6jub8+fPXfL7ql1+C\nj08+BQU/M3r06CpMqWoqLe5KOUlCArRpA5s3tweu3u9++jT8738QEjKXli2bkJCgT15SZdPirpST\nGGM9el+zxouWLW+8anH/7DPrSJBHjrzI6NGj9cYlZRMt7ko50QMPgJ8f+Pq+wIIFC1ixYsUl83Nz\n4a23oEWLQ4gka5eMspkWd6WcKDQUHnwQUlI60bjxQG699VY2bdpUPH/iRDhyBDIzn6FXr146fruy\nmRERp2w4Pj5ekpL0ikmlTp2C1q0hIiKbEydiOXcugzfffJN9+yy8/voYRJYBN7Bw4UK9BFJhjEkW\nkWvf8YY+Zk8ppwsJgffegzFjfLn11i2kpPRh3LiXgdm4uwt//GMqiYlfc9111zk7qqpBtLgrVQ2M\nHg27d8OLLwYQGJiMh4fg7Q0zZxoGDHjc2fFUDaTFXalq4oUXoFcv+OYbQ926hnHjoEULZ6dSNZUW\nd6Wqkf79rV9KVZZeLaOUUi5Ii7tSSrkgLe5KKeWCtLgrpZQL0uKulFIuSIu7Ukq5IC3uSinlgrS4\nK6WUC3LawGHGmBPAQads/OpCgZPODuEEut+1T23dd1fY7wgRqV9WI6cV9+rIGJNky2hrrkb3u/ap\nrftem/Zbu2WUUsoFaXFXSikXpMX9UpOdHcBJdL9rn9q677Vmv7XPXSmlXJAeuSullAty+eJujDlg\njNlqjNlkjEkqmna7MWa7McZijIkv0dbTGPN5UfsUY8zzJeb9wxiTZow574z9KC977Lcxxs8Y87Mx\nZmfRcq85a39sZcf3e44xZnPRch8aY9ydsT/lYa99L9FmljFmW1XuQ0XY8T1fYoxJLVrPJmNMA2fs\nj92IiEt/AQeA0MumxQBtgCVAfInpo4EZRf/2K1o2suh1NyAMOO/sfaqq/S7693VF072A5cAQZ+9b\nFb3fgUXfDfAtMMrZ+1ZV+1407RbgS2Cbs/erCt/zS9rW9K9a+SQmEUkBMMaUmgX4G2M8AF8gD8gs\nWmbNVZapMcq73yKSBSwuWjbPGLMBCK+6xPZRwfc7s6iNB9Y/bDXy5FRF9t0YEwA8AYwDvqqysHZU\nkf12NS7fLYP1zVxgjEk2xowro+03wAXgKPAr8JaIZDg6oIPYdb+NMUHAjcBCR4S1I7vttzFmLpAO\nnCtqW93Za99fBt4GshyW1L7s+bP+eVGXzF9NTT6So3Y8Q7WXiBwu6j+bb4zZKSLLrtI2ASgEGgPB\nwHJjzAIR2VdVYe3IbvtddJQzHXi3Bvxf2G2/RWSQMcYH+C/QD5hfBfkro9L7DgQCLUXkcWNMZFWE\ntgN7vedjitZTB2tX3F3AF1WxA47g8kfuInK46Hs68D3WN/dqRgNzRCS/qP1KoEbeqmzn/Z4M7BaR\ndxyV117s/X6LSA7wAzDCMYntx0773h2IN8YcAFYArY0xSxyZu7Ls9Z6XWM85rOcbrrWeas+li7sx\nxr/orzDGGH9gIHCts/+/Yj1Cu9i+G7DT0TntzZ77bYx5BagL/MmRme3BXvttjAkwxoQVTfcAbqCa\n/xzYa99FZJKINBaRSKAXsEtE+joye2XY8T33MMaEFk33BIaVsZ7qz9lndB35BbQANhd9bQf+XDT9\nZuAQkAscB+YWTQ8Avi5quwN4usS63ihaxlL0/UVn75+j9xvryVMBUoBNRV/3O3v/qmC/GwLrgS1Y\nf8HfAzycvX9V9bNeYp2RVPOrZez4nvsDyUXv+Xbg34C7s/evMl96h6pSSrkgl+6WUUqp2kqLu1JK\nuSAt7kop5YK0uCullAvS4q6UUi5Ii7tSSrkgLe5KKeWCtLgrpZQL+v9DhL6htLJ0pwAAAABJRU5E\nrkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#let's check the resulting fit\n", "plt.plot(spec[0], spec[1], 'k-')\n", "plt.plot(spec[0], result.best_fit, 'b-')\n", "\n", "#and extract the fitted lambdas\n", "for p in result.params:\n", " if 'center' in p:\n", " print('%s %.2f'%(p.split('_')[0], result.params[p].value))" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 2", "language": "python", "name": "python2" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.13" } }, "nbformat": 4, "nbformat_minor": 2 }