docs/PopSyCLE_example.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "from popsycle import synthetic\n",
    "from popsycle import utils\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from astropy.table import Table"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "In additional to the usual python/astropy things, you'll need to install Galaxia, SPISEA, PopSyCLE.\n",
    "For SPISEA, it goes faster once you've made isochrones.\n",
    "Same goes from PopSyCLE with the initial-final cluster mass ratio file."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "# Galaxia\n",
    "\n",
    "This is where you generate the stellar model and set the survey direction and area.\n",
    "\n",
    "Galaxia requires a text file containing the parameters you want for the model,\n",
    "following the format given in the Galaxia example. That file can either be \n",
    "generated by hand, or with the `write_galaxia_params` function. \n",
    "\n",
    "To create the file `example_galaxia_params.txt`, run:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "metadata": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "synthetic.write_galaxia_params(output_root = 'example',\n",
    "                               longitude = 1.25,\n",
    "                               latitude = -2.65,\n",
    "                               area = 0.001)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "metadata": false
    }
   },
   "source": [
    "\n",
    "Then create the stellar model, \n",
    "by running `galaxia -r example_galaxia_params.txt` on the command line. \n",
    "PopSyCLE can also combine these two steps by executing:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "metadata": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "synthetic.run_galaxia(output_root = 'example',\n",
    "                      longitude = 1.25,\n",
    "                      latitude = -2.65,\n",
    "                      area = 0.001,\n",
    "                      galaxia_galaxy_model_filename= '/g/lu/code/galaxia/docs/galaxyModelParams_PopSyCLEv3.txt')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "metadata": false
    }
   },
   "source": [
    " \n",
    "The results should look like this:\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "\n",
    "CODEDATAPATH=/u/casey_lam/scratch/GalaxiaData/\n",
    "\n",
    "Reading Parameter file-             example_galaxia_params.txt\n",
    "___\n",
    "\n",
    "outputFile               example                 \n",
    "outputDir                ./                      \n",
    "photoSys                 UBV                     \n",
    "magcolorNames            V,B-V                   \n",
    "appMagLimits[0]          -1000                   \n",
    "appMagLimits[1]          1000                    \n",
    "absMagLimits[0]          -1000                   \n",
    "absMagLimits[1]          1000                    \n",
    "colorLimits[0]           -1000                   \n",
    "colorLimits[1]           1000                    \n",
    "geometryOption           1                       \n",
    "longitude                1.25                    \n",
    "latitude                 -2.65                   \n",
    "surveyArea               0.001                   \n",
    "fSample                  1                       \n",
    "popID                    -1                      \n",
    "warpFlareOn              1                       \n",
    "seed                     12                      \n",
    "r_max                    20                      \n",
    "starType                 0                       \n",
    "photoError               0                       \n",
    "___\n",
    "Reading tabulated values from file- /u/casey_lam/scratch/GalaxiaData/Model/vcirc.dat\n",
    "Using geometry:                     Patch at l ,b : (1.25 -2.65) d_theta=0.0178412\n",
    "Reading Isochrones from dir-        /u/casey_lam/scratch/GalaxiaData/Isochrones/padova/UBV\n",
    "Isochrone Grid Size:                (Age bins=182,Feh bins=34,Alpha bins=1)\n",
    "Time Isocrhone Reading              2.73553     \n",
    "Generating populations................\n",
    "___\n",
    "Thin Disk,ID=0:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_0_E1.ebf\n",
    "Time Tree generation/reading =      1.34589     \n",
    "Completed % <0..10..20..30..40..50..60..70..80..90..>\n",
    "Stars spawned =                     15          \n",
    "Time Spawning=                      0.519746    \n",
    "___\n",
    "Thin Disk,ID=1:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_1_E1.ebf\n",
    "Time Tree generation/reading =      0.894184    \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     785         \n",
    "Time Spawning=                      3.01123     \n",
    "___\n",
    "Thin Disk,ID=2:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_2_E1.ebf\n",
    "Time Tree generation/reading =      1.21969     \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     2304        \n",
    "Time Spawning=                      2.31239     \n",
    "___\n",
    "Thin Disk,ID=3:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_3_E1.ebf\n",
    "Time Tree generation/reading =      1.12221     \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     3299        \n",
    "Time Spawning=                      2.63636     \n",
    "___\n",
    "Thin Disk,ID=4:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_4_E1.ebf\n",
    "Time Tree generation/reading =      1.12869     \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     9040        \n",
    "Time Spawning=                      5.16435     \n",
    "___\n",
    "Thin Disk,ID=5:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_5_E1.ebf\n",
    "Time Tree generation/reading =      0.924357    \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     9983        \n",
    "Time Spawning=                      5.99943     \n",
    "___\n",
    "Thin Disk,ID=6:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_6_E1.ebf\n",
    "Time Tree generation/reading =      1.03473     \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     11927       \n",
    "Time Spawning=                      5.80737     \n",
    "___\n",
    "ThickDisk:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_7_E0.ebf\n",
    "Time Tree generation/reading =      1.00016     \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     6561        \n",
    "Time Spawning=                      1.83702     \n",
    "___\n",
    "Spheroid:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_8_E0.ebf\n",
    "Time Tree generation/reading =      0.320754    \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     589         \n",
    "Time Spawning=                      0.693446    \n",
    "___\n",
    "Bulge:\n",
    "Reading tree from file-             /u/casey_lam/scratch/GalaxiaData/BHTree-2.2/bhtree_with_wf/bhtree_9_E0.ebf\n",
    "Time Tree generation/reading =      1.0575      \n",
    "Completed % <0..9..19..29..39..49..59..69..79..89..99..>\n",
    "Stars spawned =                     620676      \n",
    "Time Spawning=                      28.1113     \n",
    "___\n",
    "Total stars written                 665179                  \n",
    "File written-                       .//example.ebf\n",
    "Calulating magnitudes................\n",
    "Reading Isochrones from dir-        /u/casey_lam/scratch/GalaxiaData/Isochrones/padova/UBV\n",
    "Isochrone Grid Size:                (Age bins=182,Feh bins=34,Alpha bins=1)\n",
    "Time Isocrhone Reading              2.67512     \n",
    "Calulating Extinction................\n",
    "Time for extinction calculation     0.876306    \n",
    "Total Time=                         73.89 \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "# Population Synthesis\n",
    "\n",
    "The stellar model now needs compact objects, which we inject with SPISEA. \n",
    "This produces an .h5 file with both stars and compact objects, sorted by latitude and longitude.\n",
    "Here you can set the black hole and neutron star kick speed means, as well as the SPISEA IFMR object used to produce the compact objects."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "*********************** Starting popid 0\n",
      "Starting age bin  6.995522499084473\n",
      "Starting age bin  6.3\n",
      "Starting sub-bin  0\n",
      "Starting age bin  7.0\n",
      "Starting sub-bin  0\n",
      "Starting age bin  7.7\n",
      "Starting sub-bin  0\n",
      "Starting age bin  8.0\n",
      "Starting sub-bin  0\n",
      "*********************** Starting popid 1\n",
      "Starting age bin  8.178046226501465\n",
      "Starting sub-bin  0\n",
      "Starting age bin  8.38562906551361\n",
      "Starting sub-bin  0\n",
      "Starting age bin  8.593211904525756\n",
      "Starting sub-bin  0\n",
      "Found 2 stars out of mass range\n",
      "Starting age bin  8.800794743537903\n",
      "Starting sub-bin  0\n",
      "Found 14 stars out of mass range\n",
      "*********************** Starting popid 2\n",
      "Starting age bin  9.000550270080566\n",
      "Starting sub-bin  0\n",
      "Found 7 stars out of mass range\n",
      "Starting age bin  9.077987862825394\n",
      "Starting sub-bin  0\n",
      "Found 12 stars out of mass range\n",
      "Starting age bin  9.15542545557022\n",
      "Starting sub-bin  0\n",
      "Found 1 stars out of mass range\n",
      "Starting age bin  9.232863048315048\n",
      "Starting sub-bin  0\n",
      "Found 12 stars out of mass range\n",
      "*********************** Starting popid 3\n",
      "Starting age bin  9.301044464111328\n",
      "Starting sub-bin  0\n",
      "Found 12 stars out of mass range\n",
      "Starting age bin  9.347432012081146\n",
      "Starting sub-bin  0\n",
      "Found 25 stars out of mass range\n",
      "Starting age bin  9.393819560050964\n",
      "Starting sub-bin  0\n",
      "Found 23 stars out of mass range\n",
      "Starting age bin  9.44020710802078\n",
      "Starting sub-bin  0\n",
      "Found 11 stars out of mass range\n",
      "*********************** Starting popid 4\n",
      "Starting age bin  9.477165222167969\n",
      "Starting sub-bin  0\n",
      "Found 37 stars out of mass range\n",
      "Starting age bin  9.535034675598144\n",
      "Starting sub-bin  0\n",
      "Found 52 stars out of mass range\n",
      "Starting age bin  9.59290412902832\n",
      "Starting sub-bin  0\n",
      "Found 257 stars out of mass range\n",
      "Starting age bin  9.650773582458495\n",
      "Starting sub-bin  0\n",
      "Found 108 stars out of mass range\n",
      "*********************** Starting popid 5\n",
      "Starting age bin  9.69898796081543\n",
      "Starting sub-bin  0\n",
      "Found 150 stars out of mass range\n",
      "Starting age bin  9.7379761531353\n",
      "Starting sub-bin  0\n",
      "Found 110 stars out of mass range\n",
      "Starting age bin  9.77696434545517\n",
      "Starting sub-bin  0\n",
      "Found 127 stars out of mass range\n",
      "Starting age bin  9.815952537775038\n",
      "Starting sub-bin  0\n",
      "Found 19 stars out of mass range\n",
      "*********************** Starting popid 6\n",
      "Starting age bin  9.845110893249512\n",
      "Starting sub-bin  0\n",
      "Found 78 stars out of mass range\n",
      "Starting age bin  9.886327680826188\n",
      "Starting sub-bin  0\n",
      "Found 120 stars out of mass range\n",
      "Starting age bin  9.927544468402862\n",
      "Starting sub-bin  0\n",
      "Found 102 stars out of mass range\n",
      "Starting age bin  9.968761255979537\n",
      "Starting sub-bin  0\n",
      "Found 4 stars out of mass range\n",
      "*********************** Starting popid 7\n",
      "Starting age bin  9.94097840309143\n",
      "Starting sub-bin  0\n",
      "Found 172 stars out of mass range\n",
      "*********************** Starting popid 8\n",
      "Starting age bin  10.038501205444335\n",
      "Starting sub-bin  0\n",
      "Found 14 stars out of mass range\n",
      "*********************** Starting popid 9\n",
      "Starting age bin  9.9\n",
      "Starting sub-bin  0\n",
      "Found 9219 stars out of mass range\n",
      "perform_pop_syn runtime : 17.248314 s\n",
      "******************** INFO **********************\n",
      "Total number of stars from Galaxia: 367004\n",
      "Total number of compact objects made: 24775\n",
      "Total number of things binned: 391779\n"
     ]
    }
   ],
   "source": [
    "synthetic.perform_pop_syn(ebf_file = 'example.ebf', \n",
    "                          output_root = 'example',\n",
    "                          iso_dir = '/g/lu/models/PopSyCLE_isochrones',\n",
    "                          IFMR = 'SukhboldN20',\n",
    "                          bin_edges_number = None, \n",
    "                          BH_kick_speed_mean = 100, \n",
    "                          NS_kick_speed_mean = 350,\n",
    "                          overwrite=True);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "# Finding microlensing events\n",
    "\n",
    "Here is when all microlensing events are found.\n",
    "You can set the survey duration, cadence, maximum impact parameter, blend radius.\n",
    "This also has the option of being parallelizable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Working on loop ll, bb =  0 0\n",
      "Candidate events detected:  21\n",
      "calc_events runtime : 16.382421 s\n"
     ]
    }
   ],
   "source": [
    "synthetic.calc_events(hdf5_file = 'example.h5', \n",
    "                      output_root2 = 'example', \n",
    "                      radius_cut = 2, \n",
    "                      obs_time = 1000, \n",
    "                      n_obs = 11, \n",
    "                      theta_frac = 2, \n",
    "                      blend_rad = 0.65, \n",
    "                      overwrite = False, \n",
    "                      n_proc = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "# Photometry\n",
    "\n",
    "This is the last thing, where you choose the photometric band for the observations and the reddening law.\n",
    "The final file of microlensing events is produced here as a .fits file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Original candidate events:  21\n",
      "Candidate events in survey window:  19\n",
      "refine_events runtime : 0.160328 s\n"
     ]
    }
   ],
   "source": [
    "synthetic.refine_events(input_root = 'example', \n",
    "                        filter_name = 'I',\n",
    "                        photometric_system = 'ubv',\n",
    "                        red_law = 'Damineli16', \n",
    "                        overwrite = False, \n",
    "                        output_file = 'default')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "# Things you can do with PopSyCLE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "First, read in the table, and print out a list of all the possible keys:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['zams_mass_L', 'mass_L', 'px_L', 'py_L', 'pz_L', 'vx_L', 'vy_L', 'vz_L', 'age_L', 'popid_L', 'exbv_L', 'glat_L', 'glon_L', 'mbol_L', 'grav_L', 'teff_L', 'feh_L', 'rad_L', 'rem_id_L', 'obj_id_L', 'ubv_J_L', 'ubv_H_L', 'ubv_K_L', 'ubv_U_L', 'ubv_I_L', 'ubv_B_L', 'ubv_V_L', 'ubv_R_L', 'vr_L', 'mu_b_L', 'mu_lcosb_L', 'zams_mass_S', 'mass_S', 'px_S', 'py_S', 'pz_S', 'vx_S', 'vy_S', 'vz_S', 'age_S', 'popid_S', 'exbv_S', 'glat_S', 'glon_S', 'mbol_S', 'grav_S', 'teff_S', 'feh_S', 'rad_S', 'rem_id_S', 'obj_id_S', 'ubv_J_S', 'ubv_H_S', 'ubv_K_S', 'ubv_U_S', 'ubv_I_S', 'ubv_B_S', 'ubv_V_S', 'ubv_R_S', 'vr_S', 'mu_b_S', 'mu_lcosb_S', 'theta_E', 'u0', 'mu_rel', 't0', 't_E', 'ubv_I_app_S', 'ubv_I_app_L', 'cent_glon_I_N', 'cent_glat_I_N', 'ubv_I_app_N', 'ubv_I_app_LSN', 'delta_m_I', 'f_blend_I', 'pi_rel', 'pi_E', 'pps_seed', 'gal_seed']\n"
     ]
    }
   ],
   "source": [
    "t = Table.read('example_refined_events_ubv_I_Damineli16.fits')\n",
    "print(t.colnames)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "## Figure out how many events due to BH, NS, WD, star"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of BHs: 0\n",
      "Number of NSs: 0\n",
      "Number of WDs: 4\n",
      "Number of stars: 20\n"
     ]
    }
   ],
   "source": [
    "bh_idx = np.where(t['rem_id_L'] == 103)[0]\n",
    "ns_idx = np.where(t['rem_id_L'] == 102)[0]\n",
    "wd_idx = np.where(t['rem_id_L'] == 101)[0]\n",
    "st_idx = np.where(t['rem_id_L'] == 0)[0]\n",
    "print('Number of BHs: ' + str(len(bh_idx)))\n",
    "print('Number of NSs: ' + str(len(ns_idx)))\n",
    "print('Number of WDs: ' + str(len(wd_idx)))\n",
    "print('Number of stars: ' + str(len(st_idx)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "## Histogram of tE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, '$t_E$ (days)')"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEOCAYAAACTqoDjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAASRElEQVR4nO3df6zddX3H8efLYocwVoNA3FpcqyXVhkWriBplOOdMgbVMGI6GGSykRDL8kc0sNVuC29iGzJiNDGpAGUskdNDxq9KBc6NjI+gAQdZamYXAuEGDitSBY+XHe3+cw3o9fm7vObfn3sPtfT6Sm/T7+X6+n+/79NPe1/18v+eeb6oKSZJ6vWzUBUiSXpoMCElSkwEhSWoyICRJTQaEJKnJgJAkNR0w6gKG5bDDDqvFixePugxJmlXuueee71fV4a19+01ALF68mLvvvnvUZUjSrJLkkYn2eYlJktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1zfrfg0iyCli1dOnSUZciNS1ef/O0jPvwhSdNy7jSi2b9CqKqNlfVOQsWLBh1KZK0X5n1ASFJmh4GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1zfqASLIqyWW7du0adSmStF+Z9QHhh/VJ0vSY9QEhSZoeBoQkqcmAkCQ1GRCSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQeMuoC9SXIwcCmwG9haVVeNuCRJmjP6WkEkeWWSTUm+lWRHkndM5WRJrkjyeJJtjX0rkzyQZGeS9d3mU4BNVbUOWD2Vc0qSpqbfS0x/BdxSVa8H3gjsGL8zyRFJDulpW9oY50pgZW9jknnAJcAJwHJgTZLlwCLg0W635/usVZI0BJMGRJKfA34Z+AJAVe2uqid7uh0P3JjkwO4x64CLe8eqqtuBJxqnORbYWVUPVdVuYCNwMjBGJyT6qlWSNDz9fNN9LfA94G+S3Jvk8917A/+vqq4FbgE2JjkDOAv4wAB1LGTPSgE6wbAQuA44NckGYHPrQD/uW5KmRz8BcQDwZmBDVa0AngbW93aqqouAZ4ANwOqqemqAOtJoq6p6uqrWVtW5E92g9uO+JWl69BMQY8BYVX2tu72JTmD8hCTHAUcD1wPnD1jHGHDkuO1FwGMDjiFJGqJJA6Kqvgs8mmRZt+lXgW+O75NkBXA5nfsGa4FDk1wwQB13AUclWZJkPnA6cNMAx0uShqzfG78fAa5Kcj/wJuDPevYfBJxWVQ9W1QvAmcAjvYMkuRq4E1iWZCzJ2QBV9RxwHnArnXdIXVNV26fygiRJw9HXL8pV1X3AMXvZf0fP9rN0VhS9/dbsZYwtwJZ+6pEkTT/fOipJajIgJElNBoQkqcmAkCQ1GRCSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQaEJKmpr89ieilLsgpYtXRp6wmn0mAWr7951CVILxmzfgXhA4MkaXrM+oCQJE0PA0KS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1vaSfKJfkYOBSYDewtaquGnFJkjRn9L2CSDIvyb1JvjTVkyW5IsnjSbY19q1M8kCSnUnWd5tPATZV1Tpg9VTPK0ka3CCXmD4G7GjtSHJEkkN62loPib4SWNk4fh5wCXACsBxYk2Q5sAh4tNvt+QFqlSTto74CIski4CTg8xN0OR64McmB3f7rgIt7O1XV7cATjeOPBXZW1UNVtRvYCJwMjNEJiQlrTbIqyWW7du3q56VIkvrU7wriL4HfB15o7ayqa4FbgI1JzgDOAj4wQB0L2bNSgE4wLASuA05NsgHYPMG5N1fVOQsWLBjgdJKkyUx6kzrJrwOPV9U9Sd49Ub+quijJRmAD8LqqemqAOtIesp4G1g4wjiRpSPpZQbwTWJ3kYTqXft6T5Iu9nZIcBxwNXA+cP2AdY8CR47YXAY8NOIYkaYgmDYiq+mRVLaqqxcDpwD9X1W+P75NkBXA5nfsGa4FDk1wwQB13AUclWZJkfvc8Nw1wvCRpyIb1i3IHAadV1YNV9QJwJvBIb6ckVwN3AsuSjCU5G6CqngPOA26l806pa6pq+5BqkyRNwUC/KFdVW4GtjfY7erafpbOi6O23Zi9jbwG2DFKPJGn6+FEbkqQmA0KS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1GRCSpKaBnig305IcDFwK7Aa2VtVVIy5JkuaMSVcQSQ5M8u9JvpFke5I/murJklyR5PEk2xr7ViZ5IMnOJOu7zacAm6pqHbB6queVJA2un0tM/wu8p6reCLwJWJnk7eM7JDkiySE9bUsbY10JrOxtTDIPuAQ4AVgOrEmyHFgEPNrt9nwftUqShmTSgKiOp7qbL+9+VU+344EbkxwIkGQdcHFjrNuBJxqnORbYWVUPVdVuYCNwMjBGJyT6qlWSNDx93YPo/oR/D7AUuKSqvjZ+f1Vdm2QJsDHJtcBZwK8NUMdC9qwUoBMMb6MTMn+d5CRg8wS1rQJWLV3aWrBI+6/F62+elnEfvvCkaRlXs09fP5VX1fNV9SY6P80fm+ToRp+LgGeADcDqcauOfqR92nq6qtZW1bkT3aCuqs1Vdc6CBQsGOJ0kaTIDXbapqieBrbTvIxwHHA1cD5w/YB1jwJHjthcBjw04hiRpiPp5F9PhSV7Z/fMrgPcC3+rpswK4nM59g7XAoUkuGKCOu4CjkixJMh84HbhpgOMlSUPWzwri54HbktxP5xv5P1bVl3r6HAScVlUPVtULwJnAI70DJbkauBNYlmQsydkAVfUccB5wK7ADuKaqtk/1RUmS9t2kN6mr6n5gxSR97ujZfpbOiqK335q9jLEF2DJZPZKkmeFbRyVJTQaEJKnJgJAkNRkQkqQmA0KS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1GRCSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQaEJKnJgJAkNRkQkqQmA0KS1HTAqAvYmyQHA5cCu4GtVXXViEuSpDlj0hVEkiOT3JZkR5LtST421ZMluSLJ40m2NfatTPJAkp1J1nebTwE2VdU6YPVUzytJGlw/l5ieA36vqt4AvB34nSTLx3dIckSSQ3raljbGuhJY2duYZB5wCXACsBxY0z3HIuDRbrfn+6hVkjQkk15iqqrvAN/p/vm/k+wAFgLfHNfteODcJCdW1TNJ1gHvB07sGev2JIsbpzkW2FlVDwEk2QicDIzRCYn7mCDMkqwCVi1d2soj7a8Wr7951CXst2bT3+3DF5406hL2awPdpO5+c18BfG18e1VdC9wCbExyBnAW8IEBhl7InpUCdIJhIXAdcGqSDcDm1oFVtbmqzlmwYMEAp5MkTabvm9RJfhb4e+DjVfWj3v1VdVH3J/8NwOuq6qkB6kijrarqaWDtAONIkoakrxVEkpfTCYerquq6CfocBxwNXA+cP2AdY8CR47YXAY8NOIYkaYj6eRdTgC8AO6rqsxP0WQFcTue+wVrg0CQXDFDHXcBRSZYkmQ+cDtw0wPGSpCHrZwXxTuCDwHuS3Nf9OrGnz0HAaVX1YFW9AJwJPNI7UJKrgTuBZUnGkpwNUFXPAecBtwI7gGuqavuUX5UkaZ/18y6mf6N9j2B8nzt6tp+ls6Lo7bdmL2NsAbZMVo8kaWb4URuSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQaEJKnJgJAkNRkQkqQmA0KS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1GRCSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQaEJKnJgJAkNRkQkqQmA0KS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkJgNCktRkQEiSmgwISVKTASFJajIgJElNBoQkqcmAkCQ1GRCSpCYDQpLUZEBIkpoMCElSkwEhSWoyICRJTQaEJKnJgJAkNR0w6gJakhwMXArsBrZW1VUjLkmS5pwZW0EkuSLJ40m29bSvTPJAkp1J1nebTwE2VdU6YPVM1ShJ2mMmLzFdCawc35BkHnAJcAKwHFiTZDmwCHi02+35GaxRktQ1Y5eYqur2JIt7mo8FdlbVQwBJNgInA2N0QuI+9hJiSc4BzgF4zWteM/yitc8Wr7951CVIA5uuf7cPX3jStIw7XUZ9k3ohe1YK0AmGhcB1wKlJNgCbJzq4qi6rqmOq6pjDDz98eiuVpDlm1Dep02irqnoaWDvTxUiS9hj1CmIMOHLc9iLgsRHVIkkaZ9QBcRdwVJIlSeYDpwM3jbgmSRIz+zbXq4E7gWVJxpKcXVXPAecBtwI7gGuqavtM1SRJmthMvotpzQTtW4AtM1WHJKk/o77EJEl6iTIgJElNBoQkqSlVNeoahiLJ94BHupsLgF0TdG3tOwz4/jSVNlV7ew2jGnOQ4/vtO1m/Qedyb+3O8/CPdZ4nNx1zvK/jjj/2F6uq/ZvGVbXffQGXDbIPuHvUNQ/yGkY15iDH99t3sn6DzuUk7c7zkI91nkczxzM1z/vrJaYJP55jkn0vJdNR576OOcjx/fadrN9U5nK2zDG89OZ50GOd58lNV53TPs/7zSWmfZHk7qo6ZtR1aHo5z3OD8zw8++sKYlCXjboAzQjneW5wnofEFYQkqckVhCSpyYCQJDUZEJKkJgOiIcnBSf42yeVJzhh1PZoeSV6b5AtJNo26Fk2fJL/R/b98Y5L3jbqe2WTOBESSK5I8nmRbT/vKJA8k2Zlkfbf5FGBTVa0DVs94sZqyQea5qh6qqrNHU6n2xYDzfEP3//KHgN8aQbmz1pwJCOBKYOX4hiTzgEuAE4DlwJoky+k82e7FZ2U/P4M1at9dSf/zrNnrSgaf5z/s7lef5kxAVNXtwBM9zccCO7s/Se4GNgIn03kU6qJunznzd7Q/GHCeNUsNMs/p+DTwD1X19ZmudTab69/8FrJnpQCdYFgIXAecmmQDs+fX+TWx5jwneVWSzwErknxyNKVpiCb6//wR4L3Abyb58CgKm61m7IlyL1FptFVVPQ2sneliNG0mmucfAH7D2H9MNM8XAxfPdDH7g7m+ghgDjhy3vQh4bES1aPo4z3OD8zxkcz0g7gKOSrIkyXzgdOCmEdek4XOe5wbnecjmTEAkuRq4E1iWZCzJ2VX1HHAecCuwA7imqraPsk7tG+d5bnCeZ4Yf1idJapozKwhJ0mAMCElSkwEhSWoyICRJTQaEJKnJgJAkNRkQkqQmA0KS1GRASD2SLErS14Nlkrwiyb90n0XQu+9TST4xxLrmJ7k9yVz/kE3NEANC+mm/Cry5z75nAddV1bQ/WKr7jIN/wqeiaYYYENI4Sd4FfJbOswPuS7JkkkPOAG4cd/wfdB95+RVgWc/YNyS5J8n2JOd02/4kycfG9fnTJB/tPhf95iTfSLJt3Irmhu45pWnnZzFJPZLcAnyiqrZN0m8+8F9V9eru9lvoPArzbXSetfJ14HNV9Znu/kOr6okkr6DzyaPHA4fQWYG8OcnLgG/TeTLau4GV3Wcpk2RBVe3qXsr6blUdPuzXLfVyBSH9tGXAAwDdn/g/1/16a0+/w4Anx20fB1xfVT+uqh/x0x81/dEk3wC+Sue5BUdV1cPAD5KsAN4H3Nt9kNF/AO9N8ukkx1XVLoDupazdSQ4Z6iuWGrzZJY2T5FXArqp6NsmRwF1VNdFT5/4HOLCnrbkkT/JuOo+9fEdV/TjJ1nHHfh74EPBq4AqAqvrP7orkRODPk3y5qv642/9ngGem8PKkgbiCkH7SEvY8hewtwBu6q4dP9Xasqh8C85K8+I3+duD93Xc2HQKsGtd9AfDDbji8Hnj7uH3XAyuBt9J5lgFJfgH4cVV9EfgM3Zvm3QD7XlU9O5RXK+2FKwjpJ30LOCzJNuA7wMer6t699P8y8C7gK1X19SR/B9wHPAL867h+twAfTnI/nctXX31xR1XtTnIb8OS4d0P9EvAXSV4AngXO7bb/CrBlX1+k1A9vUksTSLKFzmriOeDWqrq+0WcF8LtV9cF9OM/L6NzQPq2qvj1J3+uAT1bVA1M9n9QvVxDSBKrqxD763JvktiTzpvK7EEmWA1+ic3N7snCYD9xgOGimuIKQJDV5k1qS1GRASJKaDAhJUpMBIUlqMiAkSU0GhCSpyYCQJDUZEJKkpv8Dgp5kam0FyeAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(t['t_E'], bins = np.logspace(0, 2.5, 15))\n",
    "plt.xscale('log')\n",
    "plt.yscale('log')\n",
    "plt.xlabel('$t_E$ (days)')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "## piE vs tE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, '$\\\\pi_E$')"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEOCAYAAACaQSCZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAQiklEQVR4nO3df2zc913H8derbjpOrLJFE6nKtSUZqbxWipBbrzCyiQ3GUgZWQwAxQEjTokVFGgOhRko0NE0gcCZPCE0MVYFVFULamKbUW+k2T1sGRWgTdeptaVW8lrJuuYC6drP5o56ahjd/3Nm52E78ufN973P3/T4fktX6e9+7e1+cb17+fn46IgQAwFauy10AAGA4EBgAgCQEBgAgCYEBAEhCYAAAkhAYAIAk1+cuoCg7d+6MPXv25C4DAIbK2bNnX4qIXZs9VtrA2LNnj+bn53OXAQBDxfYLV3uMJikAQBICAwCQpHSBYXvK9qnl5eXcpQBAqZQuMCLi0Yg4Ojo6mrsUACiV0gUGAKAYBAYAIAmBAQBIQmAAAJKUduIeBsvsQkMzc4u6sLSi3WM1HTs4rkMT9dxlAegAgYHCzS40dOL0Oa1cvCRJaiyt6MTpc5JEaABDhCYpFG5mbnEtLFatXLykmbnFTBUB6AaBgcJdWFrp6DiAwURgoHC7x2odHQcwmAgMFO7YwXHVdoxccay2Y0THDo5nqghAN+j0RuFWO7YZJQUMNwIDfXFook5AAEOOJikAQBICAwCQhMAAACQhMAAASQgMAEASAgMAkITAAAAkYR4GKo+l14E0BAYqrail1wkhlBFNUqi0IpZeXw2hxtKKQpdDaHahsc1qgbwIDFRaEUuvs/8Hyqp0gWF7yvap5eXl3KVgCBSx9Dr7f6CsShcYEfFoRBwdHR3NXQqGQBFLr7P/B8qqdIEBdOLQRF3Th/erPlaTJdXHapo+vH9bHdTs/4GyYpQUKq/XS6+z/wfKisAACsD+HygjmqQAAEm4wwAKxAQ+lAmBARSkqFnkQC40SQEFYQIfyobAAArCBD6UDYEBFIQJfCgbAgMoCBP4UDZ0egMFYQIfyobAAArEBD6UCU1SAIAkBAYAIAmBAQBIQmAAAJIQGACAJAQGACAJgQEASMI8DGyKZbkBrEdgYIMcy3ITUMDgo0kKG/R7We7VgGosrSh0OaBmFxqFvB+A7hAY2KDfy3KzbwQwHAgMbNDvZbnZNwIYDgQGNuj3stzsGwEMBwIDGxyaqGv68H7Vx2qypPpYTdOH9xfWCc2+EcBwYJQUNtXPZbnZNwIYDgQGBgL7RgCDj8DoAnMGAFQRgdGhHJPaAGAQ0OndIeYMAKgqAqNDzBkAUFUERoeYMwCgqgiMDjFnAEBVla7T2/aUpKl9+/YV8vrMGQBQVY6I3DUUYnJyMubn53OXAQBDxfbZiJjc7DGapAAASQgMAEASAgMAkITAAAAkITAAAEkIDABAEgIDAJCEwAAAJCEwAABJCAwAQBICAwCQhMAAACQhMAAASQgMAEASAgMAkITAAAAkITAAAEkIDABAEgIDAJCEwAAAJCEwAABJCAwAQBICAwCQ5PrcBQAYLrMLDc3MLerC0op2j9V07OC4Dk3Uc5eFPiAwACSbXWjoxOlzWrl4SZLUWFrRidPnJInQqACapAAkm5lbXAuLVSsXL2lmbjFTRegnAgNAsgtLKx0dR7kQGACS7R6rdXQc5UJgINnsQkMHTp7R3uOP6cDJM5pdaOQuCX127OC4ajtGrjhW2zGiYwfHM1WEfqLTG0no7IR0+WfNKKlqIjCQ5FqdnfxjUS2HJur8zCuKJikkobMTAIGBJHR2AiAwkITOTgD0YSAJnZ0ACAwko7MTqDYCA0AWLGI4fJL7MGzfXWQhAKpjdV5PY2lFocvzepgMOtg66fQ+Ikm2P9T6732FVASg9FjEcDh1M0rq5tZ/D/ayEADVwbye4dRJYNxj+28k3WH7p0X/B4AuMa9nOG0ZGLZfL0kRMSnpLyT9laRfl0TvFICuMK9nOKXcJSy0+iu+HRHnJZ23/UxEfKjg2oDCdTJSh1E9vcO8nuHkiLj2CfbLkr4l6Y2SXpL0jKS7ImJf8eV1b3JyMubn53OXgQG2fgVeqflb7vTh/Rv+4erkXGCY2T7balHaIOUO47sR8fbWC90i6Q5Jle6Z4jfNcuhkBV5W6wXSAmPU9s9Jenq1SargmgYa+0KURycjdRjVA6SNknq9pAck/bvt79j+gu2ZgusaWIwfL49ORuowqgdIC4w3R8ThiBhXsx/jg2r2aVQSv2mWRycjdRjVAyQ0SUXEf7b9/48kPdn6qqTdYzU1NgkHftMcPp2M1GFUz+Cjb7F4W46SGlZFjZJitAwweLgue+dao6TYQKlDhybqmj68X/WxmiypPlbjLyWQGX2L/cHyHl1gXwhgc7mahehb7A/uMAD0RM4lyxnF1h8EBoCeyNksxCi2/qBJCkBP5GwWYhRbfxAYAHoi95Bz+haLNzRNUrbfYPsTtj+TuxYAG9EslN/sQkMHTp7R3uOP6cDJMz3vP+pLYNh+yPaLtp9ad/xe24u2n7N9/FqvERHPR8SRYisF0C2GnOfVj0EH/WqSeljSX0v6+9UDtkckfVzSL6m5oOETtj8naUTS9LrnvzciXuxPqVhVtZmzVfu8RaBZKJ9+rKjcl8CIiMdt71l3+B5Jz0XE85Jk+1OS7ouIaUm/2s372D4q6agk3XbbbV3Xi+qtylu1z4vy6cegg5x9GHVJ32v7/ryuse2r7ZtsPyhpwvaJzc6JiFMRMRkRk7t27epttRVTtZmzVfu8KJ9+zEXJGRje5NhVF7aKiJcj4v6I+KnWXQgKVLWZs1X7vCiffgw6yDms9rykW9u+v0XShUy1YJ3cQySLcK0+ijJ+XlRLP+ai5AyMJyTdbnuvpIakd0v6nYz1oM2xg+Obrv45rEMkt+qjKNvnRTUVPeigX8NqPynpa5LGbZ+3fSQiXpP0fklzkp6R9OmIeLof9WBrZRsiuVUfRdk+L1AE9sNAJew9/timHWSW9F8nf6Xf5QADi/0wUHmsZgpsH4GBSmDZCmD7WHwQlcBqpsD2lS4wbE9Jmtq3b1/uUjBgWLYC2J7SNUlFxKMRcXR0dDR3KQBQKqULDABAMQgMAEASAgMAkITAAAAkITAAAEkIDABAEgIDAJCEwAAAJCldYNiesn1qeXk5dykAUCqlCwxmegNAMUoXGACAYhAYAIAkpVutFtszu9BgCXAAmyIwsGZ2oaETp8+t7X3dWFrRidPnJKkvoUFYAYONJimsmZlbXAuLVSsXL2lmbrHw914Nq8bSikKXw2p2oVH4ewNIQ2BgzYWllY6O91LOsAKQhsDAmt1jtY6O91LOsAKQhsDAmmMHx1XbMXLFsdqOER07OF74e+cMKwBpCAysOTRR1/Th/aqP1WRJ9bGapg/v70vHc6/CanahoQMnz2jv8cd04OQZ+kCAHmKUFK5waKKeZWTS6ntuZ5RU7lFeQNk5InLX0FO2pyRN7du3733PPvts7nLQRwdOnlFjkz6P+lhN/3b8FzJUBAwf22cjYnKzx0rXJMVaUtVFxzlQrNIFBqqLjnOgWAQGSiPnKC+gCuj0Rmn0ouMcwNURGCiVXKO8gCqgSQoAkITAAAAkITAAAEkIDABAEgIDAJCEwAAAJCEwAABJCAwAQJLSBYbtKdunlpeXc5cCAKVSusBgtVoAKEbpAgMAUAzWklpndqHB4nUAsAkCow1bfALA1dEk1WZmbnEtLFatXLykmbnFTBUBwODgDqMNW3wOPpoMgXy4w2jDFp+DbbXJsLG0otDlJsPZhUbu0oBKIDDasMXnYKPJEMiLJqk2bPE52GgyBPIiMNZhi8/BtXuspsYm4UCTIdAfNElhaNBkCOTFHQaGBk2GQF4EBoYKTYZAPjRJAQCSEBgAgCSlCwz2wwCAYpQuMNgPAwCKUbrAAAAUg8AAACQhMAAASQgMAEASAgMAkITAAAAkITAAAEkIDABAEgIDAJCEwAAAJCEwAABJCAwAQBICAwCQhMAAACQhMAAASQgMAEASAgMAkITAAAAkKV1gsKc3ABSjdIHBnt4AUIzSBQYAoBgEBgAgCYEBAEhCYAAAkhAYAIAkBAYAIAmBAQBI4ojIXUMhbH9f0gsdPGVUUqez/Tp9Tur5OyW91GEtZdfNz6ffctRY1Hv26nW38zrdPpfrcnt+MiJ2bfpIRPDVDM1TRT8n9XxJ87n/PAbtq5ufTxVqLOo9e/W623mdbp/LdVncF01Slz3ah+d08x5oGoY/uxw1FvWevXrd7bxOt8/luixIaZukhpnt+YiYzF0HgMu4Lun0HlSnchcAYIPKX5fcYQAAknCHAQBIQmAAAJIQGACAJATGELD9BtufsP2Z3LUAaLJ9yPbf2v6s7XfmrqcfCIxMbD9k+0XbT607fq/tRdvP2T4uSRHxfEQcyVMpUB0dXpezEfE+Se+R9FsZyu07AiOfhyXd237A9oikj0v6ZUl3Svpt23f2vzSgsh5W59fln7QeLz0CI5OIeFzSD9YdvkfSc607ilclfUrSfX0vDqioTq5LN31E0hci4sl+15oDgTFY6pK+1/b9eUl12zfZflDShO0TeUoDKmvT61LSH0h6h6TfsH1/jsL67frcBeAK3uRYRMTLkirxFxIYQFe7Lj8m6WP9LiYn7jAGy3lJt7Z9f4ukC5lqAdDEddlCYAyWJyTdbnuv7RskvVvS5zLXBFQd12ULgZGJ7U9K+pqkcdvnbR+JiNckvV/SnKRnJH06Ip7OWSdQJVyX18bigwCAJNxhAACSEBgAgCQEBgAgCYEBAEhCYAAAkhAYAIAkBAYAIAmBAQBIQmAACWzfYjtpkxzbNdv/0tpHYf1jH7b9QA/rusH247ZZSBSFIzCANL8o6a7Ec98r6XREXCqwHklSa3+Gr6giO74hLwID2ILtt0j6SzX3PfiG7b1bPOV3JX227fkfbG3v+WVJ4+tee9b2WdtP2z7aOvZntv+w7Zw/t/0B2z9u+zHb37T9VNsdz2zrPYFCsZYUkMD2FyU9EBFPbXHeDZK+GxE3t76/W81tP39Gzf1nnpT0YER8tPX4T0TED2zX1FwV9ecl3ajmHcpdtq+T9Kyau769TdK9rX2kZXs0IpZbTV//ExG7ev25gXbcYQBpxiUtSlLrjuDB1teb1p23U9JS2/dvlfRIRLwSEf+rjctif8D2NyV9Xc09F26PiO9Ietn2hKR3SlpobaJ1TtI7bH/E9lsjYlmSWk1fr9q+saefGFiHjjJgC7ZvkrQcERdt3yrpiYi42g6IK5J+bN2xTW/jbb9NzS0+3xwRr9j+57bn/p2k90i6WdJDkhQR327dsbxL0rTtL0XEn7bOf52kH3Xx8YBk3GEAW9uryzus3S3pjtbdxYfXnxgRP5Q0Ynv1H/7HJf1aa+TUjZKm2k4flfTDVli8UdLPtj32iKR7Jb1JzX0YZHu3pFci4h8kfVStTvhWoH0/Ii725NMCV8EdBrC1/5C00/ZTkv5b0h9FxMI1zv+SpLdI+nJEPGn7HyV9Q9ILkv617bwvSrrf9rfUbO76+uoDEfGq7a9KWmobbbVf0ozt/5N0UdLvt46/XdLnt/shga3Q6Q10wPbn1bzbeE3SXEQ8ssk5E5L+OCJ+bxvvc52aHeS/GRHPbnHuaUknImKx2/cDUnCHAXQgIt6VcM6C7a/aHulmLobtOyX9k5qd5VuFxQ2SZgkL9AN3GACAJHR6AwCSEBgAgCQEBgAgCYEBAEhCYAAAkhAYAIAkBAYAIAmBAQBI8v+i68dhZRxjkgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(t['t_E'], t['pi_E'])\n",
    "plt.xscale('log')\n",
    "plt.yscale('log')\n",
    "plt.xlabel('$t_E$ (days)')\n",
    "plt.ylabel('$\\pi_E$')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {}
   },
   "source": [
    "## dL vs dS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "pycharm": {}
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, '$d_S$ (kpc)')"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAU/klEQVR4nO3df5BdZ33f8fcHSQ4rTxjZtYBKxJXJwE4DNghE+eGEBgzIFA9WnUyKiwMFOu40GdqSVgkaOuUvgqciTTvtQEcDjh0gzoARCk1JZAYoLoSYyl4b+ZfCNPjXiuCljmBKl1gW3/6xd81qvZL2ru4959573q8Zj3Sfe/ae7x2v7uee53nO86SqkCR119PaLkCS1C6DQJI6ziCQpI4zCCSp4wwCSeq49W0XsBYXXHBBbdu2re0yJGms3H777d+rqs3L28cyCLZt28ahQ4faLkOSxkqSB1dqt2tIkjrOIJCkjjMIJKnjDAJJ6jiDQJI6ziCQpI4by+mjkjTJDszMsvfgEY4em2fLpil275xm1/atQzufQSBJI+TAzCx79h9m/vgJAGaPzbNn/2GAoYVBY11DSa5P8miSu5e1vzvJkST3JPn3TdUjSaNo78EjT4bAovnjJ9h78MjQztnkGMENwOVLG5K8BrgSuKSqXgB8qMF6JGnkHD0231f7IDQWBFV1K/DYsuZ/DlxXVX/TO+bRpuqRpFG0ZdNUX+2D0PasoecDv5DktiRfSfKyUx2Y5Nokh5Icmpuba7BESWrO7p3TTG1Yd1Lb1IZ17N45PbRzth0E64HzgFcAu4FPJclKB1bVvqraUVU7Nm9+yuJ5kjQRdm3fygevupitm6YIsHXTFB+86uKJnjX0CLC/qgr4RpIfAxcAfuWX1Fm7tm8d6gf/cm1fERwAXguQ5PnAOcD3Wq1IkjqmsSuCJDcBvwhckOQR4P3A9cD1vSmljwNv710dSJIa0lgQVNXVp3jqmqZqkCQ9VdtdQ5KklhkEktRxBoEkdZxBIEkdZxBIUscZBJLUcQaBJHWcQSBJHWcQSFLHGQSS1HEGgSR1nEEgSR3X9n4EkqQzODAzy96DRzh6bJ4tm6bYvXN6oPsVGASSNMIOzMyyZ/9h5o+fAGD22Dx79h8GGFgY2DUkSSNs78EjT4bAovnjJ9h78MjAzmEQSNIIO3psvq/2tTAIJGmEbdk01Vf7WhgEkjTCdu+cZmrDupPapjasY/fO6YGdw8FiSRphiwPCzhqSpA7btX3rQD/4l7NrSJI6ziCQpI4zCCSp4wwCSeo4g0CSOs4gkKSOMwgkqeMMAknqOINAkjqusSBIcn2SR5PcvcJz/yZJJbmgqXokSQuavCK4Abh8eWOSnwFeDzzUYC2SpJ7GgqCqbgUeW+Gp3wV+E6imapEk/USrYwRJ3gzMVtVdqzj22iSHkhyam5troDpJ6obWgiDJRuB9wL9bzfFVta+qdlTVjs2bNw+3OEnqkDavCH4WuAi4K8kDwHOAO5I8u8WaJKlzWtuPoKoOA89cfNwLgx1V9b22apKkLmpy+uhNwNeB6SSPJHlXU+eWJJ1aY1cEVXX1GZ7f1lApkqQlvLNYkjrOIJCkjjMIJKnjDAJJ6jiDQJI6ziCQpI4zCCSp4wwCSeo4g0CSOs4gkKSOMwgkqeNaW31UksbNgZlZ9h48wtFj82zZNMXundPs2r617bLOmkEgSatwYGaWPfsPM3/8BACzx+bZs/8wwNiHgV1DkrQKew8eeTIEFs0fP8Heg0daqmhwDAJJWoWjx+b7ah8nBoEkrcKWTVN9tY8Tg0CSVmH3zmmmNqw7qW1qwzp275xuqaLBcbBYklZhcUDYWUOS1GG7tm+diA/+5ewakqSOMwgkqeMMAknqOINAkjrOIJCkjjMIJKnjDAJJ6jiDQJI6ziCQpI4zCCSp4xoLgiTXJ3k0yd1L2vYmuT/JN5N8NsmmpuqRJC1o8orgBuDyZW1fAF5YVZcAfwHsabAeSRINBkFV3Qo8tqztlqp6ovfwz4HnNFWPJGnBKI0RvBP4k1M9meTaJIeSHJqbm2uwLEmabCMRBEneBzwBfPJUx1TVvqraUVU7Nm/e3FxxkjTh+t6PIMm5wI+q6sQZD17d670duAK4rKpqEK+pbjswMzuRm4dIw3LGIEjyNOAtwFuBlwF/A/xUkjng88C+qvrWWk6e5HLgt4C/X1X/by2vIS11YGaWPfsPM3984XvK7LF59uw/DGAYSKewmq6hLwM/y8KMnmdX1c9U1TOBX2BhgPe6JNec6UWS3AR8HZhO8kiSdwH/Bfhp4AtJ7kzyX9f6RiRY2EZwMQQWzR8/wd6DR1qqSBp9q+kael1VHV/eWFWPAZ8BPpNkw5lepKquXqH5Y6s4v7RqR4/N99UuaRVXBIshkOTGpTd8JTkvyfVLj5HatmXTVF/tkvqbNXRJVR1bfFBVfw1sH3xJ0trt3jnN1IZ1J7VNbVjH7p3TLVUkjb5+Zg09Lcl5vQAgyfl9/rw0dIsDws4aGk/O+GpHPx/kvwP8WZKbgQJ+BfjAUKqSzsKu7Vv98BhDzvhqz6q7hqrq94FfAr4LzAFXVdXHh1WYpG5xxld7+uraqap7gXuHVIukDnPGV3tWHQRJng78GvDzLHQNfRX4SFX9aEi1SeqQLZummF3hQ98ZX8PXz6yh3wdeAPxnFm4E+7uAXUOSBsIZX+3pp2touqpetOTxl5PcNeiCJHWTM77a008QzCR5RVX9OUCSlwNfG05ZkrrIGV/t6CcIXg68LclDvccXAvcmOQxUb5cxSdKY6ScIfpWFaaNLvRn43ODKkSQ1rZ/B4n3AM6rqwap6EHglcM2Sx5KkMdTPFcEvAzcneSsLU0jfBrxhKFVJkhqz6iCoqr9M8hbgAPAw8Iaq8k4PSRpzq9mh7DALN5AtOh9YB9yWBAeJJWm8reaK4IqhVyFJas1qguChM20qnyRuPC9J42lVexYneXeSC5c2JjknyWuT3Ai8fTjlSZKGbTVXBJcD7wRuSnIRcAx4OgvjBLcAv1tVdw6vREnSMJ0xCHqri34Y+HBvk/oLgPml21ZKksZXv/sRHAe+M6RaJEkt6OfOYknSBDIIJKnj+g6CJH97GIVIktqxliuCDwAkeWuSryV504BrkiQ1aC1BsDhb6A0sLD73DwdXjiSpaWsJgvVJ/i3wcO9u4h8OuCZJUoNWs+jcG4G7qupor+lfs3Al8GerfQ1J0uhazRXBLwGfS/Jwki8Avw1sBZ6XZF1V/fpqTpTk+iSPJrl7Sdv5Sb6Q5Fu9P89b07uQJK3ZGYOgqv5pVe0Afgf4C+DbwGuAbwD97Ex2AwvLVSz1XuCLVfU84Iu9x5KkBvXTrfOOqnrR4oMkHwZ2r/aHq+rWJNuWNV8J/GLv7zcC/wP4rT5qkiSdpX4Gi3+Q5KWLD6rqduD5Z3n+Z1XVd3qv9x3gmac6MMm1SQ4lOTQ3N3eWp5UkLerniuCdwCeS3AvcDlwMHB9KVSuoqn3APoAdO3a494EkDciqrwiq6lvAq4DPA88C7gP+wVme/7uLdyr3/nz0LF9PktSnflcfPQF8uvffIHyOhU1truv9+UcDel1J0io1tuhckpuArwPTSR5J8i4WAuD1Sb4FvL73WJLUoMZuBquqq0/x1GVN1TBODszMsvfgEY4em2fLpil275xm1/atbZclaQJ5V/AIOjAzy579h5k/fgKA2WPz7Nl/GMAwkDRw7kcwgvYePPJkCCyaP36CvQePtFSRpElmEIygo8fm+2qXpLNh19AI2rJpitkVPvS3bJpqoRqp27owXucVwQjavXOaqQ3rTmqb2rCO3TunB36uAzOzXHrdl7jovf+dS6/7EgdmZgd+DmlcLY7XzR6bp/jJeN2k/TsxCEbQru1b+eBVF7N10xQBtm6a4oNXXTzwbyFd+SWX1qor43V2DY2oXdu3Dv3y83S/5JN26SutRVfG6wyCDuvKL7km37D68bsyXmcQdFhXfsnHwSQPSA77vQ3zvpvdO6dPem0Y3nhdmxwj6LAmB6V1apM8VtPEextmP35T43Vt84qgwxZ/mSf1m+i4mOSxmibe2yC6OE931dLEeF3bDIKO68IveRv66Q6Z5LGaJt7b2XZxuqSLXUPSwPXbHXKqD6zF9nG+1+NM720QzraLsytTRE/HIJAGrN8PltN9kI37+EET41Bn248/yVdkq2XXkDRg/X6wnG6s5tLrvjTW4wdNjUOdTRens+cMAmng1vLBcqoPskn4tjrq41BdmSJ6OnYNSQM2yO6QJvrYB2GcxzG6MkX0dLwikAZskN0h4/BtdRJm3Yz6VcuwGQTSEAzqg2Uc7vUYxr0Ck3yn9SgyCKQRN+rfVgc9jjEJVxjjxjECSWdl08YNK7avdRzDef3NMwgkrdmBmVn+74+eeEr7hnVZ8zjGJMyUGjcGgaQ123vwCMd/XE9pP/ec9Wc1r7+fdp09g0DSmp3qW/r354+f9udON93UVXGb52CxJGBtM3XWcvPcmQaDx2Gm1KQxCBrgVDiNurXO1FnLfQ6rmW466jOlJo1dQ0M27ouGqRvWOlNnLXflOhg8erwiGLJJ3nREk+NsPpz7/fbuIm+jxyuCIfPbj8ZBkzN1HAwePSMRBEnek+SeJHcnuSnJ09uuaVCcCqdx0OSHs4u8jZ7Wu4aSbAX+BfBzVTWf5FPAW4AbWi1sQMZh0TCp6Zk6DgaPltaDoGc9MJXkOLARONpyPQPjVDiNCz+cu6v1IKiq2SQfAh4C5oFbquqW5ccluRa4FuDCCy9stsiz5D8wSaOs9TGCJOcBVwIXAVuAc5Ncs/y4qtpXVTuqasfmzZubLlOSJlbrQQC8Dvh2Vc1V1XFgP/CqlmuSpM4YhSB4CHhFko1JAlwG3NdyTZLUGaMwRnBbkpuBO4AngBlgX7tVSd3msijd0noQAFTV+4H3t12HJHcI66JR6BqSNELcIax7DAJJJ3FZlO4xCCSdxGVRuscgkHQSF4XrnpEYLNbkcdbJ+HJZlO4xCDRwzjoZfy6L0i12DWngnHUijReDQAPnrBNpvBgEGjhnnUjjxSDQwDnrRBovDhZr4Jx1Io0Xg0BD4awTaXzYNSRJHWcQSFLHGQSS1HEGgSR1nEEgSR1nEEhSxzl9VI1yVVJp9BgEaoyrkkqjya4hNcZVSaXRZBCoMa5KKo0mg0CNcVVSaTQZBGqMq5JKo8nBYjXGVUml0WQQqFGuSiqNHruGJKnjDAJJ6jiDQJI6biSCIMmmJDcnuT/JfUle2XZNktQVozJY/J+AP62qX05yDrCx7YIkqStaD4IkzwBeDfwTgKp6HHi8zZokqUtGoWvoucAc8HtJZpJ8NMm5yw9Kcm2SQ0kOzc3NNV+lJE2oUQiC9cBLgI9U1Xbgh8B7lx9UVfuqakdV7di8eXPTNUrSxGq9awh4BHikqm7rPb6ZFYJgEFwLX5KeqvUrgqr6K+DhJIsLzlwG3Dvo8yyuhT97bJ7iJ2vhH5iZHfSpJGmstB4EPe8GPpnkm8CLgd8e9AlcC1+SVjYKXUNU1Z3AjmGew7XwJWllo3JFMHSuhS9JK+tMELgWviStbCS6hprgWviStLLOBAG4Fr4kraQzXUOSpJUZBJLUcQaBJHWcQSBJHWcQSFLHpararqFvSeaAB9uuY8AuAL7XdhFD5nucDJP+Hif5/f2dqnrK8s1jGQSTKMmhqhrqMhtt8z1Ohkl/j5P+/lZi15AkdZxBIEkdZxCMjn1tF9AA3+NkmPT3OOnv7ykcI5CkjvOKQJI6ziCQpI4zCEZEknVJZpL8cdu1DEOSB5IcTnJnkkNt1zNoSTYluTnJ/UnuS/LKtmsapCTTvf93i//9IMm/aruuQUvyniT3JLk7yU1Jnt52TU1wjGBEJPkNFrbrfEZVXdF2PYOW5AFgR1VN5I06SW4E/mdVfTTJOcDGqjrWdl3DkGQdMAu8vKom5sbOJFuBrwI/V1XzST4FfL6qbmi3suHzimAEJHkO8Cbgo23Xov4leQbwauBjAFX1+KSGQM9lwP+epBBYYj0wlWQ9sBE42nI9jTAIRsN/BH4T+HHbhQxRAbckuT3JtW0XM2DPBeaA3+t17300ybltFzVEbwFuaruIQauqWeBDwEPAd4DvV9Ut7VbVDIOgZUmuAB6tqtvbrmXILq2qlwBvBH49yavbLmiA1gMvAT5SVduBHwLvbbek4eh1e70Z+HTbtQxakvOAK4GLgC3AuUmuabeqZhgE7bsUeHOvD/0Pgdcm+US7JQ1eVR3t/fko8Fng77Vb0UA9AjxSVbf1Ht/MQjBMojcCd1TVd9suZAheB3y7quaq6jiwH3hVyzU1wiBoWVXtqarnVNU2Fi65v1RVE/UtJMm5SX568e/AG4C7261qcKrqr4CHk0z3mi4D7m2xpGG6mgnsFup5CHhFko1JwsL/x/tarqkRndq8Xq15FvDZhX9brAf+oKr+tN2SBu7dwCd7XSd/Cbyj5XoGLslG4PXAP2u7lmGoqtuS3AzcATwBzNCR5SacPipJHWfXkCR1nEEgSR1nEEhSxxkEktRxBoEkdZxBIEkdZxBIUscZBNIKknwsyZvOcMxUkq/09pLYluSs75ZOck6SW3urX0qNMAiklb0YuOsMx7wT2F9VJwZ10qp6HPgi8I8G9ZrSmRgEEpDk+Um+2ttF7T3As6vqkTP82FuBP1rhtZ7bW476Zb0rhfuT3Jjkm71dzDb2jntbr+2uJB9f8hIHeq8tNcIgUOf1umE+AfxGVV0MPA+4/ww/cw7w3Kp6YFn7NPAZ4B1V9b96zdPAvqq6BPgB8GtJXgC8D3htVb0I+JdLXuZu4GVn/cakVTIIJLgKuK+qvtF7fA9wJ0CS/3aKn7kAWL4L2WYWrhCuqao7l7Q/XFVf6/39E8DPA68Fbl7curOqHls8uNfV9Pjiiq3SsBkEElwCLN0Y6KXAXUm2AQ+c4mfmgeUbm38feJiFPSaWWr6yYwFZoX2pnwJ+dJrnpYExCCT4P8ALAZK8lIU19+9kYXOZO1b6gar6a2BdkqVh8DiwC3hbkn+8pP3CJK/s/f1qFjZI/yLwK0n+Vu+85y8e3Gtb3BxFGjqDQIKPAy9OcicLe0cfY2FDkpdy8pXCcrew0M3zpKr6IXAF8J4kV/aa7wPenuSbwPksbGl5D/AB4CtJ7gL+w5KXeQ3w+bN+V9IquR+BdApJ/gR4kIVNSj5dVV9Z9vx2FgaYf/U0r7EN+OOqemEf590P7KmqI2upW+qXN61Ip1BVbzzD8zNJvpxk3aDuJejNRjpgCKhJXhFIUsc5RiBJHWcQSFLHGQSS1HEGgSR1nEEgSR1nEEhSxxkEktRx/x9l4HrJlS//JwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(t['rad_L'], t['rad_S'])\n",
    "plt.xlabel('$d_L$ (kpc)')\n",
    "plt.ylabel('$d_S$ (kpc)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {}
   },
   "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.0"
  },
  "pycharm": {
   "stem_cell": {
    "cell_type": "raw",
    "metadata": {
     "collapsed": false
    },
    "source": []
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}