HajimeKawahara/exojax

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Optimization of a Voigt profile using JAXopt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:04.542522Z",
     "iopub.status.busy": "2022-10-20T05:48:04.539874Z",
     "iopub.status.idle": "2022-10-20T05:48:06.557642Z",
     "shell.execute_reply": "2022-10-20T05:48:06.557320Z"
    }
   },
   "outputs": [],
   "source": [
    "from exojax.spec.lpf import voigt\n",
    "import jax.numpy as jnp\n",
    "import matplotlib.pyplot as plt\n",
    "import jaxopt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's optimize the Voigt function $V(\\nu, \\beta, \\gamma_L)$ using exojax!\n",
    "$V(\\nu, \\beta, \\gamma_L)$ is a convolution of a Gaussian with a STD of $\\beta$ and a Lorentian with a gamma parameter of $\\gamma_L$. \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:06.568609Z",
     "iopub.status.busy": "2022-10-20T05:48:06.568321Z",
     "iopub.status.idle": "2022-10-20T05:48:07.102258Z",
     "shell.execute_reply": "2022-10-20T05:48:07.101928Z"
    },
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fbcef7569a0>]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAu50lEQVR4nO3de3xU9Z3/8dcn93vIjQQSQhIugQRQIARRFBAFtCrVahfsxVtLu1237bbd/dltf7a17W/X7rZ2W123drFrtYrWakVFQbl5B8KdEBJCCCEhN0JCboTcvr8/ZnDjmJAJmZkzl8/z8cgjk3O+Z+aTyeQ9Z77ne75HjDEopZTyX0FWF6CUUsq9NOiVUsrPadArpZSf06BXSik/p0GvlFJ+ToNeKaX8nFNBLyIrRKRURMpF5IFB1l8jIntEpFdEbh9kfZyIVIvIo64oWimllPOGDXoRCQYeA24A8oDVIpLn0KwKuBt4doi7+SnwzqWXqZRS6lI5s0dfCJQbYyqMMd3AOmDlwAbGmEpjzAGg33FjEZkLpAKbXFCvUkqpEQpxok06cHLAz9XAfGfuXESCgF8CXwSuc2ab5ORkk5WV5UxTpZRSdrt37z5tjEkZbJ0zQT8a3wA2GGOqRWTIRiKyBlgDkJmZSVFRkZvLUkop/yIiJ4Za50zQ1wATBvycYV/mjAXA1SLyDSAGCBORdmPMJw7oGmOeAJ4AKCgo0Ml3lFLKhZwJ+l3AFBHJxhbwq4A7nblzY8wXLtwWkbuBAseQV0op5V7DHow1xvQC9wMbgRLgBWNMsYg8JCK3AIjIPBGpBu4Aficixe4sWimllPPE26YpLigoMNpHr5RSIyMiu40xBYOt0zNjlVLKz2nQK6WUn9OgV0opP+fucfRK+bzjpzt4ZV8N/f2241lxkaHcOT+TqDD991G+QV+pSg3BGMMzH53g5xtK6Orp58I5f8bAMx+d4Fd/czlzMhOsLVIpJ2jQKzWIxrbzfPfP+3mnrJFFU1P4xe2zSI2LAOCjiia++8J+bn/8A/5uyWS+fd1UgoOGPvNbKatp0CvloKevnzVPF1FS28rPPjuDL8zPZOAUHlfkJPHmt6/mJ68e5rdbygkJCuJb102xsGKlLk6DXikHv9xUxt6qFn67ejY3XzZ+0DaxEaH82+2z6O3r5z82lzE/J5ErcpI8XKlSztFRN0oNsL2skf/afozVhROGDPkLRISf3TqTiUnRfGvdXs50dHuoSqVGRoNeKbuG1i6+8/w+pqbG8OBN+U5tExMewm9Xz6a5o4fv/Xk/3namuVKgQa/Ux37y2mE6unt59M45RIYFO73djPR4vn/jNLYcaeCv+5yd2FUpz9GgVwo4fKqV1w/U8pWFOUxNjR3x9nctyGJaWiz/8fZRevs+daE1pSylQa8U8Ou3y4iNCOGrV+dc0vZBQcJ3rp9KZVMnL+3VvXrlXTToVcA7WH2WTYfr+crCHOKjQi/5fq7PS2VWRjy/2XyU7l7dq1feQ4NeBbxfvVXKmKhQ7l2YNar7ERH+4fqpVDef48+7Tw6/gVIeokGvAtruE81sLW1kzTU5xEZc+t78BYunpjAncwyPbimnq6fPBRUqNXoa9Cqg/efWchKjw7hrQZZL7k9E+O6yXGrPdvGy9tUrL6FBrwJWTcs5tpY2cGdhJtHhrjtJ/MpJSUxLi+VPO0647D6VGg0NehWwnt9ZhQFWFU5w6f2KCHfOz+RQTSsHqltcet9KXQoNehWQevr6WbfrJIunppCREOXy+//s7HQiQ4P500dVLr9vpUZKg14FpM0lDTS0nefO+RPdcv9xEaHcctl41u8/RWtXj1seQylnadCrgPTszirS4iJYkpvitse4c34m53r6eEUPyiqLORX0IrJCREpFpFxEHhhk/TUiskdEekXk9gHLLxeRD0WkWEQOiMjfuLJ4pS5FVVMn75Q1sqpwAiHB7tvXmZURz4z0OP60o0onO1OWGvZVLiLBwGPADUAesFpE8hyaVQF3A886LO8EvmyMyQdWAL8WkTGjrFmpUXluVxVBAn8zz7UHYR2JCHcWTuRIXRt7qlrc+lhKXYwzuzOFQLkxpsIY0w2sA1YObGCMqTTGHAD6HZaXGWOO2m+fAhoA931WVmoY/f2Gl/fUsCR3LOPiI93+eCsvH09kaDAv7al2+2MpNRRngj4dGHg+d7V92YiISCEQBhwb6bZKucquyjPUtXZxy+UXv6iIq0SHh3BdXiobDtbSo7NaKot45GCsiIwDngbuMcZ86tUuImtEpEhEihobGz1RkgpQrx44RURoENdNT/XYY948axzNnT28X37aY4+p1EDOBH0NMLAzM8O+zCkiEge8DvzAGPPRYG2MMU8YYwqMMQUpKdqzo9yjp6+fDQfruG56qkvPhB3OotwUYiNCWL//lMceU6mBnAn6XcAUEckWkTBgFbDemTu3t38Z+KMx5sVLL1Op0fvgWBNnOrqHvRasq4WHBLMiP41NxfU60ZmyxLBBb4zpBe4HNgIlwAvGmGIReUhEbgEQkXkiUg3cAfxORIrtm38euAa4W0T22b8ud8cvotRwXt1/itjwEBZN9fynxpsvG0/7+V62lWrXpPI8pz6/GmM2ABsclj044PYubF06jts9AzwzyhqVGrWunj42Hqpj+Yw0IkKdvx6sq1w5KYmk6DBe3X+KFTPSPP74KrDpmbEqIGwva6TtfK/Hu20uCAkO4saZ49h8pJ72872W1KAClwa9Cgiv7j9FYnQYV01KsqyGWy4fT1dPP28frresBhWYNOiV3+vq6WPLkQaW56e5dcqD4czNTCA1Lpw3D9VZVoMKTBr0yu99cOw0nd19LM/33Nj5wQQFCdfnpbK9rFFH3yiP0qBXfu+tw/XEhIewwMJumwuW5aVxrqdPT55SHqVBr/xaX7/hrcP1LMpNITzE86NtHF2Rk0RseAibirWfXnmOBr3ya/tONnO6vZtledZ221wQFhLE4mljebuknr5+nbpYeYYGvfJrmw7XExIkLM4da3UpH1uWl0pTRzd7q5qtLkUFCA165dfeKq5nwaQk4iNDrS7lY4tzUwgNFjbpMEvlIRr0ym+VN7RTcbrDa7ptLoiNCGXBpGQ2FdfplaeUR2jQK7+16bBtvPp1Xhb0YOu+qWzqpLyh3epSVADQoFd+a1NxPbMy4j1yJamRut7+5qPdN8oTNOiVXzrdfp791S0sneZ9e/MAqXERXJYRz+YSDXrlfhr0yi9tK23EGFg63XtG2zhanDuWvSdbONPRbXUpys9p0Cu/tPVIA2Njw8kfH2d1KUO6dtpYjIHtZQ1Wl6L8nAa98js9ff28U9bI4twURMTqcoY0Mz2e5JhwthzRi5Eo99KgV36nqLKZtvO9XDvNe7ttwDbJ2eLcFLaXNtDb1291OcqPadArv7O1tIHQYGHhFO+/0Py108bS2tXLnqoWq0tRfkyDXvmdrUcaKMxOJCbcqStlWmrhlGRCgoStpdpPr9xHg175lZNnOjna0M4SL5rb5mLiIkKZl5XI1iMa9Mp9NOiVX7mwZ+zt/fMDLZmWwpG6NmpazlldivJTGvTKr2w50kBWUhQ5KTFWl+K0C29Kulev3MWpoBeRFSJSKiLlIvLAIOuvEZE9ItIrIrc7rLtLRI7av+5yVeFKOerq6ePDY00s8aG9eYBJKTFMSIzUoFduM2zQi0gw8BhwA5AHrBaRPIdmVcDdwLMO2yYCPwLmA4XAj0QkYfRlK/VpH1Y0cb6336vmnneGiLB46lh7/XotWeV6zuzRFwLlxpgKY0w3sA5YObCBMabSGHMAcBwMvBx4yxhzxhjTDLwFrHBB3Up9yvbSRiJCg5ifnWh1KSO2ODeFzu4+dh3Xi5Eo13Mm6NOBkwN+rrYvc8ZotlVqRLaVNrAgJ4mIUOuvDTtSCyYlERYcxDYdZqncwCsOxorIGhEpEpGixkY9HVyNXOXpDiqbOn2u2+aCqLAQ5ucksr1MX//K9ZwJ+hpgwoCfM+zLnOHUtsaYJ4wxBcaYgpQU7z+bUXmfCwG5aKrvvn4WTU3haEO7DrNULudM0O8CpohItoiEAauA9U7e/0ZgmYgk2A/CLrMvU8qltpXahlVmJUdbXcolW5xre5PS7hvlasMGvTGmF7gfW0CXAC8YY4pF5CERuQVAROaJSDVwB/A7ESm2b3sG+Cm2N4tdwEP2ZUq5TFdPHx9WNPlst80Fk1JiSB8TyfZS7b5RruXUZCDGmA3ABodlDw64vQtbt8xg2z4JPDmKGpW6qJ3Hz9DV08+iXN/ttgH7MMvcFP66t4bu3n7CQrziEJryA/pKUj5vW2kjYSFBXJGdZHUpo7Zoagod3X0UndAPvsp1NOiVz9tW1sAVOUlEhvnesEpHV05OJjRYtPtGuZQGvfJpJ890UtHY4dOjbQaKCQ9hXlYi2zTolQtp0Cufts0+rHKxj/fPD7Roagql9W3UntVhlso1NOiVT9te2siExEhyfHhYpaMLo4fe0ZOnlIto0Cufdb63jw+OnWbRVO++CPhITU2NIS0uQrtvlMto0Cuftbuymc7uPhZP9e3x844uDLN87+hpevSi4coFNOiVz9pW1khYcBALJvn+sEpHi6am0Ha+l7160XDlAhr0ymdtL21kXnYC0T5wEfCRumpKMsFBwvYynQ5BjZ4GvfJJp1rOUVrf5jfDKh3FRYQyNzNB++mVS2jQK5/0zsfDKv2rf36gRbkpFJ9qpaGty+pSlI/ToFc+aVtpI+PiI5gy1ncuAj5SFz6tvFN22uJKlK/ToFc+p6evn/fLT7M417+GVTrKHx9HSmy4TlusRk2DXvmcPSeaaTvf67f98xeICIumpvDu0dP09Rury1E+TINe+ZytpY2EBgtXTU62uhS3W5I7lrPneth3Ui8ari6dBr3yOdtKGyiYmEhsRKjVpbjdQvswy61HdPSNunQa9MqnnGo5x5G6NpZM8+9umwviI23DLLdqP70aBQ165VMujCtf4sfDKh0tnmYbZlnfqsMs1aXRoFc+ZWtpA+ljIpnsx8MqHV14U9OLkahLpUGvfMb53j7eLz/Nkmn+PazS0bS0WNLiIrT7Rl0yDXrlM3Ydt81WGUjdNmAbZrlkmm2Ypc5mqS6FBr3yGVtLGwgL8c/ZKoezOHcs7ed7KarUYZZq5JwKehFZISKlIlIuIg8Msj5cRJ63r98hIln25aEi8pSIHBSREhH5vovrVwFka6ntIuBRYf43W+VwrrJfNFzPklWXYtigF5Fg4DHgBiAPWC0ieQ7N7gOajTGTgUeAh+3L7wDCjTEzgbnA1y68CSg1EieaOqho7GCxn58NO5SY8BAKsxPZckSDXo2cM3v0hUC5MabCGNMNrANWOrRZCTxlv/0isFRsR8sMEC0iIUAk0A20uqRyFVA2l9gC7tppgdU/P9CS3LEcbWjn5JlOq0tRPsaZoE8HTg74udq+bNA2xphe4CyQhC30O4BaoAr4d2PMmVHWrALQ5iP1TEqJJsuPLgI+UkunpwLwdkm9xZUoX+Pug7GFQB8wHsgGvisiOY6NRGSNiBSJSFFjo44VVp/U1tXDjoozXGcPukCVnRxNTkq0dt+oEXMm6GuACQN+zrAvG7SNvZsmHmgC7gTeNMb0GGMagPeBAscHMMY8YYwpMMYUpKQEZh+sGto7Zafp7Tcf79EGsuump/JRRRNtXT1Wl6J8iDNBvwuYIiLZIhIGrALWO7RZD9xlv307sMUYY7B111wLICLRwBXAEVcUrgLH5pJ6xkSFMidzjNWlWG7ptLH09BnePaoXI1HOGzbo7X3u9wMbgRLgBWNMsYg8JCK32JutBZJEpBz4DnBhCOZjQIyIFGN7w/iDMeaAq38J5b/6+g1bSxtYkjuWkGA97WPuxATiI0O1n16NiFMDko0xG4ANDsseHHC7C9tQSsft2gdbrpSz9lQ109zZw9LpgTvaZqCQ4CCW5KawrbSRvn5DcFDgTAWhLp3uIimvtrmkgZAg4ZoAHT8/mGunp3Kmo1svRqKcpkGvvNrmknoKsxOJC4CLjDhr0dQUQoKEt0t09I1yjga98lpVTZ0cbWjX0TYO4iNDmZeVyGbtp1dO0qBXXmvT4ToArteg/5Tr8lIpq2+n8nSH1aUoH6BBr7zWxuI6pqXFkpkUZXUpXmdZnu3Nb2NxncWVKF+gQa+8UmPbeYpONLM8P83qUrzShMQo8sfHadArp2jQK6/0dkk9xqBBfxHL89PYU9VCg15LVg1Dg155pY3FdUxIjGT6uFirS/FaF94ENx3Wg7Lq4jTolddp6+rhg/ImVuSnBdS1YUdqamoM2cnR2n2jhqVBr7zO1tJGuvv6tdtmGCLCsvxUPjzWxNlzOsmZGpoGvfI6G4vrSI4JZ05mgtWleL3l+Wn09hu26tTF6iI06JVX6erpY9uRBq7PSyVI53EZ1uUZYxgbG86bh7T7Rg1Ng155lffLT9PR3cfyfD1JyhlBQbbum+1ljZzr7rO6HOWlNOiVV3n9QC3xkaFcOSnZ6lJ8xo0zxnGup4+tpdp9owanQa+8RldPH5sO17M8P5WwEH1pOmt+ThLJMeG8duCU1aUoL6X/TcprbC9rpP18LzfNGm91KT4lOEi4cWYaW4400HG+1+pylBfSoFde47UDtSREhbJgUpLVpficm2aNp6unX688pQalQa+8wrnuPjaX1LNixjhC9ZKBI1YwMYHUuHBeO1BrdSnKC+l/lPIKW0sb6Ozu4+ZZ46wuxScFBQk3zhzH9tJG2rr05Cn1SRr0yiu8duAUyTFhzM/RbptLddOs8XT39fOWzn2jHGjQK8t1nO9ly5EGbpgxTi92PQpzMseQPiZSu2/Up2jQK8u9XVJPV08/N2m3zaiICJ+ZNY53jzbS0tltdTnKizgV9CKyQkRKRaRcRB4YZH24iDxvX79DRLIGrJslIh+KSLGIHBSRCBfWr/zAy3trGBcfwbysRKtL8Xm3XDaenj6je/XqE4YNehEJBh4DbgDygNUikufQ7D6g2RgzGXgEeNi+bQjwDPB1Y0w+sBjQI0XqYw1tXbxT1sits9N1bhsXyB8fR25qLC/tqba6FOVFnNmjLwTKjTEVxphuYB2w0qHNSuAp++0XgaVim0h8GXDAGLMfwBjTZIzRCTnUx9bvO0W/gdvmpFtdil8QEW6bk86eqhaO64XDlZ0zQZ8OnBzwc7V92aBtjDG9wFkgCZgKGBHZKCJ7ROSfRl+y8id/2VPDZRnxTB6rV5JylZWXpyMCL+tevbJz98HYEGAh8AX791tFZKljIxFZIyJFIlLU2Njo5pKUtyipbaWktpXb5mRYXYpfSYuPYOHkZF7aW0N/v7G6HOUFnAn6GmDCgJ8z7MsGbWPvl48HmrDt/b9jjDltjOkENgBzHB/AGPOEMabAGFOQkpIy8t9C+aSX99YQEiTcfJnObeNqt81Jp7r5HEUnmq0uRXkBZ4J+FzBFRLJFJAxYBax3aLMeuMt++3ZgizHGABuBmSISZX8DWAQcdk3pypf19vXz8t4aFueOJTE6zOpy/M7y/DSiwoL1oKwCnAh6e5/7/dhCuwR4wRhTLCIPicgt9mZrgSQRKQe+Azxg37YZ+BW2N4t9wB5jzOsu/y2Uz3n/WBONbef5nB6EdYuosBBWzEjj9QO1dPXo+IdAF+JMI2PMBmzdLgOXPTjgdhdwxxDbPoNtiKVSH3th10nGRIVy7fSxVpfit26fk8FLe2p441Att87W4yCBTM+MVR7X2HaejcV1fG5OBuEhwVaX47euyEkiKymK53acHL6x8msa9MrjXtxdTW+/YXXhhOEbq0sWFCSsKsxkZ+UZjta3WV2OspAGvfKo/n7Dul1VFGYn6th5D7h9bgahwcJzO3WvPpBp0CuP+uBYEyeaOrmzMNPqUgJCckw4y/LT+Mueaj0oG8A06JVHPbezijFRoayYkWZ1KQHjC4WZnD3XwxuHdKKzQKVBrzxm4EHYiFA9COspCybpQdlAp0GvPObPu0/qQVgLiAir7Qdly/SgbEDSoFce0dPXz9MfnmBBTpIehLXAHQUTCA8J4g/vV1pdirKABr3yiDcP1VF7tov7FmZbXUpASowO47Y56by0p5rmDr36VKDRoFcesfa942QlRXHtND0T1ir3XJXN+d5+nt1ZZXUpysM06JXb7alqZt/JFu65KluvImWhqamxXD0lmT9+WEl3b7/V5SgP0qBXbrf2vePERoRw+1ydb8Vq9y7Mpr71PBsO6lDLQKJBr9yqpuUcbx6qY3VhJtHhTs2hp9xo0ZQUclKiWfvecWwziatAoEGv3OqpDyoBuOvKLEvrUDZBQcK9V2VzsOYsO4+fsboc5SEa9Mptmju6eeajE3xm5jjSx0RaXY6y+9ycDJKiw3hs2zGrS1EeokGv3OYP7x+ns7uP+6+dbHUpaoDIsGC+cnUO75Q1sv9ki9XlKA/QoFdu0drVwx8+qGRFfhpTU/UEKW/zpQUTiY8M5bdbyq0uRXmABr1yiz9+UElbV6/uzXupmPAQ7r0qm7dL6jl8qtXqcpSbadArl+s438va945z7bSxzEiPt7ocNYS7r8wiJjyEx7bqXr2/06BXLvenHSdo7uzRvXkvFx8VypcXTGTDoVrKG3SyM3+mQa9cqrWrh8e3HePqKcnMyUywuhw1jPsWZhMVGsyv3iqzuhTlRhr0yqV+/04FzZ09/NPyaVaXopyQFBPOV6/JYcPBOvZWNVtdjnITp4JeRFaISKmIlIvIA4OsDxeR5+3rd4hIlsP6TBFpF5Hvuahu5YUaWrv473ePc/Nl45mZoX3zvuIrV+eQHBPGv75xRM+W9VPDBr2IBAOPATcAecBqEclzaHYf0GyMmQw8AjzssP5XwBujL1d5s//YfJSevn6+t2yq1aWoEYgJD+GbS6ew4/gZtpU1Wl2OcgNn9ugLgXJjTIUxphtYB6x0aLMSeMp++0VgqYgIgIh8FjgOFLukYuWVKhrbWbfrJF+Yn8nEpGiry1EjtGpeJhOTonj4jSP09etevb9xJujTgYEXm6y2Lxu0jTGmFzgLJIlIDPB/gJ+MvlTlzX7xZikRIUH8/dIpVpeiLkFYSBDfW5bLkbo2XtpTbXU5ysXcfTD2x8Ajxpj2izUSkTUiUiQiRY2N+tHR17x7tJE3i+v4+qJJJMeEW12OukSfmTmOyyeM4eE3S2nt6rG6HOVCzgR9DTDwas4Z9mWDthGRECAeaALmA78QkUrg28A/i8j9jg9gjHnCGFNgjClISUkZ6e+gLHS+t48fvVJMVlIUX70mx+py1CgEBQk/XTmDpo7zPKLDLf2KM0G/C5giItkiEgasAtY7tFkP3GW/fTuwxdhcbYzJMsZkAb8G/p8x5lHXlK68wdr3jlNxuoMf3ZJPRGiw1eWoUZqZEc8X5mfy1AeVlNTq1Aj+Ytigt/e53w9sBEqAF4wxxSLykIjcYm+2FluffDnwHeBTQzCV/znVco7fbi7n+rxUluTqtWD9xfeW5RIfGcqDrxzS4ZZ+wqlL/hhjNgAbHJY9OOB2F3DHMPfx40uoT3mxn752mH5jePAmx9G2ypeNiQrjgRum8X/+cpC/7KnRS0D6AT0zVl2SDQdreeNQHd9cOoUJiVFWl6Nc7I65E5g7MYGfvnaYhtYuq8tRo6RBr0asqf08//evh5iZHs8aPQDrl4KChF/cPouunj7++WXtwvF1GvRqxB58pZi2rl7+/Y7LCA3Wl5C/mpQSw/eW5fJ2ST1/3ec40E75Ev0vVSPy+oFaXj9Yy7eum0Juml45yt/duzCbuRMT+PF67cLxZRr0yml1Z7v44V8PMisjnq9pl01ACA4S/s3ehfOPLx6gX6dH8Eka9MopvX39fPO5vZzv7eeRv7mcEO2yCRg5KTH88KY8tpc18sS7FVaXoy6B/rcqp/zH5qPsrDzDz2+dwaSUGKvLUR72xfmZfGbmOP5tYym7T5yxuhw1Qhr0aljvHT3No1vLuWNuBrfO1jHVgUhE+JfPzSR9TCR//+xeWjq7rS5JjYAGvbqoUy3n+Pbze5mcEsNPVuZbXY6yUFxEKI/eOZvG9vN8+/l9Op2xD9GgV0Pq7O7lK08Vcb6nn8e/OIeoMKdOpFZ+bFbGGH50cz7bSht5+M0jVpejnKT/uWpQ/f2G776wnyN1ray9ex6Tx+pQSmXzxSsmUlbfxhPvVDBlbAx3FEwYfiNlKd2jV4P69dtlvHGojn++cbpOWKY+5cGb8lg4OZl/fvkguyr14Ky306BXn/Lczip+s6WczxdkcN/CbKvLUV4oJDiIx+6cQ0ZCFF/9YxFH69usLkldhAa9+oQ3Dtbyg5cPsjg3hZ99dib2S/8q9SnxUaE8dU8hocFBfGntTqqbO60uSQ1Bg1597L2jp/nWun3Mzkzg8S/MJSxEXx7q4jKTovjjvYV0dvfypbU7Od1+3uqS1CD0P1kBsPP4GdY8XUROSjRP3jWPyDC9WpRyzvRxcTx59zxqz57jS2t3cqZDx9h7Gw16xQflp7nryZ2Mi4/gj/cWEh8VanVJyscUZCXyuy8VUNHYzp2//0j37L2MBn2A217WyD3/s4vMxCjWrVnA2LgIq0tSPmrR1BSevHselU0drHriI53t0oto0AewDQdr+epTRUxKieG5NVeQEhtudUnKx101OZmn7imktuUcd/zuQypPd1hdkkKDPiAZY/jvdyv4u2f3MCsjnme/Op/E6DCry1J+Yn5OEk9/ZT6t53q47fEP2FvVbHVJAU+DPsD09Rt+8uphfvZ6CTfMSOOZr8xnTJSGvHKtOZkJvPSNq4iNCGH17z9iY3Gd1SUFNA36ANLc0c3df9jJ/3xQyVcWZvPo6jlEhOroGuUe2cnR/OVvr2RaWhxff2Y3v9l8VC9cYhGngl5EVohIqYiUi8gDg6wPF5Hn7et3iEiWffn1IrJbRA7av1/r4vqVkw7VnOXmR99jR8UZ/uW2mfzwpjyCgvRkKOVeyTHhrFtzBbdens6v3ipjzdNFtHb1WF1WwBk26EUkGHgMuAHIA1aLSJ5Ds/uAZmPMZOAR4GH78tPAzcaYmcBdwNOuKlw5xxjDup1VfO7xD+jtMzz/tStYXZhpdVkqgESEBvPLz1/Gj2/OY1tpI7f89j0OVp+1uqyA4swefSFQboypMMZ0A+uAlQ5tVgJP2W+/CCwVETHG7DXGnLIvLwYiRUSHdnhIS2c3f/vMHh546SBzJybw6t8vZHZmgtVlqQAkItx9VTbPfvUKunr6ue3x9/nd9mPaleMhzgR9OnBywM/V9mWDtjHG9AJngSSHNp8D9hhj9EwKD3j3aCMrfv0um4/U8/0bpvHMffN1+KSyXGF2Im9++2qWTkvlX944wpee3KFz5HiARw7Gikg+tu6crw2xfo2IFIlIUWNjoydK8ltnO3v4xz/v50trdxIVHsxLf3sVX1s0SfvjldcYExXG41+cw7/eNpO9VS0sf+Qd/vhhpe7du5EzQV8DDLyyQIZ92aBtRCQEiAea7D9nAC8DXzbGHBvsAYwxTxhjCowxBSkpKSP7DRRg64tfv/8U1z2ynZf21vCNxZPY8M2rmZkRb3VpSn2KiLCqMJNN/3ANc7MSefCVYj7/uw8pqW21ujS/5EzQ7wKmiEi2iIQBq4D1Dm3WYzvYCnA7sMUYY0RkDPA68IAx5n0X1awclNS2suqJj/jmc3tJjQvnlb+7in9aMU2HTiqvl5EQxVP3zOOXd1xGeWM7n/nNu/zolUOc7dSROa4kxgz/cUlEbgR+DQQDTxpjfi4iDwFFxpj1IhKBbUTNbOAMsMoYUyEiPwS+DxwdcHfLjDENQz1WQUGBKSoquuRfKJA0tHXxm81HeXZHFXGRofzj8lxWzcskWLtplA9q7ujmV2+V8acdJxgTFcY3r53MnfMn6nTZThKR3caYgkHXORP0nqRBP7zWrh6e2F7B2veO09PXz53zM/nO9VP1DFflF4pPneWnrx3mo4ozZCZG8d1lU7l51ng9zjQMDXo/cbazhyffP84f3j9Oa1cvN80ax/eW5ZKVHG11aUq5lDGG7WWNPPxmKSW1rUwZG8P9107mplnj9RPrEDTofVx9axf/80ElT394gvbzvVyfl8q3lk5hRroeaFX+rb/f8NrBWh7dcpSy+nZykqNZc00On52drsegHGjQ+6jiU2dZ+95xXt1/it5+ww0z0rh/yRTyxsdZXZpSHtXfb9hYXMdvt5RzuLaV5JgwvnRFFl+4IpPkGD0/BDTofUpXTx+vH6jlmR0n2FvVQlRYMJ8vmMA9V2UxMUm7aFRgM8bw4bEmfv9uBVtLGwkNFlbMGMcX52dSmJ0Y0Bezv1jQh3i6GPVpxhgO1bTy590neWXfKc6e6yEnOZoffmY6d8ydoJf2U8pORLhycjJXTk6mvKGdP+04wYu7q3l1/ykmpURz+9wJ3DYnnVS9Uton6B69hU6e6WT9/lOs33eK0vo2wkKCWJ6fxqp5E7hyUlJA750o5azO7l5e21/LC0UnKTrRTJDAwikp3HLZeJbnpxIbERg7Stp140UqT3ewsbiONw7Vse9kCwBzJyZw6+x0br5sPPGRgfGiVModKhrb+cueal7Zd4rq5nOEhQSxeGoKN8xM49ppqX79/6VBb6G+fsO+ky1sOVLP5pIGjtS1ATAjPY7PzBzPzZeNIyMhyuIqlfIvxhj2VLXw6v5TvHmojrrWLkKChCtyklg6fSzXThvrd8e8NOg9rO5sF+8ebeSdo6d572gjzZ09BAcJBRMTWJafxvL8VA13pTykv9+wv7qFN4vrePtwPccabRcsz0mJ5popKVwzNZn52UlEh/v2IUsNejdraOti1/FmPjh2mg8rmqiwv5CSY8K5ekoyS6aNZdGUFD2oqpQXqDzdwZYjDWwva2TH8Sa6evoJCRIunzCGBZOSWJCTxOzMBCLDfGucvga9C/X1G442tLG3qoXdJ5opqjxDZZNtPu3osGDmZSeyICeJq6ekMC0tVk/bVsqLdfX0UVTZzHvltp20g9Ut9BsICRLy0+OZNzGBORMTmJ05hnHxkVaXe1Ea9JfIGEPVmU4O1pzlYPVZ9le3cKimlfbzvQAkRIVSkJVIYVYiBVkJzEiPJzRYJ2BSyle1dfVQVNlM0Ykz7DrezL7qFrp7+wFIi4tgVka8/WsM+ePjSPKik7V0HL0TOrt7Katvp7SulZLaNg7XtlJyqpU2e6iHBQcxfVwst85OZ3bmGOZkJjAxKUqHQCrlR2IjQlkybSxLpo0FoLu3n5LaVvZUNbO3qoWDNWfZdLj+4/bj4iPIHx/H9HFx5KbFMi0tjqykKEK8bIcv4Pboz3b2UN7YzjH7V3l9O2UNbZw8c+7jNlFhweSmxZI/Po788fHkj49jWlqcTpeqlOLsuR4O1Zzl8KlWik+dpfhUKxWnO+izXyErLDiInJRopqbGMnlsDJNSYpg8NoaJSVFunZ8noLpujDGc6ejmZPM5TjR1cKKpkxNNnVQ2dXD8dAdnOro/bhsWHER2cjRTUmOYmhrL1NRYpo+LZUJClPatK6Wc1tXTR3lDO0fq2jha30ZZfRtl9e3UtPzvDqQIjI+PJCclmolJUWQlRZOZGMXEpGgmJEYSFTa6DpaA6LppaOviy2t3cvJMJx3dfZ9YlxYXQVZyFMvz08hOjiIn2fYOm5EQ6XUfsZRSviciNJgZ6fGfmlG2s7uXisYOjjW2c/x0x8df6/edorWr9xNtk2PCWDApmd+unu3y+vwm6MdEhpGREMmCSUlMSIhiQmIUE5OiyEx078clpZQaSlRYyKBvAAAtnd2caOqk6oztq7q5kwQ3XTzIb4I+LCSI/75rntVlKKWUU8ZEhTEmKozLJoxx+2Npv4VSSvk5DXqllPJzGvRKKeXnNOiVUsrPORX0IrJCREpFpFxEHhhkfbiIPG9fv0NEsgas+759eamILHdh7UoppZwwbNCLSDDwGHADkAesFpE8h2b3Ac3GmMnAI8DD9m3zgFVAPrAC+E/7/SmllPIQZ/boC4FyY0yFMaYbWAesdGizEnjKfvtFYKnYJoFZCawzxpw3xhwHyu33p5RSykOcCfp04OSAn6vtywZtY4zpBc4CSU5uq5RSyo284oQpEVkDrLH/2C4ipaO4u2Tg9Oircjmta2S0rpHRukbGH+uaONQKZ4K+Bpgw4OcM+7LB2lSLSAgQDzQ5uS3GmCeAJ5yoZVgiUjTUxD5W0rpGRusaGa1rZAKtLme6bnYBU0QkW0TCsB1cXe/QZj1wl/327cAWY5sWcz2wyj4qJxuYAux0TelKKaWcMewevTGmV0TuBzYCwcCTxphiEXkIKDLGrAfWAk+LSDlwBtubAfZ2LwCHgV7g74wxfYM+kFJKKbdwqo/eGLMB2OCw7MEBt7uAO4bY9ufAz0dR40i5pAvIDbSukdG6RkbrGpmAqsvrLjyilFLKtXQKBKWU8nM+GfQicoeIFItIv4gUOKwbdsoF+4HlHfZ2z9sPMru6xudFZJ/9q1JE9g3RrlJEDtrbue9iuf/7eD8WkZoBtd04RLuLTnvhhrr+TUSOiMgBEXlZRMYM0c4jz9dopv1wY00TRGSriBy2v/6/NUibxSJydsDf98HB7ssNtV307yI2v7E/XwdEZI4Hasod8DzsE5FWEfm2QxuPPF8i8qSINIjIoQHLEkXkLRE5av+eMMS2d9nbHBWRuwZrMyxjjM99AdOBXGAbUDBgeR6wHwgHsoFjQPAg278ArLLf/i/gb91c7y+BB4dYVwkke/C5+zHwvWHaBNufuxwgzP6c5rm5rmVAiP32w8DDVj1fzvz+wDeA/7LfXgU874G/3Thgjv12LFA2SF2Lgdc89Xpy9u8C3Ai8AQhwBbDDw/UFA3XARCueL+AaYA5waMCyXwAP2G8/MNhrHkgEKuzfE+y3E0b6+D65R2+MKTHGDHZS1bBTLtinZrgW21QNYJu64bPuqtX+eJ8HnnPXY7iBM9NeuJQxZpOxnVUN8BG2cy6sMpppP9zGGFNrjNljv90GlOA7Z5qvBP5obD4CxojIOA8+/lLgmDHmhAcf82PGmHewjUgcaOBraKgcWg68ZYw5Y4xpBt7CNm/YiPhk0F+EM1MuJAEtA0LF3dMyXA3UG2OODrHeAJtEZLf9DGFPuN/+8fnJIT4uWj11xb3Y9v4G44nnazTTfniEvatoNrBjkNULRGS/iLwhIvkeKmm4v4vVr6lVDL2zZcXzBZBqjKm1364DUgdp45LnzSumQBiMiLwNpA2y6gfGmFc8Xc9gnKxxNRffm19ojKkRkbHAWyJyxP7u75a6gMeBn2L7x/wptm6le0fzeK6o68LzJSI/wHbOxZ+GuBuXP1++RkRigL8A3zbGtDqs3oOte6Ldfvzlr9hOVHQ3r/272I/B3QJ8f5DVVj1fn2CMMSLitiGQXhv0xpjrLmEzZ6ZcaML2sTHEvic26LQMrqhRbNNB3AbMvch91Ni/N4jIy9i6DUb1D+LscycivwdeG2SVU1NXuLouEbkbuAlYauwdlIPch8ufr0GMZtoPtxKRUGwh/ydjzEuO6wcGvzFmg4j8p4gkG2PcOq+LE38Xt7ymnHQDsMcYU++4wqrny65eRMYZY2rt3VgNg7SpwXYc4YIMbMcmR8Tfum6GnXLBHiBbsU3VALapG9z1CeE64IgxpnqwlSISLSKxF25jOyB5aLC2ruLQL3rrEI/nzLQXrq5rBfBPwC3GmM4h2njq+RrNtB9uYz8GsBYoMcb8aog2aReOFYhIIbb/cbe+ATn5d1kPfNk++uYK4OyAbgt3G/JTtRXP1wADX0ND5dBGYJmIJNi7WZfZl42Mu482u+MLW0BVA+eBemDjgHU/wDZiohS4YcDyDcB4++0cbG8A5cCfgXA31fk/wNcdlo0HNgyoY7/9qxhbF4a7n7ungYPAAfsLbZxjXfafb8Q2quOYh+oqx9YXuc/+9V+OdXny+Rrs9wcewvZGBBBhf+2U219LOR54jhZi63I7MOB5uhH4+oXXGXC//bnZj+2g9pUeqGvQv4tDXYLtAkbH7K+/AnfXZX/caGzBHT9gmcefL2xvNLVAjz277sN2TGczcBR4G0i0ty0A/nvAtvfaX2flwD2X8vh6ZqxSSvk5f+u6UUop5UCDXiml/JwGvVJK+TkNeqWU8nMa9Eop5ec06JVSys9p0CullJ/ToFdKKT/3/wEgrxSoFWbWQQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "nu=jnp.linspace(-10,10,100)\n",
    "plt.plot(nu, voigt(nu,1.0,2.0)) #beta=1.0, gamma_L=2.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## optimization of a simple absorption model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we try to fit a simple absorption model to mock data.\n",
    "The absorption model is \n",
    "\n",
    "$ f= e^{-a V(\\nu,\\beta,\\gamma_L)}$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:07.104653Z",
     "iopub.status.busy": "2022-10-20T05:48:07.104342Z",
     "iopub.status.idle": "2022-10-20T05:48:07.105779Z",
     "shell.execute_reply": "2022-10-20T05:48:07.106056Z"
    }
   },
   "outputs": [],
   "source": [
    "def absmodel(nu,a,beta,gamma_L):\n",
    "    return jnp.exp(-a*voigt(nu,beta,gamma_L))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Adding a noise...\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:07.108361Z",
     "iopub.status.busy": "2022-10-20T05:48:07.108032Z",
     "iopub.status.idle": "2022-10-20T05:48:07.310843Z",
     "shell.execute_reply": "2022-10-20T05:48:07.310515Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fb90d0398b0>]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAbeklEQVR4nO3dfbAddZ3n8ffn3hhYRtRAWEECSVgYVnC3FE4hs/MEwwiBoowP+xDc2YHRKcpdsHZ0rC0orWjF3XWe3KearEycTSmWC7IwoykKC1HD+scazb0ISNDIJWPkBkYy5OLsVihCcr/7x+nLdE763Nvnnu7TD+fzqrqVc7v75PzuOX2+/evv79u/VkRgZmbtNVF1A8zMrFwO9GZmLedAb2bWcg70ZmYt50BvZtZyK6puQK/Vq1fHunXrqm6GmVmjTE9P/01EnJG1rnaBft26dUxNTVXdDDOzRpG0v986p27MzFrOgd7MrOUc6M3MWs6B3sys5RzozcxazoHezKzlxjrQT++fY+vOGab3z1XdFDOz0tSujn5UpvfP8S//fBdHjs6zcsUEX/rdy7l07aqqm2VmVrix7dHv2vcCR47OMx/wytF5du17oeommZmVYmwD/eXnnc7KFRNMCl6zYoLLzzu96iaZmZVibFM3l65dxZd+93J27XuBy8873WkbM2utJQO9pO3A9cDzEfGWjPUC/itwHXAYuCkiHknW3Qh8PNn030fEF4pqeBEuXbvKAd7MSjG9f+7VjiRQaacyT4/+88CfAnf2WX8tcEHy83bgs8DbJZ0GfALoAAFMS9oRES5xMbNWSxd7rJgQSBw9Vl3hx5I5+oj4NnBokU02AndG1y7gDZLOAq4BHoqIQ0lwfwjYUESjq+JyTDPL47hij2PBKxUXfhSRoz8beCb1+2yyrN/yE0i6GbgZ4Nxzzy2gScVzOaaZ5bVQ7PHK0Xkmkx79sWPzlRV+1GIwNiK2AdsAOp1OVNycTFnlmA70Zs2XzqUv5zud9fzeYg+of45+KQeAc1K/r0mWHQCu6Fn+cAGvV4n0EdrlmGbtMOyZ+mLP7y32qLJjWEQd/Q7gt9V1OfDziHgOeBC4WtIqSauAq5NljbRwhP7I1Rc6bWPWEsNeODns80c17penvPIuuj3z1ZJm6VbSvAYgIu4AHqBbWjlDt7zyd5J1hyR9Ctid/FdbImKxQd1K5Tl9czmmWT0tN/0y7Jn6MM8f5bjfkoE+Im5YYn0At/RZtx3YvrymjY4HWs2aa5jv77AXTg7z/FGO+9ViMLYMgxzhPdBq1lzDfn+HPVNf7vNHOe7XykDfe4TffP3FzB0+0jfoe6DVrLma+v3NqszZunOmlMocdTMv9dHpdGJqamqo/2Przhk+8/W9zEd3tHliQsxHLHpaN2yJlZlVpy7f3+W2o4j0saTpiOhkrWtljz59hJe6QX6p07pBT7/qsmOZWf/v7yi/p8ME67LTx60M9OlTolWnrGTL/XuWfVqXtaN48Nas/kb9PR0mWJedfmploIfjj/AXnnlqoadToxq89VmD2fKNushimGBd9rTprQ30acsdFe+3o4xi8MdnDWbDGfZ7OmhHq4hSTZdXVqDfjjKKm5b0HmTue2TWvXuzAQzzPV1uR6uuF1U60C9isR0l/YGWcYOB3tnv7p2erXQ+a7MmWu4gbduurXGgX8JSR+jl3mBgqR0tfZB59sWXuOt7P23NTmdWpTy99abW5vczdoG+6AHO3hsMQBAsHpDznhYuHGSm989x3yOzmTudB2zNBpOnt15meraK7+xYBfoyBjiXc4OBQU8L++10HrA1O9FSgTRvb72MfHtV39mxCvRl5N2Wc4OB5ZwWZu10bcsj2ngpo2ebJ5Au1lsvu7dd1Xd2rAJ93gC7nLKqQW4wUNRpYdvyiDY+iujZZn1P8wbSrI7TKHrbVX1nxyrQ5wmwozq1KuK0cBRlnmZlGLZn2+97OkwgHUVvu6rvbKsCfRE3Dyn6wy77VLCudbtmixm2Z9vvezpMIB1Vb7uK72xrAn1RPfEiP+ymDZa6gsdGZZCAnLVfLvY9XW4gbfMZcmsCfVE98SI/7OW0qapg27SDkjVfnoDcb78sKyjnvcCqaZ2i1gT6InviRZ1aDdqmsoJtnp3SFTxWpX776GL75ahSIFk3Mtpy/55GdYpaE+jreNo1aJvKCLZ5Dx6u4LGqZO2jwKvTjFe9X/Z+L7/2xHON6xS1JtBDPQcm+82JU3awXXitZ198KXe5Wd0OlDYesibw+4tHZnPfCjSv5aZber+X177lLHb/5FCjOkWtCvR1NuyFHMt9rRUTYsXkxJJX6y68vgO8jVpvIBUcF/jnDh/hlivPH+o1hkmLZn0vl3uPi6o40I/IMBdy5JXViz82H/yLy87h7Df8vcbslDZesq4u7ze303INmxbNuiiySd8lB/oRKTsHvlgv/r2XrGnUTmnN1y9N0m95b+AseoqCcR+DcqAfkbLn10j3WNyLtyr1S5MMkj4peoqCcR+DcqAfoTLn1+jtsbgXb1XplyYZNn1SdPplnDjQV6yOF3qZDaNfmmTY9Mm4p1+GoYioug3H6XQ6MTU1VXUzRmahR7+w8zbh4guzpfS7vWb68XL286ZdkTpKkqYjopO5zoG+et55ra08tcboLBbonbqpgXHOHVq7eWqNepjIs5GkDZL2SpqRdFvG+rWSvinpcUkPS1qTWndM0qPJz44iG29m1ZreP8fWnTNM75/LXL+QV58UzqtXaMkevaRJYCvwDmAW2C1pR0Q8mdrsT4A7I+ILkn4D+DTwr5J1L0XEW4tttplVZSHVuOqUlUtO7uUigXrIk7q5DJiJiH0Aku4GNgLpQH8R8JHk8U7gKwW20cxqIp1zn5CYjyj1am8rRp7UzdnAM6nfZ5NlaY8B70kevxs4VdLCOdrJkqYk7ZL0rqwXkHRzss3UwYMH87fezEZiIUWzMNnYfMD8fDAhOS3TAEUNxn4U+FNJNwHfBg4Ax5J1ayPigKTzgG9J+kFEPJ1+ckRsA7ZBt+qmoDaZWQEWm16jqJklrVx5Av0B4JzU72uSZa+KiGdJevSSXgu8NyJeTNYdSP7dJ+lh4G3AcYHeBjNMOWbT75Rjo+fpNZovT6DfDVwgaT3dAL8JeF96A0mrgUMRMQ/cDmxPlq8CDkfEy8k2vwz8UYHtHzvD1CW34U45NnqeXqP5lgz0EXFU0q3Ag8AksD0i9kjaAkxFxA7gCuDTkoJu6uaW5OlvBv5M0jzd8YA/6KnWsQENU5fchjvl2Oi5cqb5cuXoI+IB4IGeZZtTj+8F7s143v8B/tGQbbSUYeb7aMOdcqwarpxpNk+B0EDO0ZtZL891Y4XxgcGsnjzXjRXCE1SZNZMDvS0p6160HrxtD5+ltZ8DvS1qsYtl0oO3DhbNlFVyu9wLoLwP1JcDvS0qz8UyTuk0V/rzPfLKPJu/+gTzEUNfo+F9oF5yTVNs46t3mtn3XrKGW648/7gvcVZtvzVD+vOdmDhxkrK8vA/Um3v0tqg8F8v4Xp7Nlf58F6YdLuIaDe8D9eLySss0aL7V+dl26HevV+8D9ec6ehuI863mfaB5Fgv0ztHbCZxvNe8D7eJAbyfwfT7bz/d6HS9O3VimPPlW52SbKW9axp9vs3gKBBvYUrMVOofbXHmnuvaMle3h1I0ti3O4zeW0zPhxj96WxXXTzeUbiYwf5+ht2ZzDNasP5+itFM7hmjWDc/RmZi3nQG+lWqpe28zK59SNlcYlmGb14B69lcYlmGb14EBvpXG9tlk9OHVjpXG9tlk9ONBbqVyCaVY9p26scK60MasX9+itUK60Masf9+itUK60MasfB3orlCttzOrHqRsrlCttzOonV49e0gZJeyXNSLotY/1aSd+U9LikhyWtSa27UdJTyc+NRTbe6unStau45crzHeTNamLJQC9pEtgKXAtcBNwg6aKezf4EuDMi/jGwBfh08tzTgE8AbwcuAz4hyd/+MeVqHLNq5EndXAbMRMQ+AEl3AxuBJ1PbXAR8JHm8E/hK8vga4KGIOJQ89yFgA3DX0C23RnE1jll18qRuzgaeSf0+myxLewx4T/L43cCpkk7P+VwbA67GMatOUVU3HwV+XdL3gV8HDgDH8j5Z0s2SpiRNHTx4sKAmWZ24GsesOnlSNweAc1K/r0mWvSoiniXp0Ut6LfDeiHhR0gHgip7nPtz7AhGxDdgG3VsJ5m++NYWrccyqkyfQ7wYukLSeboDfBLwvvYGk1cChiJgHbge2J6seBP5jagD26mS9jSHPe2NWjSVTNxFxFLiVbtD+IXBPROyRtEXSO5PNrgD2Svox8EbgPyTPPQR8iu7BYjewZWFg1szK5SonW6CIemVKOp1OTE1NVd0Ms0ZzldP4kTQdEZ2sdZ4CwayFXOVkaQ70Zi3kKidL81w3Vonp/XOuwCmRq5wszYHeRs7549HorXLywXV8OdDbyGXljx14yuWD63hzjt5Gzvnj0fPg7Hhzj95GzvnjYuVJySwcXF85Ou+D6xhyHb1Zg/WmZDZffzFzh49kBn3n6NttsTp69+jNGiydkjnyyjybv/oE8xGZeXhPQTG+nKM3a7D0eMfEhJiPcB7eTuAevdWK0wuDSY93rDplJVvu3+M8vJ3Agd5qwyWAy5NOyVx45qk+UNoJHOitNlxfPzzn4S2Lc/RWG66vNyuHe/RWuXRe3vX1ZsVzoLdKZeXlb7ny/KqbZdYqTt1YpXxp/vL47lE2CPforVK+NH9wrk6yQTnQW6U8783gXJ1kg3Kgt8q5JHAwPguyQTnQmzWMz4JsUA70Zg3ksyAbhKturLZcWWJWDPforZZcWWJWHPforZZcX29WHAd6qyXPe2NWHKdurJZcWWJWHAd6qy1XlpgVw6kbM7OWc6A3M2s5B3ozs5bLFeglbZC0V9KMpNsy1p8raaek70t6XNJ1yfJ1kl6S9Gjyc0fRf4CZmS1uycFYSZPAVuAdwCywW9KOiHgytdnHgXsi4rOSLgIeANYl656OiLcW2mozM8stT4/+MmAmIvZFxBHgbmBjzzYBvC55/Hrg2eKaaGZmw8gT6M8Gnkn9PpssS/sk8FuSZun25j+UWrc+Sen8b0m/mvUCkm6WNCVp6uDBg/lbb2ZmSypqMPYG4PMRsQa4DviipAngOeDciHgb8BHgf0p6Xe+TI2JbRHQionPGGWcU1CRrM094ZpZfngumDgDnpH5fkyxL+wCwASAiviPpZGB1RDwPvJwsn5b0NPCLwNSwDbfx5QnPzAaTp0e/G7hA0npJK4FNwI6ebX4KXAUg6c3AycBBSWckg7lIOg+4ANhXVONtPHnCM7PBLNmjj4ijkm4FHgQmge0RsUfSFmAqInYAvw98TtKH6Q7M3hQRIenXgC2SXgHmgQ9GxKHS/hprten9c+za9wKrTlnpW+mZDUARUXUbjtPpdGJqypkdO15vumbz9Rczd/iIJzwzS0iajohO1jpPamaN0JuumTt8hFuuPL/qZo3UwhmND242KAd6a4SF+enHLV2TTldtuX+PB6BtWRzorRHGcX76dLpqQmI+4rgB6HF4D6wYDvTWGOM2P306XUUEExNCxFid0VgxHOjNaqo3XeUBaFsuB3prvLYNUqb/nnFLV1k5HOit0dp2lWzW3zNu1UVWPN94xBqtbVfJtu3vsXpwoLdGW8hjT4pWDFK27e+xevCVsdZ4bc7Rt+HvsdHwlbHWam0ru2zb32PVc+rGzKzlHOjNzFrOgd7MrOUc6M3MWs6B3sys5RzozcxazoHezKzlHOjNzFrOgd7MrOUc6M1qYHr/HFt3zjC9f67qplgLeQoEs4q1baplqx/36K1Vmtgz9tTEVjb36K01mtoz7r1loKcmtqI50FtrZPWM6xbos6YgvnTtKt8y0ErlQG+tUfee8WJnHJ6a2MrkQG+tUfeecRPOOKydHOitVercM677GYe1lwO92YjU/YzD2suB3myE6nzGYe2Vq45e0gZJeyXNSLotY/25knZK+r6kxyVdl1p3e/K8vZKuKbLxZk3WxJp/a6Yle/SSJoGtwDuAWWC3pB0R8WRqs48D90TEZyVdBDwArEsebwIuBt4EfEPSL0bEsaL/ELMmaWrNvzVTnh79ZcBMROyLiCPA3cDGnm0CeF3y+PXAs8njjcDdEfFyRPwVMJP8f2ZjzVfD2ijlCfRnA8+kfp9NlqV9EvgtSbN0e/MfGuC5SLpZ0pSkqYMHD+Zsutni6pwaWajAmRSuwLHSFTUYewPw+Yj4jKRfAr4o6S15nxwR24BtAJ1OJwpqk42xuqdGXIFjo5Qn0B8Azkn9viZZlvYBYANARHxH0snA6pzPNStcEy5OcgWOjUqe1M1u4AJJ6yWtpDu4uqNnm58CVwFIejNwMnAw2W6TpJMkrQcuAL5XVOPN+nFqxOzvLNmjj4ijkm4FHgQmge0RsUfSFmAqInYAvw98TtKH6Q7M3hQRAeyRdA/wJHAUuMUVNzYKTo2Y/R1143F9dDqdmJqaqroZZmaNImk6IjpZ63zjETOzlnOgt7FT57JLszJ4rhsbK3UvuzQrg3v0NlZ8RaqNIwd6Gysuu7Rx5NSNjZUyyy6z7gdrVgcO9DZ2yrgi1bl/qzOnbmwslF1ps1ju31U+VjX36K31RtHb7nc/WPf0rQ4c6K31RjHBWb/cfxMmV7P2c6C31uvX2y5aVu5/VK9tthjPdWNjocqKGFfj2CgsNteNe/Q2Fqqc+93zzlvVXHVjZtZyDvRmZi3nQG9jzTXuNg6co7exNWyNuwdZrSkc6G1s9da43/fIbO7A7QuhrEkc6G1spWvcJyfEvdOzHD12YuDO6rn7QihrEgd6G1vpq1mfffEl7vreT08I3L09983XX8zc4SOsOmWlL4SyxnCgt7G2UOM+vX+O+x6ZPSFwp3vuR16ZZ/NXn2A+4rig7xy91Z0DvRn956pJp3ckMR/xaq9/7vARbrny/IpbbrY0B3qzRPoK1nRefuEAsOqUlWy5f4/TNdY4DvRmPbIqahZ67heeeapLKq1xHOjNeixWUeN5a6yJfGWsWQ/fQNzaxj16sx5l3kDcrAoO9GYZhk3ReHoEqxMHerOCeXoEqxvn6M0KljWYa1alXIFe0gZJeyXNSLotY/1/lvRo8vNjSS+m1h1LrdtRYNvNasmDuVY3S6ZuJE0CW4F3ALPAbkk7IuLJhW0i4sOp7T8EvC31X7wUEW8trMVmNefBXKubPDn6y4CZiNgHIOluYCPwZJ/tbwA+UUzzzJrJ9fZWJ3lSN2cDz6R+n02WnUDSWmA98K3U4pMlTUnaJeldfZ53c7LN1MGDB/O13MzMcil6MHYTcG9EHEstWxsRHeB9wH+R9A96nxQR2yKiExGdM844o+AmmZmNtzyB/gBwTur3NcmyLJuAu9ILIuJA8u8+4GGOz9+bmVnJ8gT63cAFktZLWkk3mJ9QPSPpHwKrgO+klq2SdFLyeDXwy/TP7ZuZWQmWHIyNiKOSbgUeBCaB7RGxR9IWYCoiFoL+JuDuiIjU098M/JmkeboHlT9IV+uYmVn5dHxcrl6n04mpqamqm2Fm1iiSppPx0BPX1S3QSzoI7B/iv1gN/E1BzSmS2zUYt2swbtdg2tiutRGRWc1Su0A/LElT/Y5qVXK7BuN2DcbtGsy4tctz3ZiZtZwDvZlZy7Ux0G+rugF9uF2DcbsG43YNZqza1bocvZmZHa+NPXozM0txoDcza7lGBnpJ/0zSHknzkjo9625PbpCyV9I1fZ6/XtJ3k+2+nEztUHQbv5y64cpPJD3aZ7ufSPpBsl3pV4pJ+qSkA6m2Xddnu0VvNlNCu/5Y0o8kPS7pLyW9oc92I3m/ctxs56TkM55J9qV1ZbUl9ZrnSNop6clk//+3GdtcIennqc93c9ntSl530c9FXf8teb8el3TJCNp0Yep9eFTS30r6vZ5tRvJ+Sdou6XlJT6SWnSbpIUlPJf9mzmst6cZkm6ck3bisBkRE437oTq1wId1J0jqp5RcBjwEn0Z0u+WlgMuP59wCbksd3AP+65PZ+BtjcZ91PgNUjfO8+CXx0iW0mk/fuPGBl8p5eVHK7rgZWJI//EPjDqt6vPH8/8G+AO5LHm4Avj+CzOwu4JHl8KvDjjHZdAdw/qv0p7+cCXAd8DRBwOfDdEbdvEvhruhcVjfz9An4NuAR4IrXsj4Dbkse3Ze3zwGnAvuTfVcnjVYO+fiN79BHxw4jYm7FqI935dl6OiL8CZujeOOVVkgT8BnBvsugLwLvKamvyev+cnlk9a+7Vm81ExBFg4WYzpYmIr0fE0eTXXXRnSa1Knr9/I919B7r70lXJZ12aiHguIh5JHv9f4If0uTdEDW0E7oyuXcAbJJ01wte/Cng6Ioa56n7ZIuLbwKGexel9qF8cugZ4KCIORcQc8BCwYdDXb2SgX0Sem6ScDryYCip9b6RSkF8FfhYRT/VZH8DXJU1LurnEdqTdmpw+b+9zupj7ZjMleT/d3l+WUbxfef7+V7dJ9qWf0923RiJJFb0N+G7G6l+S9Jikr0m6eERNWupzqXqfOmEK9ZQq3i+AN0bEc8njvwbemLFNIe9bnlsJVkLSN4AzM1Z9LCK+Our2ZMnZxhtYvDf/KxFxQNLfBx6S9KPk6F9Ku4DPAp+i+8X8FN200vuHeb0i2rXwfkn6GHAU+FKf/6bw96tpJL0WuA/4vYj4257Vj9BNT/y/ZPzlK8AFI2hWbT+XZAzuncDtGaurer+OExEhqbRa99oG+oj4zWU8Lc9NUl6ge9q4IumJLXYjlaHaKGkF8B7g0kX+j4Ubszwv6S/ppg2G+oLkfe8kfQ64P2PVIDebKaxdkm4CrgeuiiRBmfF/FP5+Zcjz9y9sM5t8zq+nu2+VStJr6Ab5L0XEX/SuTwf+iHhA0n+XtDoiSp3AK8fnUso+ldO1wCMR8bPeFVW9X4mfSTorIp5L0ljPZ2xzgO44woI1dMcmB9K21M0OYFNSEbGe7pH5e+kNkgCyE/inyaIbgbLOEH4T+FFEzGatlPQLkk5deEx3QPKJrG2L0pMXfXef18t1s5mC27UB+HfAOyPicJ9tRvV+5fn7d9Ddd6C7L32r38GpKMkYwP8AfhgR/6nPNmcujBVIuozud7zUA1DOz2UH8NtJ9c3lwM9TaYuy9T2rruL9SknvQ/3i0IPA1erexGkV3ff2wYFfqezR5jJ+6AaoWeBl4GfAg6l1H6NbMbEXuDa1/AHgTcnj8+geAGaA/wWcVFI7Pw98sGfZm4AHUu14LPnZQzeFUfZ790XgB8DjyY52Vm+7kt+vo1vV8fSI2jVDNxf5aPJzR2+7Rvl+Zf39wBa6ByKAk5N9ZybZl84bwXv0K3RTbo+n3qfrgA8u7GfArcl78xjdQe1/MoJ2ZX4uPe0SsDV5P39Aqlqu5Lb9At3A/frUspG/X3QPNM8BrySx6wN0x3S+CTwFfAM4Ldm2A/x56rnvT/azGeB3lvP6ngLBzKzl2pa6MTOzHg70ZmYt50BvZtZyDvRmZi3nQG9m1nIO9GZmLedAb2bWcv8fq2pmMcNZKTkAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from numpy.random import normal\n",
    "data=absmodel(nu,2.0,1.0,2.0)+normal(0.0,0.01,len(nu))\n",
    "plt.plot(nu,data,\".\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's optimize the multiple parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We define the objective function as $obj = |d - f|^2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:07.316921Z",
     "iopub.status.busy": "2022-10-20T05:48:07.316618Z",
     "iopub.status.idle": "2022-10-20T05:48:07.317971Z",
     "shell.execute_reply": "2022-10-20T05:48:07.318232Z"
    }
   },
   "outputs": [],
   "source": [
    "# loss or objective function\n",
    "def objective(params):\n",
    "    a,beta,gamma_L=params\n",
    "    f=data-absmodel(nu,a,beta,gamma_L)\n",
    "    g=jnp.dot(f,f)\n",
    "    return g\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:07.320127Z",
     "iopub.status.busy": "2022-10-20T05:48:07.319790Z",
     "iopub.status.idle": "2022-10-20T05:48:07.320959Z",
     "shell.execute_reply": "2022-10-20T05:48:07.321203Z"
    }
   },
   "outputs": [],
   "source": [
    "# Gradient Descent"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:07.323961Z",
     "iopub.status.busy": "2022-10-20T05:48:07.323658Z",
     "iopub.status.idle": "2022-10-20T05:48:11.061847Z",
     "shell.execute_reply": "2022-10-20T05:48:11.062106Z"
    }
   },
   "outputs": [],
   "source": [
    "gd = jaxopt.GradientDescent(fun=objective, maxiter=10)\n",
    "res = gd.run(init_params=(1.5,0.7,1.5))\n",
    "params, state = res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:11.064618Z",
     "iopub.status.busy": "2022-10-20T05:48:11.064301Z",
     "iopub.status.idle": "2022-10-20T05:48:11.066350Z",
     "shell.execute_reply": "2022-10-20T05:48:11.066576Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(DeviceArray(1.9579332, dtype=float32, weak_type=True),\n",
       " DeviceArray(1.0382165, dtype=float32, weak_type=True),\n",
       " DeviceArray(1.8850585, dtype=float32, weak_type=True))"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:11.071839Z",
     "iopub.status.busy": "2022-10-20T05:48:11.071551Z",
     "iopub.status.idle": "2022-10-20T05:48:11.180374Z",
     "shell.execute_reply": "2022-10-20T05:48:11.180644Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fb90cf3d490>]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA0AklEQVR4nO3deXhTZdr48e+dpC2UpSwtUGihrLKvpaiIoKACCiiOCjrjMjOv4zo6Or95cZzF19F3Vt9ZXQZnXEdxVxBQRAFFFsu+bwUKbVla9r1tkuf3R041xBTSNunJcn+uq1fTsyR3k5M7T+7znOcRYwxKKaXil8PuAJRSSkWWJnqllIpzmuiVUirOaaJXSqk4p4leKaXinMvuAAKlp6ebnJwcu8NQSqmYsmLFigPGmIxg66Iu0efk5LB8+XK7w1BKqZgiIruqW6elG6WUinOa6JVSKs5poldKqTiniV4ppeKcJnqllIpzmuiVUirOJXaiL8qHhU/5fiulVJyKun709aYoH14eD54KcCbDbTMgO8/uqJRSKuwSt0VfuNCX5I3H97twod0RKaVURCRuos8Z5mvJi9P3O2eY3REppVREJG7pJjvPV64pXOhL8lq2UUrFqfO26EXkBREpFZH11awXEfmbiBSIyFoRGei37jYR2Wb93BbOwMMiOw+GPaxJXikVfv6dPWzu+BFKi/4l4B/AK9WsHwN0tX6GAM8CQ0SkBfBrIBcwwAoRmWGMOVzXoJVSKqr5d/ZwOAEBr9u2jh/nbdEbY74ADp1jkwnAK8ZnKdBMRDKBq4C5xphDVnKfC4wOR9C20e6YSqlQnNXZo9L2jh/hqNG3A4r8/i62llW3/FtE5E7gToD27duHIaQI0O6YSqlQVXX2CNait6HjR1ScjDXGTAWmAuTm5hqbwwkuWHdMTfRKxb6i/Lp1ygi2f2BnD7C140c4En0JkO33d5a1rAQYEbB8QRgezx7+n9DaHVOp+FDXb+rn2j877+z7srFhGI5+9DOAW63eNxcCR40xe4E5wJUi0lxEmgNXWstiU9Un9OWPatlGqXhR1wsn67p/PZ33O2+LXkSm4WuZp4tIMb6eNEkAxpjngNnAWKAAOAXcYa07JCK/AZZZd/W4MeZcJ3XtFcrXt8BPaKVUdKht+aWu39Trsn89nvc7b6I3xkw+z3oD3FvNuheAF2oXWj3SE61Kxa66vH/reuFkXfavx/N+UXEyNiJq8gmvJ1qVil11ff/W9Zt6bfevx/N+8ZnoAz/hR/8OTh+sPunriValYlesvn+D9cxZ+FREeuaIr/ISPXJzc83y5ctrvJ/Xa5i2bDdtmzWk745/0yL/j4jxAA5wOMCYc3+tq2sXK6WUfaLl/VvDOIwxHD5VyeHNX5IzezLiqcDhSqlV+VhEVhhjcoOti5sWfenxch593zccz0BpwGvJTpLEAILD68GBweuuoGTVJ6Sk9SWjSQoi8s0d1PTrV7QcWEqp6t+/9fk+reZcgdvjpfjwaQoPnmT3oVMUHjjF7kOnKD7s+32qwsM9zuk85KrAKd6IlI/jJtG3apLCkkcup+TwaUqO9GdWYQca71tC8ZmGfPfIs7iMm0qcPLCkESsXf0ajZCcdMxrRMb0xnTMa0TmjMV1aNaZTRiNSXM5v7jjYgaInb5WKfvX8Pj2z7XNS3OUIXjzuCj54dxpPu33J3O39pnLSMMlJ+xapZLdI5aLOLclqnkovjwPHF9Mx3kokAuWnuEn0DoeQmdaQzLSG5AL0nwhM9K0sugbPji840nIwDyT3pPDASXZaP6uLDjNz7R6qKlhOh9ChZSoXtG7CiNRCvrPhHhzeSnAmI1UHSn2dvNVvDUrVXoTepxVuL9tKj7N573E27zvG5n3H2br/OO2Ou3gt2UUSvkbll+4edM9swujebchJb0TH9EZ0aJlKRuOAagIAHaHThxF7v8dNoj+n7Dyc2XlkApnA8G4ZZ60+U+lhR9lJCspOsG2/70XbvO84HY/MwTgrEPHidpfz9luvUXBBY4anducSZxLiISKfvoB+a1Cqrup6krYoH/eOLyhsPIAllV1YV3yE9SXH2FZ6nEqPr2WY7HLQrXVjhnZJp1vr8Wxw9qDLydU06X4Zf+4wpGaPF8HrdBIj0Z9HgyQnPds2pWfbpmctP7MjBcdr0/F6KvGKi9XOPnywdBf/dsNAmcLwlK0cbTmEpluaMfB0GQM7NKdxSpie0sDWyJrXtXWvVE3Uoo/73qOnWV54mLKNC7lly304jZt2uHi/4ucUpvamV9umXNqtky9fZDYlp2UqLqf/AAOdicZBejXRn0ODThfB7b6vU8k5w/h9dh5Perxs3X+CtcV9WFN8hJW7jtBo3iwqZBN/NT3omN6YMY0LaNL9Mi4YPJK0hkm1e/DA0e9WvW7reNZKxaTznKTd32Iwn5/qyNKdB8nfeYjiw6cB+HHyPFwON07x4sDDi5eV0/SKUUFKLrEhbrpX2qYoH/PyOHBX4BEnXgMO46ESF9+t/DnutoO5pGs6l3bNYGCH5iRVffqHUn+v2uZoMax42de6F6dvvJ1hD9ff/6hUnDh6upKNX81l0Oe34zCVVBoXt1T8nF2pvRmc04LBHVswOKc5PT2bcb16bUyVThOie2XIwn2Cs3Ah4qkEvLiqPjTF4BQPD3Yt5a9nHDz3+Q6enr+dJikuhnVL54ZWexnx1Q98+53rIKpqjRTlw+ppwWuNesJWqWoZY9hedoJPN5Uyb1MpK3Yf5kfyAYNdlbjEi0M8/PPSM6SPDmytXxi5OaVteM8mVqKPxAnOaiYYEGcyw0Zdx7DsPI6ermTJ9gMs2FLGvM2l5GycwTBXBS7x4nVXcHLzfJqcK47qao16wlapbzG7v2L/2k/59HQ3/r0rg50HTgLQM7Mpdw/vzNjmN+L8ZAZ4KnA4k8noPRKClWQicXLUpvdsYiX6SHS3CmGCgbSGSYzuncno3pl4vYadqxyYWR/g9lZSaZzcNi8Z5/bFTOjfjqv7ZNK8UXLwxwmMVcfoUbEszC3bgtLj5C/8mInr7iHduLkeF9tb/4HvX3s5I7u3om2zhtaWF0BmNa31SLe2bXrPJlaiD7W7VU1f7BpMMOBwCJ0HXQ6tZkLhQvY1G8TlZe2YvnoPv/hgPY/N2MBl3VtxU242Iy7ICDijX8v/R6loE46WbVE+Z7Z9zrzybjxT0IL1Jce41zWXJJfvJKpTPPy6zyG4sMO39w3WcKqP1rZN79nESvShdLeqr69W1oGWBdwH3HtZFzbuPcYHq0p4f1UJczfuJ6NJCjfmZnHLkA5+rZEa/j9KRaM6tGyNMWxZ/hmdZt+My1vJZbj4tOmTXHfNFVyX0RTn276yTI2vcamP1rZN79n4SvThmDwk3C92iN8ORIRebdPo1TaNn43uzvzNpby5rIhnF2znuc93cEWP1tw+NIchHVvUbYwepaJBLVq25W4P01fv4ZUlhVy67y0esk6oOsXD/+Udh0s6Ah1rn0jrq7Vtw3s2fhJ9uFri4XyxaxlTktPBlb3acGWvNhQdOsVrX+3mzWW7+XjDPvplN+OuSztxZa82OB1h7tOrPXhUfalBy/ZEwWLWL5rF87vb8tnJHC5o3YR+w67BubyalnttE2kcf0OOn0QfrpZ4OF/s2sQUkGyzW6QyZUx3HhzVlXdWFPP8wh3c/dpKOqU34scjuzKuX9vwJHztwaPq23kS8tFTlcz66AOuW3sPubgZ4Ehi0/hX6XfRMN+32l5twp+UQx0FM8YaRfGT6MPZEg/XV6uaxnSOZNsgycl3L+zA5Lz2fLx+H3+ft40H31zt+z2qG9f0zaz+qr1QDkrtwaPs5HeMnmg1kH8t3MG/F+7ke+65JCe5ceLFhZv+nvUgV/n2qa8SSLCJjD6eElONovhJ9NH4taumMYWQbJ0O4eq+mYzp3YaPN+zjr59u4/5pq/jXlzt5dGwP8jq2OPs+Q22paw8eZRfrGDWeCjySxP3yS46equTxzCIu7DUA55IZ9h6Xge/LTdNjrlEUP4keovPEpH9M52tZ1yDZOhzC2D6ZjO7VhvdWlfCnOVu48Z9LuKpXa355TU+yTqz/ZviEUA7KaPygVImhcCHGGscdr2Fy6iJGMQ/H4UpYknz+qUBDVdtyS+D7sscE2LUkphpF8ZXoo1koLetaJFuHQ/jOoCyu7pPJvxbu4JkF2/np/z3Pq0lP4jJuxOEEhwu8nP+gjMYPShXX9h87w3+2teYe4yJJ3IgrmSt6tkZWVn7TODl9sO5jO9XlHFSw92XrnjHVKNJEX19CrYHXMtk2THZyf7fD3OpZw+r1a5GjlYh4MV6QQbdCWnbMHJQq/hljeD1/N7+dvZkKT0uyBj7LxBaFuDpf6ttgzRvhbTHX9RxUsIsiY+i9pIm+vkS6Bm61WNI8FQx3OPE4Xbi9biq9Tt47M5TvjJ549hSJSkVSdWWSonyObZrPnwta8eLuVgzt0pInr+1DTnqjs/ev7pttuMovMVBuCScdprg+nePgr/PXwIVPwbwnvxnKeNCtnGnUjn8VteVPG5vRvU0T/j55AF1bNwnP/6JUdaorkxTl435xHHgqqMTFggv/zejR40If472uXYBjrEtkTekwxdEikuNrBLZY+t1Mg+w87gN6bt7Pz95Zy7h/fMnjE3pzw6CsmJ1AQcWAIGWSM20GsWDm24zyVHx9NeuYxtuCjxpZg/utU/klgZxjxCxVL4IdvLVRdcLo8ke/9WFxeffWzP7xMAZkN+dn76zl4bfWcLrCE6Z/QKkAVY0OcYIzmeJmg5jwj0VM3d0W40zCiLN2cy0H3G+ilV/qQks3dqvHK1I9XsPf523jr59to1fbpjx/ay6ZaUEGS1OqrqwyySpHb56au5UhspFLRl3LgPbN61Y+ifPyS12cq3SjiT4a1PPB+9mm/TzwxmoaJjuZ+r1BvjefUmFkjOHFRYXMmv0BryX/Lyni8bXiY+Aq0lh1rkSvpZtokJ3n6ydcT2+AkT1a8949F9MwyclNU5fy0bq99fK4KjF4vIZfTl/P4zM38t02Rb4kX9fSpKqTkBK9iIwWkS0iUiAiU4Ks7yAin4nIWhFZICJZfus8IrLa+pkRzuBV7XVr3YTp9w6lT7s07nl9Ja9/tdvukFQsKsr39fgqygd8QwnfP20l/1m6mx8N78SECTf5WvJaV7fVeUs3IuIEtgJXAMXAMmCyMWaj3zZvAzONMS+LyOXAHcaY71nrThhjGocaUEKWbmx0usLDPa+tYP6WMn56ZTfuvayL9shR51ZVamzY8qzBvU5Nfp//mi8sKjjIL67uwQ+HdTp7e62rR1Rdu1fmAQXGmB3Wnb0BTAA2+m3TE3jIuj0f+KDW0ap61TDZydRbc/nZO2v50ydbOVXh4f9ddYEmexWcf+cBETBeMF6Mp4J335vG0sOjeeqGflw/KOubfRK4W2O0CCXRtwOK/P4uBoYEbLMGmAj8FbgOaCIiLY0xB4EGIrIccAO/M8Z8EPgAInIncCdA+/bta/o/qDpKcjp46oZ+NEhy8syC7bicDh66opvdYaloUtUq9x8kzzjA4cAgVBgX04904umbBzC6d6bd0aoA4bpg6qfAP0TkduALoASo6qjdwRhTIiKdgHkiss4Ys91/Z2PMVGAq+Eo3YYpJ1YDDITx5bW88Xi9/+2wbSQ7h/pFd7Q5LRQP/VnzAIHnlVzzJe1+u5Z2DOdwx6UZN8lEqlERfAmT7/Z1lLfuaMWYPvhY9ItIYuN4Yc8RaV2L93iEiC4ABwFmJXtVQXWqe55gpx5Gdx28n9sXtMTw1dyupKS5+cEnHyPwPKnb4X9TnBaxB8tzth/LDT4VFB9rw55v6c03ftnZHqqoRSqJfBnQVkY74Evwk4Gb/DUQkHThkjPECjwAvWMubA6eMMeXWNkOBP4Qx/sRTlwusQpgpx5mdxx9v6MepCg9PzNpIZloDxvbRVlpCCzK8hskazH+/vZaF24r5/fV9mNC/nd1RqnM4b/dKY4wbuA+YA2wC3jLGbBCRx0VkvLXZCGCLiGwFWgNPWst7AMtFZA2+k7S/8++to2qhLkMmhDJTDr5ZrP4yqT8D2zfnwTdXs6zwUIT+GRUTggyv8edPt/HuymIeGNmVmwbrebVoF1KN3hgzG5gdsOxXfrffAd4Jst9ioE8dY1T+6jLcag1mymmQ5OT5W3P5zrOL+a9XlvPu3RfTOSPkXrIq3vj1nHlrWRF/+2wbNwzK4sFReh4nFugQCLEoQjX6YPe1++AprntmEWmpSUy/dyhNGiSF6Z9QsWjFrsNMmrqECzu15IXbB5Pk1Ivro4WOdaPqZOmOg9zyr68Y2b0Vzw334Nj9pV78koBKj59h3N+/pEGSkxn3XkJaqn7oRxMdj17VyYWdWvLzsT2YNesDPIW/xWHcER9pU0WXSo+Xe19bybHTbl66I0+TfIzR713q/Iry+b73PR5qvRLxVNZ97HwVXQLGqwnmyVmbWFZ4mN9d34cemU3rMTgVDtqiV+dmdckUTwVDHU4qxYkbcDqTEP8TwTqeSWwK1uX29MGzXsfPNu3npcWF3DE059zdKPUYiFqa6NW5+XXJFC+c7HEzL6yvpLzNxTySNRiBep08RYWZf5dbdznMfhiM+fp1LGvWj5+9s5YemU2ZMqZ79fejx0BU09KNOreA6duaX3wrLUZPYWphOv+pGto4XNMhqvrn//o6HNYgZb7X0excyP97Zw0nyt38bVJ/UlzO6u9Hj4Gopi16dW5VF8v4fSW/PcuwYEsZT8zcyEWdWtClLn37lb38X9+AYYc/OtGFBVvKeHxCL7q2bnLu+9FjIKpp90oV3HnqraXHzzD6LwvJbt6Q9+4ZirNkmdZn44H1uu9tlsuDb67mhvRCrp84CWkfOGBt9fvqMWAP7UevaibEeuuMNXv48bRV/PKanjr4WRwxxvDYMy8ypfS/aeDQuV5jhc4Zq2omxHrruL6ZXHZBBk99soWiQ6fqOUgVKe+sKKbRnqWkiFvneo0TmujVtwWcgK2u3ioiPHGdbyijX3ywnmj7dqjOoZq+82XHy3li1iaOth6CuFJ0rtc4oSdj1bcFOQEbVFE+7QoX8ochnblvYRkz1uzR4WpjwTlKc4/P3MjpCg93TLoJKe+jNfc4oYleBXe+eT79ksXVzmS+aPMET8xKYWSP1jRO0cMqqgUrzWXnsbjgAB+u2cNPRnWjS6vGgM71Gi+0dKNqx/9CKk8FD3bZT9nxcp6eX2B3ZOp8gpTm3B4vj8/cSFbzhvxoeCe7I1Rhpk0vVTsB/abb9ruCiaeS+ffCnUwe3J72LVPtjlBVJ0hp7o2lu9i87zjP3jKQBknnuDBKxSRN9Kp2giSL/047w8fr9/Hk7I3883tBe3mpaOFXmjt6qpKnPtnCkI4tGN27jc2BqUjQ0o2qvew8GPbw1wmjddMG3DOiM3M27Gfx9gM2B6dC9dfPtnHkdCW/GtcTEbE7HBUBmuhVWP1wWCfaNWvIEzM34fVqd8toV3jgJK8sKWTS4Gx6tU2zOxwVIZroVVg1SHLy06u6sXHvMWav3xvSWOfKPn/5dCsup/CTUd3sDkVFkCZ6FXbj+7Wja6vGfPzxDMzL42Hek76umJrso8qWfceZvmYPt12cQ6umDewOR0WQJnoVdk6H8PCV3Wh/dCXGXa5D10ap/5u7hcbJLu66tLPdoagI00SvIuKqXm0obTmYClwYvYw+6qwtPsKcDfv54bBONG+UbHc4KsK0e6WKCBFh3DXXcvOLJ/hl70MMuHScXmUZRf70yVaapybx/Uty7A5F1QNt0auIubRrOq4OF3LXruGUZw6yOxxlWbHrMF9sLeOu4Z1p0iDJ7nBUPdBEr8LP6mkjxcu4f2QX9h8r572VJXZHpSzPzC+geWoS37uog92hqHqipRsVXgEjI15y63T6tEvjn59v58bcbJwOvSDHTpv3HeOzzaX8ZFQ3UpP17Z8otEWvwitgZETZ9SX3jOhM4cFTzF631+7oEt6zC7bTKNnJbRdraz6RaKJX4RVkZMSrerWhU0YjnlmwXScnsdHug6f4cM0ebh7Snmap2tMmkWiiV+FVNdjZ5Y9+PaGFwyHcNbwzm/YeY8HWMrsjTFj//GI7LoeDHw7TYYgTTUiJXkRGi8gWESkQkSlB1ncQkc9EZK2ILBCRLL91t4nINuvntnAGr6JUwGBnANf2b0dmWgOenb/dxsASV+nxM7y9opjrB7WjtV4Fm3DOm+hFxAk8DYwBegKTRaRnwGZ/Al4xxvQFHgd+a+3bAvg1MATIA34tIs3DF76KFckuB4/0OUZu0YtsXzHP7nASzn+W7qbS4+W/tDWfkEJp0ecBBcaYHcaYCuANYELANj2BqnfvfL/1VwFzjTGHjDGHgbnA6LqHrWJOUT7jVt/FQ663yZ45Sce9qUflbg+vf7WLyy5oRaeMxnaHo2wQSqJvBxT5/V1sLfO3Bpho3b4OaCIiLUPcVyWCwoWIpxKXeHF4KzmxZb7dESWMmWv2cuBEBXcMzbE7FGWTcJ2M/SkwXERWAcOBEsAT6s4icqeILBeR5WVlerIuLlm9cYw4qcTFrGM6kFZ9MMbw4uKddGnVmEu6pNsdjrJJKIm+BMj2+zvLWvY1Y8weY8xEY8wA4FFr2ZFQ9rW2nWqMyTXG5GZkZNTsP1CxweqNI5c/yp8y/8gfNzaj3B1yW0DV0vJdh1lfcozbL87R2aMSWCiJfhnQVUQ6ikgyMAmY4b+BiKSLSNV9PQK8YN2eA1wpIs2tk7BXWstUIrJ641x6+dUcOFGuF1DVg5cWFdK0gYuJA7VimsjOm+iNMW7gPnwJehPwljFmg4g8LiLjrc1GAFtEZCvQGnjS2vcQ8Bt8HxbLgMetZSqBXdo1nc4ZjXhxUaFeQBVBZRsXkrPpOX7a46gOd5DgJNreaLm5uWb58uV2h6Ei7NUlhfxy+gZm3DeUvlnN7A4n/hTlU/niNYinEkdSMo7bPtRhouOciKwwxuQGW6dXxipbTBjQjoZJTqbl77Y7lLjk2fEF4rV6OXkqdXavBKeJXtmiaYMkxvXLZPrqPZwod9sdTtxZIb2pNC686OxeSocpVnYpyufBlE8oqGzE9NU9uGWIjqYYTs9ub4HD9RhTLz0NnS7Vsk2C00Sv6p81Zn2mp4LXU5z8YlEjbhlyp91RxY2SI6dZsLWM+y4bhXP4Bd+sKMr3lXByhmniTzCa6FX9s8asF+MhGWh1cDnrim+iT1aa3ZHFhTeX+S5GvzHX7xKWgAlhqkYWVYlBa/Sq/vmNWS+uJFY6evG6npQNC7fHy1vLiri0awbZLVK/WREwIYyenE0smuhV/fMbs15u+5CsviOYsbqEk3pStnasOXopymfBljL2HTvD5Lz2Z28TZEIYlTi0dKPskZ33dengJs8h3llRzEfr9/GdQVnn2VGdJaAks7fZ3fy04T5GNWkJtPlmu6oPV63RJyRN9Mp2uR2a06FlKu+uKNZEX1N+JRnjLmdS2V9xisHx6vvfrsP7fbiqxKKlG2U7EWHigCyW7DhI8eFTdocTW/xKMl4RHHhx4NU6vDqLJnoVFaoG3Vqy4KOv680qBH7nO55NvZtKSdI6vPoWLd2oqJDdIpXvZe1j3JopGPEg2gUwdNl5FKT04E+zvqD9JbmMT9uhdXh1Fk30KmrcmL4LV5kb8S89aLIKyTsrSnA6hIuGj4UmKXaHo6KMlm5U1Og8eAyVuPDg0NJDDXi8hvdXFTOiWwYZmuRVEJrolf2sfuCpyU6m5vyZv3MT5be8r635EC0qOMD+Y+Vcrz2WVDW0dKPsFdAPfMTlrzBhyzi6nchhrN2xxYgPVpXQtIGLkT1a2R2KilLaolf2Crg0v497HRlNUpi++ltTCyt/1reg8h1LmLNhH2P7ZJLictodlYpS2qJX9qrqB2616B0dh3HN0Ua89tVujp6uJK1hkt0RRh+/b0EuRxIXVE5hfL+gEwspBWiLXtnNrx94VXfKCf3bUeH2MmfDPruji04B34JGNtzGkE4t7Y5KRTFt0Sv7BVya3y8rjQ4tU5mxes/ZQ+0qH+tbkPFUUOF10rDbcJwOsTsqFcW0Ra+ijogwvl9bFm8/QOnxM3aHE32sb0Hrut3HLRU/Z9DQq+yOSEU5TfQqKk3o3xavgVlr99odSnTKzuP3J8dyqEV/+uqELeo8NNGrqNSlVRMmZpTgWvxnHfcmiNJjZ1i8/SDj+7dDRMs26ty0Rq+iU1E+vz/5C8RTifflN3Hc9qFeQOVn5tq9GAPj+7W1OxQVA7RFr6JT4UJcxo1LvOCu1CF3A8xcu4cemU3p0qqx3aGoGKCJXkWnnGGIMxk3Dipw6rg3fkqOnGbl7iNc0zfT7lBUjNBEr6KT1bNkVed7uLn85+xs2MvuiKLGR+t8J6iv7qOJXoVGE72KXtl5tB33C1aabsxep71vqsxcu5debZuSk97I7lBUjNBEr6Jau2YNGdi+GTO1myUAxYdPsbroCFdr2UbVgCZ6FfWu7tuWTXuPsaPshN2h2O6jdb5hIbRso2pCE72KemP7tAHQ8g0wc91e+rRLo0NLLduo0IWU6EVktIhsEZECEZkSZH17EZkvIqtEZK2IjLWW54jIaRFZbf08F+5/QMW/zLSG5HZonvDlm6JDp1ijZRtVC+dN9CLiBJ4GxgA9gcki0jNgs18AbxljBgCTgGf81m03xvS3fu4KU9wqwYztk8nmfccpKE3c8s1s7W2jaimUFn0eUGCM2WGMqQDeACYEbGOAptbtNGBP+EJUypfo4ZuuhYlotlW2yW6RancoKsaEkujbAUV+fxdby/w9BnxXRIqB2cD9fus6WiWdz0Uk6FUvInKniCwXkeVlZWWhR68SRpu0Bgxs34yP1ifmGPXFh0+xpvgoY6zzFUrVRLhOxk4GXjLGZAFjgVdFxAHsBdpbJZ2HgNdFpGngzsaYqcaYXGNMbkZGRphCUvFmbJ9MNu49xq6DJ7+eSi9RBjz72PqAG9Nbyzaq5kJJ9CWA/+wPWdYyfz8A3gIwxiwBGgDpxphyY8xBa/kKYDvQra5Bq8Q0urevNbty0RzfVHrznvT9ToBk//H6ffTIbEpHvUhK1UIoiX4Z0FVEOopIMr6TrTMCttkNjAQQkR74En2ZiGRYJ3MRkU5AV2BHuIJXiSWreSp9s9I4vnnBWVPpxfuAZ/uPnWH5rsOM6a1lG1U75030xhg3cB8wB9iEr3fNBhF5XETGW5s9DPyXiKwBpgG3G2MMcCmwVkRWA+8AdxljDkXg/1CJoCifnzf5iI1HXHidSSBO38TicT7gWdXcuWO1Pq9qKaTx6I0xs/GdZPVf9iu/2xuBoUH2exd4t44xKuUrz7w8niGeCvolOVnU5WcMa+fwJfk4H6d+9rq9dGnVmC6tmtgdiopRemWsig2FC8FTgRgPSeJmd3ExDHs47pP8gRPl5O88xNjebRLuBLQKH51hSsWGnGG+Mo2nAiMu3j2Yw6hjZ2jdtIHdkUXUyi/ncJdjJpOcA+Dl//Gdk3Amw20z4v5DToWPJnoVG6zx6SlcyP60Qax8/SQfr9/HbRfn2B1Z5BTlM/yrHzIyqRLHonfBeH0/VSegNdGrEGnpRsWO7DwY9jBZfUfQpVXjr/uWx6sz2z7H6a3EiRfxekEcCXMCWoWXtuhVTBrTuw1Pzy/g4IlyWjZOsTuciFjs6cFFuHCKB3Emw+jfwemDCXECWoWXJnoVk0b3bsPf5xUwd+N+JmXu85Uy4iUBFuVD4UKW7kjnjZT/4Z+XnIaOcfK/KVtoolcxqWdmU9q3SGXrinlwcEr8nKS0upEaTwU/8TqZ1v0fyKU32R2VinFao1cxSUQY07sNqXsWY+LpKln/bqS4ubLRNrsjUnFAE72KWaN7t2GxuwceiaOrZK1upB4cuMVF235X2B2RigNaulExq19WM/Y27csfW/6BR3ociI8afXYe5Te/z9MvvUSjC0bwow5D7I5IxQFt0auY5XAIV/Vqw0u7W3Ey74HYT/KW+ac68reK8fQeoq15FR6a6FVMG927DeVuL/O3lNodSth8vH4vzVKTyOvYwu5QVJzQRK9i2uCcFqQ3TuajdfFx8VS528Nnm0q5okdrkpz69lThoUeSimlOq3wzb3Mppys8dodTZ19uO8Dxcjdj++pMUip8NNGrmHd1n0xOV3r4fGvsl29mr9tH0wYuhnZOtzsUFUc00auYl9exBS0aJTM7xss3FW4vczfu44qebUh26VtThY8eTSrmuZwOrurVms827edMZeyWbxZtP8CxM26dSUqFnSZ6FRfG9snkZIWHL7aW2R1KrX20bi9NUlxc0lXLNiq8NNGruHBhp5Y0S03ioxgdurjS4+WTjfsZ1bM1KS6n3eGoOKOJXsWFJKeDK3u25tON+yl3x175Zsn2gxw5VcmY3lq2UeGniV7FjTF9Mjle7ubLbQfsDqXGNnw1lweSZzA8dafdoag4pIlexY2hndNp2sDFrLV77Q6lRty7lnLH9gf4seMtUl67Tif/VmGniV7FjWSXgzs7HaD9xuco37nE7nBCtnvlJ7iMGyfe+BhqWUUdTfQqfhTlc/euh7iPN3H959qYaRnPOtaZSnFh4mWoZRV1dJhiFT8KF+LwViLixVPVMo62ES2taQKrhlQud3uYuiMdOv2F+zvti4+hllXU0USv4kfOMMSZjMddToVxQduLaWh3TP6saQL9pz38/Fh7jpe76XvRldAtw+4IVZzS0o2KH9l5cNsM9gx8mFsqfs4nx9vbHdHZrGkC/ac9nLl2Ly0aJXNx55Z2R6fimCZ6FV+y82h3zaPsadKXD9dEWe8ba5rAqmkPz7S7mLkb9zOmdxsdklhFlJZuVNxxOIRr+mby8pJCjp6qJC01ye6QfKxvHFU1+k8PZ3G6chXj+rW1OzIV57QZoeLSNf3aUukxzNkYZUMiZOfBsIchO48P1+yhVZMUBufoTFIqskJK9CIyWkS2iEiBiEwJsr69iMwXkVUislZExvqte8Tab4uIXBXO4JWqTr+sNNq3SGX66hK7Qwnq6KlKDm9ZxG8z5uIsWWZ3OCrOnTfRi4gTeBoYA/QEJotIz4DNfgG8ZYwZAEwCnrH27Wn93QsYDTxj3Z9SESUiXDegHYu3H2TPkdN2h/MtSz7/iJedT3D53ud9PXFipM+/ik2htOjzgAJjzA5jTAXwBjAhYBsDNLVupwF7rNsTgDeMMeXGmJ1AgXV/SkXc9QOzMAbeXxV9rfoD6z8lWdyIXw8cpSIllETfDijy+7vYWubvMeC7IlIMzAbur8G+iMidIrJcRJaXlcXueOIqurQ/tZ7/Tf+ETcs+xRhjdzhf21F2gvcOdcQ4kr7ugaNXw6pIClevm8nAS8aYp0TkIuBVEekd6s7GmKnAVIDc3NzoeUeq2GVdnDTJXc51xsXWFV24IHeU3VEB8N7KElbTjWM3vkuLsny9GlZFXCiJvgTI9vs7y1rm7wf4avAYY5aISAMgPcR9lQo/6+IkB16ScLN7xSdRkei9XsP7q0oY1jWDFt3zoLu25FXkhVK6WQZ0FZGOIpKM7+TqjIBtdgMjAUSkB9AAKLO2myQiKSLSEegK6FknFXl+Fyd5HUm8sjc7KiYkWbrzICVHTnP9oCy7Q1EJ5LwtemOMW0TuA+YATuAFY8wGEXkcWG6MmQE8DDwvIj/Bd2L2duMrim4QkbeAjYAbuNcYY/+7TcU/v4uT1jt6s/BDN/M2lTKmT6atYb27ooQmKS6u7Nna1jhUYgmpRm+MmY3vJKv/sl/53d4IDK1m3yeBJ+sQo1K1k50H2Xn08xpaLfiMt1cU25roj5+p5KP1exnfry0NkrSXsao/emWsintOh3DT4Gzmbyml6NAp2+L4YFUJpyo8TM6LssHWVNzTRK8SwuS89gjwev5uX4+chU/V60VKxhheWbKLvllp9MtuVm+PqxTooGYqQbRt1pBRPVqzKf9TzPInEL8x4euja+NXOw+xrfQEf/hO34g/llKBtEWvEsb3LupAz/K1GHf5WWPC14dXl+4irWES4/rqSJWq/mmiVwljaOd0djcdSCWuer0itfTYGeas38eNuVk0TNaTsKr+aelGJQyHQxgw9ComzzrNM5ecok3fK8JbtgmYD7bKtPwi3F7DLUM6hO+xlKoBbdGrhPKdgVlsdHXnL2fGhT/Jvzwe5j151miUlR4v0/J3c2m3DHLSG4Xv8ZSqAU30KjFYPW3SDq5i4sAs3ltZwr6jZ8J3/0HmgwXfyJltj6/l180+1qGIlW000av4F9DafqDbYbzG8Nzn28P3GAHzwZIzDLfHy4JPZzIt5X/ptO4vOu68so0mehX/AlrbrQ8t47oB7ZiWv5vS42Fq1VcNuXD5o1932ZyxZg85x1eRhI47r+yliV7FvyCt7Xsv60Klx8vzX+wI3+P4zQfr8Rr+Mb+Avc1zEVeKjjuvbKW9blT88xvgrKpHTA5wbf92/Gfpbu4a3pmWjVPC+pCz1u1lR9lJHr55HNK8f9DeOErVF030KjFYA5z5u+eyLry/uoTnF+5kypjuYXsor9fwj3nb6NqqMWN6twFHpiZ4ZSst3aiE1aVVY8b1bcuLi3ay+2D4Bjt7Y1kRW/ef4P6RXXE4JGz3q1RtaaJXCe2Rsd1Jcjp49IN1YZlXtvT4GX730SYu7NSCcX3tHfteqSqa6FVCyzy2jpe7LuRkwWJmrNlT5/t7YuYmzlR6efK6Pohoa15FB63Rq8Rl9a8f6KlgWoqTuz5MYkS3H5GWmhT6/n4nWRdsKWXGmj08OKornTMaRzZ2pWpAE71KXFb/ejEekgVGlc9j8cu7GHPNDec/eVp1EZY13PHpye/zy+mn6ZTRiLtHdK6f+JUKkZZuVOLy618vDic3JX3BFfv+hfvFa86+gjXYRCV+F2EZTwUfzniLksOn+e11fUhx6QiVKrpoi14lLv/+9UeLca54GREvbk8lGxbPotdNed9quTP6d3D6IDRsCc5kjKeCSly8Udqe307sw5BOLe3+r5T6Fk30KrFV9a8vykdWT8N4KvDg4rG1zbm3fykjSv2GT3CXw+yHwRhfkh/9W+at3MzTO9sw6qpx3DRY54JV0UkTvVLwdeteChfibnsRJ2d6+OHLy5nSux0/cCYhHkAEjBeMF+Op4K0v1vDfpVfw/aEduXu41uVV9NJEr1QVq3XfCHhjzCIWz1vA8+szWZD8C+5svwdJbcFFW/+ImEoqvE5mHe/MY+N6cutFOdqVUkU1TfRKBSrKp+lb1zPaU8GVDZN4Iv333Lr1EgAGyiMMdW2mZe/L+fu4iaF3xVTKRprolQrk16PG4YVf9T7IQ9//Ll5jgCtJdjpokKQ9a1Ts0ESvVKCqbpdVPW1yhtE4Rd8qKnbp0atUoCDDGisVyzTRKxVMkGGNayRgeASl7KSJXqlwC7zIyppaUCm76BAISoVbwBy1Ok+ssltIiV5ERovIFhEpEJEpQdb/WURWWz9bReSI3zqP37oZYYxdqegUZI5apex03tKNiDiBp4ErgGJgmYjMMMZsrNrGGPMTv+3vBwb43cVpY0z/sEWsVLTTk7kqyoRSo88DCowxOwBE5A1gArCxmu0nA78OT3hKxai6nsxVKoxCKd20A4r8/i62ln2LiHQAOgLz/BY3EJHlIrJURK6tZr87rW2Wl5WVhRa5UkqpkIT7ZOwk4B1jjMdvWQdjTC5wM/AXEfnW6E/GmKnGmFxjTG5GRkaYQ1JKqcQWSqIvAbL9/s6ylgUzCZjmv8AYU2L93gEs4Oz6vVJKqQgLJdEvA7qKSEcRScaXzL/Ve0ZEugPNgSV+y5qLSIp1Ox0YSvW1faWUUhFw3pOxxhi3iNwHzAGcwAvGmA0i8jiw3BhTlfQnAW8YY4zf7j2Af4qIF9+Hyu/8e+sopZSKPDk7L9svNzfXLF++3O4wlFIqpojICut86LfXRVuiF5EyYFcd7iIdOBCmcMJJ46oZjatmNK6aice4OhhjgvZmibpEX1cisry6TzU7aVw1o3HVjMZVM4kWl451o5RScU4TvVJKxbl4TPRT7Q6gGhpXzWhcNaNx1UxCxRV3NXqllFJni8cWvVJKKT+a6JVSKs7FZKIXkRtEZIOIeEUkN2DdI9YEKVtE5Kpq9u8oIl9Z271pDe0Q7hjf9JtwpVBEVlezXaGIrLO2i/iVYiLymIiU+MU2tprtzjnZTATi+qOIbBaRtSLyvog0q2a7enm+QphsJ8V6jQusYyknUrH4PWa2iMwXkY3W8f9AkG1GiMhRv9f3V5GOy3rcc74u4vM36/laKyID6yGmC/yeh9UickxEHgzYpl6eLxF5QURKRWS937IWIjJXRLZZv5tXs+9t1jbbROS2WgVgjIm5H3xDK1yAb5C0XL/lPYE1QAq+4ZK3A84g+78FTLJuPwfcHeF4nwJ+Vc26QiC9Hp+7x4Cfnmcbp/XcdQKSree0Z4TjuhJwWbd/D/zerucrlP8fuAd4zro9CXizHl67TGCgdbsJsDVIXCOAmfV1PIX6ugBjgY8AAS4Evqrn+JzAPnwXFdX78wVcCgwE1vst+wMwxbo9JdgxD7QAdli/m1u3m9f08WOyRW+M2WSM2RJk1QR84+2UG2N2AgX4Jk75mogIcDnwjrXoZeDaSMVqPd6NBIzqGeW+nmzGGFMBVE02EzHGmE+MMW7rz6X4Rkm1Syj//wR8xw74jqWR1msdMcaYvcaYldbt48AmqpkbIgpNAF4xPkuBZiKSWY+PPxLYboypy1X3tWaM+QI4FLDY/xiqLg9dBcw1xhwyxhwG5gKja/r4MZnozyGUSVJaAkf8kkq1E6mEyTBgvzFmWzXrDfCJiKwQkTsjGIe/+6yvzy9U83Ux5MlmIuT7+Fp/wdTH8xXK///1NtaxdBTfsVUvrFLRAOCrIKsvEpE1IvKRiPSqp5DO97rYfUx9awh1P3Y8XwCtjTF7rdv7gNZBtgnL8xbKVIK2EJFPgTZBVj1qjJle3/EEE2KMkzl3a/4SY0yJiLQC5orIZuvTPyJxAc8Cv8H3xvwNvrLS9+vyeOGIq+r5EpFHATfwWjV3E/bnK9aISGPgXeBBY8yxgNUr8ZUnTljnXz4AutZDWFH7uljn4MYDjwRZbdfzdRZjjBGRiPV1j9pEb4wZVYvdQpkk5SC+r40uqyV2rolU6hSjiLiAicCgc9xH1cQspSLyPr6yQZ3eIKE+dyLyPDAzyKqaTDYTtrhE5HbgGmCksQqUQe4j7M9XEKH8/1XbFFuvcxq+YyuiRCQJX5J/zRjzXuB6/8RvjJktIs+ISLoxJqIDeIXwukTkmArRGGClMWZ/4Aq7ni/LfhHJNMbstcpYpUG2KcF3HqFKFr5zkzUSb6WbGcAkq0dER3yfzPn+G1gJZD7wHWvRbUCkviGMAjYbY4qDrRSRRiLSpOo2vhOS64NtGy4BddHrqnm8kCabCXNco4GfAeONMaeq2aa+nq9Q/v8Z+I4d8B1L86r7cAoX6xzAv4FNxpj/q2abNlXnCkQkD997PKIfQCG+LjOAW63eNxcCR/3KFpFW7bdqO54vP/7HUHV5aA5wpfgmcWqO77mdU+NHivTZ5kj84EtQxUA5sB+Y47fuUXw9JrYAY/yWzwbaWrc74fsAKADeBlIiFOdLwF0By9oCs/3iWGP9bMBXwoj0c/cqsA5Yax1omYFxWX+PxderY3s9xVWArxa52vp5LjCu+ny+gv3/wOP4PogAGljHToF1LHWqh+foEnwlt7V+z9NY4K6q4wy4z3pu1uA7qX1xPcQV9HUJiEuAp63ncx1+veUiHFsjfIk7zW9ZvT9f+D5o9gKVVu76Ab5zOp8B24BPgRbWtrnAv/z2/b51nBUAd9Tm8XUIBKWUinPxVrpRSikVQBO9UkrFOU30SikV5zTRK6VUnNNEr5RScU4TvVJKxTlN9EopFef+Py+qnZH2cU5PAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from numpy.random import normal\n",
    "model=absmodel(nu,params[0],params[1],params[2])\n",
    "plt.plot(nu,model)\n",
    "plt.plot(nu,data,\".\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:11.182786Z",
     "iopub.status.busy": "2022-10-20T05:48:11.182469Z",
     "iopub.status.idle": "2022-10-20T05:48:11.184079Z",
     "shell.execute_reply": "2022-10-20T05:48:11.183790Z"
    }
   },
   "outputs": [],
   "source": [
    "#NCG"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:11.187312Z",
     "iopub.status.busy": "2022-10-20T05:48:11.186821Z",
     "iopub.status.idle": "2022-10-20T05:48:18.162965Z",
     "shell.execute_reply": "2022-10-20T05:48:18.162656Z"
    }
   },
   "outputs": [],
   "source": [
    "gd = jaxopt.NonlinearCG(fun=objective, maxiter=100)\n",
    "res = gd.run(init_params=(1.5,0.7,1.5))\n",
    "params, state = res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:18.165241Z",
     "iopub.status.busy": "2022-10-20T05:48:18.164923Z",
     "iopub.status.idle": "2022-10-20T05:48:18.167385Z",
     "shell.execute_reply": "2022-10-20T05:48:18.167620Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(DeviceArray(1.9526778, dtype=float32),\n",
       " DeviceArray(1.0492882, dtype=float32),\n",
       " DeviceArray(1.8708111, dtype=float32))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "params"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-10-20T05:48:18.172661Z",
     "iopub.status.busy": "2022-10-20T05:48:18.172366Z",
     "iopub.status.idle": "2022-10-20T05:48:18.469013Z",
     "shell.execute_reply": "2022-10-20T05:48:18.469353Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fb90c0d6eb0>]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAz5klEQVR4nO3deXxU1dnA8d8zM1lYw5IAIQkk7DsCMSCIUlQ2FdS6gLZibaVW7ap969JP7Wtf17d2sdWqrb6gtbgvqCiiiKIgYZFN1hACSdjCvmcyM+f9Y25wCBOYJDO5szzfzyefzNxl5pmZe585c86554gxBqWUUvHLYXcASimlIksTvVJKxTlN9EopFec00SulVJzTRK+UUnHOZXcANaWnp5vc3Fy7w1BKqZiybNmyPcaYjGDroi7R5+bmsnTpUrvDUEqpmCIiW2tbp1U3SikV5zTRK6VUnNNEr5RScU4TvVJKxTlN9EopFec00SulVJxL7ERfWggLHvf/V0qpOBV1/egbTWkhzJgIXjc4k2HqLMgpsDsqpZQKu8Qt0Zcs8Cd54/X/L1lgd0RKKRURiZvoc0f6S/Li9P/PHWl3REopFRGJW3WTU+CvrilZ4E/yWm2jlIpTZy3Ri8jzIrJbRNbUsl5E5AkRKRKRVSIyOGDdVBHZZP1NDWfgYZFTACPv1CSvlAq/wM4eNnf8CKVEPx34O/BCLevHA92tv6HAP4ChItIGuB/IBwywTERmGWP2NzRopZSKaoGdPRxOQMDnsa3jx1lL9MaYz4F9Z9hkEvCC8fsKaCUimcBYYK4xZp+V3OcC48IRtG20O6ZSKhSndPaosr3jRzjq6LOA0oD7Zday2pafRkSmAdMAOnXqFIaQIkC7YyqlQlXd2SNYid6Gjh9R0RhrjHkWeBYgPz/f2BxOcMG6Y2qiVyr2lRY2rFNGsP1rdvYAWzt+hCPRlwM5AfezrWXlwKgay+eH4fnsEfgNrd0xlYoPDf2lfqb9cwpOfSwbC4bh6Ec/C7jR6n0zDDhojNkBzAHGiEhrEWkNjLGWxabqb+jR92m1jVLxoqEXTjZ0/0Zq9ztriV5EZuIvmaeLSBn+njRJAMaYp4HZwASgCDgG/MBat09E/gAssR7qAWPMmRp17RXKz7ea39BKqehQ3+qXhv5Sb8j+jdjud9ZEb4yZcpb1Bri9lnXPA8/XL7RGpA2tSsWuhpy/Db1wsiH7N2K7X1Q0xkZEXb7htaFVqdjV0PO3ob/U67t/I7b7xWeir/kNP+4ROL639qSvDa1Kxa5YPX+D9cxZ8HhEeuaIv+YleuTn55ulS5fWeT+vz/DUp0V0bNWEc8v+j5wVf0aMF3CAwwHGnPlnXUO7WCml7BMt528d4zju9lJ+4BiHNi1kwCc34vC5cbhS6lV9LCLLjDH5wdbFTYl+9+ETPD53IwCDpTkvJTtJwmBEcPi8ODH4PG42fjUbj6MneenNaJYS8PLr+vMrWg4spVTt529jnqdB2gpM9rnsPepmy56jlOw5yrZ9x9i69xil+49Ruu84e45UAnCb8x0GuNw4xBeR6uO4SfSZaU1Y/4dxbD9wnPIDBSwo7kpq+ULKK5ty1a6/YfBQZZzcuzyN5cu+AKB9yxTy0pvRNaM5XTKa061dc7q3a05mWioi4n/gYAeKNt4qFf0a8Tyt8vo4tOYTWnsqceDD66nkpZn/5n+P7udwpefkdg6Bjq2a0KlNUy7q1Y6cNk3Ibt2UHlUpOD56B+OtQiJQ/RQ3iR4gNclJFytp030iMNG/onQclCzA3XE4f2jSh617j7Flz9GTf++t2sHB41UnH6dFiovu7ZtzSYtt3LLl5zh9VeBKRqa+6z9QGqvxVn81KFV/EThPjTHsOlTJ2h0HWbfjMBt2+v+K9xyhvy+Vl5JdJOHBIy62NB/ElT2zyEtvRm56M3LbNiOrVROSXcEuX8qCzHcjdr7HVaKvlfWzrinQF+jbMe2U1cYY9hxxs7niCJt2HWbjriNs3HUY9+bPwVQh4sNT5eb5F2awMtfFJS3yuNyRhMNHRL59Af3VoFRDNbCR1mxbzOEN81mbPIAFJ/JYVXaQtdsPsfeo++Q2Wa2a0KtDC0b3bkf3dgPYbgaSc2gZqd0u5P76dNXU7pWRIyJktEgho0UKw7q0PbncbHPAC2/i81ZhHC72txvKyrIDvL8/hRfkboY511HedAhNClMZtHMb+blt6JLe7Ntqn4aoWRpZ+R8t3StVF3Xs4+72+FhdfpBlW/exf8OX/Lz8TpoaDwNx8UfPfRxtN4TRvdrRt2NL+nRMo1dmC1qmJtV4lGxgdMReUn1poj8D6TQUpr6LlCzAkTuS3+QU8Btg31E3q8sLWFl6gEPb9rN/zWe0Wb6aV3y9adkkiUmtiknuegHdhlxEj/bN65f4a45+9/V/bB3PWqmYdIZG2qrNn7M+dSBzD3dmcfFeVpQeoNLjA+CeFotJwoNTfDjFy8xLPCSNipFum0HETfdK25QWYqwqFn/fHoMYL1XGxQ3ue9nWrB8juqUzsnsGF3RPp13L1JP7nbWkUb3NwTJYNsNfuhenf7ydkXc23mtUKg4YY9iw6zDrl3zM+OW34vRVUYWL71XdizvzXM7NbUNBXmuGdG5DxoGVMVd1mhDdK0MW7gbOkgWIVcXiwmctNDgdXn4/cD/Pk84XRXt4Z8V2APp2bMn3snZx3drbEZ/Vwl7bQVRdGikthBUzg9c1aoOtUrU6UeVl4eY9fLxuN/PW7WbnoRPc5vyQy5OqcIoPh3h56WIPqaPPP3XHFhGcU9qGczaxEn0kGjhrmWBAnMkMGHEZf8kZhM9nWLfzEJ9trGDeut2Uff0RPqcbl/i7YW3/+iM6Zp2L01FLFU9tdY3aYKvUaY5vXkTJ8jm8f7Arz2/L4JjbS7NkJyO7ZzC6VztGN2+J841Z4HXjcCaT2v3C4A8UicZRm87ZxEr0kegWGcIEAw6H0LdjGn07pnHbqG4c2mjg5Xfw+ty4jYufL2rG9jXzuHxgJpPOyaJvx5an1+sHO+h0jB4Vy8JYsnV7fHy2sYJViz7itm2/ojse7hAXzXo8QZ+CixnWpQ0pLqe1dU7tpfVIl7ZtOmcTK9GH2t2qrh92HScYaNljBPzA32fWdBzOTUc68+7K7UxfWMI/F2yhe7vmXJufw5WDs0hvntLw16NUtAlHyba0kN2rP+bdg134+6Y27D9WxZ1NFpIsHpz4G1F/0nk79Mg4fd9gBafGKG3bdM4mVqIPpbtVY/20CujbPxGYOLAjB465mb16J68tK+XB2et49MP1XNKnPTeel8uwLm2Cl/IjVY+oVCQ1oGR7osrLl/Nnc/6XN9PGeLgeF/s7/YkhI8YysklbnC++BV533a9xaYzStk3nbHwl+nBMHhLuD7sOvw5aNU3m+qGduH5oJzbtOswrS0p5fXkZH6zZSY/2zZk6PJfvDs4mNcn57U46GYqKRfUo2e44eJzpC0t4ZUkp11e+y6ikb7s/3tVjN/RqB7SrfyJtrNK2Deds/HSvDFdJPJwl+jA81okqL7NWbueFRSWsKT9E22bJ3DQ8l++f15lWTZPrF9eZ4tVfB6qxhHi8bVs5n9Vfvs/07Vks83ZnbN8O3NZtH/0+/r6/x1s4f3nH8DmQGN0rw1USD+dPq/rEVONAS01ycm1+DtcMyWbxln0889lmHp+7kac/28xNI3K5ZWSX8CR87cGjGttZSrbrdx5i1rtv89OyOxmLhzFJSeyb8jrt+w7xb5AVgSqQUEfBjLEvhPhJ9OH82RWun1Z1jekMyVZEGNalLcO6tGX9zkP8fV4RT83fzIyFW7l5RC7TLuxK85RaPs5QDkrtwaPsFHCMFqf24fGPNvL+6h38ImUhKeLBgQ/w0H7fEuAC/z6NVQUSbCKjD++OqUJR/CT6aGyYrGtMISbbXh1a8vfrB/PTnYd54pNNPDGviP8UbuMXF/dg8rk5uJwBo+OFWlLXHjzKLtYxarxuPLj4jftekpwOXuyxg8G98nHMe9ve47LmebnunZgrFMVPoofobJgMjOlsJes6JtueHVrw5A2DmVZ6gAdnr+O3b69h+sISHpjUl+HJxd8OnxDKQRmNX5QqIXiLP0escdzFGO7qsJyCg3OQ0irYnnz2qUBDVd/qlprnZe9JsHVRTBWK4ivRR7NQStb1TLYDc1rxyrRhfLR2F//z/lr++K8XeTn1YZLwIA4nOFzg4+wHZTR+Uaq4tmzrfv69tCUPGRfJ4sGRlMzQvLawrOrbwsnxvQ0f26khbVDBzsv2fWKqUKSJvrGEWgdez2QrIoxtuY3RBctZu+4bHDv94+j7vCBDbkRa5cTMQani3zG3h0c/WM+MRVvJTOvKytEzGCrrkDyrIFLb2E711dA2qGAXRcbQuaSJvrFEug7cKrEked0MdDjxuVx4vR7cxsmfywfxwwuvo331yJlKRVpt1SSlhZQun8NDa9P58FAnfjAil7vG9LTmb57w7XbhHqIgwdug4qcffSw4w8Hf4J+BCx6HeQ9+O5TxkBvxtczh/UNdueurFFKTnPzv1QMY07dDeF6LUrWppZrEs/UrzIyJiLcKj7jYPP4/9B16SYMft077x1B1S10lRj/6WBDJ8TVqllgGXo8jp4DLgT5Dj/CLl1cw7cVl/GBELveM713LvJVKhUGQapLyFv359LWZTPZW4bKuZu3rXgXUIdGHu/olgejZbrdgB299VDcYjb7vtC+LrhnNef0n53HT8Fz+78sSrn56IeUHjofpBShVQ3WhQ5zgTGa5ox+XPrGA2Ye6Ii7/8nrNtVzjcROt+qUhtOrGbo18ReqHa3by69dWkpLk4JnvD2FI5zYRey6VwEoLMVsWMOtgF15YVMJlLYsZe+nVdGzVpGHVJ3Fe/dIQZ6q60UQfDRr54C3afZgfzVjK9gMneOiq/lw9JDviz6kSS6XHy31vraF4+bxvu/rGyFWksepMiV6rbqJBToG/n3AjnQDd2rXg7dtHcG5ea+56bSV/mruRaPvCV7Hr0IkqbnyukNeXlfGrHrv9Sb6hVZOqQUJK9CIyTkQ2iEiRiNwdZH1nEflERFaJyHwRyQ5Y5xWRFdbfrHAGr+qvVdNkpv+ggGuGZPPEJ5v43Tvf4PVpsld1VFro7/FVWghAxeFKJj/zFcu27uevk8/h/Iuv9JfktV7dVmftdSMiTuBJ/M3jZcASEZlljFkbsNkfgReMMTNEZDTwMPB9a91xY8w54Q1bhUOS08FjVw+gTfNknvmsmH3H3Pz52nO0R446s+qqxiZtTxnca9eVr3LdbC+7DlXyr6n5jOrZDsjSoTWiQCjdKwuAImNMMYCIvAxMAgITfR/gV9btT4G3wxijiiAR4Z7xvWnbLJmHZq/H7fHx1A2DSXJqsldBBHYeEAHjA+PDeN28+ebL7Pddwb9/NJQhnVt/u08Cd2uMFqGczVlAacD9MmtZoJXAVdbtK4EWItLWup8qIktF5CsRuSLYE4jINGubpRUVFaFHr8Jm2gVdeWBSX+au3cXPZn6Nx+uzOyQVTaqraFbO/LY7sM8H4sCIk0rjZKG3N/+5pUaSV1EhXBdM3QX8XURuAj4HygGvta6zMaZcRLoA80RktTFmc+DOxphngWfB3+smTDGpOrrxvFzcHh//8/46fvnqSv587cBThzxWiSmwFF9jkLwDFz7AK5+vZEFVL37zoxvp2zHN7mhVEKEk+nIgJ+B+trXsJGPMdqwSvYg0B75rjDlgrSu3/heLyHxgEHBKold11JDumGeZKedHI7vg8Rke+WA9zVOcPHRl/9MnJVeJJfCiPh8w5EZIy+Fwh2Fc/a6Hne7OvPjDAvpna5KPVqEk+iVAdxHJw5/gJwPXB24gIunAPmOMD7gHeN5a3ho4ZoyptLYZATwWxvgTT0MusApxppxbL+zK4RNVPPnpZrJaNeGO0d0j+5pUdAsyvMaJDkP44XOFbNt7jBd+WMCgTlpdE83OmuiNMR4RuQOYAziB540x34jIA8BSY8wsYBTwsIgY/FU3t1u79waeEREf/vaAR2r01lF11ZDxPuowU85dY3qy/cAJ/vjRRjLTmvBdvagqcdUYj92XdS53vvw1hSX7eGLKIIZ1aXv2x1C2CqmO3hgzG5hdY9nvAm6/DrweZL+FQP8GxqgCNWS41TrMlCMiPPrdAew+fILfvLGKzLRUhndLj8ALUjEhoOfMo7PX8f6qHdwzvhcTB3a0OTAVCh0CIRZFsI6+pkMnqrj6HwupOFzJuz89n+zWTcP0IlQsemdFOT9/eQXfG9aJP0zqp+03UUTHulENsmXPUSb+7Qs6pzfljctcpJQt1ItfEtC6HYe48qkv6Z+Vxks/GqYX1kUZHetGNUheejP+fN05JG9firwwCTPvQX+jrnXZu4p/B49Vceu/l9EyNYknrx+sST7G6MQj6uxKC7l47wIyOq05ORdtvSZ+UNHpLNV3xhh+9eoKth84zsvThtFOp6SMOZro1ZkFdMkc4HDiESceAw5XEo7AhmAdJzw2Betye3zvKZ/jC4u28sn63dx/eZ8zz1+gx0DU0kSvziygS6b4oGrADTy32kNJk0E8nJlPMjT65CkqjAK73HoqYfadYMzJz3Fjcm8emr2O7/TM4KbhubU/jh4DUU0r2tSZ1Zi+rem536fn1ffzRkUWj8/d4N8mXNMhqsYX+Pk6HNYgZf7P0VP8OT+b+TXNU1w8dvXAM/ew0WMgqmmJXp1ZjYtlyClgDDCloBPPfl7Mhd0zGN6Qvv3KXoGfb41hh/+9M4f1Ow/z3NR8MlqknPlx9BiIatq9UgV3lvrWY24Plz7xBVVeHx/98gKa7lqu9bPxwPrc16cO5N63VvPjTjsYe+nVoX2mWkdvK+1Hr+omxPrWxcV7ue7Zr7hlZB73XdrHhkBVJLg9Pn7953/y6NHfkiJenes1Rmg/elU3Ida3Du3SlikFOTz3xRZWlx1s5CBVpDzz2WayDiwjWed6jRua6NXpajTAnqm+9e7xvWnbPIW731ylk5XEkhpzvVbbXHGEv80rwtllJA5Xis71Gie0MVadLkgDbFClhaSVLOCJEb2Z8mElz32xhR9f2LVxY1V1V0vVnDGGe95cTZNkJzdeey0c6Kl17nFCE70K7mzzfAYki2HOZKblPcpfP9nElYOy9MrJaFfLUNezVm6ncMs+Hrmqv7+XTQud6zVeaNWNqp/AC6m8bn7SeQdVXh+Pzdlgd2TqbIJUzR13e3nkg/X0y2rJtfk5Z38MFVO0RK/qp0a/6dZ9R3OzrwXPfF7Mjed1ZkB2K7sjVLUJUjX3zMcb2XHwBH+dPAiHQ4cejjea6FX9BEkWd2RU8cbyMh54dy2v3XqejlUezQKq5rYfOM7Tn23m0gGZFOSdYSwbFbO06kbVX04BjLzzZMJokZrEXWN6snTrft5btcPm4FSoHv1wPcbAPeN72R2KihBN9CqsrsnPoU9mSx75YD1uj3a3jHaryg7wzort3DKyi84eFsc00auwcjqE34zvRfmB47yyZFut/bVVdPjjRxtp3TSJH1/Yxe5QVARpHb0Kuwu6p3Nubmvmf/w+35M/IDp0bVQq3LKPzzdWcM/4XrRITbI7HBVBWqJXYSci3DWmJz1PrMR4KnXo2ihkjOGPH20go0UKN56Xa3c4KsI00auIGNqlLcezzsONC6OX0UedBZv2ULhlHz8d3Y0myU67w1ERpoleRcyky6/k+sp7WdjpVq22iSLVpfmsVk247ly9OCoRaKJXEXNOTiva9h7JrVsv5HDGILvDUZb5GypYVXaQn13UjRSXluYTgSZ6FX4BPW3u+E43Dp/w8NLibXZHpSxPzS+iY1oqVw3OtjsU1Ug00avwqh7sbN6DMGMiA9nIiG5tee6LLZyo8todXcJbUrKPJSX7ueWCLiQ59fRPFPpJq/AKMjLibaO6UXG4kjeWl9kdXcL7x/zNtGmWzORzO9kdimpEmuhVeAUZGXF417YMzE7jmc+KdXISG63bcYh563fzg+G52tMmwWiiV+FVPdjZ6PtO9rQREX4yqivb9h1j9pqddkeYsJ7+bDPNkp3abz4BhZToRWSciGwQkSIRuTvI+s4i8omIrBKR+SKSHbBuqohssv6mhjN4FaVqDHYGMKZPB7pmNOMf8zcTbRPSJ4LSfcd4d+V2bhjWmbSmehVsojlrohcRJ/AkMB7oA0wRkT41Nvsj8IIxZgDwAPCwtW8b4H5gKFAA3C8ircMXvooVDodwT/9DfGf3i6xZ/LHd4SSc6QtLcIhw84g8u0NRNgilRF8AFBljio0xbuBlYFKNbfoA86zbnwasHwvMNcbsM8bsB+YC4xoetoo5pYVcVDiNXyW9Rs85N+ggZ43oaKWHV5eUMqF/Jh3SdJrHRBRKos8CSgPul1nLAq0ErrJuXwm0EJG2Ie6rEkHJAsTrxoUPh6+KA2vnnX0fFRZvLC/jcKWHm0bk2h2Kskm4GmPvAi4Uka+BC4FyIORO0yIyTUSWisjSioqKMIWkoorVG8eIkypcvHNAh8VtDD6fYfrCEgbmtGJwJ601TVShJPpyIHBAjGxr2UnGmO3GmKuMMYOA+6xlB0LZ19r2WWNMvjEmPyMjo26vQMUGqzeOjL6PJzv9icfXpXHM7bE7qri3oGgPxRVH+cHwXLtDUTYKJdEvAbqLSJ6IJAOTgVmBG4hIuohUP9Y9wPPW7TnAGBFpbTXCjrGWqURk9cYZddGlHDrh4c3lp33nqzD7vy+3kNEihQn9M+0ORdnorIneGOMB7sCfoNcBrxpjvhGRB0RkorXZKGCDiGwE2gMPWvvuA/6A/8tiCfCAtUwlsCGdW9M/K43pC0u0q2UEla+eT5+if3FX74Mku/SSmUQm0Xai5efnm6VLl9odhoqw15aW8uvXV/HytGEM69LW7nDiT2khVc9fhviqcCYlI1Pf1WGi45yILDPG5Adbp1/zyhaXDehIi1QXLxfqqJaR4Cn+HDFVuMSHeKt0dq8Ep4le2aJJspMrB2Uxe81O9h912x1O3Fnk7U2VceHT2b0UOjm4sktpIbe7PmaNtwlvft2LH56vV2yG01NFbWmZ8t/84/zjkDdSq20SnCZ61fisMevbe93MTHFy18KmmBHTEBG7I4sLxRVHWFS8l1+PHYPjgm7frigt9Ffh5GriTzSa6FXjCxizPgnIObiMpVv3c25uG7sjiwuvLCnF6RCuGRIwg1T1hDBet78qR+fwTShaR68aX8CY9eJKZqWrPzO1UTYs3B4fry8r46Je7WjXMmBcmyATwqjEoSV61fiqx6wvWYDkjiR3aRNeX1bG/Zf11SF06yOgSmbu/hz2HnUzZWiNGaSqv1yrS/TaOJtQNNEre+QUnKw6uE4O8NLibby3ejs3DO1sc2AxpkaVzO6WP+a/mlVwQWo60O7b7QK+XLWOPvFoole265+VRvd2zXlzebkm+roKqJIxnkq+t/dvOMXgePHN0+vhA75cVWLROnplOxHhu0OyWbZ1P1v2HLU7nNgS0N7hE8GBDwc+rYdXp9BEr6LCFedk4RBY9NlsWPC4TkwSqoA5ev+WeiseSTplYnalQKtuVJTokJbK1JzdXLX6vzDiRbQLYOhyCvjG2ZO/vP8FPS4cyoTmRVoPr06hiV5FjWvSS3Dt8iCBVQ+arELyxrJykpzCeReMh2bJdoejooxW3aio0TV/HFW48OLQqoc6qPL6mLWynIt6tae1JnkVhCZ6Zb/SQljwOClJTp7v+lee8F3Lievf0tJ8iD7fWMGeI26+G3glrFIBtOpG2atGP/BRl7zIpWsn0uVQJybZHVuMeHvFdlo3TWJUT52GUwWnJXplrxqX5veuXEVmWiqzVmy3O7LoZv0KOr55EXPX7uTSAZkkOfV0VsFpiV7Zq8al+Y68kUw80oLnvtjC/qNurXMOJuBXUJIk0cdzN5POGWp3VCqKaRFA2SugH3h1d8rLB3bE4zN8sGan3dFFp4BfQeJzc0nTTQzp1NruqFQU00Sv7JdTACPvPNn42rdjS7pmNGPWynKbA4tS1q8gI07cxkWznqNwOHQsf1U7TfQq6ogIEwdmsXjLPnYePGF3ONHH+hW0svvt3OC+lyHnj7U7IhXlNNGrqDTxnI4YA++t0kbZoHIKePjwBA5lDKZPZku7o1FRThO9ikp56c24pv12khb9Rce9CWLHweMUluxj4sCOOgWjOivtdaOiU2khDx++D7xV+Ga8gmPqu3oBVYD3Vu7AGJg4sKPdoagYoCV6FZ1KFuA0HlziA0+VDrlbw3urttM/K43c9GZ2h6JigCZ6FZ1yRyLOZDw4cOPUcW8ClO47xsqyg1w6INPuUFSM0ESvopPVs2RFt9u5vvJeilL62B1R1Ji9egcAl/bXRK9Co4leRa+cArIv/y1f0+NkclPw/uodDMxOI6dNU7tDUTFCE72Kah3SUsnv3Jr3V2miB9i29xiryg4yQUvzqg400auod2n/TDbsOkzR7sN2h2K72Wv8X3ia6FVdaKJXUW98/0xE4P1VOvbN+6t2MDCnlVbbqDoJKdGLyDgR2SAiRSJyd5D1nUTkUxH5WkRWicgEa3muiBwXkRXW39PhfgEq/rVvmcq5ndvw/urEvkp2295jrC4/yKX9O9gdiooxZ030IuIEngTGA32AKSJSswvEb4FXjTGDgMnAUwHrNhtjzrH+bg1T3CrBXDogk427jrBpV+JW37y/WqttVP2EUqIvAIqMMcXGGDfwMpw2+Y8BqgfcSAMSu+ilwm58vw6IwOzViVt9M3u1v9omu7VW26i6CSXRZwGlAffLrGWBfg98T0TKgNnATwPW5VlVOp+JSNCrXkRkmogsFZGlFRUVoUevEka7lv7eNx+sSczeN6X7/NU2E/pptY2qu3A1xk4BphtjsoEJwIsi4gB2AJ2sKp1fAf8RkdOG2jPGPGuMyTfG5Gdk6LyXKrhx/TJZv/MwW/YcPTmVXqIMePahNQnL+H5abaPqLpREXw7kBNzPtpYF+iHwKoAxZhGQCqQbYyqNMXut5cuAzUCPhgatEtM4qzS7/Ms5/qn05j3o/58Ayf6DNTvo27ElndpqtY2qu1AS/RKgu4jkiUgy/sbWWTW22QZcBCAivfEn+goRybAacxGRLkB3oDhcwavEktWqCQNzWnFkw/xTJhSP9wHPdhw8zvJtBxiv1Taqns6a6I0xHuAOYA6wDn/vmm9E5AERmWhtdidwi4isBGYCNxljDHABsEpEVgCvA7caY/ZF4HWoRFBayL0tZrPugAufMwnE6Z9YPM4HPJtTXW2jvW1UPYU0Hr0xZjb+RtbAZb8LuL0WGBFkvzeANxoYo1L+6pkZEynwuhmQ5OTLbv/FyCyHP8nH+Tj1H6zZSY/2zema0dzuUFSM0itjVWwoWQBeN2K8JImHraWlp0woHq8qDldSWLKPcf0yE64BWoWPzjClYkPuSH81jdeNERdv7svj4oMn6JCWandkEbX8yzn8xPEeU5yDYMZ/+9sknMkwdVbcf8mp8NFEr2KDNT49JQvYlTaE5f85yodrdnDTiDy7I4uc0kJGLf4RFydV4fjiDTA+/191A7QmehUirbpRsSOnAEbeSfaAUXRv15wP1sT3VbLHN87H6avCiQ/x+UAcCdMArcJLS/QqJo3v14G/f1rEniOVpDdPsTuciFjo681wXDjFiziTYdwjcHxvQjRAq/DSRK9i0rh+mTwxr4i5a3cxJXOnvyojXhJgaSGULGDx5nReSflvnjn/OOTFyWtTttBEr2JS78wWdG7blE3L5sHeu+OnkdLqRmq8bn7pczKz99+RC66zOyoV47SOXsUkEWFcvw40LV+IiaerZAO7keJhTNNNdkek4oAmehWzxvfLZKG3N16Jo6tkrW6kXhx4xEXHgZfYHZGKA1p1o2LWgKw0drQcwGNtHuPePnvio44+p4DK69/iyenTadZzFD/uPNTuiFQc0BK9ilkOhzC2bweml7bjSMHPYz/JWz49lssT7on0G6qleRUemuhVTBvfrwNuj49P1++2O5Sw+WDNTlo3TWJoXhu7Q1FxQhO9imn5uW1Ib55ycmKOWFfp8fLJut1c0qc9Lqeenio89EhSMc3pEMb1a8+89bs57vbaHU6DLdi4hyOVHp0AXIWVJnoV8yb0y+R4lZf5G2K/+mb2mh20THUxvGu63aGoOKKJXsW8grw2tG2WzOwYr75xe3zMXbuLMX07kOzSU1OFjx5NKua5nA7G9O3AvHW7OFEVu9U3X27ew+ETHib01ykDVXhpoldxYUL/Dhx1e/lsY4XdodTb7FU7aJHiYkQ3rbZR4aWJXsWFYV3a0qppEh+s3mF3KPVS5fXx0dpdXNynPSkup93hqDijiV7FhSSng7F9OvDxut1UemKv+mbR5r0cPF6lvW1URGiiV3FjfP8OHKn0sGDjHrtDqbM1X83l58mzuKBJsd2hqDikiV7FjeFd00lrksT7MVZ949n6FTcX/5yfOV4l5aUrdfJvFXaa6FXcSHY5mJa3h85rn6ZyyyK7wwnZtmUf4TIenPjiY6hlFXU00av4UVrIrVt/yR28guvfV8RMyfi9Q12pEhcmXoZaVlFHhylW8aNkAQ5fFSI+vNUl42gb0dKaJrB6SOUTVV7+uSUdR5e/cEeXnfEx1LKKOproVfzIHYk4k/F6KnEbF6bjcJraHVMga5rAwGkPPz/UicOVHvqfNwZ6ZNgdoYpTWnWj4kdOAUydxfbBd3KD+17mHu5kd0SnsqYJDJz28N1VO2jTLJnhXdvaHZ2KY5roVXzJKSDrsvvY3mIA767cbnc0p7KmCaye9vBE1nA+XruL8f06kKRDEqsI0qobFXccDuGyAZnMWFTCwWNVpDVNsjskP+sXR3Ud/dx92Ryv+prLB3a0OzIV57QYoeLS5QM7UuU1zFkbZSNa5hTAyDshp4B3V26nfcsUzs3VmaRUZIWU6EVknIhsEJEiEbk7yPpOIvKpiHwtIqtEZELAunus/TaIyNhwBq9UbQZkp9G5bVPe/rrc7lCCOnDMzYGNX/JQ+lyc5UvsDkfFubMmehFxAk8C44E+wBQR6VNjs98CrxpjBgGTgaesfftY9/sC44CnrMdTKqJEhCsHZbGoeC/lB47bHc5pvvrsQ2Y4/4fRO/7p74kTI33+VWwKpURfABQZY4qNMW7gZWBSjW0M0NK6nQZUt4JNAl42xlQaY7YARdbjKRVxVw3KxhiislRfseZjksWDBPTAUSpSQkn0WUBpwP0ya1mg3wPfE5EyYDbw0zrsi4hME5GlIrK0oiJ2xxNX0aXTsTU8lP4R6wo/xhhjdzgnba44wlv78jCOpJM9cPRqWBVJ4ep1MwWYbox5XETOA14UkX6h7myMeRZ4FiA/Pz96zkgVu6yLkyZ7KrnSuNi4rBs98y+2OyoA3lxexgp6cOjaN2hTUahXw6qICyXRlwM5AfezrWWBfoi/Dh5jzCIRSQXSQ9xXqfCzLk5y4CMJD1uXfRQVid7nM7y1vJwLemTQplcB9NKSvIq8UKpulgDdRSRPRJLxN67OqrHNNuAiABHpDaQCFdZ2k0UkRUTygO6AtjqpyAu4OMnnSOLFHTlRMSHJV8V72X7wBN8dnG13KCqBnLVEb4zxiMgdwBzACTxvjPlGRB4AlhpjZgF3Av8UkV/ib5i9yfgrRb8RkVeBtYAHuN0YY//ZpuJfwMVJ3zj6seBdD5+s2237DE6vLy+jRaqLS/q0tzUOlVhCqqM3xszG38gauOx3AbfXAiNq2fdB4MEGxKhU/eQUQE4BA3yG9p99wmtLS21N9IdPVPHhmp1MOqcjqUnay1g1Hr0yVsU9p0O4Lj+H+RsrKN13zLY43lxezjG3lykFUTbYmop7muhVQpgytBMOEf69eKu/R86Cxxv1IiVjDC9+tZWB2WkMyG7VaM+rFOigZipBZKY14ZLe7Vlf+DFm2YNIwJjwjdG1cVHxXop2H+GP1wyM+HMpVZOW6FXC+P55nenrXo3xVJ4yJnxj+PdXW2nVNInLBtjbGKwSkyZ6lTCGd23LtrTBVOFq1CtSdx06wZxvdnFdfo42wipbaNWNShgiwpARY5ny3nGeGnGMDgMvCW+1TY35YKv9Z/E2fMZw/VBthFX20BK9SihXDc5mnas3jx+/LPxJfsZEmPfgKaNRuj0+ZhZu48IeGXRu2yx8z6dUHWiiV4nB6mmTtudrrh6SzdsrytkezuGLg8wHC/5xbbKPrOZ3aR/qUMTKNproVfyrUdr+WY/9ADz92ebwPUeN+WDJHUmV18f8T95jZspD5K3+i447r2yjiV7Fvxql7Yy9hVw9JJuXl5Sy69CJ8DxH9ZALo+872WXznRXb6XJkBUnouPPKXproVfwLUtr+yYXd8PoMz3xWHL7nCZgP1uP18eSnRexsk4+4UnTceWUr7XWj4l/AAGfVPWI6AVcOyuKlxVu5dVQX2rVIDetTvrdqB1v2HOU337scSTsnaG8cpRqLJnqVGKwBzgLd/p1uvLm8jH8t2MK9E3qH7am8PsPf5m2iZ/sWjOnTARyZmuCVrbTqRiWsvPRmXHFOFtMXllBccSRsj/vS4q1srjjKzy7qjsMhYXtcpepLE71KaHeP70WKy8Fv314Tlnlldx06wWMfbmBk93Qm9O8QhgiVajhN9CqhtTu4ihe6f8GJ4kW8ubzhs1z+ftY3VHl9/M8V/RDR0ryKDlpHrxKX1b/+HK+bmSlOpr3n4ju9fkybZsmh7x/QyPrx2l18sGYnvx7bU6+CVVFFE71KXFb/ejFekgXGVH3KoulbmXD5NUinoWfet/oiLGu446NT3uT+WSfo0b45t4zs0jjxKxUirbpRiSugf704nFyX9Dljdz+HZ/rlp17BGmyikoCLsIzXzay3X2HXoRM8fNUAkl16WqnooiV6lbgC+9cfLMO5bAYiPjzeKlYueI+B1xecVnJn3CNwfC80aQvOZIzXTRUuXtuTy+PXDWRI59Z2vyqlTqOJXiW26v71pYXIipkYrxuvuHhgTWtuWbOTcfsDhk/wVMLsO8EYf5If9zBzlq7l2a0duezSK5h0Tpbdr0apoDTRKwUnS/dSsgCTNRzvBz5ue2kZd/bO5DZnEuIFRMD4wPgwXjcvffo1v907ltu/05Wbz8+z+xUoVStN9EpVs0r3qcBL4xayeN7n/GNjJl847uNH2eXQtA0XFD+Ow1Th9jmZd6Inj313ANfkZ9sduVJnpIleqZpKC2k28ypGe92MSk3i8Q7/y482d8UYGCz3MDJpPZkDL+Fvl11BsxQ9hVT006NUqZoCetQ4vPDrnhX86ubvW1fOjschokMbqJiiiV6pmqq7XVb3tMkdidMhgCZ3FZs00StVU5BhjZWKZZrolQomyLDGdVJjeASl7KSJXqlwq3mRlTW1oFJ20Wu1lQq3GnPU6jyxym4hJXoRGSciG0SkSETuDrL+zyKywvrbKCIHAtZ5A9bNCmPsSkWnIHPUKmWns1bdiIgTeBK4BCgDlojILGPM2uptjDG/DNj+p8CggIc4bow5J2wRKxXttDFXRZlQ6ugLgCJjTDGAiLwMTALW1rL9FOD+8ISnVIxqaGOuUmEUStVNFlAacL/MWnYaEekM5AHzAhanishSEflKRK6oZb9p1jZLKyoqQotcKaVUSMLdGDsZeN0Y4w1Y1tkYkw9cD/xFRLrW3MkY86wxJt8Yk5+RkRHmkJRSKrGFkujLgZyA+9nWsmAmAzMDFxhjyq3/xcB8Tq2/V0opFWGhJPolQHcRyRORZPzJ/LTeMyLSC2gNLApY1lpEUqzb6cAIaq/bV0opFQFnbYw1xnhE5A5gDuAEnjfGfCMiDwBLjTHVSX8y8LLxj/xUrTfwjIj48H+pPBLYW0cppVTkyal52X75+flm6dKldoehlFIxRUSWWe2hp6+LtkQvIhXA1gY8RDqwJ0zhhJPGVTcaV91oXHUTj3F1NsYE7c0SdYm+oURkaW3fanbSuOpG46objatuEi0uHetGKaXinCZ6pZSKc/GY6J+1O4BaaFx1o3HVjcZVNwkVV9zV0SullDpVPJbolVJKBdBEr5RScS4mE72IXCMi34iIT0Tya6y7x5ogZYOIjK1l/zwRWWxt94o1tEO4Y3wlYMKVEhFZUct2JSKy2tou4leKicjvRaQ8ILYJtWx3xslmIhDX/4rIehFZJSJviUirWrZrlPcrhMl2UqzPuMg6lnIjFUvAc+aIyKcistY6/n8eZJtRInIw4PP9XaTjsp73jJ+L+D1hvV+rRGRwI8TUM+B9WCEih0TkFzW2aZT3S0SeF5HdIrImYFkbEZkrIpus/61r2Xeqtc0mEZlarwCMMTH3h39ohZ74B0nLD1jeB1gJpOAfLnkz4Ayy/6vAZOv208BPIhzv48DvallXAqQ34nv3e+Cus2zjtN67LkCy9Z72iXBcYwCXdftR4FG73q9QXj9wG/C0dXsy8EojfHaZwGDrdgtgY5C4RgHvNdbxFOrnAkwAPgAEGAYsbuT4nMBO/BcVNfr7BVwADAbWBCx7DLjbun13sGMeaAMUW/9bW7db1/X5Y7JEb4xZZ4zZEGTVJPzj7VQaY7YARfgnTjlJRAQYDbxuLZoBXBGpWK3nu5Yao3pGuZOTzRhj3ED1ZDMRY4z5yBjjse5+hX+UVLuE8von4T92wH8sXWR91hFjjNlhjFlu3T4MrKOWuSGi0CTgBeP3FdBKRDIb8fkvAjYbYxpy1X29GWM+B/bVWBx4DNWWh8YCc40x+4wx+4G5wLi6Pn9MJvozCGWSlLbAgYCkUutEKmEyEthljNlUy3oDfCQiy0RkWgTjCHSH9fP5+Vp+LoY82UyE3Iy/9BdMY7xfobz+k9tYx9JB/MdWo7CqigYBi4OsPk9EVorIByLSt5FCOtvnYvcxddoQ6gHseL8A2htjdli3dwLtg2wTlvctlKkEbSEiHwMdgqy6zxjzTmPHE0yIMU7hzKX5840x5SLSDpgrIuutb/+IxAX8A/gD/hPzD/irlW5uyPOFI67q90tE7gM8wEu1PEzY369YIyLNgTeAXxhjDtVYvRx/9cQRq/3lbaB7I4QVtZ+L1QY3EbgnyGq73q9TGGOMiESsr3vUJnpjzMX12C2USVL24v/Z6LJKYmeaSKVBMYqIC7gKGHKGx6iemGW3iLyFv9qgQSdIqO+diPwTeC/IqrpMNhO2uETkJuAy4CJjVVAGeYywv19BhPL6q7cpsz7nNPzHVkSJSBL+JP+SMebNmusDE78xZraIPCUi6caYiA7gFcLnEpFjKkTjgeXGmF01V9j1fll2iUimMWaHVY21O8g25fjbEapl42+brJN4q7qZBUy2ekTk4f9mLgzcwEognwJXW4umApH6hXAxsN4YUxZspYg0E5EW1bfxN0iuCbZtuNSoF72ylucLabKZMMc1DvgvYKIx5lgt2zTW+xXK65+F/9gB/7E0r7Yvp3Cx2gCeA9YZY/5UyzYdqtsKRKQA/zke0S+gED+XWcCNVu+bYcDBgGqLSKv1V7Ud71eAwGOotjw0Bxgj/kmcWuN/b+fU+Zki3dociT/8CaoMqAR2AXMC1t2Hv8fEBmB8wPLZQEfrdhf8XwBFwGtASoTinA7cWmNZR2B2QBwrrb9v8FdhRPq9exFYDayyDrTMmnFZ9yfg79WxuZHiKsJfF7nC+nu6ZlyN+X4Fe/3AA/i/iABSrWOnyDqWujTCe3Q+/iq3VQHv0wTg1urjDLjDem9W4m/UHt4IcQX9XGrEJcCT1vu5moDechGOrRn+xJ0WsKzR3y/8XzQ7gCord/0Qf5vOJ8Am4GOgjbVtPvCvgH1vto6zIuAH9Xl+HQJBKaXiXLxV3SillKpBE71SSsU5TfRKKRXnNNErpVSc00SvlFJxThO9UkrFOU30SikV5/4foF69h8jW3RUAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from numpy.random import normal\n",
    "model=absmodel(nu,params[0],params[1],params[2])\n",
    "plt.plot(nu,model)\n",
    "plt.plot(nu,data,\".\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.8.8 ('base')",
   "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.9"
  },
  "vscode": {
   "interpreter": {
    "hash": "72bc7f8b1808a6f5ada3c6a20601509b8b1843160436d276d47f2ba819b3753b"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}