{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib\n",
    "matplotlib.rcParams['figure.figsize'] = [10, 6]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<BarContainer object of 5 artists>"
      ]
     },
     "execution_count": 69,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmIAAAFlCAYAAABIu4TDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOzdd3gV1drG4d9KISEEAoTQQgkldJASukiVInjAigVFRVGxoGLDhkc/saEoFmyoKAoqglQFBBSQIqF3EnonEDokpMz3x2wPIOnJ3rOTPPd1cZHMrJU8HvTwZs2adxnLshARERERz/NxOoCIiIhIYaVCTERERMQhKsREREREHKJCTERERMQhKsREREREHKJCTERERMQhfk4HyKkyZcpYERERTscQERERydSKFSuOWJYV9u/r+bYQi4iIIDo62ukYIiIiIpkyxuxK67oeTYqIiIg4RIWYiIiIiENUiImIiIg4RIWYiIiIiENUiImIiIg4RIWYiIiIiENUiImIiIg4RIWYiIiIiENUiImIiIg4RIWYiIiIiENUiImIiIg4JNNCzBjzpTHmsDFm/UXXShtj5hhjYly/l3JdN8aYUcaYWGPMWmNM04vm9HeNjzHG9L/oejNjzDrXnFHGGJPX/5A5cngT7F3hdAoREREpwLKyIvY10P1f154F5lqWFQnMdX0O0AOIdP0aCIwGu3ADhgEtgRbAsH+KN9eY+y6a9+/v5XmWBRNug1nPOZ1ERERECrBMCzHLshYA8f+63BsY6/p4LNDnouvfWLalQEljTAWgGzDHsqx4y7KOAXOA7q57JSzLWmpZlgV8c9HXco4xEHUP7FkKB9c5nUZEREQKqJzuEStnWdYB18cHgXKuj8OBPReN2+u6ltH1vWlcd17j28EvEJaPcTqJiIiIFFC53qzvWsmy8iBLpowxA40x0caY6Li4OPd+s6DS0OBGWPsjJJxw7/cSERGRQimnhdgh12NFXL8fdl3fB1S+aFwl17WMrldK43qaLMv6zLKsKMuyosLCwnIYPRuaD4CkM7Bmgvu/l4iIiBQ6OS3EpgL/vPnYH5hy0fU7XW9PtgJOuB5hzgK6GmNKuTbpdwVmue6dNMa0cr0teedFX8t54U0hvBks/8LewC8iIiKSh7LSvmI8sASobYzZa4wZALwBXG2MiQG6uD4HmAlsB2KBz4FBAJZlxQOvAstdv15xXcM15gvXnG3Ar3nzj5ZHmt8LR7bCjgVOJxEREZECxlj5dKUnKirKio6Odv83SkqAd+tARDvo+637v5+IiIgUOMaYFZZlRf37ujrrZ8Y/EJrcAZtnwIl0t6+JiIiIZJsKsayIugesVFg5NvOxIiIiIlmkQiwrSleDyKthxdeQkuR0GhERESkgVIhlVfN74fQh2DTN6SQiIiJSQKgQy6qaXaBkVXXaFxERkTyjQiyrfHztvWK7FsHhTU6nERERkQJAhVh2NLkDfAPsBq8iIiIiuaRCLDuKhUKD6+0jjxJPOZ1GRERE8jkVYtnV/F44fxrW/uB0EhEREcnnVIhlV3gzqNAY/tb5kyIiIpI7KsSyyxh7VSxuE+xa7HQaERERycdUiOVEgxsgsCQs/9zpJCIiIpKPqRDLiSJB0KSf3dz11EGn04iIiEg+pUIsp6LugdRkWPmN00lEREQkn1IhllOhNaBGZ4j+ClKSnU4jIiIi+ZAKsdxofi+c2g9bZjqdRERERPIhFWK5UasbhFRWp30RERHJERViueHjC1F3w44/df6kiIiIZJsKsdxqehf4BcKyT5xOIiIiIvmMCrHcKhYKjfrCmh/gbLzTaURERCQfUSGWF1o+AMnnYMXXTicRERGRfESFWF4oVw+qd4C/P4eUJKfTiIiISD6hQiyvtBpkt7LYOMXpJCIiIpJP+DkdwBulplpMXLmXoCK+9GpUMWuTal4NpWvYm/Yb3ujegCLiNhHPznA6goh40PPX1OWO1lUJ9Pd15PtrRSwNxsD4v3fz+szNJKWkZm2Sj4+9V2zvctiz3L0BRUREJE98sWi7o99fhVgajDE82imSfcfPMXnlvqxPbHwbBITAstHuCyciIiJ55v6raji2GgZ6NJmuDrXDaBgewofzY7m+aTh+vlmoWQOCoekd9uPJE/sgJNz9QUUkT+18o6fTEUTEA/p9sYzNB09ya4sqjubQilg6jDE83Kkmu+PPMnXN/qxPbDEQrFQdeyQiIuKlVu4+xqLYI9zXrjpFizi3GgYqxDJ0dd1y1ClfnA/nx5KSamVtUqmqUKcnrPgKzp91b0ARERHJtg/mxlAqyJ9+rao6HUWFWEZ8fAyPdIpke9wZZq47kPWJrQbBuWOw9gf3hRMREZFsW7f3BPO3xHFvu+oUC3B+h5YKsUz0aFCemmWD+XBeLKlZXRWr0hrKN7L3illZnCMiIiJu98G8GEoE+nFna+dXw0CFWKZ8fAwPd6zJlkOnmL3xUNYmGWOvisVthu3z3RtQREREsmTTgZPM3niIu9tWo3igv9NxABViWdKrUQUiQoP4YF4MVlZXuBpcD8XKwlK1shAREfEGH86PJTjAj3vaVnM6yv+oEMsCP18fHupYkw37TzJv8+EsTgqA5vdCzGw4EuvegCIiIpKh2MOnmLnuAHe2rkpIkHeshoEKsSzr0yScSqWKMmpebNZXxaLuBt8i9l4xERERccyH82IJ9PNlwJXesxoGKsSyzN/Xh0EdarJmz3EWxhzJ2qTgstDwJlj9nf0WpYiIiHjcjiNnmLpmP/1aVSE0OMDpOJdQIZYNNzQLp0JIYPb2irV8AJLOwspv3RtORERE0vTx/Fj8fX2476rqTke5jAqxbAjw8+WB9jVYvvMYS7fHZ21ShUYQ0Q7+/hxSkt0bUERERC6xJ/4sk1ft49YWVShbPNDpOJdRIZZNfZtXpmzxAD6YF5P1SS0fgBO7YcsM9wUTERGRy4z+cxs+xvBA+xpOR0mTCrFsCvT3ZeBV1Vm87SjRO7O4Kla7B5SsqlYWIiIiHrT/+Dl+it7DTVGVKB/ifathoEIsR25vWZXQYkUYNS+LbSl8fO1Vsd1LYP8q94YTERERAD79cxuWBQ928M7VMFAhliNFi/hyb7vqLNgax+o9x7M2qUk/KFIcFn/o3nAiIiLC4ZMJjF++h+ubhlOpVJDTcdKlQiyH7mhdlZJB/nyY1b1igSUg6i7YMBmO7XJrNhERkcLuswXbSU5JZVCHmk5HyZAKsRz654iE3zcdZsP+E1mb1PJBMD6w9GP3hhMRESnEjp5O5Ltlu+nTOJyIMsWcjpMhFWK50L9NBMUD/fgwq3vFQsLtBq8rv4GzWdzoLyIiItnyxaIdJCSnMKijd6+GgQqxXAkp6s/dbSL4df1BNh04mbVJbR6xG7wuH+PecCIiIoVQ/JnzfLN4Jz0bVqBm2WCn42RKhVguDbiyOsUD/Xjv961Zm1CuHkR2tc+fTDrn3nAiIiKFzKcLtnE2KYXBnSOdjpIlfk4H8Govh2Q6JAS4N/k6Rm64ifUvXUEDn51Z//qvlc9xNBFxk5ezuOdTRLxO3KlEvlm8i95XVCSyXHGn42SJVsTywN2+vxHCaUYm3+h0FBERkULrkz+3kZicwqP5ZDUMVIjliRLmHAP9pjM3tSmrUr23aZyIiEhBdehkAuOW7uL6ppWoHub9e8P+katHk8aYx4F7AQtYB9wNVAAmAKHACuAOy7LOG2MCgG+AZsBRoK9lWTtdX2coMABIAR61LGtWbnLlmWw8orgrMZkxb81nZPgnfHNPi8wnpKbAB80gqDTcOxeMyUVQERGRwu3j+bGkpFo82in/rIZBLlbEjDHhwKNAlGVZDQBf4BbgTWCkZVk1gWPYBRau34+5ro90jcMYU881rz7QHfjYGOOb01xOKRbgxwPt7W77WTqD0scX2jwM+1bArsXuDygiIlJA7Tt+jvF/22dKVgn13i76acnto0k/oKgxxg8IAg4AnYCJrvtjgT6uj3u7Psd1v7MxxriuT7AsK9GyrB1ALJCFJSXvc0erCMoEB/DunCy+Qdn4dggKhcWj3BtMRESkAPtofiwWFg/ns9UwyEUhZlnWPmAEsBu7ADuB/SjyuGVZya5he4Fw18fhwB7X3GTX+NCLr6cxJ18pWsSXQR1qsHjbUZZsO5r5BP+i0GIgbP0NDm92f0AREZECZk/8WX5cvodbmlchvGRRp+NkW24eTZbCXs2qBlQEimE/WnQbY8xAY0y0MSY6Li7Ond8qx25rWYVyJQJ4d84WLMvKfELz+8CvKCz+wP3hRERECphRc2Pw8TE8lA+66KclN48muwA7LMuKsywrCZgEtAVKuh5VAlQC9rk+3gdUBnDdD8HetP+/62nMuYRlWZ9ZlhVlWVZUWFhYLqK7T6C/Lw93rMnyncdYFHsk8wnFQqFJP1j7A5zc7/6AIiIiBcSOI2eYtGoft7esQvmQwOx/gXUT4Yd+kOBc/8DcFGK7gVbGmCDXXq/OwEZgPvBPQ63+wBTXx1Ndn+O6P8+yl4ymArcYYwKMMdWASODvXORy3M3NKxNesijvzN6atVWx1g+BlWJ32xcREZEsGTU3Bn9fw4MdctA6KvEUzH4BTuyFIs61u8jNHrFl2JvuV2K3rvABPgOeAZ4wxsRi7wH751DFMUCo6/oTwLOur7MB+BG7iPsNeMiyrJSc5vIGAX6+PNKpJqv3HOePLVl4hFq6GtTrDdFfQUIWz6wUEREpxGIPn2LK6n30bx1B2eI5WA1b8DacOgDXjLA7GTgkV29NWpY1zLKsOpZlNbAs6w7Xm4/bLctqYVlWTcuybrIsK9E1NsH1eU3X/e0XfZ3XLMuqYVlWbcuyfs3tP5Q3uKFZJaqUDuLdOVlcFWvzKCSehBVfuz2biIhIfvfe7zEE+vsy8Krq2Z8ctxWWfGxvDaoUlffhskGd9d3E39eHRztHsm7fCWZvPJT5hPCmENEOlo6G5PPuDygiIpJPbT54kulrD3B32whCgwOyN9my4NenwT8IOr/slnzZoULMjfo0rkj1MsUYOWcrqalZWBVr+xic2g/rJ2Y+VkREpJB6b04MxQP8uK9dDlbDNk+H7fOh0/MQ7PyLfyrE3MjP14fBXSLZfPAUv64/mPmEmp2hbH34a5RdsYuIiMgl1u87wW8bDnLPldUoGVQke5PPn4XfnrP/ro0akPl4D1Ah5ma9GlUksmww7/2+lZTMVsWMgbaPQtwmiJnjmYAiIiL5yHu/b6VEoB/3XFkt+5P/eg9O7IZr3gbfXB23nWdUiLmZr4/hsS61iDl8mulrs9AnrMENUCIcFo10fzgREZF8ZPWe4/y+6TADr6pOSFH/7E2O3wGL3oOGN0FEW/cEzAEVYh7Qo0F56pQvznu/x5CckprxYF9/aPMI7F4MO//yTEAREZF84N05WykV5M9dbXOwGjbrOfDxg6tfyftguaBCzAN8fAyPX13rfx2AM9W0PwSVgYUj3B9OREQkH1i+M54FW+O4v30NggOy+Vhx62zYMhPaPw0lKronYA6pEPOQrvXKcUWlEN7/PYaEpEz61RYJgjYPw7Z5sG+FZwKKiIh4KcuyePPXzYQVD+DO1lWzNzk5EX57BkJrQqtB7gmYCyrEPMQYwzPd67Dv+DnGLd2V+YSoARBYEha84/5wIiIiXmze5sNE7zrG4M6RBBXJ5mrYkg8hfjv0eAv8svmWpQeoEPOgNjXL0C6yDB/Nj+VkQlLGgwNLQKsHYcsMOLjeMwFFRES8TEqqxVu/bSEiNIi+zStnb/KJvbBgBNTpZbeI8kIqxDzs6W51OHY2iS8WbM98cIuBUKQ4LNSqmIiIFE5TVu9jy6FTDOlaG3/fbJYts18AKxW6DXdPuDygQszDGlYKoWejCnyxaAdxpxIzHhxUGlrcCxsmw5EYzwQUERHxEonJKbw7ZysNwkvQs2GF7E3e/qf99+eVT0CpbO4r8yAVYg54smttEpNT+XBeFoqrVg+BXyAsfNf9wURERLzI98t2s/fYOZ7uVgcfH5P1iSlJ9nmSJavajdK9mAoxB1QrU4y+zSvz/d+72X30bMaDg8Og2V2w9gc4ttMT8URERBx3OjGZD+fF0qZGKO0iy2Rv8t+fQdxm6P4G+Bd1T8A8okLMIYM7R+LrY3h3zpbMB7d5BHx87Y7AIiIihcAXC7dz9Mx5nu5eB2OysRp28gDMfx1qdoHaPdwXMI+oEHNIuRKB3N22GlPW7Gfj/pMZDw4Jh8a3w+rv4GQWjkkSERHJx46eTuTzBdvp0aA8jSuXzN7kWUMh5bzdriI7BZxDVIg56IH2NSgR6M/bszZnPvjKxyA1BRZ/4P5gIiIiDvpwfiznklIY0rV29ibGzrU36LcbAqE13BMuj6kQc1BIUX8e7FCD+VviWLb9aMaDS0VAo74Q/RWcjvNIPhEREU/bE3+W75bu5uaoytQsG5z1iUkJMPNJKF3DXrzIJ1SIOeyuNhGUKxHAG79txrKsjAe3ewKSE2DpR54JJyIi4mEjf9+KMTC4S2T2Ji4aaXfQ7/kO+AW4J5wbqBBzWKC/L491qcWq3ceZs/FQxoPLREL96+Dvz+FsvGcCioiIeMiWg6eYvGofd7WJoEJINt52PLoNFr0LDW6AGh3dF9ANVIh5gZuaVaJ6mWK8PWsLKamZrYoNgfOn7VdzRURECpC3Z20mOMCPBztkY3+XZcGMIXbPTS/uoJ8eFWJewM/Xhye71Sbm8Gkmrdyb8eDyDaB2T1g6GhIyedtSREQkn4jeGc/vmw7zQPsalAzKxuHcGybB9vnQ6UUoXt59Ad1EhZiX6NGgPI0qhTByzlYSklIyHnzVEEg4DtFjPBNORETEjSzL4s3fNhNWPIC720ZkfWLCCfjtOajQGJoPcFs+dzKZbhD3UlFRUVZ0dLTTMfLUX7FHuP2LZbzQsy73tque8eBvr4MDa+GxdVAkyDMBRQqBiGdnOB1BRLJomN9Y+vvOpvf5V1lnZfL3Zjp2vtEzj1OlzRizwrKsqH9f14qYF2lbswxX1izDR/NjOZmQlPHgq56Cs0dg5VjPhBMREfEi9c0O7vSdzbiULjkuwryBCjEv83T32hw7m8TnC7ZnPLBqG6jaFv4aBcmJngknIiLiBXxI5TX/McRTghHJNzsdJ1f8nA4gl2pUqSQ9G1Xg84Xbub1lVcqHBKY/+Kqn4Ns+sGIstBzouZAiBZinHlOICCQkpdBpxB+UDi7C1IeuxMcni0cSLf8CZmyH6z9nbaP8XYhpRcwLPdOtDqmpMGJ2JgeCV+8AVdrAwncg6ZwnoomIiOSZMYt2sP9EAs9fUy/rRdjpw/D7K1DtKmh4k3sDeoAKMS9UJTSI/m2q8vPKvWzYfyL9gcZAp+fh9EFYrjcoRUQk/zhyOpHRf2yjS91ytK4RmvWJs1+A5HPQ8918cah3ZlSIeamHO0YSUtSf12Zsyvjoo4groVp7+2iHxNOeCygiIpIL/7RrGnpNnaxP2rEA1v4AbQfbp80UACrEvFRIkD+DO0eyeNtR5m0+nPHgTi/Yb1Cq276IiOQDMYdOMf7v3dzesgo1wrJ4sHfyebuDfqkI+5SZAkKFmBe7vWVVqpUpxvCZm0hKSU1/YOUWUPNqWDxK3fZFRMTrDZ+5iWIBfgzuUivrkxaNhCNb4ZoR4J+Ncyi9nAoxL1bEz4dne9RhW9wZJvy9O+PBHZ+Dc8fso49ERES81KKYI8zfEsfDHWtSulgWjzKK2woLR9iHekde7d6AHqZCzMt1rVeOFtVKM/L3mIybvIY3tc+gXPKRXZCJiIh4mZRUi/+bsZFKpYrSv01E1ialpsK0weAfBN3fcGs+J6gQ83LGGF7oWZf4M+f5eP62jAd3fA4ST9jFmIiIiJf5ecVeNh88xTPd6xDo75u1SSvHwu7F0PX/ILisewM6QIVYPtCoUkmubxLOl3/tYE/82fQHlm8A9frYjyfPHPVcQBERkUycSUxmxOwtNKlSkl6NKmRt0qmDMGcYRLSDJv3cG9AhKsTyiSe71cYAb8/KpMlrh6Fw/gwsft8juURERLLiswXbOXwqkRd61sVktf/XzKcgOQGufb9A9AxLi444ysjLIU4n+J+KwH3WTXy45jru2TSAxj6ZPKb86337l4hkz8sZNFEWkRw5dDKBzxZsp2fDCjSrWjprkzbPgE1TofNLEFrDvQEdpBWxfOQBv2mU4Tj/l9SPjHq8ioiIeJMRs7aQkmrxTPcsNm9NOAkznoRyDaDNo+4N5zAVYvlIsEngCb+JRFu1+S21udNxREREMrVh/wkmrtxL/zZVqRIalLVJc1+BUwfg2lHg6+/egA7To8mMeOEjiptTUvl61ELeSH6ezo+3p4hfOrX0sZ3wQTNodhf0fMeTEUVERACwLIvhMzcRUtSfhztm8UiiPX/D8i+g5QNQqZl7A3oBrYjlM36+Pjx3TV12HT3LN0t2pj+wVAQ0uQNWjIXjmTSDFRERcYP5Ww7zV+xRBneOJCQoCytbyedh6qNQIhw6Pe/+gF5AhVg+1KF2WdpFluGDebEcP3s+/YFXPWm/ZbJghOfCiYiIAMkpqQyfuZlqZYpxe8uqWZv013sQtwl6vQsBxd0b0EuoEMunnu9Zl1MJSYyaG5v+oJBK0OxuWP0dxO/wXDgRESn0xi/fQ+zh0zzbo07622guFrcVFrwN9a+HWt3cH9BLqBDLp+qUL8HNUZX5dulOtsedTn9guyfAxw/+fMtz4UREpFA7cS6J9+ZspUW10nStVy7zCRcfY9TjTfcH9CIqxPKxIV1rE+jny6vTN6Y/qHh5aH4vrJ0AR2I8F05ERAqt93+PIf7seV7qVS9rzVsL+DFGGVEhlo+FFQ9gcJdI5m+JY97mQ+kPvPJx+6eMef/nuXAiIlIoxRw6xdglO7m1RRUahGehMXohOMYoIyrE8rk7W0dQI6wYr0zbSGJyStqDipWBNo/Axl9g7wrPBhQRkULDsixenraBYkV8ebJr7axMgJlPFvhjjDKiQiyfK+Lnw7Br67Pz6Fm+XLQz/YGtH4ZiYTDnJdSWX0RE3GHWhoP8FXuUIV1rU7pYkcwnbJgMm6ZBx6EF+hijjKgQKwCuqhXG1fXK8cG8GA6eSEh7UEAwtH8Gdi2CmDmeDSgiIgVeQlIKr07fRJ3yxbm9ZZXMJ5w5Yq+GVWwKrR9xf0AvlatCzBhT0hgz0Riz2RizyRjT2hhT2hgzxxgT4/q9lGusMcaMMsbEGmPWGmOaXvR1+rvGxxhj+uf2H6owerFnPZJTLd74dVP6g5rdBaVrwO/DIDWdx5giIiI58Omf29l3/BzDrq2Pn28WyouZT9lnSvb+CHwL70E/uV0Rex/4zbKsOsAVwCbgWWCuZVmRwFzX5wA9gEjXr4HAaABjTGlgGNASaAEM+6d4k6yrEhrE/VdV55fV+4neGZ/2IF9/+xT7wxthzQTPBhQRkQJr77GzfPxHLD0bVaB1jdDMJ2ycChsmQYdnoFw99wf0YjkuxIwxIcBVwBgAy7LOW5Z1HOgNjHUNGwv0cX3cG/jGsi0FShpjKgDdgDmWZcVblnUMmAN0z2muwuzBDjWoEBLIsKkbSElNZx9Yvd4Q3gzmvwZJ5zwbUERECqThMzdhDDx3Td3MB5+NhxlPQPlG0PYx94fzcrlZEasGxAFfGWNWGWO+MMYUA8pZlnXANeYg8E8nt3Bgz0Xz97qupXddsimoiB/PXVOXDftP8sPyPWkPMgaufgVO7oO/P/NsQBERKXAWxx5h5rqDDOpQk/CSRTOf8OszcO4Y9PnYflJTyOWmEPMDmgKjLctqApzhwmNIACzLsoA8e0XPGDPQGBNtjImOi4vLqy9boPRqVIEW1Urz9qzNnDiblPagiCshshssfMf+yURERCQHklNSeXnaBiqXLsrAq6pnPmHLr7DuR2j3JJRv6P6A+UBuCrG9wF7Lspa5Pp+IXZgdcj1yxPX7Ydf9fUDli+ZXcl1L7/plLMv6zLKsKMuyosLCwnIRveAyxvDytfU5cS6Jd+dsSX9gl2H2JslF73ounIiIFCjjlu5i66HTvNCzHoH+vhkPPncMpj0G5RpAuyGeCZgP5LgQsyzrILDHGPNPx7bOwEZgKvDPm4/9gSmuj6cCd7renmwFnHA9wpwFdDXGlHJt0u/quiY5VK9iCW5vWZVvl+5i88GTaQ8qVx8a3wbLPoPj6TzGFBERScfR04m8O2cr7SLLZO08yVnPw5k4+y1Jvyz0GCskcvvW5CPAd8aYtUBjYDjwBnC1MSYG6OL6HGAmsB2IBT4HBgFYlhUPvAosd/16xXVNcmFI11qUKOrPy1M3YKXXwLXjc/bv84d7LpiIiBQII2Zv4ez5FIZdm4XzJGPmwOrv4MrHoGJjzwTMJ3LVuMOyrNVAVBq3Oqcx1gIeSufrfAl8mZsscqmSQUV4smttXvhlPTPXHaRnowqXDwqpBC3vh8UfQOuHoHwDzwcVEZF8Z93eE0xYvocBbatRs2zxjAcnnIBpgyGsjt1YXC6hzvoF2K0tqlC3Qglem7GRc+fTaeDa7gkIDIHfX/ZoNhERyZ8sy2LY1PWEFivCo10iM58w+0U4dQB6fwx+Ae4PmM+oECvAfH0M//1PffafSGD0n9vSHlS0lL1pMnYO7Fjg2YAiIpLv/LJ6Hyt3H+fp7nUoEZhJ+4lt82HlWPu840rNPBMwn1EhVsC1qFaa/1xRkU/+3Mae+LPpDBoIJSrBnGE6EFxERNJ1OjGZ12du5opKIdzYtFLGgxNPwdRHITTywp5kuYwKsUJg6DV18PMx6W/c9w+ETs/D/pWwYbLnA4qISL7w3pytxJ1O5OX/1MfHJ5MN+nOGwYk99luS/llo9FpIqRArBCqEFOXxLrWYu/kwszceSntQo75Qtj7MfQVS0mkEKyIihdaG/Sf4avFObm1RhSZVMjkSets8iB4DrQZBlZaeCZhPqRArJO5qG0Gd8sV5eeoGTicmXz7Axxe6vAzHdsCKrz2cTkREvFlqqsXzk9dTsqg/z3Srk/Hgc8fgl4egTC3o/KJnAuZjKsQKCX9fH167riEHTiTw3pytaQ+KvBoi2sEfr8O5454NKCIiXmv88t2s3nOcF3rVJSQokw36v092jz8AACAASURBVD4Dpw/BdZ/qkWQWqBArRJpVLcWtLarw1eKdbNh/4vIBxkC31+zzJxe87fmAIiLideJOJfLmr5tpXT2UPo3DMx684RdY+wO0fxrCm3omYD6nQqyQeaZ7bUoW9ef5yetJTU1j436FK6DpHbDsEzgS4/mAIiLiVYbP3ERCUir/d12DjDvonzoE0x+Hik10lmQ2qBArZEoGFeGFXnVZvec445fvTntQpxfBP8g+F0xERAqtxbFHmLxqHw+0r06NsOD0B1oWTHsUks7ajyR9M3l8Kf+jQqwQ6tM4nNbVQ3nz183EnUq8fEBwWbjqKYiZBTG/ez6giIg4LjE5hRd+WU/V0CAGdayZ8eBV38LW3+yXvsJqeyJegaFCrBAyxvB/1zUgISmV4TM3pT2o5QNQujrMGqp2FiIihdCnf25n+5EzvNK7AYH+vukPPLYTfhtqv+zV4n6P5SsoVIgVUjXCgnmgfXUmr9rH4tgjlw/wKwLdhsORrbB8jOcDioiIY3YcOcOH82Pp1agC7WuFpT8wNQUmPwjGB/qMBh+VFdml/8UKsUEda1I1NIgXfllPYnIah4LX6g7VO8Ifw+HMUc8HFBERj7Msi5emrCfA14cXe9XLePDSj2H3YujxJpSs7JmABYwKsUIs0N+XV3o3YPuRM3z65/bLBxgD3V+HxNN2MSYiIgXetLUHWBhzhCe71aZcicD0Bx7aaJ/GUqcXXHGr5wIWMCrECrn2tcLo1agCH86PZeeRM5cPKFsXmg+A6C/h0AbPBxQREY85cS6JV6dvpGF4CP1aVU1/YPJ5mHw/BJSAXu/ZP7hLjqgQE17sVc9egp6yPu1DwTsMtf9j+22o/YqyiIgUSO/M3sLR04kMv64hvhkd6r3gLTi4Fv4zCoIz2EMmmVIhJpQrEciT3WqzMOYI09YeuHxAUGno+Dzs+BO2zPR8QBERcbs1e47z7dJd3Nk6goaVQtIfuGc5LHwHGt8OdXp6LmABpUJMAOjXqioNw0N4dfpGTpxLo11F1D0QVsdu8pqcRu8xERHJt5JTUnlu8jrCggMY0rVW+gPPn7UfSZYIt/cQS66pEBMAfH0Mw69ryNHTiYyYtSWNAX52O4tjO2DpaM8HFBERt/lmyS427D/JsGvrUzwwg674s56D+O3Q52MIzGDVTLJMhZj8T8NKIfRvE8G4ZbtYvjP+8gE1O0OtHrBghH2mmIiI5Ht74s8yYvYW2tcK45qG5dMfuGk6rPgK2j4K1a7yXMACToWYXOLJrrUJL1mUZ35eS0JSGr3Fur0GyQkw7xXPhxMRkTxlWRbPTV6HAYZf3zD9Q71P7oepD0OFxtDxBY9mLOhMmm/J5QNRUVFWdHS00zEKpIUxcdwx5m8e7FCDZ7rXuXzArOdhyUcwcD5UbOL5gCJuFPHsDKcjiHgVH1IZ5z+cxj7b6Hl+ODusCk5HylM73/DMCwfGmBWWZUX9+7pWxOQy7SLDuDmqEp8t2M76fScuH9D+aQgKhV+fVTsLEZECbqDvdNr4bmRYcv8CV4R5AxVikqbne9YjtFgRnpq4lqSU1EtvBoZA5xdhz1JY+6MzAUVExO0amW0M8fuJ6Skt+SmlvdNxCiQ9mpR0zd5wkIHfruDJrrV4uFPkpTdTU2FMFzi+Gx5eDkVLORNSRESybfra/Tz8/SqG9qjD/e1rpD0o8TR82s7uov/gIv3/fC7p0aRkW9f65enVqAKj5sYSc+jUpTd9fKDXSDh7FOa+6kxAERHJtvgz5xk2ZQONKoUw4Mpq6Q/89RmI3wHXf6YizI1UiEmGXv5PfYoF+PL0z2tJSf3X6mmFK6DFQPscyr0rnAkoIiLZ8k/j7jdvaISfbzplwPqfYfU4uOpJiGjr2YCFjAoxyVCZ4ACGXVufVbuP8/XinZcP6Pg8BJeDGY9DahrtLkRExGvM33yYyav2MahjTepWKJH2oOO7YdrjUKk5tH/GswELIRVikqnejSvSqU5ZRszawu6jZy+9GVgCug+HA2tg+RfOBBQRkUydSkjiucnrqFUumIc71kx7UGoKTBoIVipc/zn4ZtBlX/KEn9MBvNrLOr4BwACvWaXpmvQWQ9/5gHH+w0mz59+vT9u/RPKzl9No2SJSALz+62YOnUxgdL+2FPFLZx1m4Tuwewlc9xmUzmD/mOQZrYhJllQw8Qz1+56/UhvwQ0oHp+OIiEg2LNl2lO+X7WbAldVoXLlk2oN2L4M/3oCGN8MVfT0bsBBTISZZdqvvPFr7bOC15H4ctPQGjYhIfnDufArPTlpLRGgQT1xdO+1BCSdg0r0QUgl6jvBswEJOjyYzokcUlzDAG0fP0O29BbxQ/Sc+vzPq0nPJkhJgdGt75IOLwT/QqagiIuLy7pwt7Dp6lgkDW1G0iO/lAywLpj8BJ/bBPb/ZTbvFY7QiJtlSNbQYT3atze+bDjNt7YFLb/oHQs93IH4b/PW+MwFFROR/Vu85zphFO7i9ZRVaVQ9Ne9CKr2H9ROg4FCq38Gg+USEmOXB3W3uPwctTN3DkdOKlN2t0gvrX2xs+j25zJqCIiJCQlMLTE9dQrkQgz/aok/agA2vtxq01OsGVQzwbUAAVYpIDvj6Gt29sxOnEZIZOWsdlx2R1Gw6+RWDmUzoUXETEISPnbGXrodMMv74hxQPTaEORcBJ+uguCSttvSfqoJHCC/leXHIksV5ynu9VmzsZDTFyx99KbJSpApxdg21zY+IszAUVECrFl24/y2cLt3NayCh1rl718gGXBtMFwbAfcMAaCwzwfUgAVYpIL97StRstqpfnvtI3sif9Xo9fm90L5RvDbUPunLhER8YjTickM+WkNVUoH8fw1ddMeFP0lbJhk/9CsI4wcpUJMcszHx/DOzVcA8ORPa0i9+CxKXz/o9R6cOgh/vO5QQhGRwufVaRvZf/wc7958BcUC0miOcGCN/UNyjc7Q9nHPB5RLqBCTXKlUKohh19Zj2Y54vvxrx79uNoOou2HZJ/Z/+CIi4lZzNh7ih+g9PNC+Bs2qlr58wMX7wq7XvjBvoD8BybUbm1Wia71yvDVrC1sPnbr0ZueXICgUpj4CKcnOBBQRKQSOnk5k6KS11K1Qgse61Lp8wP/2he2EG7+EYmU8nlEup0JMcs0Yw/DrG1Ii0I/HJqzmfHLqhZtFS9m9xQ6sgcWjnAspIlKAWZbFc5PXcfJcMu/1bZz2WZLRYy7sC6vaxvMhJU0qxCRPlAkO4PXrG7HxwElGzY259Ga93vavP96AuC3OBBQRKcB+XrmPWRsO8WS3WtQuX/zyAf/sC6vZRfvCvIwKMckzV9crx81Rlfj4j1hW7Dp26c1rRkCRIJjyMKSmOBNQRKQA2nvsLP+duoEW1Uoz4Mrqlw9IOAk/9oegMnDdp9oX5mX0pyF56sVe9ahYsihDflzN2fMX7QkLLgs93oK9f8OyT50LKCJSgKSmWjz101pSLYt3broCXx9z6QDLsvfoHt+tfWFeSoWY5Knigf68c9MV7Io/y/CZmy692fAmqNUd5r4C8dudCSgiUoB8+dcOlmw/yrBr61O5dNDlA5Z/YTfW7vQCVG3t+YCSKRVikudaVg/lvnbVGbd0N39sOXzhhjHQayT4+sPURyE1Nf0vIiIiGYo5dIq3Zm2hS91y3BRV6fIBe1fArOeg5tXQ9jHPB5QsUSEmbvHE1bWoXa44T09cy7Ez5y/cKFERur0GOxfCiq+cCygiko+dT07l8R9XUzzAj9evb4gx/3okeToOfrwDgsurX5iX05+MuEWgvy/v9r2CY2fP88KU9ZceDN7kDqjeAea8BMf3OBVRRCTf+mBeDOv3neS16xoSVjzg0pspyTDxbjh7FPp+azdvFa+V60LMGONrjFlljJnu+ryaMWaZMSbWGPODMaaI63qA6/NY1/2Ii77GUNf1LcaYbrnNJN6hfsUQHutSixlrDzBl9f4LN4yBa0ddaC54cZEmIiIZWrHrGB/Nj+WGppXo3qD85QN+H2Y/deg1Eio29nxAyZa8WBEbDFy8K/tNYKRlWTWBY8AA1/UBwDHX9ZGucRhj6gG3APWB7sDHxhjfPMglXuD+q6oTVbUUL/yynl1Hz1y4UaoqXP1f2DYXVn/vXEARkXzkxLkkHh2/ioolizLsP/UuH7D+Z1jyITS/Dxrf5vmAkm25KsSMMZWAnsAXrs8N0AmY6BoyFujj+ri363Nc9zu7xvcGJliWlWhZ1g4gFmiRm1ziPfx8fXj/1ib4GHhk/KpLu+5HDYAqbWDWUDh5wLmQIiL5gGVZPDdpHYdOJjDq1iaUCPS/dMChjTDlEajcEroNdyakZFtuV8TeA54G/vnbNRQ4blnWPw2k9gLhro/DgT0ArvsnXOP/dz2NOZcwxgw0xkQbY6Lj4uJyGV08JbxkUd668QrW7j3B27M2X7jh4wO9P4TkRJjxhB5RiohkYPzfe5ix7gBPdqtN0yqlLr157jj80A8CguGmseBXxJmQkm05LsSMMb2Aw5ZlrcjDPBmyLOszy7KiLMuKCgsL89S3lTzQvUF57mhVlc8X7mD+5otaWoTWsPvbbJlpL6mLiMhlthw8xX+nbaBdZBkGtvtX9/zUVJj8ABzfZRdhJSo4E1JyJDcrYm2B/xhjdgITsB9Jvg+UNMb4ucZUAva5Pt4HVAZw3Q8Bjl58PY05UoA837MudcoXZ8hPazh0MuHCjVaDILwZ/Po0nDniXEARES907nwKD3+/kuKB/rx7c2N8/t09f+EI2Pqr/ThSTVvznRwXYpZlDbUsq5JlWRHYm+3nWZZ1OzAfuNE1rD8wxfXxVNfnuO7Ps+yeBlOBW1xvVVYDIoG/c5pLvFegvy8f3taEc+dTePyH1aSkuh5F+vhC748g8RTMfNLZkCIiXuaV6RuJOXyakX2vuLxVxdbZMH84NOoLLQY6E1ByxR19xJ4BnjDGxGLvARvjuj4GCHVdfwJ4FsCyrA3Aj8BG4DfgIcuydCp0AVWzbHH++5/6LN52lNF/xF64UbYutH8GNkyGtT86F1BExItMX7uf8X/v5sEONWgX+a8tOfHbYdK9UK4B9HrPbg0k+Y6x8ukG6aioKCs6OtrpGJIDlmUxeMJqZqw7wISBrWge4Wo2mJoCX10DhzfCg39BySrOBhURcdCe+LNc8/5CapYL5sf7W+Pve9HayfmzMOZqOLEXBv4Bpas5FVOyyBizwrKsqH9fV2d98ThjDK9d14DwkkUZPH4Vx8+6jkDy8YXrP7Xfnpz8gF2YiYgUQkkpqTwyfhUYGHVLk0uLMMuCaY/CoQ1wwxgVYfmcCjFxRPFAfz64tQmHTyXyzM9rLxyBVCoCrnkbdv0Ff73vaEYREae8M3srq/cc543rG1G5dNClN5d8BOt+gk7PQ2QXZwJKnlEhJo65onJJnuleh1kbDjFu6a6LbtwC9frA/Ndg/2rnAoqIOGDB1jg++XMbt7WsQs9G/2pFsXU2zHkR6l4LVw5xJqDkKRVi4qgBV1ajQ+0wXp2xiY37T9oXjbHPSCtWFibdZ++FEBEpBA6fSuCJH1dTu1xxXur1ryOMDm+GiffYm/Ov+9Ruii35nv4UxVE+PoYRN11BSFF/Hh6/krPnXYcyBJWG60bDka32T38iIgVcaqrFEz+s4XRiMh/c1oRA/4uOXT5zFMb3Bf+icOt4KFLMuaCSp1SIiePKBAfwXt/G7DhyhpembLhwo3oHaP0wLP/CXo4XESnAPlmwjUWxRxh2bX1qlSt+4UbyefjxTvtM3lvHQ0gl50JKnlMhJl6hbc0yPNyxJhNX7GXC37sv3Oj8EpStD1MGwWmdLyoiBdPi2COMmLWFXo0qcEvziw6bsSy70fWuRfbZvJUu634g+ZwKMfEaj3WpRbvIMrw0ZQNr9hy3L/oFwA2fQ8JJmPqIDgYXkQJn//FzPDJ+FdXDgnnzhkaYixuzLvsUVo6FK5+ARjc7F1LcRoWYeA1fH8P7tzQhrHgAg75bSfwZV3+xcvWhy8v2WWorvnYwoYhI3kpMTuHB71aSmJzKJ/2aUSzA78LN2N9h1lCo0ws6aa9sQaVCTLxK6WJFGN2vKXGnEhk8YdWF8yhbPgDVO8Ks5+BIbMZfREQkn3hl2kbW7DnOiJsaUbNs8IUbcVvhp3ugbD29IVnA6U9WvE6jSiV5pXd9FsYc4d05W+yLPj7QZ7T9qHLSvZCS5GxIEZFc+il6D98t28397avTvcFF/cLOxttvSPoVsTfnBwSn/0Uk31MhJl7plhZV6BtVmY/mb2P2hoP2xRIV7INt96+CP990NqCISC6s33eCF35ZT+vqoTzVtfaFGylJ8FN/+wzJvt/pzN1CQIWYeK3/9q5Pw/AQhvy4hh1HztgX6/eBxrfDwndg5yJnA4qI5MDxs+d58LsVlC5WhA9ua4LfxedI/vYs7FgA146CKi2dCykeo0JMvFagvy+j+zXF19fwwLcrLjR77fEmlK5ud5g+fdjZkCIi2ZCaajF4wmoOnkjg49ubUiY44MLNvz+3+ya2HQyNb3UupHiUCjHxapVKBTHqliZsPXyKoZPW2YeDBxSHm7+BhBPw8wBITXE6pohIlrw/N4Y/t8Yx7Nr6NKlS6sKNzTPh16ehVg/oPMy5gOJxKsTE611VK4whV9diyur9jF28075Yrj70fMdewv/jDUfziYhkxbzNh3h/bgw3NK3E7S0v2vu1d4W9wl+hMdw4Bnx80/8iUuCoEJN8YVCHmnSpW5b/m7GJ5Tvj7YtN+kHjfrDgbbvfjoiIl9p19AyPTVhNvQoleO26BheatsZvh+9vhuCycNsPOkOyEFIhJvmCj4/hnZsbU6lUUQZ9t5LDJxPsG9e8bffZ+fk++y0jEREvc+58Cg+MW4kxhk/6NbtwmPeZozDuRrBSoN/PdjEmhY4KMck3Qor688kdzTiVkMTD368iKSUVigTBzWMh5by9tK/+YiLiRSzL4vlf1rH54Eneu6UxVUKD7BtJ52D8LfYPkLdOgDKRzgYVx6gQk3ylTvkSvHlDI/7eGc+wqRvszftlIuE/o2DPMvj9Zacjioj8zxcLdzBp5T4Gd46kY23XildqCky6D/Yut8/SrdLK2ZDiKL/Mh4h4l96Nw9l04BSf/LmNyLLB3N22GjS4AXYtgSUf2v+nVvdap2OKSCE3d9Mhhv+6iWsalufRTheteM16HjZNg26vQ73ezgUUr6AVMcmXnu5Wm671yvHq9I3M3+LqJdbtNajYBH55COJ3OBtQRAq1TQdO8uj4VTSoGMI7NzXGx8e1OX/JR7BsNLQaBK0HORtSvIKxLMvpDDkSFRVlRUdHOx1DHHQmMZkbP1nCnvizTBrUhlrlisOxnfDpVVAqAu6ZDf6BTseUfCbi2RlOR5AC6hqfpXxcZBQzU1rwUNKjWFoL8Qo73+jpke9jjFlhWVbUv6/r3wLJt4oF+DGmfxSB/r4MGLuco6cT7QKszydwYA3MGup0RBERAKLMZkb6jyY6tRaPJw1SESb/o38TJF+rWLIon9/ZjEMnE3lg3AoSk1OgzjXQ5lGI/hLW/uR0RBEp5GqYfXxR5B32WmW49/wQEinidCTxIno0mZGXQ9z79SXPTE1pzaNJj3Cj75+87fcp//RKFMm2l084nUDyqY/mx/L2rC0MuboWj3R2bc4/eQC+7Gq3qxgwB0pXczakOCa9R5N6a1IKhP/4LmFbakXeT7mBmmYfD/hNdzqSiBQiv60/wNuzttC7cUUe7lTTvng2Hr69zv69/1QVYZImFWJSYAz2m0SsVZE3k2+hujlAV98VTkcSkUJg/b4TPP7DGppUKcmbNzSyjy9KPAXjbrCPMOo3EcKbOR1TvJQKsYzoEUW+4gO8k5TC3k+X8Njhp/lpYGvqVwyBFV/DtMHQ6iHoPtzpmCJSgBw6mcC9Y6MpFeTPZ3fYLw+RlADjb7VfGuo7Dqpd5XRM8WLarC8FSqC/L5/fGUWJQH/uGxvN4VMJ0OwuaHE/LP0IVn7rdEQRKSDOnU/hvm+iOZmQxJi7mhNWPABSku3j1nYuhD6j7ZeHRDKgQkwKnLIlAvmifxTHziYx8JsVJCSlQLfhUL0DTH8cdi91OqKI5HOpqRZPTlzDun0nGHVLE+pWKAGpqTDlIdgyA3q8DVf0dTqm5AMqxKRAahAewsi+jVm95zhP/LiaFOMLN34FJSvDD/3g+B6nI4pIPjZi9hZmrD3A0B516FKvHFgW/PYsrJ0AHZ+HlgOdjij5hAoxKbC6NyjPCz3rMnPdQV6ZtgGraCm49QdITrT3b5w/43REEcmHvv5rBx//sY3bWlbhvnbV7Yt/vA5/f2rvRb3qKWcDSr6iQkwKtHvbVee+dtUYu2QXH/+xDcJqwY1fwuENMPkB+1GCiEgWTV+7n/9O30i3+uV4tXcD+w3JJR/Dn29C4372mbdqZCjZoEJMCryhPerSp3FF3p61hR+X74HIq+HqV2DTVFjwltPxRCSfWBx7hCd+WEPzqqV5/5Ym+PoYWDXOPk6t7rVw7fsqwiTb1L5CCjwfH8NbN17B0TPnGTp5HaHBRejc+mE4tNF+nFC2LtTr7XRMEfFi6/edYOC3K6hWphif3+lqU7FxKkx9BKp3hBvGgK/+SpXs04qYFApF/HwY3a8Z9SqU4KHvV7Ji93HoNRIqNbcfUR5Y63REEfFSu4+e5a6vllMi0I+v72lOSJA/xP4OPw+wG7X2HQd+AU7HlHxKhZgUGsEBfnx1d3PKlwhkwNjlxB5Lgr7fQdFS9ub904edjigiXubI6UTu/HIZyampfDOgBRVCikLsXBh/G5SpDbf/BAHBTseUfEyFmBQqZYID+Oaelvj5+HDnmL85kFoCbvkezh6FCbfbB/OKiABnEpO55+vlHDyZwJj+zalZtrhdhE24DcrUgjun2D/IieSCCjEpdKqEBvH13c05mZDMXV8u50TJ+nDdJ7B3OUwcYHfGFpFC7XxyKg+MW8GG/Sf56LamNKtaCrbNs4uw0Jp2EVYs1OmYUgCoEJNCqUF4CJ/d0YztR05z3zfRJNS6Fnq8aXfEnjnEbs4oIoVSaqrF0xPXsDDmCK9f35DOdcvBtvn2FobQmnDnVBVhkmdUiEmh1aZmGd69uTHLd8UzeMIqUpoPhCsftw8J/1NtLUQKq9d/3cQvq/fzVLfa3BxV2VWE3QKla6gIkzynd22lULv2ioocOZ3If6dt5IVf1jO8z0uYUwfhj+FQvDw06+90RBHxoM8WbOPzhTvo37oqgzrUuLQI6z9NRZjkORViUujd3bYacacS+fiPbQT6+/DStaMwZ+Jg+mNQLAzqXON0RBHxgLGLdzJ85mZ6NqrAS9fWx+z486IiTCth4h56NCkCPNWtNve0rcZXf+3k9dnbsG76Gio0hol3w+5lTscTETf7btkuhk3dwNX1yvFe38b47vwTvr+4CCvjdEQpoFSIiQDGGF7sVZc7WlXlswXbGfHHPqzbfoQSFWF8X4jb4nREEXGTH5fv4fnJ6+lUpywf3tYE/10LXUVYNRVh4nYqxERcjDH89z/1ubVFZT6av433lx6DfpPAxx/G3QAn9zsdUUTy2ORVe3lm0lraRZbh49ubErB7EXzf11WETVMRJm6nQkzkIj4+htf6NOTGZpV47/cYPlqTAv0mwrnjMO5G+3cRKRCmrtnPkB/X0Lp6qH1+5PbZ8N1NUCrC9XakijBxPxViIv/i42N484ZG9GlckbdnbeGzmGDo+y0c2erqvp/gdEQRyaVf1x3g8R9WExVRmi/6RxG46Wf7v+9y9eGuGRAc5nREKSRUiImkwdfHMOKmK+jVqALDZ27mywMRdvf9XYtg8kBITXE6oojk0JyNh3hk/CoaVy7Jl3c1J2j1VzDpPqjaRm9HisfluBAzxlQ2xsw3xmw0xmwwxgx2XS9tjJljjIlx/V7Kdd0YY0YZY2KNMWuNMU0v+lr9XeNjjDFq3CRewc/Xh5F9G9O9fnlemb6Rb09HQbfhsHEKTH0EUlOdjigi2TR/82EGfbeCBuEhfH1XFMFLR8LMJ6H2NXD7RAgo7nREKWRysyKWDAyxLKse0Ap4yBhTD3gWmGtZViQw1/U5QA8g0vVrIDAa7MINGAa0BFoAw/4p3kSc5u/rw6hbm9ClbllenLKBCb7XQoehsPo7mPaoijGRfGRhTBz3j1tB7fLFGXt3c4ov+C/M/z9odAvc/C34BzodUQqhHBdilmUdsCxrpevjU8AmIBzoDYx1DRsL9HF93Bv4xrItBUoaYyoA3YA5lmXFW5Z1DJgDdM9pLpG8VsTPh49ub0qH2mEMnbyOicX7wVVPw6pvYcbjKsZE8oHF245w79hoaoQFM+7uZoTMfhyWfAgt7oc+o8FX/c3FGXmyR8wYEwE0AZYB5SzLOuC6dRAo5/o4HNhz0bS9rmvpXU/r+ww0xkQbY6Lj4uLyIrpIlgT4+fJJv2ZcWbMMT01cw4/Bd8CVT9jnUs58UoeEi3ixRTFHGPB1NBGhxfjursaUnDEQVo+D9s9CjzfBR9ulxTm5/rfPGBMM/Aw8ZlnWyYvvWZZlAXn2N5RlWZ9ZlhVlWVZUWJjeaBHPCvT35bM7omgXGcbTk9YxJuAOaDsYosfAr0+rGBPxQr+tP8g9Xy+namgQ4+6sT+kpd8CmafD/7d13dFVl3vbx7++ck05CAgk1FGkWQOlgwQEVAcXCw6AoY0F9FUVl3sdxGMtSR2d8GQd1fB2n0CwMCCigjoIoZWiK9N5EaqgJIZ3Ucz9/JONyfEQQkuxzwvVZK4uczU7Ote6VrH1l3/veu99o6P0EmHkdUc5xZ1XEzCyC8hI22Tk3s2LzkYopRyr+PVqx/QDQ5Dtfnlqx7WTbRUJOTKSfcXd2pn+7BrzwyVZedUNxPUbAirHw6RMqYyIhZMbqNEZMWUO7xglMv+MCUmbdCrsXlU9F9njQ63giwNmtmjRgArDVOffKd/7r7ZJ4nwAAFVJJREFUI+DfKx/vAj78zvY7K1ZP9gCyK6Yw5wLXmllSxUX611ZsEwlJUQE/r9/WkcGdU3ltwU6eL74d1304fPVX+OxplTGREPDWst089l75zVr/MTiVhGk3w6H15Rfld7jd63gi3zqbqxMvB+4ANprZuoptTwKjgelmdi+wF7il4v9mA9cBO4ECYBiAcy7TzF4AVlbs97xzLvMscolUuYDfxx8GXUx8dAQTl+0mt9NtvNS1DN+XfwafH675raY8RDzgnOPPC3by8uc76Nu2Pq/39hP5Tl8oyoOh70GLXl5HFPkPZ1zEnHNLgZMdaa7+gf0dMOIk32siMPFMs4h4wecrf1B47ZgIXp23g7yLhvBG5zL8y14D88PVz6iMiVQj5xwvzt7KuCW7GdQplZcuPoj/7XshJgnunVt+13yREKP1uiJnwcwYeU1r4qMDPP/xFoa1uoWJHUoJLH0FfAG46imvI4qcE8qCjidnbmTaqv3cfVlznqm/DN+0UdCgPdw2DRIaeh1R5AepiIlUgnuuOI/46ACjZmxgSPGtTGlfRuTil8AF4aqndWZMpAoVlwb5v9PW8cnGQ4y8qgW/DL6DzfkLtOkPP58AkXFeRxQ5KRUxkUoyuEsT4qMDPPLuWm4uvpUZ7SBmyRg4kQnXjSm/dkxEKtWJ4jKG/2M1i3ak81y/5tx9+Hew7WPoPrz8kWT6vZMQp7vYiVSifu0aMuGuruw+doLrdg8mp/PDsGoivD8MSou8jidSo+QUlnDnxK9Y8nU6r13fkLt3jIDts6HfHypu1KoSJqFPRUykkl3ZJoV/3NeNjPxirt7Qm4Pdnip/UPjkwVCU63U8kRphf2YBP//rF6zbn8Wb18dz06o7IX07DJkCPYZ7HU/ktKmIiVSBzs3qMOPBy4j0+7jqy/Zs7Doa9iyFtwZAfobX8UTC2tp9xxn4l2Uczi7kg34l/GzJUCgrhmGz4fz+XscT+UlUxESqSJv68Xww4nIuaJDAjUub8mn7l3Hp22BiXzi+1+t4ImHp4w0HGTJ2ObERfuZd+Q1tF94DtZvAffOhUUev44n8ZCpiIlUoJT6Kqff34Lp2DRm+oh5/b/YyLj+9vIwd2eJ1PJGw4ZzjjYU7eXjKWjo1iuHzVu9Rb/ET0KI33PMpJDY59TcRCUEqYiJVLDqi/JFID/duxejNSTxR+yWCLghv9oN9y72OJxLyikuDPP7+Bv44dzt3t/Uz2f8cURunwJW/htunQXSC1xFFzpiKmEg18PmMX/U9nzGDL2FGWgJ38DtKourAOzfDDj1aVeRksgqKuWPCV7y/Oo2Xu2Tz7KGH8GV+A0PeLb9hslZGSphTEROpRj/vnMo793RnU34i/XOfJr92S3j3Nlg/1etoIiFnd0Y+A//yBWv3HefjzmsYtHkEFpcC/2cBXHCd1/FEKoWKmEg1u7RlXWY9dBmlMXW54sh/k163C8x6AOa/AMGg1/FEQsKK3ZkM/MsyigtyWN56Eu02j4ELb4D75kFya6/jiVQaFTERD7RIqcWshy6ndWojLk97iE31b4IlY2DaUCjM8TqeiKdmrE5j6PjlXBKTzsLEF6iz91Po8zwMfhui4r2OJ1KpVMREPJIUF8mk+7oxoFNzBuy9hcl1RuB2zIUJfSBzl9fxRKpdYUkZT83ayGPvreeB+tt5q2QUkScy4I5ZcPlIPbNVaiQ9a1LEQ1EBPy8PvoRLUhN57hMfK+Lq8XLOqwTG9oZb3oYWvbyOKFIt9mcW8NDkNWw5kMnkFgu4/OCb0LAD3DoJEpt6HU+kyuiMmIjHzIy7LmvO9AcuZaW1p2/ec2QF6uIm/Rcs/xs453VEkSo1f+sRBry+lJJju1nd+JXyEtbhFxX3B1MJk5pNRUwkRHRsmsQnj/YktWVbLs94kk1xPeDTUfDRI3pguNRIpWVBXvp0G/e+vYpfxC5ndsRvSMz7BgZNgJvfgIgYryOKVDlNTYqEkKS4SN68uyt/XpjETfMe5IWExgxdOwkyvi6foqlVz+uIIpUiPbeIR99dy6Zd+5nVYBodsz6HJj3gv8ZCUjOv44lUGxUxkRDj8xmPXt2ajk0TGTk1mtXBRvzx4N/wj+0NQyZDow5eRxQ5Kyt2Z/LwlDW0KNzM8qS/E5d9GHo9CT0fA78OS3Ju0dSkSIjq2TqFTx69gj0NruXGgmfIPlGCm9gP1kzSdWMSlpxzjF38DUPHfcFwm8G7gd8SF+WHYXOg1yiVMDkn6adeJIQ1rB3D1PsvZfScJK5alsSb8X/l4o8ehm/mw4A/QUyi1xFFTsvx/GJGzdjA5i2bmJM4nlaFG6H9LXD9GIiu7XU8Ec+oiImEuMiAj2duuIguzZP4xftJDHMfMnLLe1jaSmzQBGjaw+uIIj9qwbYjjJqxkctO/IsFtd4kMmgwcCxccqvX0UQ8Zy5Mpzi6dOniVq1aVbVv8pz+SpPQkuaSebzkAQpcFOMiXyGFLN3jspI1L5zidYQaJ4XjPB/xFv39K1kTbMXIkhHsd/W9jiUCwJ7R11fL+5jZaudcl+9v1zViImEk1TKYHPEiA/1L6V/0Ip8EuwO6ZExCleNW/0LmRz3OVb51vFRyK4OLn1UJE/kOTU2KhBmfOe4OfMaVvg08VvIg88o6MzpiHJGuBJ/OjkmIaGaHGR0Yz6X+LSwPXsgTJfex2zX0OpZIyFER+zHPZXudQOSkWgDvlQX5++JdXD+vDa9F/Jl27ISOd0D/P0BknNcRw9IerwOEsTX7jvPraWu4Jms6j0XOJBAZDde+Ro+Od7LQpwkYkR+iIiYSxgJ+HyN6t+KqC+rx+LTm9M+YyINr/4Hb+yX+n4+HRh29jijngOLSIK/N38HiRfN4I2o850fshvMHwHVjIEFnwUR+jP5EEakBLmyYwMyHf8aJnk9zR/GTZGZm4sZdBXOfgqI8r+NJDbbpQDaDX59H/JIX+CDyGVrH5sMtk8pvPqwSJnJKWjUpUsOs3XecZ6cvY0jWBG4PLKA0vjGBG16FNn29jiY1SFZBMWM+287+Ff/kd5Fv04TD0Oku6PO87m8n8gNOtmpSRUykBiosKWP8kl18sfBjfusbT2tLo+zCG/H3f0lnKeSsBIOOaav28/6cz3mk7G16+dZTltQC/42vwXlXeh1PJGSpiImcg9KOFzD64w002zaBkRGz8EVEEejzHHS5F3TxtPxE6/dnMWbWMvocfZOhgfkQWQt/r1HQ7X4IRHodTySkqYiJnMMW70hn7IfzGJ7zOlf4N1PUoDNRA1+H+m29jiZhIDO/mFfmbCB23QQeCXxAnBVhXe7Bej0BcXW9jicSFlTERM5xxaVBJizZxZ6FExll71DbV4Dr8TCB3r+ByFiv40kIKgs6pny1l7VzJzHSTaKZHaG05TUE+r0IKed7HU8krKiIiQgAh7JP8KePltN5+6vcElhEQWxjYvo+i7UfrOlK+dbqvZm8M/NDbj/+N7r7tlGUdD5R1/8/aHW119FEwpKKmIj8hy92ZjBj5lTuzRvLRb695CVdSFz/57HWfdADLM9dG9OymTx7Pl33T2CgfxklkUlE9nka63QX+HXrSZEzpSImIv9LSVmQ91buY/u8N7mneArNfEfJqd+NhAEvQpOuXseTarT1UA5TP/mMTnsnMMD/JUFfJK7r/UT2fhyia3sdTyTsqYiJyEkVlZYxffkuDi34G8PKppNiORxvei1JN/xO1wLVcDuP5jL9kzlcsmsc/X0rKfVH47reR1TPR6FWPa/jidQYKmIickqFJWVMXbqVE4v/P78IfkScFZHVZhB1rn8Oaqd6HU8q0Z6MfGZ+/E/afzOWPv7VFPnjcN3uJ/qKR7QSUqQKqIiJyGnLLypl+qJ1BL54hVvcXHw+I6f9MOr2+RXE1/c6npyF/ZkFfPjxLNrtHEsv3zpO+BNw3YcT23OE7ogvUoVUxETkJ8stLOH9+V9SZ+UYbnCLKbMAGefdQP1r/xtfw/Zex5PT5JxjzZ5jrJ43jbb7p3C5bxP5gURcjxHUumI4RCd4HVGkxlMRE5Ezll1QwgfzFxG7dhzXly0k1oo4mNSNxKtHEnvRdbrtRYgqLCnj85UbyVoynt4Fc0i1DHIikuHSESRc8QBExnkdUeScoSImImetuDTIvDVbObZoHNfkfUhDyyQjqgnB7sOpd8UwHdhDxKGsAv712QfU2TqJ3sGviLQyDtXtTtLPHiS67QDwR3gdUeScoyImIpVqw750Nsx9i0v2T6G9bxd5Fk/6+bfRtN8v8Sc29jreOcc5x5ode/lm3jg6HJlFG98B8n3xZLUZTKNrHsKSW3sdUeScpiImIlUiI7eQRfP+SZ0N47ky+BXOfOxJ7kXtbkNJ6TQAAlFeR6zRjuYUsGrpZ9i6yfysaBGxVsSBuIuIvvR+6nYfAhExXkcUEVTERKSKlZQFWbJiNQVL36BH3nySLYc8iyOtQR+SLx1Kcrurwef3OmaNkJlXyOoln1K2aRaX5C2moWVSSBRpTQaQes0Iopt19jqiiHyPipiIVJu0jGw2LP6QqO2z6F74BbWskExfHQ6l9qdhzzuo06qHHqP0E2XnF7Jm6RxKN87i4txF1LcsiolgX53LiOs4iIZdB2r1o0gIUxETEU/sOZzB1kXvkbDzA7oUryLKSjnkb0T6eTfQ9Mo7SWzSVqXsJPJOFLFu6WyKN8ykXc5i6lkWRUSyJ+kyYjsMIrX7zZjKl0hYUBETEc99sy+NnYunkrz7QzqWbsRnjnRfMofqdCPQsjfNuvYjLrmp1zE9U1hcwvYNKzi2eT6xB5fTpnA9dSyPQiLZlXQ5sR0G0azHzVhUvNdRReQnUhETkZDhnOPrnV+zf/lMYg8s5YIT60iyXADS/KmkJ3cnus3VNO/Sl5jayR6nrTolpaXs2LiCjI3ziTn4Ja1PbPh2HI746nOkThdi2van5aUD8UXX8jitiJwNFTERCVmFxSVsW/clWZvnEX9oGRcUbSTOigg6Y3dkKzLrXUp0s84kt+hA/eZt8QXC7z5YzjnSjx7i0Ndryd2zlqgDX9KqYP23xeuwrz6Hk7oQ0bInTTtdS3yDlh4nFpHKFPJFzMz6Aa8BfmC8c270j+2vIiZSc+UVFLBjzb/I3zqfpCNf0qZkG5FWBkCxC5AWaEpWrZaUJl9IbGp76rXqSErjlliI3OE/M+MIB79eS+6+jbj0bcTnfE3D4j0kk/3tPoesHoeSuhBo0ZNmna6ldqNWHiYWkaoW0kXMzPzADqAPkAasBG5zzm052deoiImcO3Jzc9j/9Xqy96wneHgLcTk7aFC4iwYc+3afPGI4ENGc3NimlMXUhbhk/PH1iKhdn7jEBsTXbUjtlEZEx5z53f9LiovIyjhM9rGDnDh+hKLsI5TkHMXlZ+A/kUFswQEaFO0hhePffk0+0RyMaEZOrZYEUy6gVmp7GrTpRFKD5mczJCISZk5WxAJehPkB3YCdzrldAGY2FbgJOGkRE5FzR3x8Ahd16gmdev7H9uPH0jnw9Rpy923Ajm4hIXsHjbPXkJSVRbSV/OD3ynMxZPlqU+BPwHE6Z9AcsWW5JLhsapNPCpDyvT1KnY8sS+B4IIW9id3ZlXwBsY3b0aB1B5IbtaR1iJypE5HQEypFrDGw/zuv04Du39/JzO4H7gdo2vTcXVklIuWS6qaQVLcv9Oj7H9tdMEhefg7Z6QfJPXaIE1mHKck5SjA3HQoyCJw4RmRJFnB6MwJ50Q05GF0HF5uC1UohMiGF6KQG1EpqQO3kRsQn1iXZ56fmLisQkaoSKkXstDjnxgJjoXxq0uM4IhKizOejVnwiteITocVFXscRETmpUDlffgBo8p3XqRXbRERERGqsUCliK4HWZnaemUUCQ4CPPM4kIiIiUqVCYmrSOVdqZg8Dcym/fcVE59xmj2OJiIiIVKmQKGIAzrnZwGyvc4iIiIhUl1CZmhQRERE556iIiYiIiHhERUxERETEIypiIiIiIh5RERMRERHxiIqYiIiIiEdUxEREREQ8oiImIiIi4hEVMRERERGPmHPO6wxnxMzSgb1V/DbJQEYVv8e5RONZ+TSmlU9jWrk0npVPY1q5qms8mznnUr6/MWyLWHUws1XOuS5e56gpNJ6VT2Na+TSmlUvjWfk0ppXL6/HU1KSIiIiIR1TERERERDyiIvbjxnodoIbReFY+jWnl05hWLo1n5dOYVi5Px1PXiImIiIh4RGfERERERDyiInYKZvaCmW0ws3Vm9pmZNfI6Uzgzsz+a2baKMZ1lZoleZwp3ZjbYzDabWdDMtJLqDJlZPzPbbmY7zew3XucJd2Y20cyOmtkmr7PUBGbWxMwWmtmWit/3kV5nCndmFm1mK8xsfcWY/taTHJqa/HFmluCcy6n4/FHgIufccI9jhS0zuxZY4JwrNbM/ADjnRnkcK6yZ2YVAEPg78Cvn3CqPI4UdM/MDO4A+QBqwErjNObfF02BhzMyuBPKAd5xz7bzOE+7MrCHQ0Dm3xszigdXAzfoZPXNmZkCccy7PzCKApcBI59zy6syhM2Kn8O8SViEOUHM9C865z5xzpRUvlwOpXuapCZxzW51z273OEea6ATudc7ucc8XAVOAmjzOFNefcYiDT6xw1hXPukHNuTcXnucBWoLG3qcKbK5dX8TKi4qPaj/EqYqfBzH5vZvuBocAzXuepQe4B5ngdQoTyA9r+77xOQwc5CVFm1hzoCHzlbZLwZ2Z+M1sHHAU+d85V+5iqiAFmNs/MNv3Ax00AzrmnnHNNgMnAw96mDX2nGs+KfZ4CSikfUzmF0xlTEan5zKwWMAP45fdmbOQMOOfKnHMdKJ+d6WZm1T6NHqjuNwxFzrlrTnPXycBs4NkqjBP2TjWeZnY3MAC42ukixdPyE35G5cwcAJp853VqxTaRkFFxHdMMYLJzbqbXeWoS51yWmS0E+gHVusBEZ8ROwcxaf+flTcA2r7LUBGbWD/g1cKNzrsDrPCIVVgKtzew8M4sEhgAfeZxJ5FsVF5ZPALY6517xOk9NYGYp/165b2YxlC/WqfZjvFZNnoKZzQDOp3xV2l5guHNOfymfITPbCUQBxyo2Ldcq1LNjZgOB14EUIAtY55zr622q8GNm1wF/AvzAROfc7z2OFNbM7F2gF5AMHAGedc5N8DRUGDOzK4AlwEbKj0cATzrnZnuXKryZ2cXA25T/zvuA6c6556s9h4qYiIiIiDc0NSkiIiLiERUxEREREY+oiImIiIh4REVMRERExCMqYiIiIiIeURETERER8YiKmIiIiIhHVMREREREPPI/zB+WtujlFZUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "x=np.linspace(-3,3,50)\n",
    "V=1000*x*x\n",
    "Va=V*(1-0.04*x)\n",
    "n=np.arange(0,5,1)\n",
    "E_arm=2000*(n+1/2)\n",
    "plt.plot(x,V)\n",
    "plt.plot(x,Va)\n",
    "plt.barh(E_arm, 3, height=100)\n",
    "χe=0.02\n",
    "E_anharm=E_arm*(1-χe*(n+1/2))\n",
    "plt.barh(E_anharm, -3, height=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0,1 : 80\n",
      "1,2 : 160\n",
      "2,3 : 240\n",
      "3,4 : 320\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlwAAAFlCAYAAAA+gTZIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAdC0lEQVR4nO3df5BdZZ3n8feXTgBHdkMgMQmB2UaIWkGWgCEyMVuFyo+AlsguO0XW0gy4heyCKyXjgpS1tJGEgAIivxnNJixKaIOlWUSooFAspUIaEsHARLJIJDFIhkAgaTeS9rt/9MHpCel0p9NPn3s771fVrb73Oc+593uenE4+ec6PG5mJJEmSytmn7gIkSZKGOwOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFTai7gJ2ZcyYMdna2lp3GZIkSX164okn/ikzx+5sWUMHrtbWVjo6OuouQ5IkqU8Rsba3ZR5SlCRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhDf3l1fx+BbSNqrsKSZLUzNo2112BM1ySJEmlGbgkSZIKa+xDioccC20ddVchSZK0R5zhkiSpyX3+85/n3HPPrbsM7YKBS5KkJrZ161YWLlzI4sWLee211+ouR70wcEmS1MRuvvlmWltbmTRpEtdff33d5agXBi5JkprU1q1bueaaa/jwhz/M9OnTueGGG5zlalAGLkmSmtQ999zDiSeeyLhx4zjooIP4xCc+wV133VV3WdoJA5ckSU1q1qxZLFy48C+vb731Vk+eb1CNfVsISZLUq5EjRzJy5Mi/vN53331rrEa74gyXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpML6DFwRsX9EPB4Rv4qIVRHx1ap9YUT8NiJWVo8pVXtExLciYk1EPBURx/V4r9kR8Vz1mF1usyRJkhpHf75LcRvwkczcEhEjgUcj4ifVsi9l5pId+p8GTKoeHwRuAT4YEQcBlwNTgQSeiIilmfnqYGyIJElSo+pzhiu7balejqweuYtVzgDuqNb7JXBgREwATgWWZeamKmQtA2buWfmSJEmNr1/ncEVES0SsBF6mOzQ9Vi2aWx02vC4i9qvaJgIv9lh9XdXWW/uOn3VeRHRERMfGjRt3c3MkSZIaT78CV2Z2ZeYU4FBgWkS8H/gy8D7geOAg4JLBKCgzb8/MqZk5dezYsYPxlpIkSbXarasUM/M14CFgZmZuqA4bbgP+JzCt6rYeOKzHaodWbb21S5IkDWv9uUpxbEQcWD1/B3Ay8I/VeVlERACfBH5drbIU+Ex1teIJwObM3AA8AJwSEaMjYjRwStUmSZI0rPXnKsUJwKKIaKE7oLVn5r0R8bOIGAsEsBI4v+p/H3A6sAboBM4ByMxNEfE1YHnVb05mbhq8TZEkSWpMfQauzHwKOHYn7R/ppX8CF/SybAGwYDdrlCRJamreaV6SJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIK689tIerz+xXQNqruKiRJamyPbIM3E9purLuSxtS2ue4KnOGSJEkqrbFnuA45Fto66q5CkqTGNncudHZC29y6K1EvnOFqQn/+85/5zW9+U3cZkrRX8+9h7Q4DVxP6+te/zlFHHcUbb7xRdymStFdav34973vf+/jhD39YdylqEgauJtPV1cUtt9zCiBEjuOmmm+ouR5L2SldddRX77LMPV111Fd3faCftmoGrySxevJjRo0fT2trKtdde6yyXJA2x9evXc+eddzJ69GhefvllHn744bpLUhMwcDWRrq4u5syZw+c+9zn2339/TjrpJGe5JGmIXXXVVZx77rm0tLRw0UUX0dbW5iyX+mTgaiKvvfYaJ598MtOmTQOgra2Nzs7OmquSpL1LRPClL30JgDPPPJPW1la6urpqrkqNzsDVRA4++GBuvPFGIgKA97znPcyZM6fmqiRp73L99dczbtw4AEaMGMGiRYsYMaKx77Kk+hm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFjai7gF36/QpoG1V3FY1nQxe89EfHRpLqtOUN+MYkOKAB5i4e2QZvJrTdWHcljaltc90VOMMlSZJUWmPPcB1yLLR11F1F41mxAh4/F9pW1F2JJO29bh0Pf78Sxo+vuxKYOxc6O6Ftbt2VqBfOcGmPPf3002zfvr3uMiQNc52dnaxevbruMqQBMXBpj3R2dvKBD3yAyy67rO5SJA1zs2fPZtq0aXR1ddVdirTbDFzaI7fddhv77rsvixYtYtu2bXWXI2mY2rRpE/fffz8A3//+92uuRtp9Bi4NWGdnJ1dffTWTJk3i3e9+NwsWLKi7JEnD1HXXXce0adM45phjmDNnjrNcajoGLg3YbbfdxvTp03nnO9/JOeecw5VXXuksl6RBt2nTJm6++WbOPPNMxo0bx6hRo5zlUtMxcGnA/vSnP3HFFVcAMHnyZD71qU+xadOmmquSNNy89NJLXHLJJbzrXe8iIrjmmmt4/fXX6y5L2i2NfVsINbRLLrnkX7y+8sora6pE0nA2efJkJk+eTHt7OwDTp09n+vTpNVcl7R5nuCRJkgozcEmSJBVm4JIkSSrMwCVJklRYn4ErIvaPiMcj4lcRsSoivlq1Hx4Rj0XEmoi4OyL2rdr3q16vqZa39nivL1ftqyPi1FIbJUmS1Ej6M8O1DfhIZh4DTAFmRsQJwFXAdZl5JPAq8Nmq/2eBV6v266p+RMRk4GzgKGAmcHNEtAzmxkiSJDWiPgNXdttSvRxZPRL4CLCkal8EfLJ6fkb1mmr5RyMiqvbFmbktM38LrAGmDcpWSJIkNbB+3Yermol6AjgSuAn4v8Brmbm96rIOmFg9nwi8CJCZ2yNiM3Bw1f7LHm/bc52en3UecB7AX48KaBu1m5u0F9jQBS/9sXHG5ndbYcGp8KC3dZNU0Ko34Zk3G+fvvi1vwDcmwQENcDr0I9vgzYS2G+uupDG1ba67gv6dNJ+ZXZk5BTiU7lmp95UqKDNvz8ypmTl17F9FqY+RJEkaMrsVyzPzNeAh4G+AAyPirSmNQ4H11fP1wGEA1fJRwCs923eyjiRJ0rDV5zGgiBgLvJmZr0XEO4CT6T4R/iHgLGAxMBv4UbXK0ur1L6rlP8vMjIilwPci4lrgEGAS8PguP/yQY6GtYyDbNbytWAGPnwttK+qupNuDM+Dc+TBjRt2VSBrO2tthyRJoa6+7km63joe/Xwnjx9ddCcydC52d0Da37krUi/7McE0AHoqIp4DlwLLMvBe4BPhiRKyh+xyt71T9vwMcXLV/EbgUIDNXAe3AM8D9wAWZ2TWYGyPNmzePY489tu4ypGEhMzn88MP53ve+V3cpUtPrc4YrM58C3vYvWGY+z06uMszM/wf8x17eay5g/FYRXV1d3HjjjWzcuJEVK1YYvKQ9dP/99/O73/2Oq6++mlmzZtF9wbmkgWiASyukwbF48WJGjRrFUUcdxZw5c+ouR2pqmUlbWxtHHnkkr7zyCg8//HDdJUlNzcClYaGrq4s5c+Zw1llnccQRR/DYY4+xYkWDnOMmNaH777+frVu3MmbMGGbPnk1bWxuZWXdZUtMycGlYWL58OYcccgjvf//7aWlp4Stf+Qp33nln3WVJTevOO+9kzpw5RAQnnXQSr7zyCmvXrq27LKlpNfadKn+/onFucNdIvPHp25wAPPjvknuWnAPPvMl/nfwAXe9MaPt2bTVJzeyOI5KWp+7l2t9tZcQdH+NX/76FloXH1FuUNz7tnTc+3bVmufGp1Axa9oldvpbUf/4+SYOrsWe4vA/Xznkfrt412n16pGbn73fvvA+XdoMzXFJBv/jFLzzRWJJk4JJKWbt2LdOnT6e9vUH+Ny5Jqo2BSypk3rx5jBgxgquvvtpZLknayxm4pALWrl3LkiVLGD16NFu2bOGBBx6ouyRJUo0MXFIB8+bN4/zzz2efffbh4osv9qaRkrSXM3BJBYwZM4YvfvGLAHz84x/n+OOPp6vL72qXpL1VY98WQmpSc+f+86XZ++yzDzfccEON1UiS6uYMlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIKM3BJkiQVZuCSJEkqzMAlSZJUmIFLkiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKmwPgNXRBwWEQ9FxDMRsSoivlC1t0XE+ohYWT1O77HOlyNiTUSsjohTe7TPrNrWRMSlZTZJkiSpsYzoR5/twMWZ+WRE/CvgiYhYVi27LjO/0bNzREwGzgaOAg4BHoyI91SLbwJOBtYByyNiaWY+MxgbIkmS1Kj6DFyZuQHYUD1/IyKeBSbuYpUzgMWZuQ34bUSsAaZVy9Zk5vMAEbG46mvgkiRJw9puncMVEa3AscBjVdOFEfFURCyIiNFV20TgxR6rravaemvf8TPOi4iOiOjYuHHj7pQnSZLUkPoduCLiAOAe4KLMfB24BTgCmEL3DNg1g1FQZt6emVMzc+rYsWMH4y0lSZJq1Z9zuIiIkXSHre9m5g8AMvMPPZb/A3Bv9XI9cFiP1Q+t2thFuyRJ0rDVn6sUA/gO8GxmXtujfUKPbmcCv66eLwXOjoj9IuJwYBLwOLAcmBQRh0fEvnSfWL90cDZDkiSpcfVnhutDwKeBpyNiZdV2GTArIqYACbwAfA4gM1dFRDvdJ8NvBy7IzC6AiLgQeABoARZk5qpB3BZJkqSG1J+rFB8FYieL7tvFOnOBuTtpv29X60mSJA1H3mlekiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIKM3BJkiQVZuCSJEkqzMAlSZJUmIFLkiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIK6zNwRcRhEfFQRDwTEasi4gtV+0ERsSwinqt+jq7aIyK+FRFrIuKpiDiux3vNrvo/FxGzy22WJElS4+jPDNd24OLMnAycAFwQEZOBS4GfZuYk4KfVa4DTgEnV4zzgFugOaMDlwAeBacDlb4U0SZKk4azPwJWZGzLzyer5G8CzwETgDGBR1W0R8Mnq+RnAHdntl8CBETEBOBVYlpmbMvNVYBkwc1C3RpIkqQHt1jlcEdEKHAs8BozLzA3VopeAcdXzicCLPVZbV7X11i5JkjSs9TtwRcQBwD3ARZn5es9lmZlADkZBEXFeRHRERMfGjRsH4y0lSZJq1a/AFREj6Q5b383MH1TNf6gOFVL9fLlqXw8c1mP1Q6u23tr/hcy8PTOnZubUsWPH7s62SJIkNaT+XKUYwHeAZzPz2h6LlgJvXWk4G/hRj/bPVFcrngBsrg49PgCcEhGjq5PlT6naJEmShrUR/ejzIeDTwNMRsbJquwyYD7RHxGeBtcDfVsvuA04H1gCdwDkAmbkpIr4GLK/6zcnMTYOyFZIkSQ2sz8CVmY8C0cvij+6kfwIX9PJeC4AFu1OgJElSs/NO85IkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIKM3BJkiQVZuCSJEkqzMAlSZJUmIFLkiSpMAOXJElSYQauJjRhwgQ+/OEP112GJEnqpxF1F7ArT6/fTOulP667jMa070f5QYOMzUsvvMpZt/6c/e/dXHcpbH32STpXb2iY/Wbjlm0cf8WDtBwwuu5SpAHx97t3jfT7vfnnq/nz9m18t0HGptG8MP9jdZfgDJckSVJpDT3DdfTEUXQ0QCrVrs149Ermnz+dGTNm1F0K7e1bWbLkedobZL8Zv3A/ln/lJMaPH193KdKA+Pvdu0b6/Z47dyWdnZ3MndsYY6O3c4ZLkiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIKM3BJkiQVZuCSJEkqzMAlSZJUmIFLkiSpMAOXJElSYX0GrohYEBEvR8Sve7S1RcT6iFhZPU7vsezLEbEmIlZHxKk92mdWbWsi4tLB3xRJkqTG1J8ZroXAzJ20X5eZU6rHfQARMRk4GziqWufmiGiJiBbgJuA0YDIwq+orSZI07I3oq0NmPhIRrf18vzOAxZm5DfhtRKwBplXL1mTm8wARsbjq+8xuVyxJktRk9uQcrgsj4qnqkOPoqm0i8GKPPuuqtt7aJUmShr2BBq5bgCOAKcAG4JrBKigizouIjojo2Lhx42C9rSRJUm36PKS4M5n5h7eeR8Q/APdWL9cDh/XoemjVxi7ad3zv24HbAaZOnZoDqU9Da968eRx33HF1lyFJe61Zs2axffv2usvQLgwocEXEhMzcUL08E3jrCsalwPci4lrgEGAS8DgQwKSIOJzuoHU28J/6+pyn12+m9dIfD6REDbX7Hqq7AgC2Pvsknas3NMx+s3HLNo6/4kFaDhjdd2epAb30wqucdevP2f/ezXWX4u93vzxXdwEN6YX5H6u7hL4DV0TcBZwIjImIdcDlwIkRMQVI4AXgcwCZuSoi2uk+GX47cEFmdlXvcyHwANACLMjMVYO+NZIkSQ0oMhv3qN3UqVOzo6Oj7jLURNrb21myZAnt7e11lwLA+PHjWblyJePHj6+7FGlAZsyYwfz585kxY0bdpfj7rYYXEU9k5tSdLfNO85IkSYUZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTCDFySJEmFGbgkSZIKM3BJkiQVZuCSJEkqzMAlSZJUmIFLkiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqbARdRcgDaaTTz6Z9773vXWXIWkvcPfddzNu3Li6y1CTaOjA9fT6zbRe+uO6y1Azumtd3RUAsHHLNo6/4kFaDhhddynSgLz0wqucdevP2f/ezXWXwtZnn6Rz9YbG+nfhJ/fVXYH64YX5H6u7BA8pSpIkldbQM1xHTxxFRwOkUmmgxi/cj+VfOYnx48fXXYo0IDMevZL5509nxowZdZdCe/tWlix5nnb/XVATcoZLkiSpMAOXJElSYQYuSZKkwgxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEkFnXPOORx44IF1lyENC5MnT+aUU06puwxpQCIz666hV/tNmJQTZn+z7jIkaa/10p3/nQNPnM3+hx5VdynSgL0w/2ND8jkR8URmTt3ZMme4JEmSChtRdwG7cvTEUXQMUSqVJL3djEevZP7505kxY0bdpUhNzRkuSZKkwgxckiRJhfUZuCJiQUS8HBG/7tF2UEQsi4jnqp+jq/aIiG9FxJqIeCoijuuxzuyq/3MRMbvM5kiSJDWe/sxwLQRm7tB2KfDTzJwE/LR6DXAaMKl6nAfcAt0BDbgc+CAwDbj8rZAmSZI03PUZuDLzEWDTDs1nAIuq54uAT/ZovyO7/RI4MCImAKcCyzJzU2a+Cizj7SFOkiRpWBroVYrjMnND9fwlYFz1fCLwYo9+66q23tolSQ3sm9/8JkcffXTdZUhNb49vC5GZGRGDdvfUiDiP7sORtPzrsbRe+uPBemtJ0kAs+UPdFUh7ZKhufLorA71K8Q/VoUKqny9X7euBw3r0O7Rq6639bTLz9sycmplTW/5q1ADLkyRJahwDneFaCswG5lc/f9Sj/cKIWEz3CfKbM3NDRDwAzOtxovwpwJf7+hBvfCpJkoaDPgNXRNwFnAiMiYh1dF9tOB9oj4jPAmuBv6263wecDqwBOoFzADJzU0R8DVhe9ZuTmTueiC9JkjQsNfSXV0+dOjU7OjrqLkOSJKlPfnm1JElSjQxckiRJhRm4JEmSCjNwSZIkFWbgkiRJKszAJUmSVJiBS5IkqTADlyRJUmEGLkmSpMIMXJIkSYUZuCRJkgozcEmSJBXW0F9eHRFvAKvrrmOYGwP8U91FDHOOcVmOb3mOcVmOb3lDNcb/JjPH7mzBiCH48D2xurdv3dbgiIgOx7gsx7gsx7c8x7gsx7e8RhhjDylKkiQVZuCSJEkqrNED1+11F7AXcIzLc4zLcnzLc4zLcnzLq32MG/qkeUmSpOGg0We4JEmSml5DBK6ImBkRqyNiTURcupPl+0XE3dXyxyKideirbG79GOO/i4iNEbGyevznOupsVhGxICJejohf97I8IuJb1fg/FRHHDXWNzawf43tiRGzusf/+j6GusdlFxGER8VBEPBMRqyLiCzvp4348QP0cX/fjPRAR+0fE4xHxq2qMv7qTPrXlidoDV0S0ADcBpwGTgVkRMXmHbp8FXs3MI4HrgKuGtsrm1s8xBrg7M6dUj28PaZHNbyEwcxfLTwMmVY/zgFuGoKbhZCG7Hl+A/9Nj/50zBDUNN9uBizNzMnACcMFO/p5wPx64/owvuB/viW3ARzLzGGAKMDMiTtihT215ovbABUwD1mTm85n5J2AxcMYOfc4AFlXPlwAfjYgYwhqbXX/GWHsgMx8BNu2iyxnAHdntl8CBETFhaKprfv0YX+2hzNyQmU9Wz98AngUm7tDN/XiA+jm+2gPVfrmlejmyeux4onpteaIRAtdE4MUer9fx9p3wL30yczuwGTh4SKobHvozxgD/oTpMsCQiDhua0vYa/f0z0MD9TXUo4ScRcVTdxTSz6jDLscBjOyxyPx4EuxhfcD/eIxHREhErgZeBZZnZ6z481HmiEQKXGsP/Bloz898Cy/jn/wFIzeBJur9S4xjgBuCHNdfTtCLiAOAe4KLMfL3ueoabPsbX/XgPZWZXZk4BDgWmRcT7667pLY0QuNYDPWdTDq3adtonIkYAo4BXhqS64aHPMc7MVzJzW/Xy28AHhqi2vUV/9nMNUGa+/tahhMy8DxgZEWNqLqvpRMRIusPAdzPzBzvp4n68B/oaX/fjwZOZrwEP8fZzP2vLE40QuJYDkyLi8IjYFzgbWLpDn6XA7Or5WcDP0huI7Y4+x3iH8zA+Qff5BRo8S4HPVFd5nQBszswNdRc1XETE+LfOw4iIaXT/3eZ/ynZDNX7fAZ7NzGt76eZ+PED9GV/34z0TEWMj4sDq+TuAk4F/3KFbbXmi9i+vzsztEXEh8ADQAizIzFURMQfoyMyldO+k/ysi1tB94uzZ9VXcfPo5xv8tIj5B95U0m4C/q63gJhQRdwEnAmMiYh1wOd0nbJKZtwL3AacDa4BO4Jx6Km1O/Rjfs4D/EhHbgT8CZ/ufst32IeDTwNPVOTAAlwF/De7Hg6A/4+t+vGcmAIuqK/P3Adoz895GyRPeaV6SJKmwRjikKEmSNKwZuCRJkgozcEmSJBVm4JIkSSrMwCVJklSYgUuSJKkwA5ckSVJhBi5JkqTC/j+ItBpOSMGfZwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "B=40\n",
    "lev=5\n",
    "J=np.arange(0,lev,1)\n",
    "Ej=B*J*(J+1)\n",
    "#plt.plot(x,Va)\n",
    "for i in range(2):\n",
    "    plt.barh(E_anharm[i]+Ej,3,height=20)\n",
    "for i in range(lev-1):\n",
    "    print(\"{},{} : {}\".format( i, i+1, Ej[i+1]-Ej[i]))\n",
    "    \n",
    "for i in range(lev-1):\n",
    "    da=E_anharm[0]+Ej[i]\n",
    "    a=E_anharm[1]+Ej[i+1]\n",
    "    plt.annotate(\"\",(2.0+i*0.2, a),(2.0+i*0.2, da), arrowprops={'arrowstyle':'->'})\n",
    "    da=E_anharm[0]+Ej[i+1]\n",
    "    a=E_anharm[1]+Ej[i]\n",
    "    plt.annotate(\"\",(1.0-i*0.2, a),(1.0-i*0.2, da), arrowprops={'arrowstyle':'->'})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<BarContainer object of 14 artists>"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAAFlCAYAAADPim3FAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAV80lEQVR4nO3dfYxlZ30f8O8PL07bhCwvO4moXzquaqqYJhSYGtSkxREkWduV3TSQ2k0DxJBVX0zJi6psRGUsUCWD2yiNZOKuqOUQURsClG66poa2EEsNpl43jvFLDItx4nUoXl7qVo2CY/HrH3NsxuOdmWueu56Znc9HurrPec6z9/7m0ZnRd5977jnV3QEA4NvzrM0uAABgOxOmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYMCuzXrjPXv29OLi4ma9PQDAzG6//favdPfC8fZtWphaXFzM4cOHN+vtAQBmVlV/uNY+H/MBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKdqIrdy8/2BYW9x/K4v5Dm10GsAZhCgBggDAFADBAmAIAGCBMAQAM2DBMVdV1VfVwVd21zpjzquqOqrq7qn5nviUCAGxds6xMXZ9k71o7q+q5Sd6T5KLufnGS182nNACArW/DMNXdtyT52jpD/kGSj3T3H03jH55TbQAAW948zpl6UZLnVdWnqur2qnr9HF4TAGBb2DWn13h5klcn+fNJPl1Vt3b351YPrKp9SfYlyZlnnjmHtwYA2FzzWJk6muTm7v5/3f2VJLckecnxBnb3ge5e6u6lhYWFObw1AMDmmkeY+o9JfqiqdlXVX0jyiiT3zuF1AQC2vA0/5quqG5Kcl2RPVR1N8vYkz06S7r62u++tqv+c5M4k30zy3u5e8zIKAAAnkw3DVHdfOsOYq5NcPZeKAAC2EVdABwAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwB374rdy8/OK7F/YeyuP/QZpcBnGDCFADAAGEKthMrQQywUgYnhjAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBgwzBVVddV1cNVddcG4/5GVT1WVa+dX3kAAFvbLCtT1yfZu96AqjolybuSfHwONQEAbBsbhqnuviXJ1zYY9pYkH07y8DyKAgDYLobPmaqq05L8eJJfHy8HAGB7mccJ6L+a5Je6+5sbDayqfVV1uKoOHzt2bA5vDQCwuXbN4TWWktxYVUmyJ8kFVfVYd3909cDuPpDkQJIsLS31HN4bAGBTDYep7j7r8XZVXZ/kPx0vSAEAnIw2DFNVdUOS85LsqaqjSd6e5NlJ0t3XntDqAAC2uA3DVHdfOuuLdfcbh6oBANhmXAEdAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABmwYpqrquqp6uKruWmP/T1XVnVX12ar63ap6yfzLBADYmmZZmbo+yd519n8xyau6+/uTvDPJgTnUBQCwLezaaEB331JVi+vs/90Vm7cmOX28LACA7WHe50y9KcnH1tpZVfuq6nBVHT527Nic3xoA4Jk3tzBVVT+c5TD1S2uN6e4D3b3U3UsLCwvzemsAgE2z4cd8s6iqH0jy3iTnd/dX5/GaAADbwfDKVFWdmeQjSX66uz83XhIAwPax4cpUVd2Q5Lwke6rqaJK3J3l2knT3tUmuSPKCJO+pqiR5rLuXTlTBAABbySzf5rt0g/1vTvLmuVUEALCNuAI6AMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYMCGt5MB5uDK3dPzI5tbB8zB4v5DT7QfuOrCTawEtgYrUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABG4apqrquqh6uqrvW2F9V9WtVdaSq7qyql82/TACArWmWlanrk+xdZ//5Sc6eHvuS/Pp4WQAA28OGYaq7b0nytXWGXJzkfb3s1iTPraoXzqtAAICtbB7nTJ2W5MEV20envqeoqn1VdbiqDh87dmwObw0AsLme0RPQu/tAdy9199LCwsIz+dYAACfEPMLUQ0nOWLF9+tQHAHDSm0eYOpjk9dO3+l6Z5JHu/tIcXhcAYMvbtdGAqrohyXlJ9lTV0SRvT/LsJOnua5PclOSCJEeS/EmSnzlRxQIAbDUbhqnuvnSD/Z3kn86tIgCAbcQV0AEABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEK2HxX7l5+bBGL+w9lcf+hzS4D2CaEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGzBSmqmpvVd1XVUeqav9x9p9ZVZ+sqt+rqjur6oL5lwoAsPVsGKaq6pQk1yQ5P8k5SS6tqnNWDfsXST7Y3S9NckmS98y7UACArWiWlalzkxzp7vu7+9EkNya5eNWYTvLdU3t3kj+eX4kAAFvXrhnGnJbkwRXbR5O8YtWYK5N8vKrekuQ7k7xmLtUBAGxx8zoB/dIk13f36UkuSPKbVfWU166qfVV1uKoOHzt2bE5vDQCweWYJUw8lOWPF9ulT30pvSvLBJOnuTyf5c0n2rH6h7j7Q3UvdvbSwsPDtVQwAsIXMEqZuS3J2VZ1VVadm+QTzg6vG/FGSVydJVX1flsOUpScA4KS3YZjq7seSXJ7k5iT3Zvlbe3dX1Tuq6qJp2C8m+dmq+v0kNyR5Y3f3iSoaAGCrmOUE9HT3TUluWtV3xYr2PUl+cL6lAQBsfa6ADgAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMCAma4zBTvWlbun50c2tw44CS3uP/RE+4GrLtzESmCMlSkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYMFOYqqq9VXVfVR2pqv1rjPnJqrqnqu6uqn8/3zIBALamXRsNqKpTklyT5EeSHE1yW1Ud7O57Vow5O8kvJ/nB7v56VX3PiSoYAGArmWVl6twkR7r7/u5+NMmNSS5eNeZnk1zT3V9Pku5+eL5lAgBsTbOEqdOSPLhi++jUt9KLkryoqv57Vd1aVXvnVSAAwFa24cd8T+N1zk5yXpLTk9xSVd/f3f975aCq2pdkX5KceeaZc3prAIDNM8vK1ENJzlixffrUt9LRJAe7+8+6+4tJPpflcPUk3X2gu5e6e2lhYeHbrRkAYMuYJUzdluTsqjqrqk5NckmSg6vGfDTLq1Kpqj1Z/tjv/jnWCQCwJW0Yprr7sSSXJ7k5yb1JPtjdd1fVO6rqomnYzUm+WlX3JPlkkn/e3V89UUUDAGwVM50z1d03JblpVd8VK9qd5BemBwDAjuEK6AAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwIBdm10AnBBX7p6eH9ncOoATbnH/oSfaD1x14SZWwk5lZQoAYMBMYaqq9lbVfVV1pKr2rzPuJ6qqq2ppfiUCAGxdG4apqjolyTVJzk9yTpJLq+qc44x7TpK3JvnMvIsEANiqZlmZOjfJke6+v7sfTXJjkouPM+6dSd6V5E/nWB8AwJY2S5g6LcmDK7aPTn1PqKqXJTmjuw9lHVW1r6oOV9XhY8eOPe1iAQC2muET0KvqWUl+JckvbjS2uw9091J3Ly0sLIy+NQDAppslTD2U5IwV26dPfY97TpK/luRTVfVAklcmOegkdABgJ5glTN2W5OyqOquqTk1ySZKDj+/s7ke6e093L3b3YpJbk1zU3YdPSMUAAFvIhmGqux9LcnmSm5Pcm+SD3X13Vb2jqi460QUCAGxlM10BvbtvSnLTqr4r1hh73nhZAADbgyugAwAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwIBdswyqqr1J/k2SU5K8t7uvWrX/F5K8OcljSY4luay7/3DOtbITXLl7en5kc+sAdpzF/YeeaD9w1YWbWAnbzYYrU1V1SpJrkpyf5Jwkl1bVOauG/V6Spe7+gSQfSvLueRcKALAVzfIx37lJjnT3/d39aJIbk1y8ckB3f7K7/2TavDXJ6fMtEwBga5olTJ2W5MEV20envrW8KcnHjrejqvZV1eGqOnzs2LHZqwQA2KLmegJ6Vf3DJEtJrj7e/u4+0N1L3b20sLAwz7cGANgUs5yA/lCSM1Zsnz71PUlVvSbJ25K8qru/MZ/yAAC2tllWpm5LcnZVnVVVpya5JMnBlQOq6qVJ/m2Si7r74fmXCQCwNW0Yprr7sSSXJ7k5yb1JPtjdd1fVO6rqomnY1Um+K8lvVdUdVXVwjZcDADipzHSdqe6+KclNq/quWNF+zZzrAgDYFlwBHQBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAbs2uwCOElcuXt6fmRz6wDYZIv7Dz3RfuCqCzexEp4pVqYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA9xOZqdzGxiALcltabYPK1MAAAOEKQCAAcIUAMCAmcJUVe2tqvuq6khV7T/O/u+oqg9M+z9TVYvzLhQAYCva8AT0qjolyTVJfiTJ0SS3VdXB7r5nxbA3Jfl6d/+VqrokybuS/P0TUTATJ44DMAMnsp94s6xMnZvkSHff392PJrkxycWrxlyc5Dem9oeSvLqqan5lAgBsTbOEqdOSPLhi++jUd9wx3f1YkkeSvGAeBZ50rtz9rVUlANhGFvcfetJKF8uqu9cfUPXaJHu7+83T9k8neUV3X75izF3TmKPT9hemMV9Z9Vr7kuybNv9qkvvm9YNsIXuSfGXDUSTm6ukwV7MzV7MzV0+P+ZrdyThXf6m7F463Y5aLdj6U5IwV26dPfccbc7SqdiXZneSrq1+ouw8kOTBLxdtVVR3u7qXNrmM7MFezM1ezM1ezM1dPj/ma3U6bq1k+5rstydlVdVZVnZrkkiQHV405mOQNU/u1Sf5bb7TkBQBwEthwZaq7H6uqy5PcnOSUJNd1991V9Y4kh7v7YJJ/l+Q3q+pIkq9lOXABAJz0Zro3X3fflOSmVX1XrGj/aZLXzbe0beuk/hhzzszV7MzV7MzV7MzV02O+Zrej5mrDE9ABAFib28kAAAwQpjZQVddV1cPT5R8e7/tAVd0xPR6oqjtW7Pvl6bY691XVj63oX/eWPCeDNebqr1fVrdNcHa6qc6f+qqpfm+bjzqp62Yp/84aq+vz0eMPx3utksMZ8vaSqPl1Vn62q366q716xb0ceW1V1RlV9sqruqaq7q+qtU//zq+oT03Hyiap63tS/Y4+tdebqddP2N6tqadW/2ZHHVbLufF1dVX8wHT//oaqeu+Lf7Mj5Wmeu3jnN0x1V9fGq+otT/876Pexuj3UeSf52kpcluWuN/f86yRVT+5wkv5/kO5KcleQLWT5p/5Sp/ZeTnDqNOWezf7ZnYq6SfDzJ+VP7giSfWtH+WJJK8sokn5n6n5/k/un5eVP7eZv9sz2D83VbkldN7cuSvHOnH1tJXpjkZVP7OUk+N83Hu5Psn/r3J3nXTj+21pmr78vytf0+lWRpxfgde1xtMF8/mmTX1P+uFcfWjp2vdebqu1eM+WdJrp3aO+r30MrUBrr7lix/Q/EpqqqS/GSSG6aui5Pc2N3f6O4vJjmS5dvxzHJLnm1vjbnqJI+vruxO8sdT++Ik7+tltyZ5blW9MMmPJflEd3+tu7+e5BNJ9p746p95a8zXi5LcMrU/keQnpvaOPba6+0vd/T+n9v9Ncm+W77qw8jZWv5Hk707tHXtsrTVX3X1vdx/vIsk79rhK1p2vj/fy3TyS5NYsX18x2cHztc5c/Z8Vw74zy3/zkx32ezjTt/lY099K8uXu/vy0fVqWf/Eet/LWO6tvyfOKE1/elvBzSW6uqn+V5Y+V/+bUv9Ztima5fdHJ7O4s/xH6aJa/Ifv4BXMdW0mqajHJS5N8Jsn3dveXpl3/K8n3Tm3HVp4yV2txXE3Wma/LknxgapuvPHWuqupfJnl9lm8l98PTsB31e2hlasyl+daqFMf3j5P8fHefkeTns3xNMtZ2WZJ/UlW3Z3kp/dFNrmfLqKrvSvLhJD+36n/D6eXPD3w1ebLeXPFUa81XVb0tyWNJ3r9ZtW01x5ur7n7b9Df+/UkuX+/fn6yEqW9TLd825+/lW/9jSda+9c4st+Q5Wb0hyUem9m9leTk8MVfH1d1/0N0/2t0vz3JQ/8K0a0fPV1U9O8t/wN/f3Y8fT1+ePjbI9Pzw1G+unjpXa9nRc5WsPV9V9cYkfyfJT01hPdnh8zXDsfX+fOvUhJ01V5t90tZ2eCRZzKoT0LP8Ge/vrOp7cZ58cuL9WT4xcdfUPivfOjnxxZv9cz0Tc5Xlz9XPm9qvTnL71L4wTz458X9M/c9P8sUsn5j4vKn9/M3+uZ7B+fqe6flZSd6X5LKdfmxNx8j7kvzqqv6r8+QT0N+904+tteZqxf5P5cknoO/Y42qDY2tvknuSLKzq37Hztc5cnb2i/ZYkH5raO+r3cNML2OqPLK8OfCnJn2X5s903Tf3XJ/lHxxn/tiyvJtyX6VtsU/8FWf72wxeSvG2zf65naq6S/FCS26c/Lp9J8vJpbCW5ZpqPz676A39Zlk/sPJLkZzb753qG5+ut03HyuSRXZbqw7k4+tqZjqJPcmeSO6XFBkhck+a9JPp/kvzz+B3knH1vrzNWPT8fYN5J8OcnNO/242mC+jmT5vJ7H+67d6fO1zlx9OMldU/9vZ/mk9B33e+gK6AAAA5wzBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYMD/B3gfPzT/TZyJAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ν0=2000\n",
    "B=12\n",
    "J=np.arange(0,15,1)\n",
    "Ej=B*J*(J+1)\n",
    "rR=ν0+Ej[1:]-Ej[:-1]\n",
    "rP=ν0+Ej[:-1]-Ej[1:]\n",
    "Ir=J*np.exp(-J*(J+1)/20)\n",
    "plt.bar(rR,Ir[1:],width=3)\n",
    "plt.bar(rP,Ir[1:],width=3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\omega=\\sqrt {\\frac k \\mu}$\n",
    "\n",
    "$B=\\frac{\\hbar^2}{2\\mu}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<BarContainer object of 14 artists>"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAAFlCAYAAADPim3FAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAWl0lEQVR4nO3dfYxld3kf8O+DF9M2IbztJqJ+ybqqXcU0ocDEQQopTqHJ2q3spgkRLm0gkKyS1pSkURQjKlMZVQq4L1EkE3eVWg4pteMATTetqaEpYKnBjteJY/wSw2LSeB2Kl5fSlyg4Fk//uGfx7LAzc+3f3Z07O5+PNNpzf+c35z777Jk73/3dc++t7g4AAE/PM7a6AACA7UyYAgAYIEwBAAwQpgAABghTAAADhCkAgAG7tuqOd+/e3Xv37t2quwcAmNvdd9/9+e7ec6J9Wxam9u7dm0OHDm3V3QMAzK2q/sd6+zzNBwAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKWAuVZWq2uoyWId/H9g6whQAwABhCgBggDAFADBAmAIAGLBpmKqqG6rqsaq6b4M5F1fVPVV1f1V9bLElAgAsr3lWpm5Msm+9nVX13CTvTnJZd78oyWsWUxoAwPLbNEx19+1JvrjBlL+X5APd/UfT/McWVBsAwNJbxDVTFyR5XlV9tKrurqofWcAxAQC2hV0LOsbLkrwqyZ9P8vGquqO7P7l2YlXtT7I/Sc4999wF3DUAwNZaxMrUkSS3dff/6+7PJ7k9yYtPNLG7D3T3Snev7NmzZwF3DQCwtRYRpv5jkldU1a6q+gtJvivJgws4LgDA0tv0ab6quinJxUl2V9WRJG9P8swk6e7ru/vBqvovSe5N8tUkv9zd676NAgDA6WTTMNXdV8wx59ok1y6kIgCAbcQ7oAMADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADNi11QUAO0dVJUm6e4srWR7HepLoC2xXVqYAAAYIU7ADVNVxKyBwIs4TeHqEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAZuGqaq6oaoeq6r7Npn3nVX1RFX90OLKAwBYbvOsTN2YZN9GE6rqjCTvTPKhBdQEALBtbBqmuvv2JF/cZNqbk7w/yWOLKAoAYLsYvmaqqs5K8gNJfmm8HACA7WURF6D/QpKf6+6vbjaxqvZX1aGqOnT06NEF3DUAwNbatYBjrCS5uaqSZHeSS6vqie7+jbUTu/tAkgNJsrKy0gu4bwCALTUcprr7vGPbVXVjkv90oiAFAHA62jRMVdVNSS5OsruqjiR5e5JnJkl3X39SqwMAWHKbhqnuvmLeg3X3G4aqAQDYZrwDOgDAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAzYNExV1Q1V9VhV3bfO/tdV1b1V9Ymq+u2qevHiywQAWE7zrEzdmGTfBvs/k+SV3f3tSd6R5MAC6gIA2BZ2bTahu2+vqr0b7P/tVTfvSHL2eFkAANvDoq+ZelOSD663s6r2V9Whqjp09OjRBd81AMCpt7AwVVXfm1mY+rn15nT3ge5e6e6VPXv2LOquAQC2zKZP882jqr4jyS8nuaS7v7CIYwIAbAfDK1NVdW6SDyT5B939yfGSAAC2j01XpqrqpiQXJ9ldVUeSvD3JM5Oku69PcnWSFyR5d1UlyRPdvXKyCgYAWCbzvJrvik32/1iSH1tYRQAA24h3QAcAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAM2PTjZICtM33eZZKku7ewEpjPsXPW+cpOYmUKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMGDTMFVVN1TVY1V13zr7q6p+saoOV9W9VfXSxZcJALCc5lmZujHJvg32X5Lk/Olrf5JfGi8LAGB72DRMdfftSb64wZTLk7ynZ+5I8tyqeuGiCgQAWGaLuGbqrCSPrLp9ZBr7OlW1v6oOVdWho0ePLuCuAQC21im9AL27D3T3Snev7Nmz51TeNQDASbGIMPVoknNW3T57GgMAOO0tIkwdTPIj06v6Xp7ky9392QUcFwBg6e3abEJV3ZTk4iS7q+pIkrcneWaSdPf1SW5NcmmSw0n+JMmPnqxiAQCWzaZhqruv2GR/J/lHC6sIAGAb8Q7oAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAVsO1WVqtrqMpamDmBrCVMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADJgrTFXVvqp6qKoOV9VVJ9h/blV9pKp+r6rurapLF18qAMDy2TRMVdUZSa5LckmSC5NcUVUXrpn2T5Pc0t0vSfLaJO9edKEAAMtonpWpi5Ic7u6Hu/vxJDcnuXzNnE7yTdP2c5L88eJKBABYXrvmmHNWkkdW3T6S5LvWzPlnST5UVW9O8g1JXr2Q6gAAltyiLkC/IsmN3X12kkuT/GpVfd2xq2p/VR2qqkNHjx5d0F0DAGydecLUo0nOWXX77GlstTcluSVJuvvjSf5ckt1rD9TdB7p7pbtX9uzZ8/QqBgBYIvOEqbuSnF9V51XVmZldYH5wzZw/SvKqJKmqb8ssTFl6AgBOe5uGqe5+IsmVSW5L8mBmr9q7v6quqarLpmk/k+THq+r3k9yU5A3d3SeraACAZTHPBejp7luT3Lpm7OpV2w8k+e7FlgYAsPy8AzoAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA+Z6nyngqauqJIn3r4Wn5tjPTuLnh+3ByhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMmCtMVdW+qnqoqg5X1VXrzPnhqnqgqu6vqn+/2DIBAJbTrs0mVNUZSa5L8jeTHElyV1Ud7O4HVs05P8lbk3x3d3+pqr75ZBUMALBM5lmZuijJ4e5+uLsfT3JzksvXzPnxJNd195eSpLsfW2yZAADLaZ4wdVaSR1bdPjKNrXZBkguq6r9X1R1VtW9RBQIALLNNn+Z7Csc5P8nFSc5OcntVfXt3/6/Vk6pqf5L9SXLuuecu6K4BALbOPCtTjyY5Z9Xts6ex1Y4kOdjdf9bdn0nyyczC1XG6+0B3r3T3yp49e55uzQAAS2OeMHVXkvOr6ryqOjPJa5McXDPnNzJblUpV7c7sab+HF1gnAMBS2jRMdfcTSa5McluSB5Pc0t33V9U1VXXZNO22JF+oqgeSfCTJz3b3F05W0QAAy6K6e0vueGVlpQ8dOrQl9w2nQlUlSUZ+xo4dY1HHOR2OscjjjPJvfHKPMXocWKSquru7V060zzugAwAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpAIABwhQAwABhCgBggDAFADBAmAIAGCBMAQAMEKYAAAYIUwAAA4QpOIGqSlVtdRnAAD/HnCrCFADAgLnCVFXtq6qHqupwVV21wbwfrKquqpXFlQgAsLw2DVNVdUaS65JckuTCJFdU1YUnmPfsJG9JcueiiwQAWFbzrExdlORwdz/c3Y8nuTnJ5SeY944k70zypwusDwBgqc0Tps5K8siq20emsa+pqpcmOae7//NGB6qq/VV1qKoOHT169CkXCwCwbIYvQK+qZyT5V0l+ZrO53X2gu1e6e2XPnj2jdw0AsOXmCVOPJjln1e2zp7Fjnp3kryb5aFX9YZKXJznoInQAYCeYJ0zdleT8qjqvqs5M8tokB4/t7O4vd/fu7t7b3XuT3JHksu4+dFIqBgBYIpuGqe5+IsmVSW5L8mCSW7r7/qq6pqouO9kFAgAss13zTOruW5Pcumbs6nXmXjxeFgDA9uAd0AEABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABghTAAADhCkAgAHCFADAAGEKAGDAXGGqqvZV1UNVdbiqrjrB/n9SVQ9U1b1V9VtV9a2LLxXmU1Wpqq0uAzgNeDxhHpuGqao6I8l1SS5JcmGSK6rqwjXTfi/JSnd/R5L3JXnXogsFAFhG86xMXZTkcHc/3N2PJ7k5yeWrJ3T3R7r7T6abdyQ5e7FlAgAsp3nC1FlJHll1+8g0tp43JfngiXZU1f6qOlRVh44ePTp/lQAAS2qhF6BX1d9PspLk2hPt7+4D3b3S3St79uxZ5F0DAGyJXXPMeTTJOatunz2NHaeqXp3kbUle2d1fWUx5AADLbZ6VqbuSnF9V51XVmUlem+Tg6glV9ZIk/ybJZd392OLLBABYTpuGqe5+IsmVSW5L8mCSW7r7/qq6pqoum6Zdm+Qbk/x6Vd1TVQfXORwAwGllnqf50t23Jrl1zdjVq7ZfveC6AAC2Be+ADgAwQJgCABggTAEADBCmAAAGCFMAAAOEKQCAAcIUAMAAYQoAYIAwBQAwQJgCABggTAEADBCmAAAGCFMAAAOEKZZGVaWqtrqM5Vc1+2LL9fTl34ONeGw7/QlTAAADhCk41ZZhZWkZalgGS9CHr61uAduWMAUAMECYAgAYIEwBAAzYtdUFwLZy7PqadpULS2I6JzuJq+Bga1iZAgAYIEwBAAwQpoCnbkneUmCrawBIXDPFTuOaJ/h6QikMsTLFQvi4BICTy+Ps8hKmAAAGCFMAAANcM8X24XonWE6eemKHm2tlqqr2VdVDVXW4qq46wf5nVdWvTfvvrKq9iy4UAGAZbRqmquqMJNcluSTJhUmuqKoL10x7U5IvdfdfTvKvk7xz0YVy8pyyixqX4OX0wJLy+HBKuIj95JhnZeqiJIe7++HufjzJzUkuXzPn8iS/Mm2/L8mryr8WALADzBOmzkryyKrbR6axE87p7ieSfDnJCxZRIOs79j+MuXLr6P/6/K8RWGbHHqM8zp10Vre+3im9AL2q9ifZP938v1X1UJLdST5/KuvYZubqz9wn9gbz5g5lo3WMHuP4/ev2Z8PjzNmvdY/xFB5IFnGMEx5nvu8/rj9P8xjr1/E0H1AX8UD8tSOMHWvWn0X8YljE32n5jvG0H5+fdh1rvu9pHWf0GCeYv84xnlJ/FhVAlvA8Wc/p9vv9W9fbMU+YejTJOatunz2NnWjOkaraleQ5Sb6w9kDdfSDJgdVjVXWou1fmqGNH0p+N6c/G9Gdj+rMx/dmY/mxsJ/Vnnqf57kpyflWdV1VnJnltkoNr5hxM8vpp+4eS/Ldur18HAE5/m65MdfcTVXVlktuSnJHkhu6+v6quSXKouw8m+bdJfrWqDif5YmaBCwDgtDfXNVPdfWuSW9eMXb1q+0+TvOZp1nBg8yk7mv5sTH82pj8b05+N6c/G9GdjO6Y/5dk4AICnz2fzAQAMWHiYqqobquqxqrpvzfibq+oPqur+qnrXqvG3Th9D81BVff+q8Q0/wma7OlF/po/iuWf6+sOqumfVPv2p+mtVdcfUn0NVddE0XlX1i1MP7q2ql676ntdX1aemr9ef6L62o3X68+Kq+nhVfaKqfrOqvmnVvp12/pxTVR+pqgemx5q3TOPPr6oPT+fDh6vqedP4jjqHNujPa6bbX62qlTXfs2POoQ36c+30++veqvoPVfXcVd+jP1XvmHpzT1V9qKr+4jS+c36+unuhX0n+epKXJrlv1dj3JvmvSZ413f7m6c8Lk/x+kmclOS/JpzO7yP2MafsvJTlzmnPhomvdiq8T9WfN/n+Z5Gr9Oe78+VCSS6btS5N8dNX2BzN7y6GXJ7lzGn9+koenP583bT9vq/9uJ7E/dyV55bT9xiTv2MHnzwuTvHTafnaST059eFeSq6bxq5K8cyeeQxv059uS/JUkH02ysmr+jjqHNujP9yXZNY2/c9X5oz+zHnzTqjn/OMn10/aO+fla+MpUd9+e2Sv6VvvJJD/f3V+Z5jw2jV+e5Obu/kp3fybJ4cw+vmaej7DZltbpT5JZik/yw0lumob0ZxpOcmy15TlJ/njavjzJe3rmjiTPraoXJvn+JB/u7i9295eSfDjJvpNf/cm3Tn8uSHL7tP3hJD84be/E8+ez3f270/b/SfJgZp/QsPojr34lyd+ZtnfUObRef7r7we5+6ATfsqPOoQ3686GefbpHktyR2fstJvpzrD//e9W0b8jsMTvZQT9fp+qaqQuSfE9V3VlVH6uq75zG1/uomnk+wuZ09D1JPtfdn5pu68/MTyW5tqoeSfIvkrx1Gtefmfvz5AP1a/Lkm+zu6P5U1d4kL0lyZ5Jv6e7PTrv+Z5JvmbZ3bI/W9Gc9+vP1/XljZqstif58rT9V9c+nx+jXJTn2av8d059TFaZ2Zbac9/IkP5vklmkVhuNdkSdXpXjSTyb56e4+J8lPZ/a+ZjzpjUn+YVXdndnS++NbXM+Wq6pvTPL+JD+15n/N6dnzDDv6Zcwb9Yf1+1NVb0vyRJL3blVty+BE/enut02P0e9NcuVW1rcVTlWYOpLkA9NS3+8k+Wpmn9mz3kfVzPMRNqeVmn0Mz99N8murhvVn5vVJPjBt/3pmS+iJ/iRJuvsPuvv7uvtlmYXxT0+7dmR/quqZmT3Qv7e7j503n5ueXsj057FLDXZcj9bpz3r058nxNyT520leNwXyRH9OdP68N09earBz+nMyLsRKsjfHXyD7E0mumbYvyGx5r5K8KMdfvPdwZhfu7Zq2z8uTF++96GTUuhVfa/szje1L8rE1Y/ozu/1gkoun7VcluXva/ls5/uLG35nGn5/kM5ld2Pi8afv5W/33Oon9OfaCjmckeU+SN+7U82c6F96T5BfWjF+b4y9Af9dOPIfW68+q/R/N8Reg76hzaIPzZ1+SB5LsWTOuP7Px81dtvznJ+6btHfPzdTKafVOSzyb5s8xWpN40nUz/Lsl9SX43yd9YNf9tmf1P+qFMr9iaxi/N7JUCn07ytq1u1MnszzR+Y5KfOMH8Hd+fJK9Icvf0gHRnkpdNcyvJdVMPPrHml8AbM7sY9HCSH93qv9dJ7s9bpnPhk0l+PtOb8e7Q8+cVmT2Fd2+Se6avS5O8IMlvJflUZq8sfv5OPIc26M8PTOfTV5J8LsltO/Ec2qA/hzNbBDg2dr3+HNef92f2+/3eJL+Z2UXpO+rnyzugAwAM8A7oAAADhCkAgAHCFADAAGEKAGCAMAUAMECYAgAYIEwBAAwQpgAABvx/2D7ghsEJEcEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "M12=12*16/(12+16)\n",
    "M13=13*16/(13+16)\n",
    "ν0_13=ν0*np.sqrt(M12/M13)\n",
    "B_13=B*(M12/M13)\n",
    "\n",
    "Ej_13=B_13*J*(J+1)\n",
    "rR_13=ν0_13+Ej_13[1:]-Ej_13[:-1]\n",
    "rP_13=ν0_13+Ej_13[:-1]-Ej_13[1:]\n",
    "Ir_13=J*np.exp(-J*(J+1)/20)*0.06\n",
    "plt.bar(rR,Ir[1:],width=3, color='black')\n",
    "plt.bar(rP,Ir[1:],width=3, color='black')\n",
    "plt.bar(rR_13,Ir_13[1:],width=3, color='red')\n",
    "plt.bar(rP_13,Ir_13[1:],width=3, color='red')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1800.0, 2200.0)"
      ]
     },
     "execution_count": 74,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmAAAAFlCAYAAABMTlT+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAXDElEQVR4nO3dfYxld3kf8O+D16YvIbztJkV+yTqqqeI0ocDUoU1a3ELJ2mrtRgmV3bRAIFmprSPyojZGRKaCf0pcpVVSE3eVWA5RYpcklGxTU0NTiKUGU4+DMX6Jyca08RoaLy+FtjQ4Fk//mLv4ejyzc8d79zcz934+0mjPOfc35zz32XvPfvfcc8+p7g4AAOM8a6cLAABYNgIYAMBgAhgAwGACGADAYAIYAMBgAhgAwGD7dmrD+/fv74MHD+7U5gEAZnb33Xd/trsPzGt9OxbADh48mNXV1Z3aPADAzKrqf8xzfT6CBAAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYbN9OFwB7VVU9Zb67d6gS2J28R2BzjoABAAwmgAEADCaAAQAMJoABAAy2ZQCrqpuq6rGquu8UYy6tqnuq6v6q+u35lggAsFhmOQJ2c5JDmz1YVc9L8q4kV3T3tyZ57XxKAwBYTFsGsO6+I8nnTzHk7yd5b3f/4WT8Y3OqDQBgIc3jHLAXJ3l+VX24qu6uqtfNYZ0AAAtrHhdi3Zfk5UleleRPJ/lIVd3Z3Z9cP7CqDic5nCQXXHDBHDYNALD3zOMI2PEkt3f3/+3uzya5I8lLNhrY3Ue6e6W7Vw4cODCHTQMA7D3zCGC/keS7qmpfVf2ZJN+R5ME5rBcAYCFt+RFkVd2S5NIk+6vqeJK3JTk7Sbr7xu5+sKr+U5J7k3w1yc9396aXrAAAWHZbBrDuvnqGMdcnuX4uFQEALDhXwgcAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYbN9OFwA8c1X1lPnu3qFK2E28LmD3cwQMAGAwR8DYc/zvHpiV/QW7lSNgAACDCWAAAIMJYAAAgwlgAACDCWAAAIMJYAAAgwlgAACDCWAAAIMJYAAAg20ZwKrqpqp6rKru22LcX66qJ6rq++ZXHgDA4pnlCNjNSQ6dakBVnZXknUk+MIeaAAAW2pYBrLvvSPL5LYb9cJJfT/LYPIoCAFhkp30OWFWdm+R7kvzc6ZcDALD45nES/r9O8hPd/dWtBlbV4apararVEydOzGHTAAB7z745rGMlya1VlST7k1xeVU909/vWD+zuI0mOJMnKykrPYdsAAHvOaQew7r7w5HRV3ZzkNzcKXwAArNkygFXVLUkuTbK/qo4neVuSs5Oku288o9UBACygLQNYd18968q6+w2nVQ0AwBJwJXwAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwbYMYFV1U1U9VlX3bfL491fVvVX1iar6nap6yfzLBABYHLMcAbs5yaFTPP6pJK/s7m9L8o4kR+ZQFwDAwtq31YDuvqOqDp7i8d+Zmr0zyXmnXxYAwOKa9zlgb0ry/s0erKrDVbVaVasnTpyY86YBAPaGuQWwqvobWQtgP7HZmO4+0t0r3b1y4MCBeW0aAGBP2fIjyFlU1bcn+fkkl3X35+axTgCARXXaR8Cq6oIk703yD7v7k6dfEgDAYtvyCFhV3ZLk0iT7q+p4krclOTtJuvvGJNcleWGSd1VVkjzR3StnqmAAgL1ulm9BXr3F4z+Y5AfnVhEAwIJzJXwAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMG2vBURTJvc7/NrunuHKgHYO+w7Wc8RMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMEEMACAwQQwAIDBBDAAgMG2DGBVdVNVPVZV923yeFXVz1TVsaq6t6peNv8yAQAWxyxHwG5OcugUj1+W5KLJz+EkP3f6ZQEALK4tA1h335Hk86cYcmWSd/eaO5M8r6peNK8CAQAWzTzOATs3ySNT88cny56mqg5X1WpVrZ44cWIOmwYA2HuGnoTf3Ue6e6W7Vw4cODBy0wAAu8Y8AtijSc6fmj9vsgwAgA3MI4AdTfK6ybchX5Hki939mTmsFwBgIe3bakBV3ZLk0iT7q+p4krclOTtJuvvGJLcluTzJsSRfTvIDZ6pYAIBFsGUA6+6rt3i8k/yTuVUEALDgXAkfAGAwAQwAYDABDABgMAEMAGAwAQwAYDABDABgMAEMAGAwAQwAYDABDABgMAEMAGAwAQwAYDABDABgsC1vxg2wlap6ynx371AlO08vgFk4AgYAMJgABgAwmAAGADCYAAYAMJgABgAwmAAGADCYAAYAMJgABgAwmAAGADCYAAYAMJgABgAwmAAGADCYAAYAMJgABgAwmAAGADDYTAGsqg5V1UNVdayqrt3g8Quq6kNV9bGqureqLp9/qQAAi2HLAFZVZyW5IcllSS5OcnVVXbxu2E8meU93vzTJVUneNe9CAQAWxSxHwC5Jcqy7H+7ux5PcmuTKdWM6yddPpp+b5NPzKxEAYLHsm2HMuUkemZo/nuQ71o3550k+UFU/nOTPJnn1XKoDAFhA8zoJ/+okN3f3eUkuT/JLVfW0dVfV4apararVEydOzGnTAAB7yywB7NEk50/NnzdZNu1NSd6TJN39kSR/Ksn+9Svq7iPdvdLdKwcOHHhmFQMA7HGzBLC7klxUVRdW1TlZO8n+6Loxf5jkVUlSVd+StQDmEBcAwAa2DGDd/USSa5LcnuTBrH3b8f6qentVXTEZ9uNJfqiqPp7kliRv6O4+U0UDAOxls5yEn+6+Lclt65ZdNzX9QJLvnG9pAACLyZXwAQAGE8AAAAYTwAAABhPAAAAGE8AAAAYTwAAABhPAAAAGm+k6YOy8qnrKvOvcArAd/h3ZXRwBAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGGymAFZVh6rqoao6VlXXbjLm71XVA1V1f1X9ynzLBABYHPu2GlBVZyW5IcnfSnI8yV1VdbS7H5gac1GStyT5zu7+QlV9w5kqGABgr5vlCNglSY5198Pd/XiSW5NcuW7MDyW5obu/kCTd/dh8ywQAWByzBLBzkzwyNX98smzai5O8uKr+a1XdWVWH5lUgAMCi2fIjyG2s56IklyY5L8kdVfVt3f2/pgdV1eEkh5PkggsumNOmAQD2llmOgD2a5Pyp+fMmy6YdT3K0u/+kuz+V5JNZC2RP0d1Hunulu1cOHDjwTGsGANjTZglgdyW5qKourKpzklyV5Oi6Me/L2tGvVNX+rH0k+fAc6wQAWBhbBrDufiLJNUluT/Jgkvd09/1V9faqumIy7PYkn6uqB5J8KMk/7e7PnamiAQD2suruHdnwyspKr66u7si296Kqesr8Tv297YY6dryGdds/ObeUvdhldewGu6EXO1rDJu+P4XWc3P4u+PvYLXXshhr2sqq6u7tX5rU+V8IHABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhs304XADwDk1uKnLyRSG0+kmXidQF7hgDG3uEfF2BW9hfscj6CBAAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGGzfThewF1TVU+a7e4cqAYC9a/2/p8ny/pvqCBgAwGAzBbCqOlRVD1XVsaq69hTjvrequqpW5lciAMBi2TKAVdVZSW5IclmSi5NcXVUXbzDuOUnenOSj8y4SAGCRzHIE7JIkx7r74e5+PMmtSa7cYNw7krwzyR/PsT4AgIUzSwA7N8kjU/PHJ8u+pqpeluT87v6Pp1pRVR2uqtWqWj1x4sS2iwUAWASnfRJ+VT0ryU8n+fGtxnb3ke5e6e6VAwcOnO6mAQD2pFkC2KNJzp+aP2+y7KTnJPmLST5cVf89ySuSHHUi/oKpSqrSSZbzC8MA2zTZb57cd8K0WQLYXUkuqqoLq+qcJFclOXrywe7+Ynfv7+6D3X0wyZ1Jruju1TNSMQDAHrdlAOvuJ5Jck+T2JA8meU93319Vb6+qK850gQAAi2amK+F3921Jblu37LpNxl56+mUBACwuV8IHABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhMAAMAGEwAAwAYTAADABhsplsRAWyoKknSJ2d3rpKdVU8+86XvBTATR8AAAAYTwAAABhPAAAAGE8AAAAYTwAAABhPAAAAGE8AAAAYTwAAABhPAAAAGE8AAAAZzK6Ldbt2tXhK3OAFgG9wybFdyBAwAYDABDABgMAEMAGAwAQwAYDABDABgMAEMAGAwAQwAYLCZAlhVHaqqh6rqWFVdu8HjP1ZVD1TVvVX1W1X1TfMorqqe8gMAcLp2Q77YMoBV1VlJbkhyWZKLk1xdVRevG/axJCvd/e1Jfi3JT827UACARTHLEbBLkhzr7oe7+/Ektya5cnpAd3+ou788mb0zyXnzLRMAYHHMEsDOTfLI1PzxybLNvCnJ+zd6oKoOV9VqVa2eOHFi9ioBABbIXE/Cr6p/kGQlyfUbPd7dR7p7pbtXDhw4MM9NAwDsGbPcjPvRJOdPzZ83WfYUVfXqJG9N8sru/sp8ygMAWDyzHAG7K8lFVXVhVZ2T5KokR6cHVNVLk/zbJFd092PzLxMAYHFsGcC6+4kk1yS5PcmDSd7T3fdX1dur6orJsOuTfF2SX62qe6rq6CarAwBYerN8BJnuvi3JbeuWXTc1/eo51wUAsLBcCR8AYDABDABgMAEMAGAwAQwAYLCZTsJfWpMbdPbJ2Z2rBAD2Lv+ePo0jYAAAgwlgAACDCWAAAIMJYAAAgwlgAACDCWAAAIMJYAAAgwlgAACDCWAAAIMJYAAAgwlgAACDCWAAAIPt3ptxV33tpp2JG3cCAKehnkwSu+Gm4I6AAQAMJoABAAwmgAEADCaAAQAMJoABAAwmgAEADCaAAQAMJoABAAwmgAEADCaAAQAMNlMAq6pDVfVQVR2rqms3ePzZVfXvJo9/tKoOzrtQAIBFsWUAq6qzktyQ5LIkFye5uqouXjfsTUm+0N1/Psm/SvLOeRcKALAoZjkCdkmSY939cHc/nuTWJFeuG3Nlkl+cTP9akldVlftnAwBsYJYAdm6SR6bmj0+WbTimu59I8sUkL5xHgQAAi2bfyI1V1eEkhyezX6mq+7b5+/Mv6hk4w3XsT/LZHa5hZgPr2LQvS9iLk57WkyXuxbSZ3kOj7Ia/k0kNO96XXdSLaTvSl13ai5OG9mSX92LaX5jnNmcJYI8mOX9q/rzJso3GHK+qfUmem+Rz61fU3UeSHEmSqlrt7pVnUvQi05eN6cvT6cnG9GVj+rIxfXk6PdlYVa3Oc32zfAR5V5KLqurCqjonyVVJjq4bczTJ6yfT35fkv3R3z69MAIDFseURsO5+oqquSXJ7krOS3NTd91fV25OsdvfRJL+Q5Jeq6liSz2ctpAEAsIGZzgHr7tuS3LZu2XVT03+c5LXb3PaRbY5fFvqyMX15Oj3ZmL5sTF82pi9Ppycbm2tfyieFAABjuRURAMBgcw1gVXVTVT02fXmJqvpLVXVnVd1TVatVdclkeVXVz0xuX3RvVb1s6ndeX1W/P/l5/Ubb2iu22ZNLq+qLk+X3VNV1U79zyttB7TWb9OUlVfWRqvpEVf2Hqvr6qcfeMnnuD1XVd08tX9q+VNXBqvp/U6+XG6d+5+WT8ccm77Od/573M1RV51fVh6rqgaq6v6rePFn+gqr64GQ/8cGqev5k+bLsW7bbl6XYv5yiL6+dzH+1qlbW/c7C71+225dl2L+coifXV9XvTfYf/76qnjf1O/N7rXT33H6S/PUkL0ty39SyDyS5bDJ9eZIPT02/P0kleUWSj06WvyDJw5M/nz+Zfv486xz5s82eXJrkNzdYx1lJ/iDJNyc5J8nHk1y808/tDPTlriSvnEy/Mck7JtMXT57zs5NcOOnFWfqSg9Pj1q3nv03eVzV5n12208/tNHryoiQvm0w/J8knJ6+Jn0py7WT5tUneOZleln3LdvuyFPuXU/TlW7J2HacPJ1mZGr8U+5dn0JeF37+coievSbJvsvydU++hub5W5noErLvvyNq3IJ+yOMnJIxnPTfLpyfSVSd7da+5M8ryqelGS707ywe7+fHd/IckHkxyaZ50jbbMnm5nldlB7yiZ9eXGSOybTH0zyvZPpK5Pc2t1f6e5PJTmWtZ4se182NHkffX1339lre413J/m78651lO7+THf/7mT6fyd5MGt335i+Bdov5snnuCz7lu32ZTML9T7arC/d/WB3P7TBryzF/uUZ9GVDi7R/OUVPPtBrd/VJkjuzdv3TZM6vlRHngP1Ikuur6pEk/zLJWybLN7vF0Sy3PtrrNutJkvyVqvp4Vb2/qr51smwZepIk9+fJF+1r8+QFgJf5tZJs3pckubCqPlZVv11Vf22y7Nys9eKkhelLVR1M8tIkH03yjd39mclD/zPJN06ml+71MmNfkiXbv6zry2aW/fVyKkuzfzlFT96YtaN8yZxfKyMC2D9K8qPdfX6SH83aNcOW3WY9+d0k39TdL0nys0net0P17ZQ3JvnHVXV31g4HP77D9ewWm/XlM0ku6O6XJvmxJL9SU+fNLZqq+rokv57kR7r7S9OPTf4nvpRf6d5GX5Zq/3KqviyzbfRlafYvm/Wkqt6a5Ikkv3wmtjsigL0+yXsn07+atUN1yea3OJrl1kd73YY96e4vdff/mUzfluTsqtqf5ehJuvv3uvs13f3yJLdk7TP1ZLlfK5v2ZXIY/HOT6bsny1+ctR6cN7WKPd+Xqjo7azvIX+7uk++dP5p8HHLyY5HHJsuX5vWynb4s0/5lk75sZtlfLxtalv3LZj2pqjck+dtJvn/yH5lkzq+VEQHs00leOZn+m0l+fzJ9NMnras0rknxxctj89iSvqarnT76985rJskWyYU+q6s+d/DZJrX0z8llZu6fmLLeD2vOq6hsmfz4ryU8mOfmtm6NJrqqqZ1fVhUkuytpJoEvdl6o6UFVnTaa/OWt9eXjyPvpSVb1i8np6XZLf2JHi52DyHH4hyYPd/dNTD03fAu31efI5LsW+Zbt9WZb9yyn6spml2L9sty/LsH/ZrCdVdSjJP0tyRXd/eepX5vta2eos/e38ZO1/559J8idZ+wz0TUm+K8ndWftWwEeTvHwytpLckLVU/Yk89dsXb8zayW3HkvzAPGsc/bPNnlyTtfN9Pp61E//+6tR6Ls/aNzT+IMlbd/p5naG+vHnyHD+Z5F9kcqHgyfi3Tp77Q5n6xs0y9yVrJ+Pfn+SerH289Hem1rOS5L5JX/7NdC/32s/k/dJJ7p0813smf+8vTPJbWfsPzH9O8oLJ+GXZt2y3L0uxfzlFX75n8p76SpI/SnL71O8s/P5lu31Zhv3LKXpyLGvndJ1cduOZeK24Ej4AwGCuhA8AMJgABgAwmAAGADCYAAYAMJgABgAwmAAGADCYAAYAMJgABgAw2P8Hq4Qq7F/lEwgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "M35=35*1/(35+1)\n",
    "M37=37*1/(37+1)\n",
    "ν0_37=ν0*np.sqrt(M35/M37)\n",
    "B_37=B*(M35/M37)\n",
    "\n",
    "Ej_37=B_37*J*(J+1)\n",
    "rR_37=ν0_37+Ej_37[1:]-Ej_37[:-1]\n",
    "rP_37=ν0_37+Ej_37[:-1]-Ej_37[1:]\n",
    "Ir_37=J*np.exp(-J*(J+1)/20)*0.35\n",
    "plt.bar(rR,Ir[1:],width=3, color='black')\n",
    "plt.bar(rP,Ir[1:],width=3, color='black')\n",
    "plt.bar(rR_37,Ir_37[1:],width=3, color='red')\n",
    "plt.bar(rP_37,Ir_37[1:],width=3, color='red')\n",
    "plt.xlim(1800,2200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9722222222222222 0.9736842105263158\n"
     ]
    }
   ],
   "source": [
    "print(M35,M37)"
   ]
  },
  {
   "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.11.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
