SimluateInverted.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import numpy as np\n",
    "from scipy.stats import logistic\n",
    "from scipy.optimize import basinhopping\n",
    "from scipy import interpolate\n",
    "from IPython.display import clear_output, Javascript, Markdown\n",
    "import matplotlib.colors as colors\n",
    "import matplotlib.pyplot as plt\n",
    "import scipy.io\n",
    "import gc\n",
    "import datetime\n",
    "import time\n",
    "import random\n",
    "import matplotlib\n",
    "import nbimporter\n",
    "from _create_pdf import create_spins, create_spins_inhomo_g\n",
    "from _helper_functions import tsallis, cauchy, find_nearest, odeintz, find_nearest_sorted_left\n",
    "from _steady_state import a_steady\n",
    "import _spin_params as spin_pars\n",
    "from _mbes import mbes, mbes_soc\n",
    "from _do_calculation import solve_mbes_wait, solve_mbes\n",
    "#matplotlib.pyplot.switch_backend('nbAgg')\n",
    "#%qtconsole"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4.174817177486928\n"
     ]
    }
   ],
   "source": [
    "gcoll = spin_pars.gcoll\n",
    "q = spin_pars.q\n",
    "quali = spin_pars.quali\n",
    "fr = spin_pars.fr;\n",
    "kappa = np.pi*fr/quali #HWHM\n",
    "delta = spin_pars.delta #detuning cavity to spin central frequency ws-wc\n",
    "deltac = spin_pars.deltac #detuning drive to cavity wc-wp\n",
    "gperp = spin_pars.gperp #HWHM\n",
    "gpar = spin_pars.gpar#0.001 #HWHM\n",
    "gammaq = spin_pars.gammaq\n",
    "nbins = spin_pars.nbins#701#5000#20\n",
    "nbins = spin_pars.nbins\n",
    "g_number = spin_pars.g_number #how many bins for inhomogeneous coupling\n",
    "#g0 = gcoll/sqrt(nspins)\n",
    "print(kappa)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def pulse_rwa(t, args):\n",
    "    wgrund = args[0]\n",
    "    fcomps = args[1]\n",
    "    pulse = 0\n",
    "    for i in range(len(fcomps)):\n",
    "        pulse = pulse+(fcomps[i])*np.sin((i+1)*wgrund*t)\n",
    "    return pulse\n",
    "\n",
    "def pulse_rwa_gauss(t, args):\n",
    "    wgrund = args[0]\n",
    "    fcomps = args[1]\n",
    "    nfcomps = len(fcomps)\n",
    "    tend = np.pi/wgrund\n",
    "    \n",
    "    width = 2*tend/(nfcomps+5)\n",
    "    t0 = 3*width/2\n",
    "    pulse = 0\n",
    "    for i in range(len(fcomps)):\n",
    "        pulse += fcomps[i]*np.exp(-4*(t-t0)**2*np.log(2)/width**2)\n",
    "        t0 += 0.5*width\n",
    "    pulse *= logistic.cdf(t, loc=0.01, scale=.002)*(1-logistic.cdf(t, loc=tend-0.01, scale=.002))\n",
    "    pulse *= 0.46971947398573455\n",
    "    return pulse\n",
    "\n",
    "def pulse_rwa_gauss_const(t, args):\n",
    "    wgrund = args[0]\n",
    "    fcomps = args[1]\n",
    "    maxampl = args[2]\n",
    "    nfcomps = len(fcomps)\n",
    "    tend = np.pi/wgrund\n",
    "\n",
    "    #xint = np.linspace(0, tend, nfcomps)\n",
    "    \n",
    "    #spline = interpolate.BSpline(xint, fcomps, 2)\n",
    "    spline = 0\n",
    "    for i in range(len(fcomps)):\n",
    "        spline += (fcomps[i])*np.sin((i+1)*wgrund*t)\n",
    "\n",
    "    #spline = spline/nfcomps\n",
    "    pulsei = 0\n",
    "    pulseq = 0\n",
    "    pulse = 0\n",
    "    pulse = maxampl*logistic.cdf(t, loc=0.02, scale=0.004)*(1-logistic.cdf(t, loc=tend-0.02, scale=0.004))\n",
    "    \n",
    "    pulsei = pulse*np.cos(spline)\n",
    "    pulseq = pulse*np.sin(spline)\n",
    "\n",
    "    return pulsei, pulseq\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def g_inhomo_function(x, sigma):\n",
    "    return np.exp(-x**2/sigma**2)\n",
    "\n",
    "x_start = 5\n",
    "x_stop = 20 #µm\n",
    "x = np.linspace(x_start, x_stop, g_number)\n",
    "sigma = 5\n",
    "g_temp = g_inhomo_function(x, sigma)\n",
    "g_list = g_temp/np.sum(g_temp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEICAYAAABWJCMKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAjj0lEQVR4nO3deZRcd3nm8e/ba/W+qFvqRWpLsuQFZFs2QthgjCGQAWySwORMcEjADBnHJyeBLAOTTCAbhCSQMI5hnIyJPRCHw3JywIQlNvEubPBYtmVbNliWtba6tbTU+76880fdK8qllrpaqlu3ludzTp+uuvdW1auSVE/d33J/5u6IiEjpKYu7ABERiYcCQESkRCkARERKlAJARKREKQBEREpURdwFLEdbW5uvXbs27jJERArGk08+OeDu7YvtK6gAWLt2Ldu3b4+7DBGRgmFm+0+3T01AIiIlSgEgIlKiFAAiIiVKASAiUqIUACIiJUoBICJSoiILADO718yOmdnHF9mXMLOvmNm24HciqjpERGRxUZ4BfAj46Gn23Qj81N3fCLwY3BcpOJMz83z9iQP865O9TM/Nx12OyLJENhHM3XvN7HS73wR8Jrj9HeBjwD8udqCZ3QTcBNDT05PlKkXO3tTsPO/7px/z1IEhAL71dC9f+uBWKsvVsiqFIa5/qSuAweD2ENB6ugPd/XZ33+LuW9rbF53NLBKLLzywm6cODHHLr2zmL9+9iUd3H+eftu2NuyyRjMUVACeA5uB2U3BfpGAMjs/wTz/cwy9u7uKXLu/mfa87j7devJJ/eGg349NzcZcnkpG4AuBh4J3B7XcG90UKxteeOMjU7AK/de2Gk9t+680bGJma45tPH4qxMpHMRTkK6IskO4FvNLO7zWyzmYWdwl8CLjGzbcAlwX2RgvGdZ/q4vKeZCzsaTm67oqeFjSvr+c4zfTFWJpK5KDuB/9sim3cE+yaBG6J6bZEo7Tk2xgv9I3zi+ledsu/6S7u45f5dHBmZYlWjRjdLftNwBZFl+vedhwF45yUdp+y77tJO3OHe5w/nuiyRZVMAiCzTo7sHuLizkc6mmlP2bVhZz5rWGn740kAMlYksjwJAZBmmZufZvn+Q15+/4rTHvOH8Nn685zjzC57DykSWTwEgsgxP7R9kZm7hjAFw1fkrGJma4/m+4RxWJrJ8CgCRZfjxnuOUGWxdd9q5i1wVhMOPXj6eq7JEzooCQGQZdvQOc8GqBhoSlac9ZmVDgjWtNTzTO5S7wkTOggJAJEPuznO9Q1y6umnJYzevaeGZg2oCkvymABDJUO/gJIMTs1y6unnJYy9b3cShoUmOjk5FX5jIWVIAiGTo2d7kN/rLMgiAzWuSx+gsQPKZAkAkQ8/2DlFVXvaKyz+czqu7migzeE79AJLHFAAiGXq2d5iLOxuoqlj6v01NVTlr2+r46eHRHFQmcnYUACIZcHde6B/h1d1LdwCHLu5o5MUjCgDJXwoAkQwcHZ1meHKWC1ct3fwTurCjgf3HJ7Q+gOQtBYBIBnYF3+Q3rqrP+DFhX8EunQVInlIAiGTgxaAtfzlnABcFAfCi+gEkTykARDLw0pExVtRVsaK+OuPHrGmppbaqXB3BkrcUACIZePHI6LKafwDKyowLVjXw08MjEVUlcm4UACJLcHd2Hx3jgmU0/4QuWFXP7qPjEVQlcu4UACJL6BueYmx67qwCYF1bPQNj04xMzUZQmci5UQCILGH30TEgudrXcq1vrwNg7zGdBUj+UQCILGH/8eSH97q2umU/9vwwAAYUAJJ/FAAiS9g7ME5tVTkrGzIfARRa01pLmcGeY2MRVCZybhQAIkvYNzDOeSvqMLNlP7a6opw1rbW8rDMAyUMKAJEl7Ds+wbq22rN+/Lq2OvUBSF5SAIicwdz8AgdPTLB2xfLb/0Pr2+rZOzDOwoJnsTKRc6cAEDmDQ0OTzC34uQVAex2Ts/McHtHqYJJfFAAiZxCO3ll7FiOAQuHoof3HJ7JSk0i2KABEziD80F57Dn0Aa1qSjz14QgEg+UUBIHIGewfGqasqp30ZF4FL19mcoLzMOKAAkDyjABA5g33Hx1nbdnZDQEOV5WV0NSc4OKgAkPyiABA5g30D4+fU/h/qaa3VGYDknUgDwMxuNLPHzOxRM7sibd96M3vEzB4yswfNbHWUtYgs1/yCc2hokp7Ws2//D61pqVUfgOSdyALAzFqADwPXAr8G3Jp2yG8Bd7j7tcCXgd+JqhaRs3F0dIrZeWd1S805P9ea1loGxma0PrDklSjPALYC29x9xt33Ag1mltqT9jzQHNxuAY5GWIvIsvUOTgKwuuXczwDCs4jwOUXyQUWEz70CGEy5PwS0Av3B/fuAe83sQ0A1ycA4hZndBNwE0NPTE1WtIqfoDTpts3EGEAbAgRMTJxeLF4lblGcAJ/jZN3yApmBb6G+Aj7v7pcCfAZ9e7Enc/XZ33+LuW9rb2yMqVeRUvSeS39a7m7MbACL5IsoAeBy42swqzawHGHP36ZT9BgwEt4+SPDsQyRu9g5O0N1STqCw/5+dqrq2kvrpCHcGSVyJrAnL3QTO7DXgYcOAjZrYZeJu7fxb4FPB/zGwOqAR+M6paRM5G79BEVpp/AMyMNRoKKnkmyj4A3P1O4M60zTuCfc8DV0f5+iLnondwkktXN2ft+Xpaa3hZl4WWPKKJYCKLmF9w+oYms3YGAMnRRH1Dk7jrstCSHxQAIovI5hyAUFdzDRMz8wxNzGbtOUXOhQJAZBHZnAMQ6m5OAMk1BkTygQJAZBHZnAMQ6m5OhkmfAkDyhAJAZBHZnAMQ6tIZgOQZBYDIIrI5ByDUWldForJMZwCSNxQAIovI5hyAkJnR1VxD35DWBpb8oAAQWUTv4GRWO4BD3c019OoMQPKEAkAkzcKC0z80ldX2/1BXU42agCRvKABE0hwfn2FmfuFkp202dbfUcGx0mum5+aw/t8hyKQBE0vQPJ7+hdzRmPwC6grOKfvUDSB5QAIikCTtpu6JoAgrOKtQMJPlAASCS5nBwBtDZlP0zgNXBZDDNBZB8oAAQSdM/PEVVRRmtdVVZf+5VTdWYKQAkPygARNL0DU/R2ZTAzLL+3NUV5bTXV6sJSPKCAkAkTf/QZCTNP6HuFk0Gk/ygABBJ0z88RWdT9juAQ13NNWoCkrygABBJMb/gHBmZivYMoLlGC8NIXlAAiKQYGJtmbsHpjGAIaKijMcH03IIWhpHYKQBEUvQPJ9vmOyOYBBYKzy7C1xKJiwJAJEV/0DbfGcFlIEIdQQAcHlE/gMRLASCSoi/4Vt4VYSdw2MGsMwCJmwJAJMXh4UmqK8porq2M7DXaG6opLzMOKwAkZgoAkRR9w1N0NddEMgksVF5mrGyo1lwAiZ0CQCRF1JPAQp1NCfUBSOwUACIpDg9PneykjVJnU436ACR2CgCRwPyCc2R0OtIO4FBHU4LDw1OaDCaxUgCIBI6OTjG/4JEOAQ11NiWYmJlnZGou8tcSOR0FgEgg7JTNRR9Ax8nJYOoHkPgoAEQC4bDMKC8EF9JsYMkHCgCRQPhtPBd9AGHIaC6AxEkBIBLoG5qiprKcxpqKyF+rvaGaMtMZgMQr0gAwsxvN7DEze9TMrlhk//8ws/vM7CEze0uUtYgs5fDIJJ3N0awElq6yvIz2huqT6w+LxCGyrzpm1gJ8GLgS6AbuAq5O2f8OoMnd3xpVDSLL0Tc0lZPmn1CH5gJIzKI8A9gKbHP3GXffCzSYWXXK/v8CJMzsfjO7y8yaFnsSM7vJzLab2fZjx45FWK6Uuv7hyZxMAgt1NiYUABKrKANgBTCYcn8IaE253wUsuPvPAY8Df7TYk7j77e6+xd23tLe3R1WrlLjZ+QWOjk7TlcMACCeDicQlygA4ATSn3G8KtqXuvye4fQ9waYS1iJzR0dFp3Il0JbB0nU0JxqbnGJ3SymASjygD4HHgajOrNLMeYMzdp1P2PwRsCW5vAXZHWIvIGZ1cCCaXTUDNGgoq8YqsE9jdB83sNuBhwIGPmNlm4G3u/lngS8AXzexBYBZ4f1S1iCylP4eTwEKpk8E2rmrI2euKhCId8OzudwJ3pm3eEeybRh/6kifCSWC5uA5QqCNYd1hnABKXjJqAzOwjmWwTKVT9w1PUVZXTUB39JLDQqiAA+jQXQGKSaR/ABxbZdmMW6xCJVbgOQC4mgYWqKspoq6/WGYDE5oxfd8zsBuBXgXVm9m8puxp45YgekYIWLgWZa51Nmgsg8VnqfPcxoB9oA/4uZfso8GxURYnk2uHhSS5Ymft5Jh1NCQ4cn8j564rAEgHg7vuB/cBVuSlHJPfCSWC5HAIa6mpK8Pie4zl/XRHIvBP4PWb2kpkNm9mImY2a2UjUxYnkQhyTwEIdTTWMTM0xPq2VwST3Mu0E/gzwC+7e5O6N7t7g7o1RFiaSK+EVOXN5HaCQFoaROGUaAEfc/SeRViISk59NAst9AISho5FAEoelRgG9J7i53cy+DtwNnLycg7t/M7rSRHKjfyj3s4BDnVobWGK01Cigd6XcngB+PuW+AwoAKXj9w1PUVpXTmMjdJLDQKs0GlhgtNQrog7kqRCQuh0cmcz4JLJSoLKe1ror+EQWA5F5GX3nM7NZFNg8D293929ktSSS3cr0SWLpOrQsgMcm0EzgBbAZeCn4uBVYDHzKzWyKpTCRHwstAxKWzKUHfkPoAJPcybfS8FHiDu88DmNk/ANtIrvH7XES1iURubn6Bo6NTsYwACnU0JXhi3+DSB4pkWaZnAC1Afcr9OqA1CITpxR8ikv+Ojk6z4PGMAAp1NtUwPDnL5Mx8bDVIacr0DOAzwA4zewgw4Brg02ZWB9wXUW0ikYtzDkCoq/lnl4U+v71+iaNFsiejAHD3O8zs+8DWYNP/dPe+4PZHI6lMJAfCztd4+wCSZx99QwoAya0zNgGZ2UXB7yuATuBg8NMRbBMpaOEErDhHAXUH1yAKJ6SJ5MpSZwC/D9zEKy8FHXLgLVmvSCSH+oenqKksp7Em95PAQqsaE5hpZTDJvaUmgt0U/H5zbsoRya3+4Uk6Y5oEFgpXBtNQUMm1TC8HXWtmHzez24P7G83s+mhLE4le//BUTheCP52u5hpdEVRyLtNhoP8XmAFeH9w/BHwqkopEcujw8BQdjfG1/4e6NBlMYpBpAJzv7p8BZgHcfYLkcFCRgjU3v8CRkXgngYU6m2roG5rC3eMuRUpIpgEwY2Y1JDt+MbPz0QQwKXDHxoJJYHnRBJRgcnae4cnZuEuREpJpAPwpcA+wxsy+AtwPfCyyqkRyIB8mgYW6msO5AOoHkNzJdOzbB4DvAf8K7AE+4u4DkVUlkgPhuPt86AMIQ6hvaJJXdWm1VcmNTAPgDuCNwNuA84GnzewRd//7yCoTidjJSWB50AR0cjKY5gJIDmV6KYgHzewR4LXAm4GbgVcDCgApWIeHp0hUltFUUxl3KbTVV1NZbvRpKKjkUKYLwtxP8gqgPyJ5GejXuvvRKAsTiVr/8BSdTTWxTgILlZUZHRoKKjmWaSfwsyTnAWwiuTbApmBUkEjBOjQ0mRfNP6HOphpdD0hyKqMAcPffc/drgPcAx0lODBuKsC6RyPUNTcZ6Ebh0XU0JXQ9IcirTS0H8tpl9HXga+EXgTuAdGTzuRjN7zMwePd3VQ83sz81s93KKFjlX03PzHB2dprsljwKguYbDw1PML2gymORGpqOAEsDngCfdfS6TB5hZC/Bh4EqgG7iL5BKSqcesAi7IuFqRLAnXAQhH3+SDzuYa5hacgbFpVjXmT9OUFK9Mm4D+1t0fz/TDP7AV2ObuM+6+F2gws+q0Yz4B/NUynlMkKw4NJpta8ikAulLmAojkQqadwGdjBZC60vUQ0BreMbONQL27P3umJzGzm8xsu5ltP3bsWCSFSunpDT5k860JCDQbWHInygA4ATSn3G8KtoX+DPjkUk/i7re7+xZ339Le3p7VAqV09Q1NYhbvUpDpwg5pTQaTXIkyAB4HrjazSjPrAcbcPfUCcuuB/21m9wCdZnZrhLWIvMKhwUna66upriiPu5STGmsqqK0q1xmA5Exk6+C5+6CZ3QY8TPIqoh8xs83A29z9s+5+VXisme129w9HVYtIukNDk3nV/ANgZnQ116gPQHIm0oVQ3f1OkkNGU+1Y5LgNUdYhkq5vaJJN3U1xl3GKzqaEmoAkZ6JsAhLJSwsLTt/QVF6NAAp1N9foekCSMwoAKTkD49PMzC/kXRMQJEcCHRudZmp2Pu5SpAQoAKTk5OMcgNDqlnAoqJqBJHoKACk5h4bCdQDyLwDWtNYCcHBQASDRUwBIyenLw0lgofAM4OCJiZgrkVKgAJCSc2hwkoZEBY2J+BeCSbeqIUFVeRkHBxUAEj0FgJScQ0OTedn+D8mFYbpbaug9oSYgiZ4CQEpO7+Akq1tq4y7jtFa31OgMQHJCASAlxd05cGKCntb8DYA1rbXqA5CcUABISTk+PsPEzDw9rfnZBASwpqWWwYlZxqaXc/V1keVTAEhJORB8s16T12cAyXDqVTOQREwBICUlbFrJ5yagsH/ioDqCJWIKACkpB44nAyCfO4HXaC6A5IgCQErKwcEJVjZUU1OVP+sApGutq6K2qlwjgSRyCgApKfk+AgiS6wKsaalVE5BETgEgJeXgicm8DwBIdgSrE1iipgCQkjEzt0Df8GRejwAKrW6ppXdwEnePuxQpYgoAKRmHhiZxz+8hoKGe1lrGpuc4MT4TdylSxBQAUjIOFMAQ0NC6tjoA9h0fj7kSKWYKACkZP5sElr+zgENrgwDYO6B+AImOAkBKxt5j49RUlrOqIRF3KUta3VJDeZmxd2As7lKkiCkApGTsHRhjbVsdZWUWdylLqiwvo6e1ln06A5AIKQCkZOwdGGd9e13cZWRs7Ypa9g6oD0CiowCQkjAzt8DBwUnWtxVQALTVse/4uIaCSmQUAFISDpyYYH7BT46uKQTr2uqYmJnn6Oh03KVIkVIASEnYcyzZmbq+vT7mSjK37uRIIDUDSTQUAFISwg/RdSsK5wxgbVDrPgWAREQBICVh78A4K+qqaKqtjLuUjHU111BVXqYzAImMAkBKwp5j4wXV/g9QXmb0rKhljwJAIqIAkJKwp8CGgIY2tNfz8lFNBpNoKACk6I1MzTIwNs26tsLpAA5dsKqefcfHmZqdj7sUKUIKACl6u4Nv0OcX4BnABR0NLHiyCUsk2yINADO70cweM7NHzeyKtH0fM7PHg32fN7P8n58vBWnX4VEALupojLmS5btgVQMAu46MxlyJFKPIAsDMWoAPA9cCvwbcmnbIt9z9de7+BmAV8JaoapHS9uKRUWoqy1ndkv9XAU23dkUdFWWmAJBIRHkGsBXY5u4z7r4XaDCz6nCnu7+Ucuw0MLfYk5jZTWa23cy2Hzt2LMJypVjtOjLKBavqC+IicOmqKspY11bHriPqCJbsizIAVgCDKfeHgNb0g8zsTUAn8MhiT+Lut7v7Fnff0t7eHkWdUuRePDx2simlEF2wqoGXjuoMQLIvygA4ATSn3G8Ktp1kZpcCfw2813XFK4nAifEZBsamubCjsAPgwIkJJmc0EkiyK8oAeBy42swqzawHGHP3k1e1MrMNwJ0kP/wHIqxDSljYdl7YZwD1uMPLx9QMJNkVWQC4+yBwG/Aw8FXgd81ss5l9NDjkFpJnCF82s4fM7LqoapHS9eLhwg+AjUHt4Z9FJFsqonxyd7+T5Lf8VDuCfddH+doiAM/3DdNaV8WqxuqlD85Ta1fUkqgs4/m+Ef7za+KuRoqJJoJJUdt5aIRN3U0U8jSTivIyXtXZyM6+4bhLkSKjAJCiNTU7z64jo2zqKrwJYOk2dTfxQt8ICwsaKyHZowCQorXryChzC84l3U1xl3LONnU3MTY9x97juiSEZI8CQIrWc4eSTSabiiEAupJ/hp2H1Awk2aMAkKK189AITTWVBXkJiHQbV9VTVVGmAJCsUgBI0Xru0BCbuhsLugM4VFlexsUdDTzbqwCQ7FEASFEan57jJ/2jXL6mJe5Ssubynhae6R1idn4h7lKkSCgApCg9c3CI+QVny9riCYAta1uYml3ghb6RuEuRIqEAkKK0ff8gZnDFeUUUAOclr6W4ff/gEkeKZEYBIEXpiX0nuHBVA42JyrhLyZqOpgTdzTU8uf/E0geLZEABIEVnfsF5+sAQrymib/+hLWtb2L5vEF08V7JBASBF5yf9I4xNz/HatacsP1Hwtqxt5ejoNPuPT8RdihQBBYAUnUdeSq4c9/oNK2KuJPuu3tAGwLaXtDqenDsFgBSdR3Yd4+LORlY2JOIuJevWrqhlTWsND+/SEhpy7hQAUlTGp+d4cv8g12xsi7uUSJgZb9zYzo9eHtB8ADlnCgApKj/ec5zZeeeaC4p3/ehrNrYzPjPPUxoOKudIASBF5f6fHqWmsryoJoCle/2GFVSUGQ/89GjcpUiBUwBI0ZibX+DenYf5uYtXUl1RHnc5kWlMVHL1xja+91y/hoPKOVEASNF4fO8Jjo/PcP2lnXGXErnrLumkd3BSF4eTc6IAkKLx3Wf7qa0q59oLV8ZdSuR+/lUdVJYb33uuP+5SpIApAKQoTM3O8+87+3nrxatIVBZv80+oqbaSaza28+0dhzQaSM6aAkCKwnef7WdoYpYbtvbEXUrOvHdrD0dGprn/J0fiLkUKlAJAisJdP9rHhpX1XLm++C7/cDpvuWgl3c013PXj/XGXIgVKASAF78n9gzzTO8yvX3leUaz+lanyMuNXX9fDo7uP85N+rREgy6cAkIJ3y327aK2r4pdfszruUnLufa/roaG6gs/9x664S5ECpACQgvbY7gG2vTTAzW9aT111Rdzl5FxzbRU3XbOe/3jhCE8d0MxgWR4FgBSsmbkFPvHtnaxpreH9V62Nu5zYfPDqdbTVV/Pxb+3UiCBZFgWAFKzP3vtTXj42zl/8wqaSGPp5OvXVFXzqlzbxQv8IX3hgd9zlSAFRAEhB+vaOQ3xx215+/crzePNFxT/xaylv39TBuy/v5tYHXuLe5w/HXY4UCAWAFJx7dvbzB994hq3rWvnj6y6Ou5y88el3X8Jlq5v58Fef5sEXdaE4WZoCQArG9Nw8n/vBi9z8L0+xqbuJOz6wpaSbftLVVJVzxwe2sGFlPR/60hPc9tBu5tQnIGcQaQCY2Y1m9piZPWpmV6TtS5jZV8xsW/C7+JZvkqyYmJnj608c4B23bOPWB3bzniu6+dpNV9KQqIy7tLyzor6ab/zmVbx9UwefuedFrrv1h9z99CGmZufjLk3ykEV1OVkzawHuB64EuoG73P3qlP03A+3u/kkz+xPgqLv/45mec8uWLb59+/ZI6pV4uTsTM/MMT84yNDHL/uPjvHxsjCf3D/L43hNMzMxzUUcDf/iOi0riYm/ZcM/Ofv72B7vYfXSM2qpy3rChjctWN3FhRyOdTQlWNSZorq2kslwNAcXMzJ509y2L7Yty4PRWYJu7zwB7zazBzKrdfTrY/ybgM8Ht7wAfA84YAGfr+s9vY2o2eSqcHninxJ+f8e6Sj0/PU087InX/UtmbzddKf/ypr535Y8+qtjMd7zA1N8/s/KlvyPq2On75Nat512VdbDmvpaRm+p6rt2/q5Odf1cGjLw9wz87DPLp7gP944dTrBlWUGYnKchKV5VSV28n3uKwMDMMMyswwAAMD/T3kWGttFd+4+aqsP2+UAbACSJ2ZMgS0Av2L7A/3ncLMbgJuAujpObsLfW1or3/lh0vav930f8rp/7hP3X9uj0/dYGl7l37u9P1LPP4ML5711zrl8Wf+kEjdnagsp6mmkqaaSpprKlnTWsv69jpqq0pvclc2lZUl1xB+48bkEpkjU7O8fHSMIyPTHBmZYnRqlsnZeSZnFpicnWd2fgH34MuEw4Inv1Z4yu1TvzVJ1BoS0fw/iPJ/1wmgOeV+U7Btsf3p+05y99uB2yHZBHQ2hdzy3svP5mEiRacxUcnlPcW7XKYsT5SNf48DV5tZpZn1AGMpzT8ADwPvDG6/M7gvIiI5ElkAuPsgcBvJD/avAr9rZpvN7KPBIV8CLjGzbcAlwX0REcmRSBtY3f1O4M60zTuCfZPADVG+voiInJ7Gf4mIlCgFgIhIiVIAiIiUKAWAiEiJUgCIiJSoyK4FFAUzOwbsj7uOJbQBA3EXkQHVmV2qM7tUZ/ac5+7ti+0oqAAoBGa2/XQXXsonqjO7VGd2qc7cUBOQiEiJUgCIiJQoBUD23R53ARlSndmlOrNLdeaA+gBEREqUzgBEREqUAkBEpEQpALLIzH7fzB4xs0fN7J/NrDLYfkWw7TEzuzHmMjGzW8zsx8HPH6Zsf7uZ/Sj4+U9x1hjU824z+4mZTaVtz7f388aglkfN7Iq46wmZ2b1mdszMPh7cNzP7vJltM7Pvmtmiq/DlmpldHrx3j5jZA2a23swSZvaVoNavmFkiD+psDP6eHzKz/2dmP5ev72nG3F0/WfoBqlJu/zNwXXD7UWAdUAU8BbTEXOfG4HcZ8CPgfKAceJbkKm3NwDNAecx1rgASwO607XnzfgItQQ1VQU0/jPM9S6ttNXAj8PHg/tuBO4Lb7wf+Ou4ag1o6gIbg9juBu4CbgU8E2/4EuDkP6iwDKoLb64En8vU9zfRHZwBZ5O4zkPymRfIfy24zqwbq3H1vsH8bsDXGMnH3l4LfC8AcMA9sAPa6+5C7DwH7gm2xcffj7p7+7T/f3s+twDZ3n3H3vUBDUGPs3L03bdObgO8Gt78T3I+dux9299Hg7jTJf5N5V6u7L7j7XHC3keQXpryrczkUAFlmZn8M7CK5yP1Bkt9ih1IOGQr2xc7M3gfscfd9JOscTNk9RJ7UmSbf3s9Ced/glbUOkTx7yRtmVgd8Cvgsp9aaF++pmXWb2Q+BHwDfIs/f06VEuiJYMTKzcpJNEOm+5+6fdPe/NLNPA18gefp9J8kmlVATcCLuOs3srcAHgXcF20+Qh3Uusj2WOs8g3+o5k9Ram3hlcMUq6C/7OvA37v6CmaXXmhfvqbsfIrnW+VrgIeBr5Ol7mgkFwDK5+zxw5WL7zCzh7lPu7mY2DEy4+5SZjZtZD9APXA38ecx1vg74JPAOTy7NCfASsM7MGoP764DdcdZ5muNjeT/P4HHgU8EHWCcw5u7TMdZzJg8D7wbuJtnW/nCs1QTMrAz4F+Bud7872PwwyRp3kCe1mll1yt/tCDBKnr6nmVIAZNffmdmrCdr/gT8Ntn8E+CpgwG3uHve3hDuC33cnuyv4A3d/0sz+CLg32PdHwYdzbMzsjSTfwy4zu4/ke/dN8uj9dPdBM7uN5H98D2rLC2b2ReD1QLWZbQHeA1xvZttIfoC9P876UrwHuA5YZWa/BjwHfAy4M6i1l+TZatw2mdn/ItlnVgH8LvAg+fmeZkQzgUVESpQ6gUVESpQCQESkRCkARERKlAJARKREKQBEREqUAkBEpEQpAERESpQCQEqWmf2SmbmZXZS2/e/M7Bkz+3za9hozezi4fAXBY/8lZX9FcPnl76ZsG0t7jhvN7AunqacquCSyJmhKTigApJTdAGwPfgNgZucDb3D3y9z9d9KO/6/AN1NmSI+TnB1aE9x/G3DobIsJrm56P/ArZ/scIsuhAJCSZGb1wLXAbxAEgJldSPICX+eZ2dPB1SlTvQ/4dtq275O8jAHB83w1w9e/2cx2BD97zezBYNfdweuIRE4BIKXqF4H73P0ZYMzMXuPuLwJfJrkQyeXuPh4ebGZVwPrg0tmpvga8N1ix6lKSF4dLVZPyQb8D+AsAd/9Hd98MvJbktW4+Fxy/M9gmEjkFgJSqG4BvBLe/wc+agS4huRpaujZeuQ4BAO7+LLA2ePz3F3ncpLtvDn9Irm6V6u+BB9z9O8HzzQMzZtawnD+MyNlQAEjJCdZtfR1wT7DpG8CvBCu5vZrkt/B0kySXp1zMvwF/S4bNPyl13Aicx6mXs64Gpk55gEiWabSBlKJfBr4fXtvd3feYWT9wDTCbskbCScFln8vDNR/Sdt8JDLn7c2Z2bSYFmNlrgP8OvDFYmjPcvgIYcPfZs/mDiSyHAkBK0Q3AZWa2L2XbCuDPWPzbf+gHJBeguS91Y7D27q3LrOG3SS5z+GCwJsN2d/8N4M3A95b5XCJnResBiGTIzK4Afs/dfz3C1/gm8Ifuviuq1xAJqQ9AJEPu/hTJb+zlUTx/MNLobn34S67oDEBEpETpDEBEpEQpAERESpQCQESkRCkARERKlAJARKREKQBERErU/wfWnUkGtX6vTQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "spins, gs, pdf = create_spins_inhomo_g(gcoll, g_list, q, delta, gammaq, nbins, holes=False, inhomo=True, hyperfine=True)\n",
    "plt.figure()\n",
    "plt.plot(spins[:len(pdf)], pdf/max(pdf))\n",
    "plt.xlabel('$\\Delta f$ (MHz)')\n",
    "plt.ylabel('weight');\n",
    "#%qtconsole"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fb1af2020d0>]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEICAYAAACwDehOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAxsUlEQVR4nO3dd3jUZbr/8fedSa+kUgIpQAi9I0VQEOvadV3rHnXlWI6KW86e3f2t7uq66/GsetS1s+pRsbvr2isC0tXQe4AUCJBGSO8zz++PmcQQQxLCzHyTmft1XXMl8512zwWZzzzfp4kxBqWUUgogwOoClFJK9R4aCkoppVppKCillGqloaCUUqqVhoJSSqlWgVYXcLISEhJMWlqa1WUopVSfsn79+lJjTGL7430+FNLS0sjKyrK6DKWU6lNEJL+j43r6SCmlVCsNBaWUUq00FJRSSrWyrE9BRB4DZriuvmeMedB1PBdoOdf1pTHmLxaUp5RSfsnKjuanjDE/F5EAYLWIvGOM2QfYjTFzLaxLKaX8lmWnj4wxe1w/HUAzYHfdJCKyTEQ+E5GJVtWnlFL+yPIhqSJyLZBjjMlzHZpujCkVkQnAa8DYDh5zM3AzQEpKirdKVUopn+fRloKI2ERkXQeXe1y3nwncCNza8hhjTKnr52agVkRi2z+vMWaRMWaqMWZqYuIP5l4o1WOfbj3Mmn2lVpfRaldhJW9/d4DGZofVpSg/4dGWgjHGzvedyccQkenA/cB5xpg617EQQIwx9SKSDPQDyj1Zo1ItvtpZxG2vbQDg07vmMGpgtKX11DQ0c93z31Ja3cDB8jp+cdYIS+tR/sHKIakvAFHAeyKyXESmAEnAGhFZCbwD3GJ0FyDlJYvX5RMebCNA4J2sAqvLYcnOIkqrG4gMCeTN7/Zjd+ifgvI8Kzuax7ouc12X9caYA8aYycaYOcaYWcaYZVbVp/xLQ7OdtfuOcOW0IczNTGLpriKrS+KrncXERwTzp4vHUFTZwI5DlVaXpPyATl5TCti4v5yGZgenDktgenoceUdqOVLdYGlN3+QeYXZGAjOHxQPwXV6ZpfUo/6ChoBSw7WAFAJNTY5k4pB8AmwvKLaunpKqBosoGxiXHMDAmjMGxYWTlaygoz9NQUArYV1JDbHgQcRHBjBscgwhsKaiwrJ7th5yvPWZQDADjB8fo6SPlFRoKSgE5JdUMTYwEIDw4kJS4cPYUV1tWz3ZXAIwe5BwBldk/mvyyWmobmy2rSfkHDQWlgJzSGoYmRLReH54YyT4LQyG7qIrkfmHEhAUBMHJgFMZAdpF1NSn/oKGg/F5lfRMlVQ2tLQWA4UmR5JTW0Gy3ZtJYbmkNQxO/D6lRA5wthp2H9RSS8iwNBeX3ckpqAI75EB6WFEljs4MDR+u8Xo8xhtzSGtLiv69ncGwYoUEBlrZelH/QUFB+L6fE+UE7rE0oDE9ythr2WvAhXFbTSFV9M2ltTmcFBAhp8RHkltZ4vR7lXzQUlN/LKanBFiCkxLVpKSQ4QyG31PuhkHfE+cGfnhB+zPGhiRoKyvM0FJTfyymtJiUunODA7/8cYsKDiAkLIv9IrdfryS11vmZ6QuQxx9MTIthfVkuTRf0cyj9oKCi/l1Ny7MijFqnx4ewvsyIUqrEFCINjw445np4QSbPDUGBBP4fyHxoKyq85HOYHI31apMRZEwp5pbUMiQ0jyHbsn2e6K7isOKWl/IeGgvJrB8vraGh2HDMctUVqfDgHj9Z5fVhqbmnNMZ3MLVpaMy2jpZTyBA0F5ddyXB23HZ0+SokLp9lhOFxR77V6jDHkHzl2OGqL2Ihg+oUHaWez8igNBeXXWoajdtRSaBmN5M3O5tLqRmoa7aTFh3d4e3qCjkBSnqWhoPxaTkkNUaGBJEQG/+C2VNcHc36Z9z6EW/owUjtoKYAzFPT0kfIkDQXl13JKnQvhicgPbhsQHUqwLYD9Xmwp7HcF0JC447QU4iMorKynrtHutZqUf9FQUH4tp6SGYR30J4BzFvHguDCvjkDKP1KLCD8Yjtoi3TVKqmWCm1LupqGg/FZtYzOHK+o7HI7aIjUu3Kt9CvvLahkQHUpokK3D21s6oLVfQXmKhoLyW98vhPfDTuYWqfHOWcTGGK/UtP9ILSnHOXUEbecqaCgoz9BQUH6rdThqJy2FlLhwqhuaKatp9EpN+8tqWzu4OxIREkj/6BANBeUxGgrKb+WUVCNCh3MCWqQltIxA8vwppLpGO8VVDZ22FABdLVV5lIaC8ls5JTUk9ws77vl7+H5oaL4XOnZbOrRTOgkpcLZs8jQUlIdoKCi/1TIctTODY8MIEOd6RJ7WGgpdtBTSEyI4UtNIRW2Tx2tS/kdDQfklYwy5x1kdta2QQBsDY8K80lJoeY3Ubpw+AsjVYanKAzQUlF8qqmygptF+zG5rx5OWEE6eF4al7i+rJSo0kH7hQZ3er6VjXE8hKU/QUFB+6fstODs/fQTOfgVv9SmkxIV3OLu6rSFx4QTI96OnlHInDQXll/aVdj1HoUVafDhHa5s8fg5//5HOh6O2CAm0kRwbpi0F5REaCsov7SuuJiLYRv/okC7v2zoCyYML49ldO6odb82j9tITInVYqvIIDQXll/aVVDMsqeOF8Npr6dj1ZL9CYWU9jXYHqXFd93EApMeHk1ta47WZ1sp/aCgov7SvuLpb/Qnw/RDRfA9+M28dedSN00fgHJZa3dBMabV3Zlor/6GhoPxOTUMzhyrqGZ7UvVAIC7YxIDrUoy2FlnWY0rsYItsiTddAUh6ioaD8TssHcHeGo7ZIiQ9v3evAE3JKaggLcoZPdwxNcAaadjYrd7MsFETklyKyQkRWi8grIhLkOj7ZdWyNiNxgVX3Kd+07geGoLdLiPTtXIae0mvSECAICuu7jAEiODSPIJuwrrfZYTco/WdlSeNIYc5ox5lTX9bNdP58ArgPmAgtFJNaK4pTv2ltcjS1AjrvlZUdS4yMoqWqgpqHZIzXllNR0ulpre7YAIT0hgn3F2lJQ7mVZKBhjGgHEOfwjANgrIiFAhDEm13X7SuCU9o8VkZtFJEtEskpKSrxat+r79pVUkxoXTnBg9//7e3Jzm4ZmOwVHa7s1Z6KtjP5RZBdVub0e5d8s7VMQkd8D2UAccACIB8rb3KXcddsxjDGLjDFTjTFTExMTvVCp8iV7i53DUU/EsCRnKLScenKn/CO1OMyJ9XEAjEiK4sDRWt2vWbmVR0NBRGwisq6Dyz0Axpi/ACOAXOAGoAzo1+YpYlzHlHKLZruDvCM1J9SfAM5RQbYAYW+x+0OhZcmNls7j7hrRPxJj8EhNyn8FevLJjTF2YEZHt4lIqDGm3hhjRKQCqDXG1ItIjYikAIeB2cB9nqxR+Zf8slqa7OaEv5WHBNpIjQ9nT5H7P4D3tQxHPdGWwoAoAHYXVTFucIzb61L+yaOh0IVHRGQMrv4E4I+u43cBbwACPG2MOWpRfcoH7S50noMfOSD6hB+bkRRJdrH7z+HnlNTQPzqEyJAT+3NMjQsn2BbAHu1XUG5kWSgYY24/zvEs4NSOblPqZO06XEmAQEb/EztVA5CRFMWSncU0NNsJCTz+bm0nal9J9QmfOgIItAUwNDFCO5uVW+nkNeVXdhZWkZ4Q0ekWnMeT0T8Su8O4dRc2h8OQXVRFputU0Ika0T+KbA+c0lL+S0NB+ZVdhZWMHHjip46A1mUx9rjxFNKBo7XUNtoZ2cNQyBwQxcHyOqo9NH9C+R8NBeU3quqbOFBWx6gefgAPS4wkQHBrZ/Oulj6OHgZVhiuo9BSSchcNBeU3Wj44e9LJDBAaZCMlLtytLYVdh6sQcQ4v7YlRrjDZebjSbTUp/6ahoPzGzsMt38p71lIAZ6C0PI877C6qJDUunPDgno35GBwbRkxYENsOaigo99BQUH5jV2ElUSGBJPcL6/FzjE2OJre0hqp692zNuetwVY9bLgAiwphB0ew4VOGWepTSUFB+Y+fhKkYOjOrWbmvHMybZOUlsx6GT/2Ze12gn90hNj0cetRibHMPOwiqa7I6TrkkpDQXlF+wOw45DlYzuYYdui7GDnKGwzQ2hsLOwEmO+7xfoqTGDomlsduhyF8otNBSUX9hbXE1dk53xg/ud1PMkRoXQPzqE7QdP/nTNlgPlAEwccnI1jXEF1XY3BJVSGgrKL2wuKAdgwpCTXyNo7KAYtrnhHP7mggr6R4cwIKZ7u60dT3pCBOHBNra5IaiU0lBQfmFLQTmRIYE9Wk6ivTHJMc6Wx0kuWb35QPlJt1zAueHOqIHRbNfOZuUGGgrKL2wtqGBscnS3t7vszMQhMTjM962PnqioayKntOakTx21GDsomu2HKrE7jFueT/kvDQXl8xqbHew8XMUEN3wrB5iS4tz3KSuv51t9bC1wfqt3V00TU/pR22hvXQVWqZ7SUFA+b1dhJY12h1tO1QDEhAeR2T+K7/J6vqr7t3llBAiMd0MfB8DkFOdW5hv260rz6uRoKCift9n1rXy8GzeimZIWy4b8oz0+XbMu5whjk2OIDg1ySz0pceHERwRrKKiTpqGgfN7G/KMkRAYzOLbnM5nbm5YWS1VDM7sKT3wYaH2TnU37y5k5NN5t9YgIk1Ji2bS/3G3PqfyThoLyed/llzE1Ne6kZjK3Nz3d+YG+dt+RE37shvyjNNodzHBjKABMTu1HTmkNR2sa3fq8yr9oKCifVlhRz4GyOqamxbr1eQf1CyMjKZKvs0tO+LFrc44QILi9ppZ+hY0H9BSS6jkNBeXTvnONEDolPc7tz336iES+ySmjtvHENrhZtruYySmxRLmpP6HF+MEx2AKEDfnlbn1e5V80FJRPy8orIzzYdtJrHnVkbmYSjXYH63K6fwrpcEUd2w5WMn9Uf7fXEx4cyJhB0Xyb2/OhskppKCif9m3eUSal9CPQ5v7/6tPSYwkPtvHljqJuP+arncUAnDU6ye31AMwcGs/GA0dPera18l8aCspnVdY3sauwkqmp7j91BBASaOPs0f35ZGshjc3dW7b6022HSYsPZ1jiyS+30ZEZQ+NpshvW52u/guoZDQXls77NKcMYmO6B/oQWF09MpqKuiRXd6HA+VF7Hmn1HuGRSsltHQrU1LT0OW4CwNqfUI8+vfJ+GgvJZq/aWEhoUwORU947yaWt2RgKx4UH8a+PBLu/73qaDGAOXTRrssXoiQwIZlxzDuhztV1A9o6GgfNbqvaVMS4sjNMjmsdcIsgVw2eTBfLa9kEPldce9X7Pdwevf7OeU9DhS4sM9Vg/AzGHxbD5QTk3DiY2KUgo0FJSPKqyoZ09xNXMyEjz+WjeemgbAi6tyj3ufT7YVUnC0jgWz0z1ez8yh8TQ7DFnar6B6QENB+aTVe53n1E8d7vlQGBwbzoXjB/LqN/kc7KC10NBs59Evs8lIiuRMDwxFbW9qWizBtgBW9mBinVIaCsonrdpbSnxEMKMGuH9+Qkf+85xMjIE/vr8NY45dJO/Z5TnkltZw9wWj3bKfQ1fCgwOZPjSO5RoKqgc0FJTPMcawam8ppw5P8MqHMDhbC78+J5MlO4t56PPdrcHw0ZZDPP5VNhdPHMTpIxK9Ugs4J9btLa7mQFmt115T+YZAqwtQyt22H6qkpKqB07z4IQxw0+x09hZX8/TyfSzfXUJMWBBrc44wNTWWBy4d59Va5mUmcv9HsDy7hJ/OSPXqa6u+TVsKyud8uaMIEecHozeJCP992Tj+5/JxhAYFUFHXxC/OHMGrC6YTEeLd71/pCRGkxoezfFexV19X9X3aUlA+56tdRUxJiSU+MsTrry0iXDkthSunpXj9tdvXMS8ziTe/2099k92jw3KVb7GspSAivxSRFSKyWkReEZEg1/HlIrLW9fMJq+pTfZMnF5zra+ZmJlLfdGIL9ill5emjJ40xpxljTnVdP7vNbVcYY+YaY+60ojDVd3l6wbm+ZMbQ+BNesE8py0LBGNMIIM5FYAKAvS03AW+KyFIROaOjx4rIzSKSJSJZJSU67E59b8nOIlI9uOBcXxIaZGNeZhKfby/q8V7Syv9Y2tEsIr8HsoE44IDr8BXGmNnA9cAzIhLV/nHGmEXGmKnGmKmJid7tTFS9V1V9E2v2HWH+yP4eW3Curzl37ABKqxvIytO1kFT3eDQURMQmIus6uNwDYIz5CzACyAVucB0rdf08AGwGhnuyRuU7vtxRRGOzg/PHD7C6lF5j3sgkQgID+HRbodWlqD6iy9FHIjITuA6YAwwE6oBtwMfAq8aYiuM91hhjB2Yc53lDjTH1xhgjIhVAretUUpQxptLVQhgH5J/om1L+6aMthxkUE8qkIZ5bFbWviQwJ5LQRiXy2rZA/eGlGterbOm0piMinwALgc+BcnKEwGrgbCAXeF5GLevjaj7hGGK0ABgCv4QypZSKyClgC3GuM0Xav6lJFbRMr95Rw/viB+sHXznljB1BYWc/GA+VWl6L6gK5aCj9tOZ3TRjWwwXV5RER6tOKYMeb249w0pSfPp/zb5zsKabIbLhg/yOpSep35o/oTbAvgw82HmOLBvSWUb+i0pdBBICAiCdKmF6+j+yjlbR9tOUxKXDjjB8dYXUqvExMWxPxRSXy4+RBN9u5tG6r8V1enj2a4TvG8KyKTRGQbzv6EIhE51zslKtW54qp6Vu8t5YLxA3XU0XFcOimZIzWN3do2VPm3rkYfPQk8ALwBLAUWGGMGAKcB/+3h2pTqlvc2HsTuMFw+xXPbXPZ1czOTiA0P4t1ubBuq/FtXoRBojPnCGPMOUGiMWQdgjNnl+dKU6poxhneyCpic0k8nrHUiODCACycM4ssdRVTWN1ldjurFugqFticg228ppVMkleU2F1Swp7iaK6YOsbqUXu/SSck0Njv4dOthq0tRvVhXoTBBRCpFpAoY7/q95bp3F4hXqgPvZB0gNCiA88cPtLqUXm/ikH4MS4zgze8OdH1n5be6Gn1kM8ZEG2OijDGBrt9brgd5q0ilOlLfZOeDzYc4d8wAokP1v2NXRIRrpqeycX85Ow5VWl2O6qW6Gn0U19nFW0Uq1ZEPNx+iqr6Zn0zTU0fddfnkZIIDA3j9W10oQHWsq9NH64Es188SnIvX7XH9vt6zpSnVucXr8slIimTm0HirS+kz+oUHc8H4gby38RA1Dc1Wl6N6oa5OH6UbY4biXHLiQmNMgjEmHrgA+MIbBSrVkU0HytlSUMFPZ6bq3IQTdO30VKobmvlg8yGrS1G9UHdXSZ1hjPmk5Yox5lNglmdKUqprr6zNIyLYxqWTkq0upc+ZnNKPkQOieHVdPsboIEJ1rO6GwiERuVtE0lyX3wP6NUNZoqymkY+2HOayyYOJ0g7mEyYiXD8rje2HKlmXo+tNqmN1NxSuBhKBfwHvun6/2lNFKdWZ17/Jp7HZwU9nplpdSp916aRk4iOCeWFVjtWlqF6my/0UAFzLV9/l4VqU6lJ9k52X1uQxNzOREf1/sCmf6qbQIBvXzUjl8a/2kFNSzVCdDa5cuhqS+ncR6XCSmohEiMjPRORaz5Sm1A+9u+EgpdWN3HLaMKtL6fOum5FKcGAAL6zKtboU1Yt0dfroKeAeEdkpIu+IyNMi8qKIrATWAFHAPzxepVKA3WH4+8ocJgyOYcZQnSZzshKjQrh0YjL/3FBAWU2j1eWoXqKrIambjDE/AabhDIiVwAc4V0udYIx53BjT4IU6leLLHYXkltZwy+nDdBiqmyyYk059k4P/W62tBeXU3T6FamC5Z0tR6viMMTzzdQ6p8eGcM2aA1eX4jIz+UZw3dgAvrc5jwZyhxITpaC5/193RR0pZ6uvsEjYfKOfm04Zi0z2Y3eqOM4ZT1dDMS6vzrC5F9QIaCqrXM8bw6JfZJPcL44opus6Ru40ZFMOZo/rz4upcqnSvBb/X7VAQEQ0QZYmlu4rZXFDBwvnDCQ7U/4aesHD+cCrqmnhlrS6U5++6GpL6goh8KCKpwLNeqkmpVsYY/vfLbFLiwrlssm636SnjB/djXmYif1+Zozuz+bmuvnYZ4CfAPYAuRam87osdRWw/VMnC+RkE2bSV4Em/OjuT8tomnvt6n9WlKAt19VdWYIypA34N6JoCyqua7Q4e+WI36QkRXDJxkNXl+LyxyTFcNGEQL6zKpaiy3upylEW6mqdwr+vnUWPMVK9UpJTLO+sLyC6q5r/OySRQWwle8Z9nZ2J3GB7/ao/VpSiLdDpPQUT+rZvPs8kYs8UN9SgFQE1DM498kc3U1FjOHavzErwlJT6ca6ensnhdPjfNTmeYronkd7r6+pXezUs/z5Wo/NFzK3IorW7g9+eP0tnLXnbHGcMJDQzgr5/tsroUZYFOWwrGmPuOd5uIBBljdJiCcrvCinoWrdjHBeMHMikl1upy/E5CZAi3zR3Gw19ks3pvKacOT7C6JOVFJ3SiVpzmi8gLQIGHalJ+7qHPd+NwwG/OHWl1KX5rwZyhpMSFc+8H22myO6wuR3lRt0JBRGaIyN+AfOB9YAWgf7HK7b7LK+OfGwq4aU46Q+LCrS7Hb4UG2bjngtHsKa5msU5o8ytdTV57QET2AH8BtgCTgBJjzMvGmKPeKFD5j2a7g3ve20ZyvzDuPGO41eX4vTNHJXHaiEQeXZJNabUuhuwvumopLACKgGeAxcaYIzgntCnldi+tyWNXYRX3XDCa8OBuLeCrPEhE+MMFo6lrtGunsx/pKhQGAn8GLgT2ichiIExE3PYXKyIvi8iSNtfPFZG1rss57nod1bsVVdbz2JI9zM1M5Jwx/a0uR7kMT4rkpjnpvJ1VwLqcI1aXo7ygq8lrdmPMZ8aY64FhwHvAauCgiLx+si8uIuNpM5xVRGzAX4HzXJe/uo4pH/enD3fQaHdw30VjdAhqL/Pz+SNIiQvnd+9upb7JbnU5ysO6PfrIGNNgjPmnMebHQAbwmRte/x7ggTbXhwO5xphyY0w5kOc6pnzYZ9sO8/HWw9w5bzip8RFWl6PaCQu28cCl48gtreGJpTrT2dd1NaM5pZObl7e5vdwYU3kiLywic4FsnH0WLeKBth3Y5cAPNuMVkZuBmwFSUjorUfV25bWN3P3edsYMiubWucOsLkcdx+yMBC6fPJjnvs7hgvGDGDUw2uqSlId01TfwMs6O5c7a8wZ4CXil/Q2uUz+rO3jMx8CpwFUcOxu6rN31GNexY1/QmEXAIoCpU6dqx3cf9qcPd1Be28jLP5umq6D2cnefP4rlu4v57T+38M/bZul6VD6qqxnN807myY0xdmBG++MiEgVcDrwJhAFjROT3wINAuoi0fA1JB/aeTA2q91q6q4h3Nx5k4RnDGTMoxupyVBdiI4K57+Ix3PH6Rp5evo+F8zOsLkl5gCXj/owxVcBEABFJA543xvzFdf13wOeuu/7OFSzKx5TXNvL/3t3GiP6R3K5zEvqMC8YP4ovtRfztqz3My0xi3GANc19jefvPGJNnjDmzzfVPjDEzXZdPrKxNeYYxhv/3r62UVjfwyBUTCQnUAWZ9yf0XjyUhMoSfv7VRRyP5IMtDQfmfd9YX8MnWQn51dqZ+0+yDYsKDeOiK8ewrqeHBT3VSm6/RUFBelVtaw70fbGfG0DhuPm2o1eWoHpqTkcgNs9J4aU0eK7JLrC5HuZGGgvKaJruDn7+5kSBbAI9eORFbgE5S68t+c+5IMvtH8Yu3Nun2nT5EQ0F5zcNf7GZzQQUPXjaOgTFhVpejTlJYsI2nrp1EbaOdhW9spFmX2PYJGgrKK77YXshzX+dwzfQUzhs30OpylJsMT4riz5eM5ZvcMv6m+zr7BA0F5XF5pTX86p3NjB8cwx8uGG11OcrNLp8ymB9PGcwTy/ayak+p1eWok6ShoDyqvsnOba9tIECEp66ZTGiQDj/1RX+6eAzDEyO5682NHCqvs7ocdRI0FJRH3fPeNnYVVvLYVRN1JzUfFh4cyDPXTaah2cEti9fr/IU+TENBecwra/N4Z30Bd84bzrzMJKvLUR42PCmKR6+cyNaDFfzu3a0Yo8uS9UUaCsojVu4p4b4Pd3DmqCTuOnOE1eUoLzlrdH9+edYI/rXxIC+syrW6HNUDGgrK7XJKqrn9tQ0MT4zksasm6XwEP3PHvOGcO2YAD3yyUzue+yANBeVWFbVNLHg5i0BbAM9fP5XIEN1r2d8EBAiP/GQCGUlR/Mdr69lbXGV1SeoEaCgot2myO7jjjQ0cOFrLs9dN0Y5lPxYREsjz108lONDG9S9+R3GVznjuKzQUlFsYY/jdu1tZuaeUv1wyjlPSf7BhnvIzQ+LCefGGqZTVNLLg5SxqG5utLkl1g4aCcotHvsjmH+sLuGt+Bj+ZNsTqclQvMX5wP564ehLbDlaw8I2N2B06Iqm301BQJ23xunyeXLaXq6YN4edn6m5c6lhnju7PvReNYcnOYu77cLsOVe3ltBdQnZTPtxfyx/e3MX9kEn++ZCwiOtJI/dC/zUyj4Ggdi1bkEBsezC/O0mHKvZWGguqxdTlHWPjGRucpgmsm6UbuqlO/PXckZTWNPP7VHmLCgvjZ7HSrS1Id0FBQPbJh/1Fueuk7V2fiNMKD9b+S6lxAgPDgZeOorm/mTx/tICo0kCumav9Tb6Nf7dQJ23awgutf/JaEqBBeXzCduIhgq0tSfUSgLYDHr57InIwEfvPPLXy2rdDqklQ7GgrqhGQXVfFvL35LdGgQry2YTlJ0qNUlqT4mJNDGs9dNYcKQfix8YyPLdhdbXZJqQ0NBdVtuaQ3XPv8NgQHCawumMzhWJ6epnokICeSlG05hxIBIbnllPct2aTD0FhoKqlv2lVRz1aK12B2G1xZMJy0hwuqSVB8XEx7EazfNIHNAFLcsXs/SXUVWl6TQUFDdsLuwiiufW4fdYXjj32eQ0T/K6pKUj4gJD+LVm6YzcmAUty7ewFc7NRispqGgOrXtYAVXLVpLgMCbN88kc4AGgnKvmPAgFrcEw6vrWbJDg8FKGgrquDYfKOeav68jLMjG27fMZHhSpNUlKR8VE+YMhtEDo7n11fW8v+mg1SX5LQ0F1aFvc8u49vlviAkP4q1bZmofgvK4mLAgXl0wnSmpsfz8rU0sXpdvdUl+SUNB/cDn2wu57oVvSIoK4a2bZ+oS2MprokKDePlnpzB/ZBL3vLeNJ5fu0bWSvExDQR3jzW/3c9ur6xk1MJp/3DaLQf3CrC5J+ZnQIBvPXDeFSycl8/AX2TzwyU4NBi/StQkU4NwP4alle3n4i2xOH5HIM9dN1qUrlGWCbAE8csUEokMD+fvKXI5UN/Lg5eMJDtTvsZ6mf/UKu8Nw/0c7eGlNHpdOSuavPx5PkC5upywWECDce9EYEiJDeOTLbA5X1PPsT6cQExZkdWk+Tf/y/VxNQzO3LF7PS2vyWDA7nUeumKCBoHoNEeHO+Rk8euUEsvLL+PEzayg4Wmt1WT7N8r9+EXlZRJa0uZ4rIstdl99bWZuvO1xRxxXPrmXpriLuvXA0d18wmoAA3Q9B9T6XThrMyz87hcLKei59eg3bDlZYXZLPsjQURGQ80K/dYbsxZq7r8hcLyvILWwsquOSp1ewvq+WF66dxw6m6tr3q3WYNS+Dd22YRbAvgJ8+t5Uud5OYRVrcU7gEeaHdMRGSZiHwmIhMtqMnnfbatkJ88t5bAgAD+cdtM5o1Msrokpbolo38U/7p9FsOTIrl5cRZPLdurI5PczLJQEJG5QDbQPu6nG2PmAb8BXj3OY28WkSwRySopKfFonb6kZYTRba+tZ8QA5x/XyAHRVpel1AlJigrl7VtmcvGEQTz0+W7ueGMjtY3NVpflMzw6+khEbMDqDm76GDgVuIp2p4+MMaWun5tFpFZEYo0xR9vdZxGwCGDq1Kn6NaEbqhua+dXbm/h8exEXTRjEX388ntAgm9VlKdUjoUE2Hr1yIqMGRvPgZ7vILanh79dPJVnn1Zw0saLpJSJRwEqgEAgDxgCPAg+7aqoXkWRgGZBpOily6tSpJisrywtV9137Sqq5ZfF6cktr+N15I7lpdjoi2qGsfMOy3cUsfGMjwbYAnr52MtOHxltdUp8gIuuNMVPbH7fk9JExpsoYM9EYcy5wPbDJ1amcBKwRkZXAO8AtnQWC6tqXO4q45MnVHK1pZPFNp7BgzlANBOVT5mUm8d7tpxITFsQ1z3/D31fkaD/DSbCkpeBO2lLomN1heHxJNn9bupfxg2N49ropumSF8mmV9U385h9b+HRbIWeP7s9DV0zQiW6d6FUtBeVZRZX1XPv8Ov62dC9XTBnM27fM1EBQPi86NIinr53MPReMZumuYi58YpXOZ+gBDQUf83V2CT96fCWbD1Tw8BUTeOiKCdqhrPyGiHDT7HTeumUGTXYHlz2zhje+3a+nk06AhoKPaLY7+J/PdnH9i9+SEBnCB3ecyo+nDLa6LKUsMSU1jo/unM309Dh+9+5W7npzE5X1TVaX1Sfogng+4FB5HQvf2EhW/lGumjaEP144hrBgbR0o/xYfGcJLN57CM8v38uiSPWzYf5THr5rIlNQ4q0vr1bSl0Me9v+kg5z62gp2HK3n8qok8ePl4DQSlXGwBwh1nZPDOrTMRgSueXcujX2bTbHdYXVqvpaHQR1XUNrHwjY3c9eYmhiVF8vHCOVw8MdnqspTqlSanxPLJwjlcMjGZx7/aw5WL1nGgTFdb7YiGQh+0Zm8p5z6+go+3HuaXZ43gHd1DWakuRYUG8b9XTuSxKyeyu7CKHz2+kn+uL9BO6HY0FPqQ+iY793+0g2ue/4awIBvv3jaLhfMzCNT9D5TqtksmJfPpXXPIHBDFr97ZzL+/sp7iqnqry+o19NOkj9i4/ygXPrGKF1bl8tMZqXy8cA4ThvSzuiyl+qQhceG8dctM7j5/FCv2lHD2oyv4YPMhbTWgodDr1TfZeeCTnVz+zBqqG5r5vxuncf8lY7UzWamTZAsQFswZyicL55AaH8HCNzbyH69t4Eh1g9WlWUqXuejFvs0t4zf/3EJuaQ3XTE/hd+eNJCpUp+0r5W7NdgeLVubw2Jd7iAoN5L6Lx3D+uIE+vU6YLnPRh9Q0NPPH97fxk+fW0uxw8PqC6Txw6TgNBKU8JNAWwH/MHc6Hd85mUL8w7nh9IwtezuJgeZ3VpXmdthR6mWW7i7n7X9s4VFHHDbPS+PU5mYQH6xxDpbyl2e7gpTV5PPJFNiLwn2dncv2sNGw+tn/58VoKGgq9xOGKOv704Q4+3VbIsMQI/ufy8UxN05mXSlnlQFkt97y/jeW7Sxg/OIb/vmwcYwbFWF2W22go9FIt30oe/TKbZodh4fwM/n3OUIID9cyeUlYzxvDRlsPc9+F2jtY2sWB2OnedmeETrffjhULff2d92Ib9R/n9v7ax83Al8zITue+isaTEh1tdllLKRUS4cMIgTstI5MHPdvLcihw+2HyIu88fzY/GDfDJjmhtKVigrKaRhz7fzZvf7ad/VCj3XjSac8b45n8wpXxJVl4Zf3h/OzsOVzJrWDz3XTSGjP5RVpfVI3r6qBdosjt4dV0+j36ZTU2jnRtnpfHzs0YQGaINNqX6CrvD8Po3+Tz0+W5qG+3cMCuNu87M6HOjAzUULLZyTwl/+nAHe4qrmT08gT9cOJoRffQbhlIKjlQ38PAXu3nzuwMkRIbw/340kksmJveZFr+GgkXySmv488c7WbKziNT4cO4+fzRnjkrqM/9xlFKd23SgnD++v43NBRVMHNKPey4Y1Sf2bNBQ8LKq+iaeXLaXF1flEmwL4I4zMvjZ7DRCAnV5CqV8jcNh+MeGAh7+fDfFVQ38aNwAfnPuSFLje+/qxRoKXtLY7OC1b/J5YuleymoauWLKYH59TiZJ0aFWl6aU8rDaxmYWrcjhua9zaHY4uH5mGneekUFMeO/rb9BQ8DCHw/Dx1sM89Plu9pfVMmtYPL89byTjB/ezujSllJcVVdbz8Oe7+ceGAmLCgrhrfgbXzUglqBctc6+h4EFr9pXy4Ke72FJQwcgBUfz2vJGcPiJR+w2U8nPbD1XwwCc7Wb33CGnx4fzirBFcOH4QAb1gyQwNBQ/YVVjJg5/uYvnuEgbFhPLLszO5dFKyz62RopTqOWMMy3YX89fPdrOrsIqRA6L49TmZnDHS2gEnGgpulFNSzeNf7eGDzYeICgnk9nnDuX5WGqFB2omslOqYw2H4cMsh/vfLbPKP1DI1NZZfn5PJ9KHxltSjoeAG+4/U8vhXe/jXxgJCAm1cPyuNW08fSr/wYK+8vlKq72uyO3g76wB/+2oPRZUNnD4ikV+fk8nYZO8utqehcBIOltfx5NI9vJNVgC1AuG5GKreePozEqBCPvq5SynfVN9l5ZW0eTy/fR3ltE+eOGcDC+RmMHhTtldfXUOiBosp6nlq2lze/PQDAVacM4fZ5w+mvw0uVUm5SWd/E8ytz+b9VuVQ1NHP26P4snJ/h8ZaDhsIJOFhex6Kv9/HmdwewOwxXTB3CHWcMJ7lfmFtfRymlWlTUNvHi6lxeXJ1LVX0zZ47qz13zMxg32DPhoKHQDXmlNTyzfB/vbizAGLhscjJ3zMvQ5ayVUl5TUdfES6vzeGFVDpX1zcwfmcTC+RlMGNLPra+jodCJ7KIqnlq2lw83HyLQFsDV04Zw8+nDtGWglLJMZX0TL6/O4/lVuVTUNXH6iERunzecU9Lds65SrwsFEbkXuBIoch2ab4yxi8hk4AlAgEXGmJc6e56TCYUtBeU8uXQvX+woIjzYxk9npHLTnHSSorTPQCnVO1TVN/HK2nxeWJVLWU0j09Ji+Y+5w5mbeXITZHtrKOw1xrza7vhq4DrgILAOZ1gcPd7z9CQU7A7Dv7+SxdJdxUSHBnLDqencOCuN2AgdWqqU6p3qGu289d1+Fq3I4VBFPaMGRvPyjdN6vK5ab92O879E5FbgbWPM30QkBIgwxuQCiMhK4BTgc3e+qC1ASIuP4DfnjuS6GSl9bnMMpZT/CQu2ccOp6VwzPZX3Nx1kyc4iEiLdPyzeypZCPFAGhAIfAvcB+4DXjTFzXfe5D9hljHmj3WNvBm4GSElJmZKfn+/FypVSqu+zpKUgIjZgdQc3fWyMud/1e52IvAtMBb4D+rW5XwzO4DiGMWYRsAicp4/cWbNSSvkzj4aCMcYOzOjoNhHpZ4wpF2dPyVzgJWNMvYjUiEgKcBiYjbMFoZRSygus7FN4TEQycY4yWm6M+cR1/C7gDdfxpzvrZFZKKeVeloWCMeaG4xzPAk71bjVKKaUAes82QEoppSynoaCUUqqVhoJSSqlWGgpKKaVa9fkF8USkBOhrs9cSgFKri/Ayfc/+Qd9z35FqjElsf7DPh0JfJCJZHc0k9GX6nv2Dvue+T08fKaWUaqWhoJRSqpWGgjUWWV2ABfQ9+wd9z32c9ikopZRqpS0FpZRSrTQUlFJKtdJQsIiIvCwiS9pcP1dE1rou51hZm7uJyGMiss51+W2b4778nm8QkTUistq177jPEZFJrve3QkSWishQEQkVkddEZKXrp09ueC4iI0SkSURm+9p71j4FC4jIeOB+nFuPnunajGgjcJrrLl8Dk137UfR5IpJhjNkjIgE4N126DsjDR9+ziMQCX+HcSyQZWGyMmW1tVe4nIgOAGmNMlYj8CLga579vojHmfhH5A1BsjHnW0kI9QEQWAwOBe4Gx+NB71paCNe4BHmhzfTiQa4wpN8aU4/zAHG5BXR5hjNnj+ukAmgE7vv2eTwFWGmMaXfuNR7n2H/cpxphCY0yV62oDzn/b04GPXMc+dF33KSIyHSgEClyHfOo9ayh4mYjMBbKBojaH44G2mwmVA3FeK8pLRORaIMcYk4dvv2dffm8/ICIRwJ+Bhzj2vZfjm+/798CDba771Hu2cuc1n9XZ3tQ4NxC6imP3oi6jG3tT92Zd7cctImcCNwIXuo73+ffcCV9+b8cQkSDgLeB/jDE7RKTte/e59y0i5wNZxpgjzp2EgWP/vfv8e9ZQ8IDj7U0tIlHA5cCbQBgwRkRavnWki0i0667pwF4vlesWXezHPR1nH8p5xpg61+E99PH33IlvgD+7PjAHAtXGmAaLa3I7Vx/Rq8B7xpj3XIe/Bn4EbHL9/NqS4jxnIjBXRGYB44CRwKf40HvWjmaLiEga8Lwx5kzX9R/h7GsAuL/NntV9nohsc/3aspLkr4wx6338Pf8MWAAY4C7XNrM+RUR+DLwEtLy3rcB/AS8Cg3Gec7/RGFNvSYEeJiIvAc8D6/Gh96yhoJRSqpV2NCullGqloaCUUqqVhoJSSqlWGgpKKaVaaSgopZRqpaGglFKqlYaCUkqpVhoKSrUhIpeIiBGRke2OPyIim0XkiXbHw0Tka9cyH7ge+2qb2wNFpEREPmpzrLrdc9wgIk8ep55g19LUuvqA8goNBaWOdTXOGbpXtxwQkWHAqcaYCcaYO9vd/2fAu22W/K4BxopImOv6WcDBnhZjjGnEuQz3lT19DqVOhIaCUi4iEgnMxbk8xdWuY5nAciBVRDa6VgRt61rg/XbHPgHOd/1+NfBGN1//VhHZ5Lrkisgy103vuV5HKY/TUFDqexcDS4wxm4FqEZlijNkNvAzcY4yZZIypabmziAQDQ11Lgbf1JnCVaweu8TgXyGsrrM2H/ybgTwDGmGeNMROBaTjX0Plf1/23uY4p5XEaCkp972rgbdfvb/P9KaRxwOYO7p+Ac/38YxhjtgBprsd3tMhfnTFmYssF+EO72x8HlhpjPnQ9nx1odK2yq5RHaSgoBYhIHDAd+Mx16G3gSnEumj8G57f19uqA4+3H+wHwMN08ddSmjhuAVOC+djeFAH125U3Vd+iIBqWcfgx80rLvgTEmR0QO49xDuqnNPhCtjDFHRcQmIqEdLJX8IlBujNnq2m2vSyIyBfhPYI5r69KW4/FAqTGmqSdvTKkToaGglNPVwAQRyWtzLB7nxuwdtRJafAHMBpa0PWiMKQD+doI13IFzK8dlrl29sowxC4B5OHftU8rjdD8FpU6CiEwGfmGM+akHX+Nd4LfGmGxPvYZSLbRPQamTYIzZgPObvc0Tz+8a4fSeBoLyFm0pKKWUaqUtBaWUUq00FJRSSrXSUFBKKdVKQ0EppVQrDQWllFKtNBSUUkq1+v+8HBU316I8EAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "initsteady = (0.1,0.1)\n",
    "fstart = 50\n",
    "deltaarr = np.linspace(-fstart*2*np.pi,fstart*2*np.pi, 2000)\n",
    "etasteady = 0.00001\n",
    "losasteady = np.zeros(np.shape(deltaarr))\n",
    "losareal = np.zeros(np.shape(deltaarr))\n",
    "losaimag = np.zeros(np.shape(deltaarr))\n",
    "\n",
    "newinit = initsteady\n",
    "for i, delt in enumerate(deltaarr):\n",
    "    ar, ai = a_steady(etasteady, delt, gs, spins, gpar, gperp, kappa, newinit)\n",
    "    newinit = (ar, ai)\n",
    "    losasteady[i] = ai**2+ar**2\n",
    "    losareal[i] = ar\n",
    "    losaimag[i] = ai\n",
    "\n",
    "fig = plt.figure()\n",
    "plt.xlabel('$\\Delta f$ (MHz)')\n",
    "plt.ylabel('|A|² (dB)');\n",
    "plt.plot(deltaarr/2/np.pi, 10*np.log10(losasteady/0.0000000001))\n",
    "#deltaarr[find(losasteady == max(losasteady))]/2/pi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SOC pulses"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "pumptime = 20e-3#5*np.pi/gcoll\n",
    "decaytime = 12000e-3\n",
    "inittime = 2000e-3\n",
    "dt = 4e-3 #The time resolution of the AWG. Don't change that because the solution of the MBEs are written in a way that use that as time resolution\n",
    "numsteps = int(pumptime/dt)\n",
    "tlist = np.arange(0,pumptime+dt,dt)\n",
    "tlistdecay = np.arange(0, decaytime+dt, dt)\n",
    "tlistinit = np.arange(0, inittime+dt, dt)\n",
    "# test time it takes to calculate\n",
    "middle = int(nbins/2)\n",
    "\n",
    "init = np.ones(len(spins)*2+1)*(-1)\n",
    "init[0] = 0\n",
    "init[1::2] = 0\n",
    "init[2::2] = np.array([-1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -0.99999999,\n",
    "       -0.99999999, -0.99999999, -0.99999999, -0.99999999, -0.99999999,\n",
    "       -0.99999998, -0.99999998, -0.99999998, -0.99999997, -0.99999997,\n",
    "       -0.99999996, -0.99999995, -0.99999994, -0.99999993, -0.99999992,\n",
    "       -0.9999999 , -0.99999989, -0.99999986, -0.99999984, -0.99999981,\n",
    "       -0.99999977, -0.99999973, -0.99999968, -0.99999962, -0.99999954,\n",
    "       -0.99999946, -0.99999936, -0.99999924, -0.9999991 , -0.99999893,\n",
    "       -0.99999873, -0.9999985 , -0.99999823, -0.9999979 , -0.99999752,\n",
    "       -0.99999707, -0.99999654, -0.99999592, -0.99999519, -0.99999432,\n",
    "       -0.99999331, -0.99999213, -0.99999074, -0.99998912, -0.99998722,\n",
    "       -0.99998499, -0.9999824 , -0.99997938, -0.99997585, -0.99997175,\n",
    "       -0.99996698, -0.99996143, -0.999955  , -0.99994755, -0.99993892,\n",
    "       -0.99992894, -0.99991741, -0.9999041 , -0.99988877, -0.99987111,\n",
    "       -0.99985081, -0.9998275 , -0.99980075, -0.9997701 , -0.99973503,\n",
    "       -0.99969495, -0.99964919, -0.99959701, -0.99953759, -0.99947002,\n",
    "       -0.99939325, -0.99930617, -0.99920749, -0.99909582, -0.99896962,\n",
    "       -0.99882715, -0.99866655, -0.99848572, -0.99828239, -0.99805404,\n",
    "       -0.99779792, -0.99751103, -0.99719008, -0.99683149, -0.99643137,\n",
    "       -0.99598547, -0.99548921, -0.9949376 , -0.99432529, -0.99364646,\n",
    "       -0.99289488, -0.99206383, -0.99114611, -0.99013401, -0.98901929,\n",
    "       -0.98779315, -0.98644626, -0.98496866, -0.98334984, -0.98157865,\n",
    "       -0.97964335, -0.97753156, -0.97523029, -0.97272593, -0.97000425,\n",
    "       -0.96705043, -0.96384906, -0.96038416, -0.95663923, -0.95259724,\n",
    "       -0.94824074, -0.94355183, -0.93851226, -0.93310351, -0.9273068 ,\n",
    "       -0.92110321, -0.91447376, -0.90739948, -0.89986154, -0.89184133,\n",
    "       -0.88332057, -0.87428145, -0.86470673, -0.85457986, -0.84388515,\n",
    "       -0.83260784, -0.82073429, -0.80825208, -0.79515015, -0.78141894,\n",
    "       -0.7670505 , -0.7520386 , -0.73637889, -0.72006895, -0.70310842,\n",
    "       -0.6854991 , -0.66724497, -0.64835231, -0.62882971, -0.60868811,\n",
    "       -0.58794084, -0.5666036 , -0.54469444, -0.52223374, -0.49924418,\n",
    "       -0.47575063, -0.45178012, -0.42736172, -0.40252643, -0.37730708,\n",
    "       -0.35173818, -0.32585577, -0.2996973 , -0.27330143, -0.2467079 ,\n",
    "       -0.21995733, -0.19309109, -0.16615106, -0.13917954, -0.11221902,\n",
    "       -0.08531203, -0.05850101, -0.03182809, -0.00533503,  0.020937  ,\n",
    "        0.04694749,  0.07265674,  0.09802593,  0.12301726,  0.14759403,\n",
    "        0.17172067,  0.19536285,  0.21848749,  0.24106287,  0.26305856,\n",
    "        0.28444554,  0.30519615,  0.32528412,  0.34468455,  0.36337393,\n",
    "        0.3813301 ,  0.39853222,  0.41496079,  0.43059757,  0.44542558,\n",
    "        0.45942905,  0.47259341,  0.48490522,  0.49635216,  0.506923  ,\n",
    "        0.51660753,  0.52539655,  0.53328186,  0.5402562 ,  0.54631322,\n",
    "        0.55144747,  0.55565438,  0.55893025,  0.56127219,  0.56267818,\n",
    "        0.56314697,  0.56267818,  0.56127219,  0.55893025,  0.55565438,\n",
    "        0.55144747,  0.54631322,  0.5402562 ,  0.53328186,  0.52539655,\n",
    "        0.51660753,  0.506923  ,  0.49635216,  0.48490522,  0.47259341,\n",
    "        0.45942905,  0.44542558,  0.43059757,  0.41496079,  0.39853222,\n",
    "        0.3813301 ,  0.36337393,  0.34468455,  0.32528412,  0.30519615,\n",
    "        0.28444554,  0.26305856,  0.24106287,  0.21848749,  0.19536285,\n",
    "        0.17172067,  0.14759403,  0.12301726,  0.09802593,  0.07265674,\n",
    "        0.04694749,  0.020937  , -0.00533503, -0.03182809, -0.05850101,\n",
    "       -0.08531203, -0.11221902, -0.13917954, -0.16615106, -0.19309109,\n",
    "       -0.21995733, -0.2467079 , -0.27330143, -0.2996973 , -0.32585577,\n",
    "       -0.35173818, -0.37730708, -0.40252643, -0.42736172, -0.45178012,\n",
    "       -0.47575063, -0.49924418, -0.52223374, -0.54469444, -0.5666036 ,\n",
    "       -0.58794084, -0.60868811, -0.62882971, -0.64835231, -0.66724497,\n",
    "       -0.6854991 , -0.70310842, -0.72006895, -0.73637889, -0.7520386 ,\n",
    "       -0.7670505 , -0.78141894, -0.79515015, -0.80825208, -0.82073429,\n",
    "       -0.83260784, -0.84388515, -0.85457986, -0.86470673, -0.87428145,\n",
    "       -0.88332057, -0.89184133, -0.89986154, -0.90739948, -0.91447376,\n",
    "       -0.92110321, -0.9273068 , -0.93310351, -0.93851226, -0.94355183,\n",
    "       -0.94824074, -0.95259724, -0.95663923, -0.96038416, -0.96384906,\n",
    "       -0.96705043, -0.97000425, -0.97272593, -0.97523029, -0.97753156,\n",
    "       -0.97964335, -0.98157865, -0.98334984, -0.98496866, -0.98644626,\n",
    "       -0.98779315, -0.98901929, -0.99013401, -0.99114611, -0.99206383,\n",
    "       -0.99289488, -0.99364646, -0.99432529, -0.9949376 , -0.99548921,\n",
    "       -0.99598547, -0.99643137, -0.99683149, -0.99719008, -0.99751103,\n",
    "       -0.99779792, -0.99805404, -0.99828239, -0.99848572, -0.99866655,\n",
    "       -0.99882715, -0.99896962, -0.99909582, -0.99920749, -0.99930617,\n",
    "       -0.99939325, -0.99947002, -0.99953759, -0.99959701, -0.99964919,\n",
    "       -0.99969495, -0.99973503, -0.9997701 , -0.99980075, -0.9998275 ,\n",
    "       -0.99985081, -0.99987111, -0.99988877, -0.9999041 , -0.99991741,\n",
    "       -0.99992894, -0.99993892, -0.99994755, -0.999955  , -0.99996143,\n",
    "       -0.99996698, -0.99997175, -0.99997585, -0.99997938, -0.9999824 ,\n",
    "       -0.99998499, -0.99998722, -0.99998912, -0.99999074, -0.99999213,\n",
    "       -0.99999331, -0.99999432, -0.99999519, -0.99999592, -0.99999654,\n",
    "       -0.99999707, -0.99999752, -0.9999979 , -0.99999823, -0.9999985 ,\n",
    "       -0.99999873, -0.99999893, -0.9999991 , -0.99999924, -0.99999936,\n",
    "       -0.99999946, -0.99999954, -0.99999962, -0.99999968, -0.99999973,\n",
    "       -0.99999977, -0.99999981, -0.99999984, -0.99999986, -0.99999989,\n",
    "       -0.9999999 , -0.99999992, -0.99999993, -0.99999994, -0.99999995,\n",
    "       -0.99999996, -0.99999997, -0.99999997, -0.99999998, -0.99999998,\n",
    "       -0.99999998, -0.99999999, -0.99999999, -0.99999999, -0.99999999,\n",
    "       -0.99999999, -0.99999999, -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        , -1.        , -1.        , -1.        , -1.        ,\n",
    "       -1.        ])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA24AAAE8CAYAAAC4iOKyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA+X0lEQVR4nO3deZRcd3ng/e9TS6tlbZZs2ZhFls0SQsjYgCAEMIuHDGAmGxOSlwMkJBM8yWQGPHMGJhAIyQTOIZCZyQt5mcQTHHZCFkLAToCwGWEDRgabNWw2EIyNZVmyFquXqnreP+pWd3WrurskdXV13fp+zunTVffeqnquW93Xz31+v+cXmYkkSZIkaf2qDDsASZIkSdLyTNwkSZIkaZ0zcZMkSZKkdc7ETZIkSZLWORM3SZIkSVrnTNwkSZIkaZ2rDTuAjrPPPjt379497DAkSQN244033pWZO4cdx6jw+ihJ42O5a+S6Sdx2797Nvn37hh2GJGnAIuK7w45hlHh9lKTxsdw10qGSkiRJkrTOmbhJkiRJ0jpn4iZJ0imKiBdExPURcV1EPHLRvsdFxJciYioi7t+1fXdEfKx4zcvXPmpJ0igycZMk6RRExHbgRcCTgecBb1h0yFeAnwQ+s2j7a4FXZebjgUsj4qEDDlWSVAImbpIknZrHAHszcyYzbwW2RMSGzs7MvCczj/Z43cWZubd4fA3wpDWIVZI04kzcJEk6NWcBB7ueHwJ29PG67mtvz9dExOURsS8i9u3fv/90YpQklYSJmyRJp+Zu4Myu59uKbStprfSazLwyM/dk5p6dO13yTpJk4iZJ0qn6LPCEiKhHxC7gaGZO9/G6myPiccXjZwCfHFiEkqTSMHGTJOkUZOZB4E3AtcC7gSsi4uKIeAlARDwkIj4CXAS8OyJ+s3jpy4DXRMR1tOfIfW0I4UuSRkxt2AHo1H3zh0dI4CHnbhl2KJI0ljLzKuCqRZtvKvZ9A3hqj9fcAjxl4MF1+avP/QsETFQrTNQq1Ivv7efBRLVKvRZMVNv7NtQWHlerBBGxliFLkhYxcRthz3/zDdxxeIrvvPaZww5FkrSOveLvv8xMo7XygUuIoJ3EdSV8yyV6CxLDE7bFwue1+ffofEa9tvCz2se0X7dtY50tk/VV/K8jSaPBxG2E3XF4atghSJJGwKde+hSmGy1mmy1mmi1mG8lMs8lMI5lptpjp7Gu0TnxePJ7uet1s1+sWH39susF0sX12bl8uOPZ01CrBz158P171Mw9jqwmcpDFi4lYCU7NNJuvVYYchSVqnztk6OewQ5mQms83sSupaPRK9VpFo5oLkcabR4mt3HObtn/4u3zlwjHe/8LFM1JyuL2k8mLiNqO47lsdnTNwkSaMhItrz6moVNm1Y+fheLn7Ambz4L2/iPfv+hec/9vzVDVCS1ilvU42oe2cac4+PzzaHGIkkSWvrZy66L486fzt/vvcWMnPY4UjSmjBxG1HTXRW3KRM3SdIYiQh+cc/9+e6Be/nSbfcMOxxJWhMmbiNqwVBJEzdJ0pj5Nw+7DwCf/Mb+IUciSWvDxG1ELay4nV6HLkmSRs32TRM8+JzN3Pjdg8MORZLWhInbiJpuzFfZHCopSRpHjzp/O5//3iHnuUkaCyZuI2rGOW6SpDH3o+dt5Z7js+w/Mj3sUCRp4EzcRpRz3CRJ4+7CnZsA+Pb+Y0OORJIGz8RtRHXPcWs0HSIiSRo/F+7cDMAtdx0dciSSNHgmbiOqu+I227Q5iSRp/Jy3dZLJeoVbrbhJGgMDTdwi4lER8eGI+HhEvG6QnzVuZprdiZsVN0nS+KlUgvtu28jt90wNOxRJGrjaoN44IiaA1wLPyswjg/qccdXdVbLRsuImSRpP526d5I7DJm6Sym+QFbefBI4C74qIj0XEJQP8rLHTPVSy+7EkSePkPtsmucOKm6QxMLCKG3Bf4CLgYmAL8NGI+NHsWmwlIi4HLgfYtWvXAEMpn+5krdFyqKQkaTydu3WSO49M0WollUoMOxxJGphBVtzuBq7PzMOZeRtwF7Cz+4DMvDIz92Tmnp07d/Z8E/W2sKukFTdJ0ni6z9YNzDaTu++dGXYokjRQg0zcPgs8JCJqEbEFOAc4MMDPGyvdiduMzUkkSWPqnK2TAC7CLan0BjZUMjMPRcQbgU8AdeC/Z6YrRa+SzlDJaiWsuEmSxtaZZ9QBOGjFTVLJDXKOG5n5duDtg/yMcTXTbFGvBrVKxTlukqSxdebGCQDuuXd2yJFI0mANNHHT4Mw2WtSrFaqVsKukJGlsbd/UqbiZuEkqNxO3EdVoJbVKUKtWXMdNkjS2OhU3h0pKKrtBNifRADVaLWrVCvVqMNtwqKQkaTxtnKiyoVbhnuNW3CSVm4nbiGp2Km6VCrNW3CRJY+zMM+ocPGbFTVK5mbiNqEaznbhN1Co0XA5AkjTGtp8xwSErbpJKzsRtRDVaSbUa1CrBrMsBSJLG2LaNdbtKSio9E7cR1W5OUqFWrTBrxU2SNMY2b6hxdLox7DAkaaBM3EZUs9VqD5Wshl0lJUljbdOGGsdmTNwklZuJ24hqNJNqsRyAQyUlSeNs04Yax6y4SSo5E7cR1WgltWpQrQTNlkMlJUnja8ukQyUllZ+J24hqtJJqpUI1TNwkSeNt00SNqdkWDUegSCoxE7cR1Wy1qFeCWjVomLhJksbYpg1VAI7NNIcciSQNjonbiOrMcatWgpaJmyRpjG3eUANwnpukUjNxG1GdOW61ihU3SdJ422TiJmkMmLiNqM4ct4pz3CRpaCLiBRFxfURcFxGPXLRvMiLeGRF7i++Txfa3RMQXIuITEfHXw4m8XDoVNxuUSCozE7cR1T3HzcRNktZeRGwHXgQ8GXge8IZFh7wA+OfMvAT4evG84z9n5pMz89mDj7T85ituznGTVF4mbiNqfo5bxcRNkobjMcDezJzJzFuBLRGxoWv/k4Cri8cfKJ53/K+iEvdLaxRrqXWak1hxk1RmtWEHoFPjHDdJGrqzgINdzw8BO4Dbe+zv7AP4b5l5V0TsAD4aEZ/LzFu63zgiLgcuB9i1a9dAgi+TTRPt/505PmviJqm8rLiNqKZz3CRp2O4Gzux6vq3Y1mv/3L7MvKv4fjfwT8BFi984M6/MzD2ZuWfnzp2rHnjZTNbbFbepWddxk1ReJm4jqtGZ41YxcZOkIfks8ISIqEfELuBoZk537b8WuKx4fFnxnIg4s/g+ATwe+MaaRVxSk/X2/84cdx03SSXmUMkR1ezMcXMBbkkaisw8GBFvop2QJfDiiLgY+KnMfD3wFuCqiNgLfB/41eKl74mIzUAdeEdmfmXNgy+ZuYpbw8RNUnmZuI2o2WKOWzWCZsuhIZI0DJl5FXDVos03FfuOA8/p8ZqnDT6y8bKhViECpqy4SSoxh0qOqPYct6BqcxJJ0piLCDbUKkw1vJEpqbxM3EZUo9miVqlQqwQtEzdJ0pjbWK8yNWvFTVJ5mbiNqGYrqTnHTZIkoD3PzeYkksrMxG1EzbaS6twcNxM3SdJ4m6xXHSopqdRsTjKi5ipuETTTxE2SNN4mHSopqeSsuI2gzCwStwrVSoVMnOcmSRprk/WKiZukUjNxG0GdoZG1SlCrBoDz3CRJY22yZsVNUrmZuI2gTpJWrQaVaCduznOTJI2zjRNVjpu4SSqxgSZuEXE8Ij5RfP37QX7WOGl0V9wqReLmPDdJ0hhrD5W0OYmk8hp0c5LbMvPJA/6MsdNsdhK3ygnbJEkaRw6VlFR2gx4qeZ+IuDYi3hsRuwf8WWOj0WrfUaxVg2olFmyTJGkcTU6YuEkqt0FX3HZn5l0R8TTgzcC/HvDnjYW5OW6VoDNC0jlukqRxtqHmUElJ5TbQiltm3lV8/xBw/uL9EXF5ROyLiH379+8fZCil4hw3SZIWmqhVmGmauEkqr4ElbhGxOSKqxeN/Bdy1+JjMvDIz92Tmnp07dw4qlNLpnuM2N1TSOW6SpDG2oVphptEivZEpqaQGOVTyYcCfRcQRIIH/MMDPGivdc9w6NxcdKilJGmcTtfa96NlmMlGLIUcjSatvYIlbZt4APGJQ7z/Ouue4dThUUpI0zjqJ20yzNfdYksrEv2wjqNHsnuPW/hFacZMkjbOJapG4NZznJqmcBt1VUgPQbHWv49Z+7Bw3SdI4q88NlTRxk1ROJm4jqDPHrVoNOumaFTdJ0jiz4iap7EzcRlD3cgCdzM05bpKkcdaZ1zZt4iappEzcRlBnWOTCBbi9UEmSxteGmhU3SeVm4jaCOsMi69UKme0LlHPcJEnjrLurpCSVkV0lR9DcHLdKzC0J4Bw3SdI4qzvHTVLJWXEbQd3LAbQ6iZtz3CRJY6zTnMSukpLKyorbCOpegLtTcWtYcZMkjbEJ57hJKjkTtxHUPcdtbqikc9wkSWPMrpKSys7EbQT1nOPmUElJ0hjbYHMSSSVn4jaCuue41SqVBdskSRpHE9Uq4FBJSeVl4jaCmgvmuBXbrLhJksZYvdYegWJzEkllZeI2ghoL5ri1f4Qtm5NIksbYhMsBSCo5E7cR1Oye4xZ2lZQkya6SksrOxG0EzXbNcatW24mbFTdJ0jibsDmJpJIzcRtBC+a4WXGTJGluqKTLAUgqKxO3EdQ9x61icxJJkogI6tVwqKSk0jJxG0Hdc9w6ywE0HRoiSRpzE9WKXSUllZaJ2wjqzHGrxvxQSZdxk6S1FxEviIjrI+K6iHjkon2TEfHOiNhbfJ8stu+OiI8Vr3n5cCIvp4laxYqbpNIycRtBzVZSCah0NSfpVOEkSWsjIrYDLwKeDDwPeMOiQ14A/HNmXgJ8vXgO8FrgVZn5eODSiHjoWsQ7DkzcJJWZidsIarRybojkXMXN65QkrbXHAHszcyYzbwW2RMSGrv1PAq4uHn+geA5wcWbuLR5f07Vdp6lerdhVUlJpmbiNoGarRa2otFUrxXIANieRpLV2FnCw6/khYMcS+7v3dV97F78GgIi4PCL2RcS+/fv3r1K45Vd3jpukEqutdEBEPAR4CXB+9/GZeekA49IyGq2cS9g63xtOcpOktXY3cGbX823Ftl77u/e1lnkNAJl5JXAlwJ49e/wD36d6NbweSiqtFRM34K+BPwX+L9AcbDjqR7OV1IqErfjmcgCStPY+C7w6IurAecDRzJzu2n8tcBlwU/H92mL7zRHxuMy8HngGcMWaRVxytUqFhnO+JZVUP4lbIzP/z8AjUd/aFbf2SJuIoFoJm5NI0hrLzIMR8SbaCVkCL46Ii4GfyszXA28BroqIvcD3gV8tXvoy4M0RMQH8Y2Z+bc2DL6l6NeY6L0tS2fSTuH0gIv4j8HfA3J3EzDxhaIfWRrM5X3EDisRtiAFJ0pjKzKuAqxZtvqnYdxx4To/X3AI8ZeDBjaFa1YqbpPLqJ3H7leL7S7q2JXDh6oejfnTPcYN2Z0krbpKkcVerWHGTVF4rJm6ZecFaBKL+Nbq6SkL7QmXFTZI07urVCvfONIYdhiQNRD9dJevAbwJPLDZ9AvizzJwdYFxaxuKKW8U5bpIkUasGjZYVN0nl1M9Qyf8D1IE3Fc+fX2z79UEFpeUtnuNWq4RdJSVJY69WqThUUlJp9ZO4PTozL+p6/rGIuLnfDyjWgfsK8JTM/NTJBqgTNVpJrTK/fmu74uaFSpI03iZq4QLckkqrsvIhNCPigZ0nEXEhJ7ee2yuZX7tGq6DZc46biZskabzVKhUaJm6SSqqfittLgI9HxC1AAOczvxbNsiLiJ4A7cOHuVXXCHLdwTL8kSTXXcZNUYv10lfxoRDwY+JFi09czc3q513T5HdpJ3v/stTMiLgcuB9i1a1efb6lma9Ect2rQMnGTJI25esV13CSV15KJW0Rcmpkfi4hnLdr1oIggM9+73BtHxDOBfZl5ICJ6HpOZVwJXAuzZs8fMo0+91nGz4iZJGne1atCw4iappJaruD0J+Bjw0z32JbBs4gZcDDw5Ih4H/Djw0Ij4pcz87qkEqnnNVjJRr849r1aCll0lJUljrl6t2JxEUmktmbhl5quKh/8jM2/t3hcRKy7KnZmvAV5THP8W4M9N2lbHCRW3incYJUmqVRyBIqm8+ukq+bc9tv3NyXxIZr7ApQBWT7PVWjDHzYqbJElQq1a8kSmptJab4/ZQ4MeAbYvmuW0FJgcdmJbWaPaouHmHUZI05iaqwUyzRWay1Px6SRpVy81x+xHg3wJnsnCe2xHghQOMSStotnLBOm5V13GTJIlatT2QaPF1UpLKYLk5bn8P/H1E/GRmfnoNY9IK2nPc5ke5VsPETZKkTrLWaCW16goHS9KI6WcB7i9ExG/RHjY5N0QyM39tYFFpWY0ec9xM3CRJ465e3NScbbaYrJu5SSqXfpqTvB24D/A04Frg/rSHS2pIms00cZMkaZG5ipsNSiSVUD+J24My85XAscx8K/BM4CcGG5aW0+g1x82ukpKkMdeZ4zbbci03SeXTT+I2W3w/FBEPB7YB5wwuJK2k2WMdNytukrS8iHhWRHwzIu6JiMMRcSQiDg87Lq2eesWKm6Ty6meO25URsR14JfB+YDPwuwONSstqtJKazUkk6WS9DvjpzPzasAPRYHQqbiZukspoxcQtM/+8eHgtcOFgw1E/rLhJ0in5oUlbudWLaQQzTYdKSiqf5Rbg/q/LvTAz/9fqh6N+2FVSkk7Jvoh4D/A+YLqzMTPfO7SItKrqnYqbc9wkldByFbctaxaFTkrPipvNSSRpJVuBe4F/07UtARO3kqg5x01SiS23APfvr2Ug6l97jpsVN0k6GZn5q8OOQYPVqbjNOlRSUgmtOMctIv6C9h3JBVyAezharSQTqt3NSUzcJGlFEXF/4I3A44tNe4EXZ+b3hxeVVtPcOm5eEyWVUD/LAVwNXFN8fZT2UJOjgwxKS+tcjBas42ZXSUnqx1/Q7o583+LrA8U2lUSn47IVN0ll1E9Xyb/tfh4R7wY+NbCItKzOhOvuOW61qombJPVhZ2Z2J2pviYgrhhWMVl+nq6Rz3CSVUT8Vt8UejAtwD81cxa0rcatYcZOkfhyIiOdFRLX4eh5wYNhBafXU7CopqcT6meN2hPYctyi+3wH89wHHpSU0i7uICypudpWUpH78Gu05bv+b9vXsesCGJSXSuak5a8VNUgn1M1TSZQHWkfk5bvPF0kol5hI6SVJvmfld4GeGHYcGZ6LmHDdJ5bVi4gYQEc8CnkD7DuXezHzfIIPS0po9hkpacZOkpUXESzPzdRHxRnp3SX7REMLSALiOm6Qy62eo5JuABwHvLjb9RkT8VGb+1kAjU0+9mpNUKmHrY0la2teK7/uGGoUGznXcJJVZPxW3S4EfzWyXdCLircBXBhqVlrRUxa1l4iZJPWXmB4rvb+1si4gKsDkzDw8tMK0613GTVGb9dJX8FrCr6/kDim0ags7FqLviVg0rbpK0koh4V0RsjYhNwJeBr0bES4Ydl1ZPZx23hhU3SSXUT+K2BfhaRHwiIj4BfBXYGhHvj4j3DzQ6nWC+4jb/o6sWj626SdKyHlZU2H4O+EfgAuD5Q41Iq6qzjptdJSWVUT9DJX934FGob40eywF0Gkw2WslE13ZJ0gL1iKjTTtz+JDNnI+KU/w8/InYAbwO2ATcBL+pMK+g65unAq4qnv5eZH4qI3cAXgJuL7a/PzGtONQ7Ncx03SWXWz3IA10bEucCji003ZOadgw1LS+k1x22u4mZnSUlazp8C36GdMH0yIs4HTmeO20uB92Tm2yPiKuBpwAc7OyOiCrwOeGKx6dqI+Ejx+MbMfOppfLZ6mOsq6QgUSSW04lDJiPhF4Abg2cAvAp+NiF8YdGDqba6rZLV3xU2SdKKiGckPM/N+mXlZURn7HvCU03jbJwFXF48/UDzv9iDg1sw8lJmHaCeNDyr2XRQReyPibRFx1mnEoC4uByCpzPqZ4/Y7wKMz81cy85eBxwCvHGxYWspyFbemiZsk9ZSZLdoVsu5tmZmN03jbHcCh4vGh4nm3s4CDXc87x9wOXJiZlwCfol2VWyAiLo+IfRGxb//+/acR4nipziVuDpWUVD79zHGrLBoaeYD+Ej4NwGyvOW7FQxM3SVrWRyLivwHvAY51Nmbm3Uu9oBjueF2PXdfQTsq20U7ItgGL3+du4Myu59uAuzNzGpgutr0DOGFd1My8ErgSYM+ePf5x71NEUK8Gs14PJZVQP4nbByPiQ8wvwP1LwD8MLiQtp2dXyaoVN0nqwy8V37sTpQQuXOoFmdkEHttrX0RsBi4D3lV8f++iQ74JXBARW4vnFwDfiohtmXlPse1S4OsncxJaXq1SseImqZT6aU7ykoh4FvCEYtOVmfl3gw1LS+nMcatVF67jBiZukrSczLxgld/ydcDbIuI3gS8CHwaIiD8GXpOZ+yPiZcCHiuNflpnNiHhKRPwucASYAl64ynGNtVrVtU0llVM/FTeA64Em0AI+N7hwtJJec9w6j5t2lZSkJUXEGcB/BXZl5uUR8WDgRzLz6hVe2lNmHgCe2WP7FV2P/4FFo1Qy833A+07lM7WyWiVsTiKplPrpKvnrtLtK/jzwC8BnIuLX+njd1oi4vli4+4aI+NenH646dxG757hVOombFypJWs5fADPA44rntwGvHl44GoRateI6bpJKqZ+K20uARxR3FinaFl8PXLXC644CT8zMRkRcSHsy+KNXeI1W0HOOW/HQipskLeuBmflLEfEcgMy8NyJipRdptNQrMdfIS5LKpJ/E7QDtcfgdR4ptyypaL3dueW2lPf5fp6lXxW1+OQDvMErSMmYiYiPthiRExAOZ7+6okqhVbU4iqZz6Sdy+RXvR7b+nfbH7WeCLEfFfATLzfy31woi4H+1K20OAE4ZXRsTlwOUAu3btOungx1EnOVuwjttcc5KhhCRJo+L3gA8CD4iIdwKPB14wzIC0+mxOIqms+kncvl18dfx98X3LSi/MzNuAJ0TEbuATwNWL9rtOzUlq9FrHrbPgqBU3SVpSZn44Im6k3d4/gBdn5l1DDkurzOYkksqqn+UAfv9U3jgiNhSLjAIcZuFwS52iuTlu1RMTN/M2SVpaRHyA9ppr78/MYysdr9FUq9icRFI5rdhV8jQ8PCI+GREfp12lu2KAnzU2es1xczkASerLHwGXAF+NiL+JiF+IiMlhB6XVVa/anERSOfW7jttJy8wbgScO6v3HVa+uknPLAXiHUZKWlJnXAtdGRBW4lPbC11fRbqClknA5AElltWLiFhFndZYC0PDNFh1IelbcvE5J0rKKrpI/DfwS8EjgrcONSKut6hw3SSXVz1DJz0TEX0fEZa53M3zzFbeuBbjD5iSStJKI+Cvga7SrbX9Ce123/zzcqLTa6naVlFRS/QyVfAjwVNrt/N9QXPjekpnfGGhk6qnnHLeqzUkkqQ9vBp6Tmc1hB6LBqVUqNJqNYYchSauun66SCfwT8E8R8RTgHcB/jIibgd/OzE8POEZ16VTc6tWuOW5W3CRpRZn5oYh4XLFETa1r+9uGF5VWm81JJJVVX3PcgOcBzwd+CPxn4P3AxcBfAxcMMD4t0qm4dRXc5oZNtuwqKUlLioi3Aw8EbgI6VbcETNxKxOUAJJVVP0MlPw28Hfi5zPx+1/Z9EfGngwlLS2m2WlQrQfd0w7kFuL3DKEnL2QM8rBhJopKqOsdNUkn105zkFZn5B91JW0Q8GyAz/3Bgkamn2WZSry7sEVO14iZJ/fgycJ9hB6HBqttVUlJJ9VNx+23grxZtexntYZJaY7PNFvXKwnx7ruLmHUZJWs7ZtBffvgGY7mzMzJ8ZXkhabbVqhYbr40gqoSUTt4h4BnAZcL+IeEPXrq2A7ZqGpNHMuS6SHdW5ddxM3CRpGb837AA0ePVqMOv1UFIJLVdx+wGwD/gZ4Mau7UeA/zLIoLS0Rqu1oKMkQDVM3CRpJZl57bBj0OC1lwOw4iapfJZM3DLzZuDmiHhnZlphWydmGnli4mbFTZKWFBFHaHePPGEX7VVvtq5xSBqgasXmJJLKabmhkn+Vmb8IfCEiTvgLmJn/aqCRqadGq+VQSUk6CZm5ZdgxaO3UqzYnkVROyw2VfHHx/d+uRSDqT6OZc+u2dXSeN+0qKUkac7Wq67hJKqflhkreXjz8d8BfZuYP1iYkLWemeeIct4oVN0mSgPZyALPNJDMXrHkqSaOun3XctgD/FBF7I+I/RcS5gw5KS2v0SNxqJm6SJAHtiht4TZRUPismbpn5+5n5Y8BvAecB10bERwYemXpqtE5cDsCKmyRJba5tKqms+qm4ddwJ3AEcAM4ZTDhaSa8FuK24SZLUVq+auEkqpxUTt4j4jxHxCeCjwFnAC+0oOTyzzaReW1RxCy9SkiRBex03wLXcJJXOcl0lOx4AXJGZNw04FvWh0WxR27Dwx9YZFtIycZMkjblOxW3WJQEklcyKiVtmvgwgIs4BJru2f2+AcWkJs82cuyh1VK24SZIEzDcncUkASWXTz1DJn46IbwK3AtcC3wH+ccBxaQmNVu/lACKg5TpukqQxN9ecxIqbpJLppznJq4HHAt/IzAuAfw18ZqBRaUmzzZy7m9itVgkrbpKksWdzEkll1U/iNpuZB4BKRFQy8+PAngHHpSW0u0qeuKBoJcI5bpKksWdzEkll1U9zkkMRsRn4JPDOiLgTODbYsLSURvPEddygXXFzOQBJ0rizOYmksuqn4vazwL3AfwE+CHwb+OlBBqWlzTZPnOMG7XluDguRJI27TsXNm5mSymbJxC0iHhQRj8/MY5nZysxGZr4V+Dxw5ppFqAWWStxqlbA5iSRp7FU7FTe7SkoqmeUqbn8MHO6x/Z5in4ag0UpqPea4Va24SZJEfW6Om9dESeWyXOJ2bmZ+afHGYtvugUWkZTWW6CpZrdicRJLWUkTsiIirI2JvRLwxIk64qxYRV0XE7RHx54u2vzwirouIj0XE7jULegx05oHbnERS2SyXuJ25zL6NqxyH+pCZzDRbTPRoTlINK26StMZeCrwnMy8BNgFP63HMK4HndG+IiIcCl2bm44HfA1474DjHylxzEq+JkkpmucRtX0S8cPHGiPh14MbBhaSldCZa96y4Va24SdIaexJwdfH4A8XzBTLztiVed02x/5PARYMKcBxV55qTWHGTVC7LLQdwBfB3EfFc5hO1PcAE8PMrvXFEPAL4E6AJNIBfz8xbTivaMdeYS9ysuEnSOrADOFQ8PlQ878dZwA+6nlcXHxARlwOXA+zateuUAxxHnXngLgcgqWyWTNwy84fA4yLiKcDDi83XZObH+nzv24GnZ+aRiLgM+H3g+acV7ZibLcbrTywxx61pV0lJWlURUQWu67HrGuAgsI120rYNuLvPt72bhdMRmosPyMwrgSsB9uzZ4x/3k9DpvGxzEklls+IC3Jn5ceDjJ/vGmXlH19Np2lU3nYbO3cOluko2vUhJ0qrKzCbw2F77ImIzcBnwruL7e/t822tpd2f+44h4HHDz6UeqjrnmJA6VlFQy/SzAfVoiYhPwauD1PfZdHhH7ImLf/v37Bx3KyOt0yOrdVbJixU2S1tbrgOdGxF5gFvgwQET8cUTsLB6/GvifwDMi4iMRsSkzvwZ8KiKuA14DvHw44ZdTZzkAh0pKKpsVK26nIyLqwHuAP8zMry7e71CQk9PpkFXvNcetMt+8RJI0eJl5AHhmj+1XdD1+BfCKHsf8AfAHg4xvXHUW4LY5iaSyGVjFLSIqwDuA92Xm+wb1OeNkttG+CNWXqriZuEmSxlzd5iSSSmqQQyWfRftO5PMi4hMR8cYBftZY6IzX7zlUMqy4SZJUm2tOYsVNUrkMbKhkZv4N8DeDev9x1Ll7WO/RnKRmxU2SpK7mJF4TJZXLwJuTaPV0Whv3qrhVnOMmSZLNSSSVlonbCJlptpf66dWcpGZXSUmSqFZsTiKpnEzcRsh00ZxkQ616wr5KJRwWIkkae52bm1bcJJWNidsImSkSt4naiT+2WiVombhJksZcRFCthAtwSyodE7cRMjNXcesxxy2suEmSBO2bmQ0rbpJKxsRthMw0rbhJkrSSerXiUElJpWPiNkLmhkr2XIDbYSGSJEH7mmhzEkllY+I2Qpab41apBBbcJElqNyiZ9aIoqWRM3EbISkMlrbhJktReIqfR9JooqVxM3EbIshW3CMzbJEmCWtXmJJLKx8RthEwvM8fNipskSW31asWhkpJKx8RthCzbnKQaNL1ISZJkcxJJpWTiNkJmmi0mqhUqlThhX70Stj6WJIn2KBSviZLKxsRthMw0Wj3nt0F7WIgTsSVJ8pooqZxM3EbIcolbzcVGJUkCiuYkTh+QVDImbiNkptHqOb8NOmvWtMj0QiVJGm/1SoVZK26SSsbEbYRMN5rLDpXMxAYlkqSx125O4vVQUrmYuI2QmeZyQyXbDUscGiJJGne1qs1JJJWPidsIWXaoZKW93aEhkqRxV69WXNtUUumYuI2Q6WW7SrYrbt5hlCSNu1olaHg9lFQyJm4jZKWukoDtjyVJY69dcTNxk1QuJm4jZKbZYsMSiVtnCOWsFypJ0pirVsIbmZJKx8RthMw0lk7cOs1JZhteqCRJ483mJJLKyMRthByfbTJZr/bcNzdU0snYkqQxV6/0bk5y43cP8p7Pfc81TyWNpNqwA1D/pmdbSyZuEzYnkSQJaFfcFjcnmWm0eMFf3MCRqQZbJ+s848fPG1J0knRqrLiNkOOzTTYuVXFzOQBJkoCiq+SiOd9f+N5Bjkw1APjbz39/GGFJ0mkxcRshx2eabJxYaqikFTdJkqA9fWBxc5Kv3X4YgKf+6Ll85pa7adrMS9KIMXEbEZm57By3ua6SVtwkSWOuVo0Tuix/7fYjbD+jzk9fdB5HpxtziZwkjQoTtxExXXSLnKyvtI6bdxAlSeOtXjmx4vYvB+9l99mbePTuHQB87jt3DyM0STplJm4jYmq2CbD0HLfOUEm7SkqSxlytGrQSWl1Vtx8enuI+Wye575kbue+2SW787sEhRihJJ8/EbUQcXyFxmxsq6TpukqQxV6u0b2Z2Nyi58/A0526dBOAR52/nC987NIzQJOmUDSxxi4gPRcT+iHjFoD5jnByfaSduS6/jduJFSpKkcbR4bdNj0w2OTDfmErdH7trObYeOc8c9U0OLUZJO1iArbv8eeMkA33+sdCpuSyVudZuTSNKaiogdEXF1ROyNiDdGRPQ45qqIuD0i/rxr2+6IOBgRnyi+nrm2kZdfp+LW6bT8w8PtBO3crRsAeNT52wH4/PccLilpdAwscctMF0lZRVOz7YRsqeUA6nPruFlxk6Q18lLgPZl5CbAJeFqPY14JPKfH9hsz88nF1zWDDHIc1ecadrWvnXcemQaYq7g97LytbKhV+Lzz3CSNkKHOcYuIyyNiX0Ts279//zBDWfc6zUkma0t1lSyGSlpxk6S18iTg6uLxB4rnC2TmbUu89qKiUve2iDhrUAGOq8XTBw7dOwvAmWfUAZioVfjx+23jRitukkbIUBO3zLwyM/dk5p6dO3cOM5R1rzPHbcmKW2eopHPcJGmt7AAOFY8PFc/7cTtwYVGp+xTwusUHeGPz9CxuTnJkqp24bZ2szx3zyPO385XbDjPdaK59gJJ0CmrDDkD9WamrZL2zHIBdJSVp1UREFbiux65rgIPANtpJ2zagr4XBMnMamC6evgP4rR7HXAlcCbBnzx7vyJ2kWmXhUMnDUw0AtkzO/2/Po3fv4MpP3sINt97NJQ+ev3k8Ndtk7zfv4nt338v2M+pccPYmHnTOZrZ0JX2SNAwDS9wi4v8CjwM2RMSezPy5QX3WOJhaoTnJ4g5akqTTl5lN4LG99kXEZuAy4F3F9/f2854RsS0z7ymeXgp8fRVCVZe5tU2bCytumzfM/2/PJQ8+m80balx98+1ziduXb7uH//Suz/OdA/ee8J7nbt3Ajk0bmKhVqFWCSkBEEEAlgkql/b1WCWrVytz3eiWoVbsfd/YFtUqFenV+W71aoVqJ9rZKhVq197a5fZWufcu8b60SVCtBj/45kkbIwBK3zHzhoN57HB2bbt8tPGPJoZILL1KSpIF7HfC2iPhN4IvAhwEi4o+B12Tm/oh4NfAM4D4R8RHgZ4GnRMTvAkeAKcDr5SqrL7qZeWSqwaaJ6txNTmjfCH36w+/D+2/+AS9+6oP56D/fyR9c/VV2nDHBm39lD4/ctZ0Dx6a5Zf8xvr3/GN+68yj3HJ9lttmi0WqRCa1MWgnNVjLbTJqZNJpJo5U0mi0arWwf30warRazzfntjVbSXOPpDfVqO4GrF4lfdS7Bm08Ca9ViW6UrUSySzmqRXNYWvU8nuey8T626MIFdNvnsxNE5vuuz5+Ktdr2m8znFe5uQapw4VHJEHCvmuG2e7P0jm+8qacVNktZCZh4ATmjln5lXdD1+BbB4PdP3FV8akLk5bsXNzMPHZ3sOdXzRpQ/mAzf/gCe9/uPMNpMnPmQn//sXL+Ksze1lA7ZvmuBB52wZWJytIoHrTuqarWS2SPxmm1kkha0FyWCjmcy22glhs/Pa1vzxC15bHNedMM4nkwsTzE4C2mi15j53ptHi2EyTZvE+s10xNLq2dce91tPta11JY7XaI+lbVOmcT1zn9y9OPqvVWDKZ7K5kdhLduSR4URW0k/DOJ8Tz79crCe7eb0KqxUzcRsSRqQYT1Qobar0rbpVi2EbDipskacwt7ip5ZKqxYH5bx66zzuBdL3wsf/W5f2HP7u38u0fen0pl7f5nuVIJJirBxHB7xa26xQnpXDK5OPksHs8ln13JYPdrl0tcF29fPgmeT1Jnmy3unen+vPnj5hPX+bibQ6iQdiqKK1Uc64uSy6WH2C49ZLe2ILkMql3HdZLgFRPiuYS2U6ntToJNSFeDiduIODo9u2S1raNerVhxkySNvcXNSY5Mz7J1Y+/mIo86f/vcgtxaHWVPSHtWMhckk0skn835hLZXpXJBAtnqSj4779l13OyChHdh4np89sSEdHFi250Qr3VCWq0ESw+hXWmI7aIqZ9cc0oVV0UVDcBfNNV0pIV5Y/ZyfU7ogxkXzTNciITVxGxFHpxoLJlX30k7crLhJksbb4uYkh483OGvzxDBDUgl0ElKAjfQeATWKMnNB8tnsTkK7ksb+h9h2VSrnjuserttJXIvkc6Xkstg21VhcqV0+IW4MISH9Hz/7Yzz3J84f2GeYuI2Io9P9JG5hxU2SNPZObE4yy+6zNw0zJGndimhXjOrV8iakc0nd3NzQE+drLjfP9ITmQl1zRruH2D7svK0DPScTtxFxdLqx4lDJiZpDJSVJWtyc5MhUg60rXEMllUt3QkpJEtJyDT4usX4qbhtq1bn13iRJGldzc9xaSWZyeKp3V0lJGiUmbiOinzluG2oVphtW3CRJ422uq2SzxXSjPYypV1dJSRolJm4jop+hkpN1K26SJNU7zUla7WobsGRXSUkaFSZuI6K/oZJW3CRJ6l4O4PDxBoBz3CSNPBO3ETA122RqtsW2Fe4WTtarJm6SpLE3P1QyOVJU3BwqKWnUmbiNgIP3zgCwY9Pya9BsqFUcKilJGnvdzUmOTHUqbg6VlDTaTNxGwN3H2onb9jOWT9ysuEmS1FVxa7Xm5rjZVVLSqDNxGwEHj7UvOlbcJElaWb2ouM025ytuDpWUNOr8KzYC7p4bKrn83cINVtwkSVqwHMCRqfZ10a6SkkadFbcRcLDPoZIbahWmF1XcPnPLAf7g6q9y4Oj0wOKTJGk9qVY6QyWTw8cbVAI2TVSHHJUknR4rbiPg7mMzRNBXV8mprorbTKPFi979Be48Ms3h47O8/tkXDTpUSZKGbqLavi8902hxZGqWzRtqRMSQo5Kk02PFbQQcODbNto11atXlf1wbahVmGi0yE4DPf+8gdx6Z5qxNE1z9xdud/yZJGguVSrTnfTeaHJlq2JhEUimYuI2AO+6Z4j5bJ1c8brLeHgbSmed2/bcPUAn4g597OMdnm1z/7bsGGqckSevFxokqUzNNDk/NOr9NUimYuI2AHxya4r5nblzxuA219o+zU1n7+h2HuXDnZi596DnUq8ENtx4caJySJK0XG+tVjs82OTzVsKOkpFIwcRsBt99znPO2rVxx27ShXXG7d6aduN161zEuOHsTk/UqD7/fNvZ95+6BxilJ0nrRTtxaHJlqsNXETVIJmLitc8dnmhy8d7avitvmDe2hIEenGzRbyXcO3MuFZ28C4NG7d/DF79/jPDdJ0liYrFc5PtPk8PFZtjrHTVIJmLitc7cdOg7QV8Vtc3FH8chUgx8cOs5Mo8UFReK25/ztzDRbfPH79wwuWEmS1omNE1WmZpscmZp1qKSkUjBxW+du2X8UYC4BW87mDe0L09HpBrfcdQyAC3duBmDP7h0A7PuuwyUlSeU3Wa9wbKbBkenGisvpSNIoMHFb5755Zztxe/C5W1Y8dstcxW2WWxclfDs2TfCgczbzuVtN3CRJ5bexXuXOw9NkYldJSaVg4rbOffOHR7jvtsm5atpy5ipuUw1uvesYWzbUOHvzxNz+R+/ezr7vHqTVyoHFK0nSejBZr/LDw1MAVtwklYKJ2zr3lR8c5kfus3K1DebnuHWGSl6wcxMRMbd/z/k7ODLV4Bt3HhlIrJIkrRcb61UaxY1KEzdJZWDito7dfWyGb955lEdfsKOv4zdNzDcn6SwF0O3RxTy3z33H9dwkSeW2caI699jETVIZmLitY5+95QAAj9ndX+JWrQSbJqrcdXSa2w4dPyFxe8COjZy7dQOf+faBVY9VkqT1ZGN9PnE784yJZY6UpNFg4raOXf2l29mxaYKLHnBm3685e8sGbrj1bjLhgUVHyY6I4NKHnsvHv36n67lJkkqtuyGJFTdJZTDQxC0iXhAR10fEdRHxyEF+VtncetcxPvTlO/jZi+9Lvdr/j+m8bZNznSgfdM7mE/Y/88fP496ZJh/75ztXLVZJktabHZvmq2zbN5m4SRp9A1uRMiK2Ay8CHgvcD3g78IRBfV5ZfPaWA7z5U7fylR8cZuNEld940gNP6vXnbdsIQCV6r/322At3cL8zN/Jn136bZzz8Pgual0iSVBbbu4ZHbqhVlzlSkkbDwBI34DHA3sycAW6NiC0RsSEzpwfxYd+68+jcYtWLm93ngg25xPaFr+vel4veceG+7u1Lt9k/8bN6x/H6D32d2++Z4hG7zuQNz3kE526dXPI9e3nAjjMAuO+ZG5msn3ihqlUrvPipD+alf/NFXv53X+IRD9jOl267hy/8y0G+d+BeJmpVzj/rDC44exP3O3Mj9z1zkjMmatQqQaUSnGyaZ2IolcvGepUnPPjsYYchrWjnlg3DDkGSVtUgE7ezgO72hYeAHcDtnQ0RcTlwOcCuXbtO68Pef9NtvOFj3zqt91gvrnjqg7niqQ85pdc+9UfP4Q0f/SbPftQDljzm2Y+6P/98+xGuuu5W3n3Dv7BpospFDziTn3vE/ZiebXHrgWN88hv72X90+oSEU9J4O/+sM7j2JU8ZdhjSih5aLKXz/MeeP+RIJGl1xHJVotN644inA0/PzCuK5zcBP7FUxW3Pnj25b9++U/68Hx6eYv+R+bdeXOjprhV17+v7uEWft/B1seS+WLB90XE9XlOJ4H5nbqRSOfVK1V1Hp9lxxsSK73Hn4SnunWnygB1nUO1x7EyjxZ1HppiabdJoJY2mWZw07iZqFR5ybn9rSy4lIm7MzD2rFFLpne71cZztPzLN1o01h0pKGhnLXSMHWXH7LPDqiKgD5wFHBzVMEuDcrZMnPaywrM7e3N/wkHNW+O81Uatw/+1nrEZIkiStOYdLSiqTgSVumXkwIt4EXEt7KtiLB/VZkiSttYjYAbwN2AbcBLwou4axRMT9gXfS7uBcAV6cmfuiPfziDcDFwD3AL2fm3WsbvSRp1Ax0OYDMvCozH5eZj89Mx3lIksrkpcB7MvMSYBPwtEX7jwDPLva/EPjfxfanAWcU2/+qeB9JkpblAtySJJ2aJwFXF48/UDyfk5n3ZGZn0cxpoNHP6yRJ6sXETZKkU7ODdsdkmO+cfIKIqNIeGvmaYlN31+VDwPYer7k8IvZFxL79+/evXsSSpJE1yOYkkiSNtCLpuq7HrmtoJ1/baCdf24Cl5qn9GfCPmfmR4vndwJnF420sXDoHgMy8ErgS2l0lTy16SVKZmLhJkrSEzGwCj+21LyI2A5cB7yq+v7fHMX8E3J6Zf9K1+Vrg54H3Fa+7dnWjliSVkUMlJUk6Na8DnhsRe4FZ4MMAEfHHEbEzIvbQ7qh8SUR8IiL+unjdh4DZ4nXPBV4/hNglSSPGipskSacgMw8Az+yx/Yri4X6g3mN/C/itgQYnSSodK26SJEmStM6ZuEmSJEnSOheZ66NZVUTsB757mm9zNnDXKoQzbGU4D89h/SjDeXgO68dqnMf5mblzNYIZB6t0fYRy/Bv0HNaPMpyH57B+lOE8VusclrxGrpvEbTVExL7M3DPsOE5XGc7Dc1g/ynAensP6UZbzGEdl+Nl5DutHGc7Dc1g/ynAea3EODpWUJEmSpHXOxE2SJEmS1rmyJW5XDjuAVVKG8/Ac1o8ynIfnsH6U5TzGURl+dp7D+lGG8/Ac1o8ynMfAz6FUc9wkSZIkqYzKVnGTJEmSpNIpTeIWES+IiOsj4rqIeOSw41ksIh5RxPbJiPhYRFwYEZMR8c6I2Ft8nyyO3V0cc11EvLzrPZ4eEZ8uvp42xHN5SETMRsQTRvEcIuJREfHhiPh4RLwu2t5YnMPVEbGjOG5H8XxvsT+K7Y8szuv6iHjBkM4hIuJPiv+On4uI54zKeUTEhyJif0S8outcTjvutfwb0OMcfjkibih+v/8yIjYU20/q9yAiXl4c+7GI2L2W59C1/VcjYrbr+bo9B/VvLX8/TlZ4fVxv5zDS18giXq+PQ7o+LnEeXiNX6xwyc+S/gO3A54EJ4ALgU8OOqUeM9wG2FI8vA94O/AbwymLb7wK/UTz+S+CS4vFHgIcCVeCLwJnF181AdUjn8vYirieM2jkU/0b+qfOzKLY9HXhz8fiXgdcWj18LPL94fBXw9OLxdcW/s4ni3932IfwMHg58vHi8Bfj2qJwHcH/gBcArVuu//1r/DehxDhd2/i0DrwP+/cn+HhT7PlIc+0TgL9fyHIptk8DVwLe7tq3bc/Cr75/1ur5G4vVx3ZwDJbhG4vVxqNfHJc7Da+QqnUNZKm6PAfZm5kxm3gps6WTz60Vm3pGZR4qn00ADeBLtfwAAHyieA1ycmXuLx9cU2x8E3JqZhzLzEPCdYtuaioifAO4Avl9sGrVz+EngKPCu4m7HJSx9DidsL/5dbcrMWzNzBthL+9/fWvsBMBMRddoXprt7xVs8XlfnkZnfX7RpNeJe078Bi88hM2/JzGbxtPP7DSf3e/Ck4hgy85PARYOKv9c5FF4E/CnQPfl53Z6D+raur5FeH9fPOVCOa6TXxyFeH3udh9fI1TuHsiRuZwEHu54fAnYMJ5TlRcQm4NXA61kY9yHmY+7+uXS2r5dz/B3ad3o6Ru0c7kv7F+W5wPOB/0t7pfvuc9hePN5RPO9s75xDZ1v39rV2EPgm8A3gJtr/phb/LEbhPGB14l4P/7aIiIfSvkP6nmLTyfweLN5eHVScvUTEduCJmXn1ol0jcw5a0rr4/ViJ18d1cQ5luEZ6fTxx+3r4t+U1chXOoXayL1in7qZdiuzYVmxbV4q7P+8B/jAzvxoR3XF3x9zqelln+9DPMSKeCezLzAPFcGoWxbXuz6H4vOsz8zBwOCLuov2L04lrG/O/VAeL54dYX+cA8FPA/WjfwdlG+87ahxm984AT/w2dStyxxPY1ExH3B94K/D+ZOVVsPpnfg8Xbm6ytl9EewrLYKJ2DeltPv+89eX0c/jkUynCN9Po4f/y6uD6C18gltp/0OZSl4vZZ4AkRUY+IXcDRzJwedlDdIqICvAN4X2a+r9h8Le3x/BTfry0e3xwRjysePwP4JO27RxdExNaI2Ep7nPK31iL2LhcDT46ID9L+w/hHwNcYrXP4LPCQiKhFxBbgHOBv6X0OJ/x8ij82xyJiV/E/Gk8Abliz6OcFcLAYenCE9tj1jyyOt3i8ns+jZ3xLbV8m7qH+DYiIs2n/O/qNzPx2166T+T24tjiG4jU3r1X8hYcALy9+v8+LiM4d0VE6B/W2rq+RXh/XzTlAOa6RXh/X0fURvEaymudwspPi1usX8GvA9bQnZu4Zdjw94vsF2uPGP1F8vRHYCLyb9t2gdwOTxbEXAh8vzqV7UuRlwKeLr8uGfD5vof1HYeTOgfbwj0/R/mP287RvYPx/xTlcA5xVHHdW8Xxvsb9SbN9TnNf1wK8N6Ryqxc/gU8DnaI+7HonzoD305iu0/4i9b7XiXsu/AT3O4U9oz2vp/H53Jl6f1O8B8Mri2I8DF67lOSza962ux+v2HPw6qZ/3ur1G4vVxXZ0DI36NxOvjUK+PS5yH18hVOgcX4JYkSZKkda4sQyUlSZIkqbRM3CRJkiRpnTNxkyRJkqR1zsRNkiRJktY5EzdJkiRJWudM3KQVRMRZEXFT8XVHRNxWPD4aEW8a0GdeERG/fAqvm4iIT0ZEbRBxSZLU4fVRWlsuByCdhIj4PdqLV/7RAD+jBnweeGRmNk7h9a+ivcbIO1c9OEmSevD6KA2eFTfpFEXEkyPi6uLx70XEWyNib0R8NyKeFRGvi4gvRcQHI6JeHPeoiLg2Im6MiA9FxHk93vpS4POdi1JEfCIi/jAiboiIb0TEJcX2Hyu23RQRX4yIBxevfx/w3IH/B5AkqQevj9JgmLhJq+eBtC8qPwO8A/h4Zv44cBx4ZnFxeiPwC5n5KOAq4DU93ufxwI2LttUy8zHAFcCrim2/Afy/mXkxsAf4frH9y8CjV+mcJEk6XV4fpVXgOF9p9fxjZs5GxJeAKvDBYvuXgN3AjwAPB/4pIiiOub3H+5wHfG3RtvcW328s3gvg08DvRMT9gfdm5jcBMrMZETMRsSUzj6zGiUmSdBq8PkqrwMRNWj3TAJnZiojZnJ9A2qL9uxbAVzLzJ1d4n+PAZK/3BprFe5GZ74qIzwLPBP4hIv5DZn6sOG4DMHVaZyNJ0urw+iitAodKSmvn68DOiPhJgIioR8SP9Tjua8CDVnqziLgQuCUz3wD8PfCviu1nAXdl5uyqRS5J0uB4fZT6YOImrZHMnAF+AfjDiLgZuAl4XI9D/xF4Yh9v+YvAlyPiJtpDTN5WbH8KcM3pxitJ0lrw+ij1x+UApHUoIv4OeGlnXP5Jvva9wG9n5jdWPzJJkobH66PGmRU3aX36bdqTsE9KREwA7/OiJEkqKa+PGltW3CRJkiRpnbPiJkmSJEnrnImbJEmSJK1zJm6SJEmStM6ZuEmSJEnSOmfiJkmSJEnrnImbJEmSJK1z/z/Uy9Zf94VXegAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x360 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "alos = np.zeros((len(tlist)+len(tlistdecay)))\n",
    "szlos = np.zeros(np.shape(alos))\n",
    "\n",
    "asolinit, asolpump, asoldecay, infodict = solve_mbes_wait(mbes, init, pdf, tlistinit, tlist, tlistdecay,1e1, kappa, deltac, gs, gperp, spins, gpar)\n",
    "alos = np.hstack((abs(asolinit[:,0])**2, abs(asolpump[:,0])**2, abs(asoldecay[:,0])**2))\n",
    "szlos = np.hstack((np.real(np.sum(asolinit[:, 2::2], axis=1)), np.real(np.sum(asolpump[:, 2::2], axis=1)), np.real(np.sum(asoldecay[:, 2::2], axis=1))))\n",
    "                  \n",
    "tplot = np.linspace(0, 1000*(inittime+pumptime+decaytime), len(tlistinit)+len(tlist)+len(tlistdecay))\n",
    "\n",
    "fig, (ax1, ax2) = plt.subplots(1,2, figsize=(15,5))\n",
    "\n",
    "ax1.plot(tplot, (alos))\n",
    "ax1.set_xlabel('Time (ns)')\n",
    "ax1.set_ylabel('Cavity population');\n",
    "\n",
    "ax2.plot(tplot, szlos)\n",
    "ax2.set_xlabel('Time (ns)')\n",
    "ax2.set_ylabel('Inversion');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [],
   "source": [
    "%qtconsole"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}