TOIs_eRV_RUWE_distribution.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "if 'google.colab' in str(get_ipython()):\n",
    "    !pip install astroquery"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Get list of TOIs that are 1) confirmed , and 2) giant planets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "61\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>TIC ID</th>\n",
       "      <th>RA</th>\n",
       "      <th>Dec</th>\n",
       "      <th>TOI</th>\n",
       "      <th>TFOPWG Disposition</th>\n",
       "      <th>Planet Radius (R_Earth)</th>\n",
       "      <th>Period (days)</th>\n",
       "      <th>Predicted Mass (M_Earth)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>394137592</td>\n",
       "      <td>00:17:47.32</td>\n",
       "      <td>-66:21:32.14</td>\n",
       "      <td>120.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>12.215808</td>\n",
       "      <td>11.535076</td>\n",
       "      <td>90.88</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>290131778</td>\n",
       "      <td>21:18:47.9</td>\n",
       "      <td>-26:36:58.95</td>\n",
       "      <td>123.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>17.215002</td>\n",
       "      <td>3.308877</td>\n",
       "      <td>317.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>37</th>\n",
       "      <td>267263253</td>\n",
       "      <td>00:29:18.94</td>\n",
       "      <td>-76:18:14.52</td>\n",
       "      <td>135.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>16.980000</td>\n",
       "      <td>4.126910</td>\n",
       "      <td>317.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>52</th>\n",
       "      <td>271893367</td>\n",
       "      <td>07:31:51.7</td>\n",
       "      <td>-73:36:21.73</td>\n",
       "      <td>150.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>13.889276</td>\n",
       "      <td>5.857408</td>\n",
       "      <td>125.80</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>59</th>\n",
       "      <td>140691463</td>\n",
       "      <td>04:54:48.35</td>\n",
       "      <td>-76:40:50.18</td>\n",
       "      <td>157.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>13.043500</td>\n",
       "      <td>2.084540</td>\n",
       "      <td>112.45</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4218</th>\n",
       "      <td>257060897</td>\n",
       "      <td>15:10:07.72</td>\n",
       "      <td>72:42:37.13</td>\n",
       "      <td>4138.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>17.547500</td>\n",
       "      <td>3.660038</td>\n",
       "      <td>317.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4415</th>\n",
       "      <td>256722647</td>\n",
       "      <td>21:58:15.34</td>\n",
       "      <td>71:07:53.78</td>\n",
       "      <td>4329.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>21.536300</td>\n",
       "      <td>2.922492</td>\n",
       "      <td>317.00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4494</th>\n",
       "      <td>206541859</td>\n",
       "      <td>01:12:11.64</td>\n",
       "      <td>-56:55:31.39</td>\n",
       "      <td>4406.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>12.015503</td>\n",
       "      <td>30.081780</td>\n",
       "      <td>98.34</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5253</th>\n",
       "      <td>257527578</td>\n",
       "      <td>03:05:10.23</td>\n",
       "      <td>-21:56:01.13</td>\n",
       "      <td>5152.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>11.634800</td>\n",
       "      <td>54.321242</td>\n",
       "      <td>92.86</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5254</th>\n",
       "      <td>124029677</td>\n",
       "      <td>06:06:09.68</td>\n",
       "      <td>-19:57:12.4</td>\n",
       "      <td>5153.01</td>\n",
       "      <td>CP</td>\n",
       "      <td>11.354700</td>\n",
       "      <td>20.910841</td>\n",
       "      <td>89.82</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>61 rows × 8 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "         TIC ID           RA           Dec      TOI TFOPWG Disposition  \\\n",
       "20    394137592  00:17:47.32  -66:21:32.14   120.01                 CP   \n",
       "23    290131778   21:18:47.9  -26:36:58.95   123.01                 CP   \n",
       "37    267263253  00:29:18.94  -76:18:14.52   135.01                 CP   \n",
       "52    271893367   07:31:51.7  -73:36:21.73   150.01                 CP   \n",
       "59    140691463  04:54:48.35  -76:40:50.18   157.01                 CP   \n",
       "...         ...          ...           ...      ...                ...   \n",
       "4218  257060897  15:10:07.72   72:42:37.13  4138.01                 CP   \n",
       "4415  256722647  21:58:15.34   71:07:53.78  4329.01                 CP   \n",
       "4494  206541859  01:12:11.64  -56:55:31.39  4406.01                 CP   \n",
       "5253  257527578  03:05:10.23  -21:56:01.13  5152.01                 CP   \n",
       "5254  124029677  06:06:09.68   -19:57:12.4  5153.01                 CP   \n",
       "\n",
       "      Planet Radius (R_Earth)  Period (days)  Predicted Mass (M_Earth)  \n",
       "20                  12.215808      11.535076                     90.88  \n",
       "23                  17.215002       3.308877                    317.00  \n",
       "37                  16.980000       4.126910                    317.00  \n",
       "52                  13.889276       5.857408                    125.80  \n",
       "59                  13.043500       2.084540                    112.45  \n",
       "...                       ...            ...                       ...  \n",
       "4218                17.547500       3.660038                    317.00  \n",
       "4415                21.536300       2.922492                    317.00  \n",
       "4494                12.015503      30.081780                     98.34  \n",
       "5253                11.634800      54.321242                     92.86  \n",
       "5254                11.354700      20.910841                     89.82  \n",
       "\n",
       "[61 rows x 8 columns]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from astropy.utils.data import download_file\n",
    "\n",
    "# Basically getting  the list of TOIs from \n",
    "# - it can be slow, so I download it locally\n",
    "toi_url =  \"https://exofop.ipac.caltech.edu/tess/download_toi.php?sort=toi&output=csv\"\n",
    "toi_local_path = download_file(toi_url, cache=True, timeout=120)  \n",
    "df = pd.read_csv(toi_local_path)\n",
    "\n",
    "disp_codes, min_radius, plot_title = [\"CP\"], 11.2, \"TESS Confirmed Giant Planets\"\n",
    "# disp_codes, min_radius, plot_title = [\"CP\", \"PC\"], 11.2, \"TESS Confirmed or Candidate Giant Planets\"\n",
    "# disp_codes, min_radius, plot_title = [\"KP\"], 11.2, \"Known Giant Planets in TOI\"\n",
    "# disp_codes, min_radius, plot_title = [\"FP\", \"FA\",], 11.2, \"FP/FA Giant Planets in TOI\"\n",
    "# disp_codes, min_radius, plot_title = [\"APC\"], 11.2, \"APC Giant Planets in TOI (usually deemed to be EB by SG1)\"\n",
    "# disp_codes, min_radius, plot_title = [\"FP\", \"FA\", \"APC\"], 11.2, \"FP/FA/APC Giant Planets in TOI\"\n",
    "\n",
    "# disp_codes, min_radius, plot_title = [\"CP\"], 0, \"TESS Confirmed Planets\"\n",
    "# disp_codes, min_radius, plot_title = [\"CP\", \"PC\"], 0, \"TESS Confirmed or Candidate Planets\"\n",
    "# disp_codes, min_radius, plot_title = [\"KP\"], 0, \"Known Planets in TOI\"\n",
    "# disp_codes, min_radius, plot_title = [\"FP\", \"FA\",], 0, \"FP/FA Planets in TOI\"\n",
    "# disp_codes, min_radius, plot_title = [\"APC\"], 0, \"APC Planets in TOI (usually deemed to be EB by SG1)\"\n",
    "# disp_codes, min_radius, plot_title = [\"FP\", \"FA\", \"APC\"], 0, \"FP/FA/APC Planets in TOI\"\n",
    "\n",
    "\n",
    "df = df[df[\"TFOPWG Disposition\"].isin(disp_codes)]  \n",
    "df = df[df[\"Planet Radius (R_Earth)\"] > min_radius]  # ~ > 1 R_jupiter\n",
    "print(len(df))\n",
    "display(df[[\"TIC ID\", \"RA\", \"Dec\", \"TOI\", \"TFOPWG Disposition\", \"Planet Radius (R_Earth)\", \"Period (days)\", \"Predicted Mass (M_Earth)\"]])\n",
    "\n",
    "df_src = df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "61\n"
     ]
    }
   ],
   "source": [
    "# Use TIC catalog to get the RA/ DEC in J2000\n",
    "\n",
    "from astroquery.mast import Catalogs\n",
    "\n",
    "df_tics = Catalogs.query_criteria(catalog=\"Tic\", ID=df_src[\"TIC ID\"]).to_pandas()\n",
    "print(len(df_tics))\n",
    "# df_tics[[\"ID\", \"ra\", \"dec\"]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "61\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING: OverflowError converting to IntType in column SolID, reverting to String. [astropy.io.ascii.core]\n",
      "WARNING: OverflowError converting to IntType in column Source, reverting to String. [astropy.io.ascii.core]\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>angDist</th>\n",
       "      <th>TIC_ID</th>\n",
       "      <th>TIC_RA</th>\n",
       "      <th>TIC_DEC</th>\n",
       "      <th>DR3Name</th>\n",
       "      <th>RAdeg</th>\n",
       "      <th>DEdeg</th>\n",
       "      <th>errHalfMaj</th>\n",
       "      <th>errHalfMin</th>\n",
       "      <th>errPosAng</th>\n",
       "      <th>...</th>\n",
       "      <th>B_AG</th>\n",
       "      <th>E(BP-RP)</th>\n",
       "      <th>b_E(BP-RP)</th>\n",
       "      <th>B_E(BP-RP)</th>\n",
       "      <th>Lib</th>\n",
       "      <th>RAJ2000</th>\n",
       "      <th>DEJ2000</th>\n",
       "      <th>e_RAJ2000</th>\n",
       "      <th>e_DEJ2000</th>\n",
       "      <th>RADEcorJ2000</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.415850</td>\n",
       "      <td>11561667</td>\n",
       "      <td>149.823787</td>\n",
       "      <td>-27.392782</td>\n",
       "      <td>Gaia DR3 5466556141521710592</td>\n",
       "      <td>149.823663</td>\n",
       "      <td>-27.392819</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>161</td>\n",
       "      <td>...</td>\n",
       "      <td>0.4571</td>\n",
       "      <td>0.2383</td>\n",
       "      <td>0.2346</td>\n",
       "      <td>0.2425</td>\n",
       "      <td>MARCS</td>\n",
       "      <td>149.823786</td>\n",
       "      <td>-27.392781</td>\n",
       "      <td>0.698193</td>\n",
       "      <td>0.718706</td>\n",
       "      <td>-0.1825</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.423098</td>\n",
       "      <td>15692883</td>\n",
       "      <td>17.474855</td>\n",
       "      <td>25.681695</td>\n",
       "      <td>Gaia DR3 306410392895767680</td>\n",
       "      <td>17.474736</td>\n",
       "      <td>25.681648</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>87</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>17.474857</td>\n",
       "      <td>25.681696</td>\n",
       "      <td>0.339290</td>\n",
       "      <td>0.233187</td>\n",
       "      <td>0.2832</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.673149</td>\n",
       "      <td>24358417</td>\n",
       "      <td>81.344384</td>\n",
       "      <td>-34.668327</td>\n",
       "      <td>Gaia DR3 4822916912388338176</td>\n",
       "      <td>81.344589</td>\n",
       "      <td>-34.668247</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>165</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0325</td>\n",
       "      <td>0.0044</td>\n",
       "      <td>0.0007</td>\n",
       "      <td>0.0172</td>\n",
       "      <td>MARCS</td>\n",
       "      <td>81.344383</td>\n",
       "      <td>-34.668326</td>\n",
       "      <td>0.165260</td>\n",
       "      <td>0.175828</td>\n",
       "      <td>-0.1996</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.135939</td>\n",
       "      <td>48103627</td>\n",
       "      <td>147.406744</td>\n",
       "      <td>-33.218506</td>\n",
       "      <td>Gaia DR3 5438787082011390208</td>\n",
       "      <td>147.406699</td>\n",
       "      <td>-33.218502</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>9</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0052</td>\n",
       "      <td>0.0013</td>\n",
       "      <td>0.0003</td>\n",
       "      <td>0.0028</td>\n",
       "      <td>MARCS</td>\n",
       "      <td>147.406744</td>\n",
       "      <td>-33.218507</td>\n",
       "      <td>0.146376</td>\n",
       "      <td>0.169321</td>\n",
       "      <td>0.2822</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0.036756</td>\n",
       "      <td>65412605</td>\n",
       "      <td>108.123125</td>\n",
       "      <td>-24.953562</td>\n",
       "      <td>Gaia DR3 5617241426979996800</td>\n",
       "      <td>108.123114</td>\n",
       "      <td>-24.953560</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>6</td>\n",
       "      <td>...</td>\n",
       "      <td>0.3792</td>\n",
       "      <td>0.1986</td>\n",
       "      <td>0.1932</td>\n",
       "      <td>0.2014</td>\n",
       "      <td>A</td>\n",
       "      <td>108.123125</td>\n",
       "      <td>-24.953563</td>\n",
       "      <td>0.212685</td>\n",
       "      <td>0.301223</td>\n",
       "      <td>0.0497</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 141 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "    angDist    TIC_ID      TIC_RA    TIC_DEC                       DR3Name  \\\n",
       "0  0.415850  11561667  149.823787 -27.392782  Gaia DR3 5466556141521710592   \n",
       "1  0.423098  15692883   17.474855  25.681695   Gaia DR3 306410392895767680   \n",
       "2  0.673149  24358417   81.344384 -34.668327  Gaia DR3 4822916912388338176   \n",
       "3  0.135939  48103627  147.406744 -33.218506  Gaia DR3 5438787082011390208   \n",
       "4  0.036756  65412605  108.123125 -24.953562  Gaia DR3 5617241426979996800   \n",
       "\n",
       "        RAdeg      DEdeg  errHalfMaj  errHalfMin  errPosAng  ...    B_AG  \\\n",
       "0  149.823663 -27.392819         0.0         0.0        161  ...  0.4571   \n",
       "1   17.474736  25.681648         0.0         0.0         87  ...     NaN   \n",
       "2   81.344589 -34.668247         0.0         0.0        165  ...  0.0325   \n",
       "3  147.406699 -33.218502         0.0         0.0          9  ...  0.0052   \n",
       "4  108.123114 -24.953560         0.0         0.0          6  ...  0.3792   \n",
       "\n",
       "  E(BP-RP)  b_E(BP-RP)  B_E(BP-RP)    Lib     RAJ2000    DEJ2000  e_RAJ2000  \\\n",
       "0   0.2383      0.2346      0.2425  MARCS  149.823786 -27.392781   0.698193   \n",
       "1      NaN         NaN         NaN    NaN   17.474857  25.681696   0.339290   \n",
       "2   0.0044      0.0007      0.0172  MARCS   81.344383 -34.668326   0.165260   \n",
       "3   0.0013      0.0003      0.0028  MARCS  147.406744 -33.218507   0.146376   \n",
       "4   0.1986      0.1932      0.2014      A  108.123125 -24.953563   0.212685   \n",
       "\n",
       "   e_DEJ2000  RADEcorJ2000  \n",
       "0   0.718706       -0.1825  \n",
       "1   0.233187        0.2832  \n",
       "2   0.175828       -0.1996  \n",
       "3   0.169321        0.2822  \n",
       "4   0.301223        0.0497  \n",
       "\n",
       "[5 rows x 141 columns]"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Use the RA / DEC in J2000 to do crossmatch with Gaia DR3\n",
    "\n",
    "from astropy.table import Table\n",
    "from astroquery.xmatch import XMatch\n",
    "from astropy import units as u\n",
    "\n",
    "def _xmatch_vizier_of_tics(df_tics_ra_dec, catalog, max_distance):\n",
    "    \"\"\"Return the Vizier Catalog records of the given list of TICs by coordinate match using Vizier Crossmatch.\"\"\"\n",
    "\n",
    "    src_tab = Table.from_pandas(df_tics_ra_dec[[\"ID\", \"ra\", \"dec\"]])\n",
    "    src_tab.rename_column(\"ID\", \"TIC_ID\")\n",
    "    src_tab.rename_column(\"ra\", \"TIC_RA\")\n",
    "    src_tab.rename_column(\"dec\", \"TIC_DEC\")\n",
    "\n",
    "    catalog = f\"vizier:{catalog}\"\n",
    "    return XMatch.query(cat1=src_tab, cat2=catalog, max_distance=max_distance, colRA1=\"TIC_RA\", colDec1=\"TIC_DEC\")\n",
    "\n",
    "df_xmatch = _xmatch_vizier_of_tics(df_tics, \"I/355/gaiadr3\", 10 * u.arcsec).to_pandas()\n",
    "\n",
    "\n",
    "# pick the closest one by coordinate from xmatch result as a rough approximation\n",
    "df_xmatch.sort_values([\"TIC_ID\", \"angDist\"], ascending=[True, True], inplace=True, ignore_index=True)\n",
    "df_gaia = df_xmatch.groupby(\"TIC_ID\").head(1).reset_index(drop=True)\n",
    "\n",
    "print(len(df_gaia))\n",
    "df_gaia.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA2wAAAEeCAYAAAD2GZayAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABBlklEQVR4nO3deXyM5/7/8fckIREkRKhWRJFSLa1qal+LcqpqaUQsIbWcQ+2tWtI0tVWstbUEp6QNbRCpKj16WmtriaVFORQHLbFHWickkcT8/sjPfMUSCZm5Z+T1fDzmMe65r/u+Pvcgn/nkuu5rTGaz2SwAAAAAgN1xMjoAAAAAAMDdUbABAAAAgJ2iYAMAAAAAO0XBBgAAAAB2ioINAAAAAOwUBRsAAAAA2CkXowMAgEfBhAkTtGvXLknSf//7X5UrV05ubm6SpGXLlqlv375KSEhQ8eLFsx3Xv39/tW7dWhcuXNDEiRP13//+V5Lk5uamf/zjH2rRooUk3Xf/7c6fP68ZM2bo4MGDMplMcnV1zbF9bhw6dEiDBg2Sh4eH2rdvrz/++ENhYWEPfL7cunz5surVq6fffvvtrvt37dql+fPn69SpUzKZTHJzc9Obb76pdu3aSZJmzZqlChUqqH379g/U/6ZNm7Rv3z4NGTLkjn0vv/yyChUqJDc3N5lMJqWnp6tBgwYaNWqUzpw5o7Zt2+qXX355oH5zcurUKU2ZMkVz5szJ03Hnz5/XBx98oHnz5unGjRuaO3euNmzYoGvXrqlJkyYaPXq0TCaTpf2sWbP0119/KTw8XJJ04MABLVu2TOPHj8/X6wEA5MAMAMhXzZo1M+/fvz/ba927dzf/61//uucxffv2NS9evNiyffToUbO/v7/52LFjudp/q8TERHPTpk3NX331lfnGjRtms9lsPnTokLlu3brmn3766YGva86cOebQ0NAHPv5BJSYmmqtUqXLXfZs2bTI3btzYvGvXLstrp06dMr/yyivm7777Ll/6nz17tnns2LF33Xf733VaWpo5MDDQHB0dbT516pS5Zs2a+RLD7Xbs2GFu06ZNno/r06ePed++fWaz2WxetGiRuXv37uaUlBRzWlqauVOnTuY1a9aYzWaz+ezZs+ZBgwaZn3/++TuufdSoUeYNGzY8/EUAAHKFETYAsAMXL15Uamqqbty4IScnJ/n5+WnevHny8PDI1f5bffHFF6pVq1a2EaWnn35as2fPtrTfvXu3pkyZopSUFBUqVEhDhw5V48aNFRcXp++//15OTk76/fff5ebmpsmTJ+vgwYP68ssvlZmZqdTUVDVo0EDfffed5s+fr+DgYHl6eur48ePq0qWL/v3vf+vZZ5/V3r17dfnyZQUGBurSpUvauXOnUlJSNHPmTFWtWlX/+9//9OGHH+rIkSNKT09XvXr1NGLECLm4uOjf//63ZsyYoSJFiqh69er3fN+mTZum0aNHy9/f3/Kaj4+PPvzwQ127dk2SNGrUKD311FPq3bu3YmNjtWzZMqWnp+uvv/5S37591bVr13ted3JysmJiYpSZmanixYtr2LBhOf49Fi5cWC+++KKOHz+upk2bWl6/dOmSwsPDlZiYqIsXL6pcuXKaOXOmSpUqpZdfflkdOnTQ9u3bdfbsWbVr105Dhw6VJG3YsEHz5s1Tenq63NzcNHLkSD333HMKCwvT+fPn1bt3b82fP1/jx4/Xzz//rEKFCsnHx0cREREqWrRottj27dunxMREPffcc5KkVatWaeTIkZaR4Dlz5qhQoUKSpNjYWNWuXVuVK1fWX3/9le08nTt31pgxY9SsWbMc3wsAQP7gHjYAsJEpU6aoXbt22R5JSUmSpBEjRmjJkiWqV6+e+vfvr3/+858qX768Spcunav9tzpw4IBq1ap1x+svvfSSqlatqqSkJA0ePFjvvfeevvnmG02ePFnvvvuuTp06JSlriuH777+vNWvW6Pnnn9eCBQv0+uuvKygoSK+++qqmT59+x7k9PDz07bffKjg4WJKUkJCgmJgYTZ06VVOnTlXt2rUVFxenRo0aacmSJZKkiRMn6tlnn1VcXJxWrVqlpKQkLV68WJcuXVJoaKjmzJmjuLg4lStX7q7v55UrV3TkyBE1bNjwjn3+/v5q3LhxtteuXr2qFStWaMGCBVq1apVmzJihqVOnWvbf7bqff/55y3Xfr1iTsqYcbty4UXXq1Mn2+tq1a1WzZk0tW7ZM69evl5ubm77++mvL/mvXrumLL75QTEyMFi1apFOnTunkyZOaMWOGJd7x48dr0KBBSktL04QJE+Tr66tPP/1Ue/fu1c6dO7V69WrFxcWpfPnyd50+um7dumxF1smTJ3Xs2DH17NlTbdu21RdffCFPT09J0sCBA9W9e3c5Od35MaFmzZr6448/LP9eAADWxQgbANjIiBEj1Lp167vuq1evnjZt2qS9e/dq9+7d2rhxoz755BN99tlneu655+67/1Ymk0lms/mecezfv1++vr56/vnnJUlPPfWUatWqpZ07d8pkMunZZ59V2bJlJUnPPPOMvv/++/te260jXJLUsmVLSVL58uUlSY0aNZIk+fr6aufOnZKy7g379ddfFRsbK0lKTU2VJO3Zs0dVqlSRn5+fpKwRnY8++uiOPm9e4633XA0dOlQnTpxQenq6SpUqpejoaMu+okWLKjIyUps3b9bJkyd1+PBhyyicpAe6bkkaPny43NzcdOPGDRUqVEidOnVSq1atdPr0aUubnj17avfu3Vq8eLFOnjypo0ePWt5/SWrevLkk6bHHHlOpUqX0119/ad++fbpw4YJCQkIs7Uwmk/74449s/VepUkXOzs7q1KmTGjZsqFatWt3xb0KSjh8/rldffdWynZGRoX379mnhwoW6fv26+vfvr+jo6Gz93YuPj49OnDhh+fsFAFgPBRsAGCwxMVFz5szR+++/L39/f/n7+6tfv3567733tGrVKpUrVy7H/bd/OK9Zs6b27t2r7t27Z3s9JiZGKSkpqlChQrYiR8oqfjIyMiwLaNx0v+LvJnd392zbhQsXzrZ9c6rdrW7cuKFZs2apcuXKkrJGzEwmk7Zt25atTxeXu6cqT09PVa5cWTt37rSMHM2cOVOSFB8ff8fCGOfOnVPnzp0VGBioF198Ua1bt9bGjRst+x/kuqWsaZk1atTIsc3UqVO1f/9+vfHGG6pTp44yMjKynd/V1fWOvm/cuKF69epZrkmSzp49qzJlymj37t2W1zw8PPT111/r559/1o4dOzR06FD17t1b3bp1yxbD7ddUpkwZtWnTRoULF1bhwoXVunVry8I59+Pi4iJnZ+dctQUAPBymRAKAwTw9PbVt2zZ9/vnnlg/UKSkp+uOPP/TMM8/cd//tOnfubJkid7P9gQMHNHv2bFWpUkU1a9bU8ePHtX//fknS0aNHtWvXLtWuXdtGV5ylYcOGioqKktlstozwLFmyRC+99JKOHTumw4cPS5Li4uLueY5Ro0ZpwoQJ+vnnny2vJScna9OmTXdM5ztw4IC8vLz01ltvqWHDhpZiLTMzM8c4nZ2dlZGR8aCXKUn66aef1LNnT7Vv316lSpXStm3b7ttvvXr1tHXrVsvKoJs3b9brr7+u1NRUOTs7Kz09XZK0ceNGhYSE6IUXXtCgQYPUvn17HThw4I7zVaxYMdvoXKtWrbR69WrduHFD6enp2rhx430LTymruD9z5owqVqyYl7cAAPCAGGEDABuZMmWK5s2bl+21li1bauDAgfr00081depURUdHy93dXSaTSR06dFBAQIAk3Xf/rUqUKKHo6GhNnTpV8+fPl5OTk4oUKaIPP/xQDRo0kJS1XPv48eOVmpoqk8mkiIgIVaxY0SpL0N/Le++9pw8//FBt27ZVenq66tevrz59+qhQoUKaNm2ahg8frkKFCumll1665zkaN26sjz76SPPmzVNCQoLS09NlNpvVuHFjzZ8/P1vbBg0aKDY2Vq1bt5bJZFLt2rXl5eWl33//Pcc469atq+HDh2v8+PF6//33H+haBwwYoClTpmjWrFkqVKiQatWqdcfUxtv5+flp3Lhxevvtt2U2m+Xi4qJ58+apaNGi8vPzk6urqwICArRs2TJt2bJFr732mtzd3eXp6XnXZfdbtWqlDz/8UIMHD5aUNX102rRpeu2115SZman69eurZ8+e972WX3/9Vb6+vnriiSce6L0AAOSNyZzbOR8AAMCh9e7dW0OGDLnrPW65NWrUKLVu3TrbKpgAAOthSiQAAAXE2LFj9cknn+T6/rzbHThwQCaTiWINAGyIETYAAAAAsFOMsAEAAACAnaJgAwAAAAA7RcEGAAAAAHaKgg0AAAAA7BQFGwAAAADYKQo2AAAAALBTFGwAAAAAYKco2AAAAADATlGwAQAAAICdomADAAAAADtFwQYAAAAAdoqCDQAAAADsFAUbAAAAANgpCjYAAAAAsFMUbAAAAABgpyjYAAAAAMBOUbABAAAAgJ1ysXWHqampOnDggEqXLi1nZ2dbdw8AsKHMzExdvHhR1atXl5ubm9Hh2D1yJAAUDHnJjzYv2A4cOKBu3brZulsAgIGWLl0qf39/o8Owe+RIAChYcpMfbV6wlS5dWlJWcGXLlrV1944rKCjrOSbG2DgAIA/OnTunbt26WX72I2fkSAdBTgbwkPKSH21esN2c4lG2bFn5+PjYunvHNWZM1jPvGQAHxPS+3CFHOghyMoB8kpv8aPOCDQ+oRQujIwAAABI5GYBNsUqko9i7N+sBAACMRU4GYEOMsDmKoUOznjdtMjIKAABATgZgQ4ywAQAAAICdssoIW3p6ukaNGqWEhAQ5OTlp/Pjxqly5sjW6AgDAobRv317FixeXJPn4+CgiIsLgiAAA9swqBdvmzZuVkZGhmJgYbd26VTNnztScOXOs0RUAAA4jLS1NkhQdHW1wJAAAR2GVgq1ixYrKzMzUjRs3lJycLBcXbpUDAODw4cNKSUlRr169lJGRobfffls1a9Y0OiwAMER8fLyGDh0qPz8/SdLVq1fl4+OjadOm6cKFC3r77be1fPlyS/svv/xSly5dkslkkqurq/r27StJ6tOnj5ydnTV//nxJ0uzZs+Xp6akrV65ozZo1KlOmjOUc9evXV//+/R8o3qioKF26dEnDhw+/6/6UlBQFBQXpnXfeUePGjR+oj7uxSiXl7u6uhIQE/e1vf1NSUpIiIyOt0U02T45am2375KQ2D9XO7kycaHQEAO7j9OnTev311/Xss89KyhpNcXd316xZs/Thhx+qdu3aCggIsLSPiopSUlKShg0blue+fv/9dw0YMEBr1qy5Y9+///1vTZkyRY8//rgkadCgQTKbzYqJidGMGTPy3NeECRPUs2dPlS9fPs/HIjs3Nzf17t1bnTp10smTJ9W3b1+tW7eOX2zexZOj1tpvjiYnA/mmbt262XLTO++8ow0bNqh69er3PKZhw4aKjIxU3759lZKSouTkZF29elWpqalyc3PTzp07NW7cOK1du1YhISHq0qXLQ8WYmpqqsLAw7d+/X6+88so9240bN04mk+mh+robq2SIqKgoNWzYUO+8847Onj2rnj176ptvvpGrq6s1uisY6tc3OgLAPjRteudrgYHSW29J165Jr7565/6QkKzHpUvSLQWTpHxf5c3Pzy/bdLfp06crNjZWgYGBmjVrVraC7auvvtInn3yS5z5WrVqlzz//XElJSXfdf/DgQb377rtq1aqV5bX4+Pg893PT6dOnKdbyScWKFVWhQgWZTCZVrFhRJUqU0MWLFy3FNRwEORmPKoNz7PXr13XhwgV5enrm2K5GjRo6evSozGaztm/frtq1ays5OVnx8fGqW7euEhMTValSpVz1uXv3bs2aNSvbayEhIWrevLllOy0tTe3bt1f9+vV1/Pjxu57n008/1QsvvCCz2ZyrfvPCKgWbh4eHChUqJEny9PRURkaGMjMzrdFVwbFtW9YzSQKwqfT0dH3wwQf6/fffdePGDQ0dOlR16tTJ1bFms1lnz56Vr6+v/P39dfnyZSUkJKhcuXLav3+/vL295ePjk+2Yf/zjH7p27Zplu3LlyhozZky2Np6enlqyZIlatmx5134PHjyoQ4cO6bPPPtNzzz2XbepGSkqKBg4cqHbt2umxxx7TggULVKhQIZ07d05BQUHasWOHDh8+rB49eqhr1646evSoKleurLS0NA0ZMkTJyclKTU3Vu+++m+v3Af8nNjZWR44c0ZgxY3T+/HklJyerdOnSRoeFvCInA/lmx44dCg4OVmJiopycnBQYGKh69erp9OnTd21vMpnk7OysatWq6ciRI9qyZYtee+01JScna8uWLXJ1dVXt2rUt7aOiovTtt99atvv166cGDRpYtv39/e97X7Gnp6caNmyouLi4u+7fvn27fv/9d40bN04///xzXi4/V6xSsIWEhCg0NFRdu3ZVenq6hg0bJnd3d2t0VXCEhmY9850vKOhy+j/g7p7zfm/vPP8fWrFihUqWLKmJEycqKSlJ3bt319q1a+/Z/tixYwoODtaff/6ptLQ0tW3bVh06dJAkBQQEaPXq1erfv7/i4uIUFBR0x/E359/npFmzZjnub9CggVq0aCEfHx998MEHiomJ0VNPPaVr166pX79+6tGjh5o3b674+HidO3dOq1at0sGDBzVkyBB9//33On/+vAYOHKiuXbtq48aNatasmf744w9dunRJUVFRSkxM1MmTJ+8bJ+4UEBCg0aNHq0uXLjKZTJo4cSLTIR0RORmPKhvnWOn/pkQmJSWpV69ell9kurm56fr169naXrt2zTJjr379+tq9e7f27t2rsLAwpaenKzIyUiVKlFCjRo0sx9xvSmRuRtjuJzY2VgkJCQoODtbx48d18OBBlS5dWtWqVcv1OXJilSxRtGjROy4cABzRkSNHtGfPHu3fv1+SlJGRoaSkJJUsWfKu7W9OiUxNTVW/fv1UqlQpywfydu3aKSQkRL169dLOnTsVFhZ2x/G5GWG7nzfeeEMeHh6SpObNm+u7777TU089pZ07d6pq1arZEuBTTz2lQoUKqXjx4vL19VXhwoXl6elpWc1w79696t27t5ydndWtWze9/fbbysjIUHBwcJ5iQpbChQtr+vTpRocBAHanZMmSmjp1qnr06KFVq1apdOnSunr1qo4dOyY/Pz9lZmZq27ZtloVGGjRooNGjR6tChQpycXGRi4uLPDw8tGPHDr355pu57jc3I2z3c+vP9VGjRunVV1/Nt2JNslLBBgCPikqVKqls2bLq16+fUlNTNW/evPvOrZeyfjM4bdo0tW/fXrVq1dLTTz8tLy8vVa5cWXPnzlXLli3vOrKSmxG2nJjNZr3++uuKiYlR2bJltX37dssiKE2bNtV7772nbt26qVatWpKU483Rf/75p4oXLy5nZ2f99ttvunr1qhYsWKALFy4oKCjoviN9AADkhZ+fn4KDgzVhwgTNnj1bERERCg0NlZOTk9LT09W8eXPVrVtXklShQgWdP38+273h9evX1/r161WsWDHLa7dPiaxYsaLGjRuXL/FOmTJFrVu31nPPPZcv57sXCjYAyEFQUJDCwsLUvXt3JScnq2vXrnJycsrVsd7e3hoxYoTCw8MVExNjmZt/c2XA/LR9+3bt2bNHAwcO1IQJEzRw4EC5ubmpcuXKCgwMtMyp9/b21qBBgxQaGmr5LeW9/Pjjj5ZpJU8++aQ++eQTrVq1SoUKFdLgwYPzNX4AQMFTp06dO+6HvnXJfX9//2zL+t9u/fr12bZDQkIUEhJi2R40aJAGDRqUP8FK6tixY7btESNG3NFm0qRJ+dbfTSazNZYyycHp06fVvHlzrV+//o6b7R/GI7+s/81Ve5gvD8CBWOtn/qOK9+v/2PWy/uRkAA8pLz/vGWFzFDNnGh0BgP/v448/vusy+RMnTmT5e6AgICcDsCEKNkdRs6bREQD4/wYOHKiBAwcaHQYAo5CTAdhQ7m7EgPF++CHrAQAAjEVOBmBDjLA5igkTsp5btDA2DgAACjpyMgAbYoQNAAAAAOwUBRsAAAAA2CkKNgAAAACwUxRsAAAAAGCnWHTEUcyfb3QEAABAIicDsCkKNkdRtarREQAAAImcDMCmmBLpKL75JusBAACMRU4GYEOMsDmK6dOzntu2NTYOAAAKOnIyABtihA0AAAAA7BQFGwAAAADYKQo2AAAAALBTVrmHLS4uTl999ZUkKS0tTYcOHdLWrVvl4eFhje4AAAAA4JFklYKtY8eO6tixoyRp7NixeuONNyjWHlZ0tNERAAAAiZwMwKasOiXy119/1bFjx9S5c2drdlMwlC+f9QAAAMYiJwOwIasu6z9//nwNGDDAaud/ctTafG1n15Yty3qm+AUAwFjkZAA2ZLWC7cqVKzp+/Ljq1q1rrS4Klnnzsp5JDgAAGIucDMCGrDYlcteuXapfv761Tg8AAAAAjzyrFWwnTpyQj4+PtU4PAAAAAI88q02J7NOnj7VODQAAAAAFAl+cDQAAAAB2yqqrRCIfxcYaHQEAAJDIyQBsioLNUXh7Gx0BAACQyMkAbIopkY4iKirrAQAAjEVOBmBDFGyOguQAAIB9ICcDsCEKNgAAAACwUxRsAAAAAGCnKNgAALCxxMRENWnSRP/973+NDgUAYOco2AAAsKH09HSFh4fLzc3N6FAAAA6AZf0dxbffGh0BACAfTJ48WUFBQVqwYIHRoeBBkZMB2BAjbI7C3T3rAQBwWHFxcfLy8lKjRo2MDgUPg5wMwIYo2BzF3LlZDwCAw1q5cqW2bdum4OBgHTp0SCNHjtTFixeNDsvhPDlqbY7bVkdOBmBDTIl0FMuXZz2/9ZaxcQAAHtjSpUstfw4ODtaYMWNUunRpAyPCAyEnA7AhRtgAAAAAwE4xwgYAgAGio6ONDgEA4AAYYQMAAAAAO0XBBgAAAAB2iimRjmLTJqMjAAAAUr7k5Pj4eA0dOlR+fn6SpKtXr8rHx0fTpk3ThQsX9Pbbb2v5zcVNJH355Ze6dOmSTCaTXF1d1bdvX0lSnz595OzsrPnz50uSZs+eLU9PT125ckVr1qxRmTJlLOeoX7+++vfv/0DxTpw4URUrVlSXLl3uuj8xMVEdO3bUokWLVLly5Tyde/Pmzfrkk08kSc8884w++OADmUymbOcOCwvTlStXlJmZqSlTpsjX11fh4eE6fPiwunbtqvbt2+t///ufxo4dq2nTpj3QNQL2ymojbPPnz1fnzp3VsWNHrVixwlrd4AGdPn1atWrVUnBwsIKDgxUYGKiQkBD99ddfGjFihGJjY7O1j4qK0owZMx6or3379ik4OPiu+xYvXqw2bdpY4jh+/Lji4uIe+IftkCFDlJKSkufjdu3apSZNmtx134QJE9SxY0cFBwdr3759kqRDhw4pICBAb775pq5duyZJmjdvnn755ZcHihsAUPDUrVtX0dHRio6OVlxcnAoVKqQNGzbkeEzDhg21Z88eSVJKSoqSk5N15swZpaamSpJ27txp+Z6/kJAQy/mjo6PvKNZ+/fVX/fTTTzKbzffs7/Lly+rTp0+OcaWnpys8PFxubm65uu5bJScna+rUqYqMjNTy5ctVrlw5JSUlZWszdepUtW3bVkuXLtXQoUN1/PhxJSUl6dKlS4qJidHKlSslZX32/Pvf/57nGAB7Z5WCLT4+Xr/88ou+/PJLRUdH69y5c9bopmCZNi3rkY/8/PwsP8SXL1+uGjVqKDY2VoGBgfr666+ztf3qq6/UqVOnbK9FR0frypUrOfaxcOFChYWFKS0t7a77Dx48qMmTJ1viqFSp0gNfz80+ihQpkqfjzp49q0WLFikjI+OOfRs3btSJEycUGxurWbNmaezYsZKyvktp/PjxatKkibZu3apLly7p9OnTeuGFFx44fgCAg7BCTr5+/bouXLggT0/PHNvVqFFDR48eldls1vbt21W7dm299NJLio+PV1pamhITE3OdS8uUKaOtW7cqKChIkZGRunDhwh1trl69qkGDBqldu3b3PM/kyZMVFBSUbTRv6dKl6tSpkzp37qzJkydLkrZv366PP/4427G//PKLqlSposmTJ6tr167y9vaWl5dXtjY///yzzp8/r5CQEH3zzTeqXbu2XF1dlZGRobS0NBUuXFinTp1SSkqKqlSpkqtrBxyJVQq2n376SVWqVNGAAQPUr18/NW3a1BrdFCxr1mQ9cpCenq7Q0FB169ZNXbp0UXx8fK5PbzabdfbsWXl4eMjf31+XL19WQkKCJGn//v3y9vaWj49PtmNcXV01cOBAjRo1yvLbvtv5+vpqzpw59+z34MGDWrBggbp06WKZznHT5cuXFRQUpO3btysuLk6DBg1S37591b59e8XFxWnAgAF65ZVX9MMPP0jKSgR16tTR5cuX1aNHDwUHBysoKEi//fabJGnEiBE6c+ZMtj7S0tL0wQcfaMyYMXeN79ixY2rUqJGcnJzk5eUlZ2dnXbx4Ue7u7rp27ZquXbumIkWKaO7cuerXr9+932AAwKMjFzk5N3bs2KHg4GC9+uqr6tixo1q2bKl69erds73JZJKzs7OqVaumI0eOaMuWLWrcuLEaN26sLVu26JdfflHt2rUt7aOioiwzWIKDg7V169Zs53vsscc0cuRILV26VI8//rjatm17x0ha+fLl9fzzz98zpri4OHl5eVlG9W59/b333tOyZctUvnx5ZWRkqF69eho4cGC2dklJSYqPj9fw4cO1cOFCffbZZzpx4kS2NgkJCfLw8FBUVJQef/xxLVy4UO7u7nr55Zf19ttva+DAgZo7d6569OihCRMmaOLEiZbZL8CjwCr3sCUlJenMmTOKjIzU6dOn1b9/f61bty7bfGTkvxUrVqhkyZKaOHGikpKS1L17d61du/ae7Y8dO6bg4GD9+eefSktLU9u2bdWhQwdJUkBAgFavXq3+/fsrLi5OQUFBdxwfGBiowMBAHTlyROHh4SpbtqxmzpyZrU2rVq10+vTpe8bQpk0bde3aVcWKFdPAgQO1ceNGSVnz1fv376/Q0FA9//zziouL09WrV7Vo0SKtXbtWUVFRWr58ueLj4/X555+rRYsW2rRpk/7+979r//79Kl68uKZPn65jx44pOTlZkjRlypQ7+h83bpx69eqlxx577K7xVatWTYsXL1a3bt107tw5HTt2TCkpKerevbsmTZokb29veXt7y93dXf/5z3/0z3/+U40bN1bz5s3vec0AAEhZUyJnzJihpKQk9erVy/KLUTc3N12/fj1b22vXrsnV1VVS1r1ou3fv1t69exUWFqb09HRFRkaqRIkS2QqnkJCQe95zdtOBAwe0fPlynThxQm+//bbq1q2bp2tYuXKlTCaTtm/frkOHDmnkyJGaN2+eIiIitGjRIk2bNk01a9a857TLEiVKqEaNGpYvkPf399ehQ4dUsWLFbG1efvllSdLLL79suUUjKChIQUFB+vnnn+Xr66vt27fL399fkrRmzRoFBgbm6VoAe2WVEbYSJUqoYcOGKly4sCpVqiRXV1ddvnzZGl3hFjd/2xYcHKzBgwcrIyPjjnngt7o5JXLFihV64oknVKpUKbm4ZNXw7dq107/+9S+lpaVp586datas2R3HX758WYsXL9aYMWNUs2ZNDRkyJE/xms1m9ezZU15eXipcuLCaNGmi//znP5KkH3/8UdevX9eNGzcs7atVqyZJKl68uCpXriyTySRPT0/LVMhz587piSeeUOPGjfXSSy/prbfe0uzZs+XkdPd/5ufPn9fu3bv1ySefKDg4WH/99ZeGDRuWrU3Dhg3l7++vnj17avHixXr22WdVokQJlSlTRh999JFCQ0O1cOFC/f3vf1dMTIzGjh2rpUuX5ul9AAAUbCVLltTUqVMVFhamCxcuqFSpUrp69aqOHTsmScrMzNS2bdtUo0YNSVKDBg20du1aVahQQS4uLipSpIg8PDy0Y8eOPBVcGzZs0LJlyxQYGKjo6Gh17txZ7u7ueYp96dKlWrJkiaKjo1WtWjVNnjxZpUuX1vLlyzV27FgtWbJEhw4duuc93tWrV9eRI0d0+fJlZWRkaN++fZaFWG568cUXtXnzZklZ95zfvj8qKkohISFKTU2Vs7OzTCYTI2x4pFhlhO3FF1/U559/rjfffFMXLlxQSkqKSpQoYY2ucItKlSqpbNmy6tevn1JTUzVv3rz7zoWXsn6TN23aNLVv3161atXS008/LS8vL1WuXFlz585Vy5YtLYXcrYYPH64OHTooKipKhQsXznO8ycnJeu211/Ttt9/K3d1d8fHxeuONN5SYmKj27durffv2GjJkiGXRmpxGaA8fPqyqVatKyrqHskyZMlq0aJF++eUXffTRR3f9gtrHHntM3333nWW7QYMGdyyscuLECZUqVUpffPGFzp49qxEjRsjDw8Oyf/PmzapZs6Y8PDwshSNJAgCQV35+fgoODtaECRM0e/ZsRUREKDQ0VE5OTkpPT1fz5s0txViFChV0/vx5BQQEWI6vX7++1q9fr2LFillei4qK0rfffmvZrlixosaNG2fZfvnlly0jV3k1YsQIDR06VE888cRd91etWlUBAQEqWbKkHnvsMT3//PPavn279uzZk21apJeXl9555x316dNHktS6dWtVqVJFx44d05IlSzRmzBiNHDlSYWFhiomJUbFixTR9+nTL8WvXrlWzZs1UpEgRtW7dWkOHDpWTk9MDL5QG2COrFGzNmjXTrl27FBAQILPZrPDwcDk7O1ujq4IjFwtpBAUFKSwsTN27d1dycrK6du16z9Gl23l7e2vEiBEKDw9XTEyMnJycFBgYqL59+2rdunV3PWbRokV5uoSbvvnmG127dk2dO3fWsGHD1KNHDxUuXFj16tVTkyZNFBcXJykreb3++uuKiIi472IemzZtstwr+fTTT2vYsGH67LPP5OTkpAEDBki6f3K51ZQpU9S6dWtVrVpVP/74o2JjY+Xq6qrw8HBLm8zMTC1fvtwyDbRBgwYKDAx84OQHAHAQeVzc6m7q1KmjOnXqZHvt1lUc/f39sy3rf7v169dn2w4JCVFISIhle9CgQRo0aNBDx3nr+W51t9sMbv3laKdOne5YrKxevXp3vUevTZs2atOmTbbX/Pz8LPeXlytXTosXL75rXLceV7ZsWcXExOR8IYADMplzWsvVCk6fPq3mzZtr/fr1dyxikVdPjrr3/VknJ7XJczsAQP7Kz5/5BQHv1/95ctTae+bo2/fl1BYA7FFeft7zxdmPuI8//viuq0VOnDhR5cuXNyAiAAAAALlFweYoxo/Pen7//TwdNnDgwDuW0AUAAA/hAXMyADwIq6wSCStYvz7rAQAAjEVOBmBDFGwAAAAAYKco2AAAAADATlGwAQAAAICdYtERR1GqlNERAAAAiZwMwKYo2BzFypVGRwAAACRyMgCbYkokAAAAANgpCjZHMXp01gMAABiLnAzAhpgS6Si2bzc6AgAAIJGTAdgUI2wAAAAAYKco2AAAAADATlGwAQAAAICd4h42R+HjY3QEAABAIicDsCkKNkexZInREQAAAImcDMCmmBIJAAAAAHaKgs1RDB2a9QAAAMYiJwOwIatNiWzfvr2KFy8uSfLx8VFERIS1uioY9u41OgIAACCRkwHYlFUKtrS0NElSdHS0NU4PAAAAAAWCVaZEHj58WCkpKerVq5d69OihvfwmCgAAAADyzCojbG5uburdu7c6deqkkydPqm/fvlq3bp1cXOxvUconR621/PnkpDYGRgIAAAAA2VmlgqpYsaIqVKggk8mkihUrqkSJErp48aIef/xxa3RXMFSpYnQEAABAIicDsCmrFGyxsbE6cuSIxowZo/Pnzys5OVmlS5e2RlcFx4IFRkcAAHhImZmZCgsL04kTJ+Ts7KyIiAj5+voaHRbyipwMwIascg9bQECA/ve//6lLly4aNmyYJk6caJfTIQEAsKWNGzdKkmJiYjR48GBWUAYA3JdVqqjChQtr+vTp1jh1wfX3v2c981s9AHBYLVq0UNOmTSVJZ86ckbe3t7EB4cGQkwHYEF+c7SiOHMl6AAAcmouLi0aOHKnx48erVatWRodjV25dCOx+r9+rbV7bPBByMgAbomADAMDGJk+erO+++07vv/++rl27ZnQ4AAA7RsEGAICNrFq1SvPnz5ckFSlSRCaTSc7OzgZHBQCwZ6wEAgCAjbzyyisaPXq0unXrpoyMDIWGhsrV1dXosAAAdoyCzVHUrGl0BACAh+Tu7q5Zs2YZHQYeFjkZgA3lakrkihUrsm1//vnnVgkGOZg5M+sBALAb5McCipwMwIZyHGFbs2aNNmzYoPj4eO3YsUNS1pd+Hj16VD169LBJgAAA2BvyIwDAVnIs2Bo1aqTSpUvrzz//VOfOnSVJTk5OKl++vE2Cwy26d896XrLE2DgAAOTHgo6cDMCGcizYPD09VadOHdWpU0eJiYlKS0uTlPVbRNjY6dNGRwAA+P/IjwUcORmADeVq0ZGxY8dq8+bNKlOmjMxms0wmk2JiYqwdGwAAdo38CACwtlwVbPv27dMPP/wgJye+tg0AgJvIjwAAa8tVhqlQoYJlugcAAMhCfgQAWFuuRtjOnj2rZs2aqUKFCpLElA8j1KtndAQAgNuQHwsocjIAG8pVwTZ9+nRrx4H7iYgwOgIAwG3IjwUUORmADeWqYPvqq6/ueG3gwIH5HgwAAI6E/AgAsLZcFWze3t6SJLPZrP/85z+6ceOGVYPCXbzxRtbzypXGxgEAsCA/FlDkZAA2lKuCLSgoKNt2nz59rBIMcpCYaHQEAIDbkB8LKHIyABvKVcF24sQJy58vXryos2fP3veYxMREdezYUYsWLVLlypUfPEIAAOzUg+RHAADyIlcFW3h4uOXPrq6uGjFiRI7t09PTFR4eLjc3t4eLDgAAO5bX/AgAQF7lqmCLjo5WUlKSTp06JR8fH3l5eeXYfvLkyQoKCtKCBQvyJUgAAOxRXvMjAAB5lasvzv7Xv/6loKAgRUZGqnPnzvr666/v2TYuLk5eXl5q1KhRvgUJSc2bZz0AAHYjL/kRjxByMgAbytUIW1RUlOLi4lS0aFElJyerZ8+eateu3V3brly5UiaTSdu3b9ehQ4c0cuRIzZs3T6VLl87XwO/nyVFrbdqf1b3/vtERAABuk5f8iEcIORmADeWqYDOZTCpatKgkqVixYnJ1db1n26VLl1r+HBwcrDFjxti8WAMAwBbykh8BAHgQuSrYfH19NWnSJPn7+2vPnj3y9fW1dly43d/+lvX8r38ZGwcAwIL8WECRkwHYUK7uYQsMDJSnp6e2bdumuLg4devWLVcnj46OZkn//JKSkvUAANiNB82PcHDkZAA2lKuCbdKkSWrZsqXCw8MVGxurSZMmWTsuAADsHvkRAGBtuSrYXFxc5OfnJ0kqX768nJxydRgAAI808iMAwNpydQ/bE088oY8++kg1a9bU/v37VaZMGWvHBQCA3SM/AgCsLVe/CoyIiJCXl5c2b94sLy8vRUREWDsu3O6117IeAAC7QX4soMjJAGwoVyNsrq6uCgkJsXIoyNHw4UZHAAC4DfmxgCInA7AhJtsDAAAAgJ2iYHMUTZtmPQAAgLHIyQBsiIINAAAAAOwUBRsAAAAA2CkKNgAAAACwUxRsAAAAAGCncrWsP+xAYKDREQAAAImcDMCmKNgcxVtvGR0BAACQyMkAbIqCzVFcu5b17O5ubBwAgAeWnp6u0NBQJSQk6Pr16+rfv7+aN29udFjIK3IyABuiYHMUr76a9bxpk6FhAAAe3OrVq1WiRAlNnTpVSUlJ6tChAwWbIyInA7AhCjYAAGykdevWatWqlWXb2dnZwGgAAI6AVSIBALCRokWLqlixYkpOTtbgwYM1dOhQo0Oye0+OWpvt+fZ9t75++3ZOx9kbe4wJgH2gYAMAwIbOnj2rHj16qF27dmrbtq3R4QAA7JxVpkRmZmYqLCxMJ06ckLOzsyIiIuTr62uNrgAAcBiXLl1Sr169FB4ernr16hkdDgDAAVilYNu4caMkKSYmRvHx8YqIiNC8efOs0VXBERJidAQAgIcUGRmpK1euaO7cuZo7d64kaeHChXJzczM4MuQJORmADVmlYGvRooWaNm0qSTpz5oy8vb2t0U3BQnIAAIcXFhamsLAwo8PAwyInA7Ahq60S6eLiopEjR+r777/X7NmzrdVNwXHpUtYzxS8AAMYiJwOwIasu6z958mQNHz5cgYGBWrt2rdz5gsl7un11qJOT2mRvEBCQ9cx3vgAAYCxyMgAbssoqkatWrdL8+fMlSUWKFJHJZOK7ZgAAAAAgj6wywvbKK69o9OjR6tatmzIyMhQaGipXV1drdAUAAAAAjyyrFGzu7u6aNWuWNU4NAAAAAAUGX5wNAAAAAHbKqouOIB/17290BAAAQCInA7ApCjZH0bmz0REAAACJnAzAppgS6ShOncp6AAAAY5GTAdgQI2yOIjg465nvfAEAwFjkZAA2xAgbAAAAANgpCjYAAAAAsFMUbAAAAABgpyjYAAAAAMBOseiIo3jnHaMjAAAAEjkZgE1RsDmKtm2NjgAAAEjkZAA2xZRIR/Hbb1kPAABgLHIyABtihM1R/OMfWc985wsAAMYiJwOwIUbYAAAAAMBOUbABAAAAgJ2iYAMAAAAAO0XBBgAAAAB2ikVHHEVYmNERAAAAiZwMwKbyvWBLT09XaGioEhISdP36dfXv31/NmzfP724KnhYtjI4AAABI5GQANpXvBdvq1atVokQJTZ06VUlJSerQoQMFW37YuzfruWZNI6MAAADkZAA2lO8FW+vWrdWqVSvLtrOzc353UTANHZr1zHe+AABgLHIyABvK94KtaNGikqTk5GQNHjxYQ2/+UAMAAAAA5IlVFh05e/asBgwYoK5du6pt27bW6MImnhy1NlftTk5qc8/jbt8HAAAAALmV7wXbpUuX1KtXL4WHh6tevXr5fXoAAAAAKDDy/XvYIiMjdeXKFc2dO1fBwcEKDg5WampqfncDAAAAAI+8fB9hCwsLUxjfT5L/Jk40OgIAACCRkwHYFF+c7Sjq1zc6AgAAIJGTAdhUvk+JhJVs25b1AAAAxiInA7AhRtgcRWho1jPf+QIAgLHIyQBsiBE2AABsbN++fQoODjY6DACAA2CEDQAAG1q4cKFWr16tIkWKGB0KAMABMMIGAIAN+fr6as6cOUaHAQBwEBRsAADYUKtWreTiUnAnuDw5am2Or9+6//a2T45ae8/jc+onp3M+qgrKdQIFQcHNGI5m5kyjIwAAABI5GYBNUbA5ipo1jY4AAABI5GQANsWUSEfxww9ZDwAAYCxyMgAbYoTNUUyYkPXcooWxcQAAHpqPj4+WL19udBh4UORkADbECBsAAAAA2CkKNgAAAACwUxRsAAAAAGCnKNgAAAAAwE6x6IijmD/f6AgAAIBETgZgUxRsjqJqVaMjAAAAEjkZgE0xJdJRfPNN1gMAABiLnAzAhqxWsO3bt0/BwcHWOn3BM3161gMAABiLnAzAhqwyJXLhwoVavXq1ihQpYo3TAwAAAECBYJURNl9fX82ZM8capwYAAACAAsMqI2ytWrXS6dOnrXFqq3py1FqbnvPkpDa5Oi6ndka6/docIU57jRGwF/x/AQDAvrDoCAAAAADYKZb1dxTR0UZHAAAAJHIyAJuiYHMU5csbHQEAAJDIyQBsympTIn18fLR8+XJrnb7gWbYs6wEAAIxFTgZgQ4ywOYp587KeO3c2Ng4AAAo6cjIAG2LREQAAAACwUxRsAAAAAGCnKNgAAAAAwE5RsAEAAACAnWLREUcRG2t0BAAAQCInA7ApCjZH4e1tdAQAAEAiJwOwKaZEOoqoqKwHAAAwFjkZgA1RsDkKkgMAAPaBnAzAhijYAAAAAMBOUbABAAAAgJ2iYAMAAAAAO0XBBgAAAAB2imX9HcW33xodAQAAkMjJAGyKgs1RuLsbHQEAAJDIyQBsiimRjmLu3KwHAAAwFjkZgA1RsDmK5cuzHgAAwFjkZAA2RMEGAAAAAHbKKvew3bhxQ2PGjNFvv/2mwoULa8KECapQoYI1ugIAwGGQHwEAeWWVEbYffvhB169f17Jly/TOO+9o0qRJ1ugGAACHQn4EAOSVVUbY9uzZo0aNGkmSatasqQMHDlj2ZWZmSpLOnTv38B1dvfzw58gHp0+fzv5CLuPKdlwOx5w+fVrKyLi5kdfwrOe2mO94H+zFLXHabYyAvcjn/y83f9bf/Nlf0OWUH6V8zpH26urlu//buv3fXh5z/K3H3O14S5+3tnlQ1sjJ93pf7OV8APJVXvKjyWw2m/M7gPfee0+vvPKKmjRpIklq2rSpfvjhB7m4uGj37t3q1q1bfncJALBjS5culb+/v9FhGC6n/CiJHAkABUxu8qNVRtiKFSumq1evWrZv3LhhSUbVq1fX0qVLVbp0aTk7O1ujewCAncjMzNTFixdVvXp1o0OxCznlR4kcCQAFRV7yo1UKtlq1amnjxo169dVXtXfvXlWpUsWyz83Njd+yAkABwqIa/yen/CiRIwGgIMltfrTKlMibq2AdOXJEZrNZEydOVOXKlfO7GwAAHAr5EQCQV1Yp2PKDoy59vG/fPk2bNk3R0dFGh5Kj9PR0hYaGKiEhQdevX1f//v3VvHlzo8O6p8zMTIWFhenEiRNydnZWRESEfH19jQ7rvhITE9WxY0ctWrTI7j+UtW/fXsWLF5ck+fj4KCIiwuCI7m3+/PnasGGD0tPT1aVLF3Xq1MnokO4pLi5OX331lSQpLS1Nhw4d0tatW+Xh4WFwZHeXnp6uUaNGKSEhQU5OTho/frzd/9t9FDhqzjOao+Rce+Boed9ojvq5w2iO9LnHaHn53GWVKZH54dalj/fu3atJkyZp3rx5RoeVo4ULF2r16tUqUqSI0aHc1+rVq1WiRAlNnTpVSUlJ6tChg13/4N64caMkKSYmRvHx8YqIiLD7fw/p6ekKDw+Xm5ub0aHcV1pamiQ5xIee+Ph4/fLLL/ryyy+VkpKiRYsWGR1Sjjp27KiOHTtKksaOHas33njDbos1Sdq8ebMyMjIUExOjrVu3aubMmZozZ47RYT3yHDHnGc2Rcq49cLS8bzRH/NxhNEf63GO0vH7ussr3sOWH+y19bI98fX0d5oNN69atNWTIEMu2vd/c3qJFC40fP16SdObMGXl7exsc0f1NnjxZQUFBKlOmjNGh3Nfhw4eVkpKiXr16qUePHtq7d6/RId3TTz/9pCpVqmjAgAHq16+fmjZtanRIufLrr7/q2LFj6ty5s9Gh5KhixYrKzMzUjRs3lJycnG1BDFiPI+Y8ozlSzrUHjpb3jeaInzuM5kife4yW189ddpuJk5OTVaxYMcu2s7OzMjIy7PrDQ6tWrRzmO0+KFi0qKet9Hjx4sIYOHWpsQLng4uKikSNH6vvvv9fs2bONDidHcXFx8vLyUqNGjbRgwQKjw7kvNzc39e7dW506ddLJkyfVt29frVu3zi7/vyUlJenMmTOKjIzU6dOn1b9/f61bt04mk8no0HI0f/58DRgwwOgw7svd3V0JCQn629/+pqSkJEVGRhodUoHgiDnPaI6Uc+2BI+Z9oznS5w6jOdrnHqPl9XOX3Y6w3W/pYzy8s2fPqkePHmrXrp3atm1rdDi5MnnyZH333Xd6//33de3aNaPDuaeVK1dq27ZtCg4O1qFDhzRy5EhdvHjR6LDuqWLFinr99ddlMplUsWJFlShRwm7jLVGihBo2bKjChQurUqVKcnV11eXLefuCXVu7cuWKjh8/rrp16xodyn1FRUWpYcOG+u677/T1119r1KhRlqkbsB5yHmzBEfO+0Rzlc4fRHO1zj9Hy+rnLbgu2WrVqacuWLZJ016WP8XAuXbqkXr166d1331VAQIDR4dzXqlWrNH/+fElSkSJFZDKZ7Ho6x9KlS7VkyRJFR0erWrVqmjx5skqXLm10WPcUGxurSZMmSZLOnz+v5ORku433xRdf1I8//iiz2azz588rJSVFJUqUMDqsHO3atUv169c3Ooxc8fDwsNwE7enpqYyMDGVmZhoc1aOPnAdrc7S8bzRH+9xhNEf73GO0vH7usttf37Vs2VJbt25VUFCQZelj5J/IyEhduXJFc+fO1dy5cyVl3cBtrzeKvvLKKxo9erS6deumjIwMhYaGytXV1eiwHhkBAQEaPXq0unTpIpPJpIkTJ9rtb/ebNWumXbt2KSAgQGazWeHh4XafRE+cOCEfHx+jw8iVkJAQhYaGqmvXrkpPT9ewYcPk7u5udFiPPHIerM3R8r7R+NwBa8rr5y67XdYfAAAAAAo6u50SCQAAAAAFHQUbAAAAANgpCjYAAAAAsFMUbAAAAABgpyjYAAAAAMBO2ee63YAdOn36tF5//XU9++yzkqS0tDS5u7tr1qxZ+vDDD1W7du1s320TFRWlpKQkDRs2LNt5IiMj1bBhQ23cuFHe3t7q0qVLnmP5xz/+Yfl+mFulpqbqgw8+0KRJk2QymfJ8XgAAHkZ8fLyGDh0qPz8/SdLVq1fl4+OjYcOGaeTIkVq+fLml7ZdffqlLly7JZDLJ1dVVffv2lST16dNHzs7Oljw3e/ZseXp66sqVK1qzZo3KlCljOUf9+vXVv39/G14hYHuMsAF54Ofnp+joaEVHR2v58uWqUaOGYmNjFRgYqK+//jpb26+++kqdOnXK9trZs2d15MgRVa9e/YFjOHPmjMqVK3fXfW5ubnrhhRe0atWqBz4/AAAPo27dupZcGRcXp0KFCmnDhg33bN+wYUPt2bNHkpSSkqLk5GSdOXNGqampkqSdO3eqUaNGkrK+K/LmuaOjoynWUCBQsKFAS09PV2hoqLp166YuXbooPj4+18eazWadPXtWHh4e8vf31+XLl5WQkCBJ2r9/v7y9ve/4suQvv/xSrVq1yvba77//rjfeeEOHDx/WnDlzNHz4cPXu3VsBAQGKi4tTv3791KpVK+3du1eStHHjRjVp0kQnTpxQUFCQunfvrp49e+r8+fOSpL/97W/64osvHuJdAQAgf1y/fl0XLlyQh4fHPdvUqFFDR48eldls1vbt21W7dm299NJLio+PV1pamhITE1WpUiUbRg3YF6ZEokBbsWKFSpYsqYkTJyopKUndu3fX2rVr79n+2LFjCg4O1p9//qm0tDS1bdtWHTp0kJT1rfWrV69W//79FRcXp6CgoDuO37lzpzp27GjZPnHihFauXKnp06frySef1Pfffy83Nzd9+umnWrBggTZv3qzIyEitXLlSa9euVc2aNRUfH6+AgADFxsbq2Wef1ahRo7R792799ddfeuyxx+Tp6amkpCT973//U/HixfP/TQMAIAc7duxQcHCwEhMT5eTkpMDAQNWrV08rVqy4o63JZJKzs7OqVaumI0eOaMuWLXrttdeUnJysLVu2yNXVVbVr17a0j4qK0rfffmvZ7tevnxo0aGCT6wKMQsGGAu3IkSPas2eP9u/fL0nKyMhQUlKSSpYsedf2N6dEpqamql+/fipVqpRcXLL+G7Vr104hISHq1auXdu7cqbCwsDuOT0pKkre3t2V7y5YtcnFxkbOzs+W1Z555RpJUvHhxyz0Anp6eSktLU0pKipydneXq6qqAgAAtXLhQffr0UfHixbPdK+ft7a0///yTgg0AYHN169bVjBkzlJSUpF69esnHx0dubm66fv16tnbXrl2Tq6urpKx70Xbv3q29e/cqLCxM6enpioyMVIkSJSzTIaWsKZEPcu834MiYEokCrVKlSmrTpo2io6O1cOFCtW7dWp6envc9zs3NTdOmTdPcuXN1+PBhSZKXl5cqV66suXPnqmXLlpZC7lZeXl66cuWKZbtnz54KDQ3ViBEjlJmZKUk5LhaydetW1atXT5K0fv16vfjii/rss8/UunVr/fOf/7S0u3Lliry8vHL3JgAAYAUlS5bU1KlTFRYWphs3bujq1as6duyYJCkzM1Pbtm1TjRo1JEkNGjTQ2rVrVaFCBbm4uKhIkSLy8PDQjh07VLduXSMvAzAcI2wo0IKCghQWFqbu3bsrOTlZXbt2lZNT7n6P4e3trREjRig8PFwxMTGWaR99+/bVunXr7npM7dq1tW/fPj3xxBOW1+rXr69169Zp4cKF9+1zy5YtGjBggCSpevXqevfddzVnzhw5OTlp9OjRkrKKNQ8PDxUtWjRX1wEAgLX4+fkpODhYEyZMUEREhEJDQ+Xk5KT09HQ1b97cUoxVqFBB58+fz7bacv369bV+/XoVK1bM8trtUyIrVqyocePG2e6CAAOYzGaz2egggIIiISFBkydP1uzZs63Wx9KlS1WsWDG1a9fOan0AAADANhhhA27z8ccf33W1yIkTJ6p8+fIPde5y5cqpatWq+vXXXy3TQPJTamqqfv75Z02dOjXfzw0AAADbY4QNAAAAAOwUi44AAAAAgJ2iYAMAAAAAO0XBBgAAAAB2ioINAAAAAOwUBRsAAAAA2CkKNgAAAACwU/8Pp87mI1/J16wAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 1080x288 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot results\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.style.use(\"seaborn-white\")\n",
    "\n",
    "df_f = df_gaia\n",
    "\n",
    "fig, axs = plt.subplots(1, 2, figsize=(15, 4))\n",
    "fig.suptitle(f\"{plot_title} ({len(df_f)})\");\n",
    "\n",
    "ax = axs[0]\n",
    "ax.hist(df_f[\"e_RV\"], bins=100);\n",
    "ax.axvline(1.5, linestyle=\"--\", c=\"red\", label=\"e_RV = 1.5km/s\");\n",
    "ratio_large_e_rv = len(df_f[df_f[\"e_RV\"] > 1.5]) / len(df_f)\n",
    "ax.text(1.7, ax.get_ylim()[1] * 0.7, f\"e_RV > 1.5km/s: {100* ratio_large_e_rv:.1f}%\")\n",
    "ax.legend();\n",
    "ax.set_xlabel(\"e_RV (km/s)\");\n",
    "ax.set_ylabel(\"count\");\n",
    "\n",
    "ax = axs[1]\n",
    "ax.hist(df_f[\"RUWE\"], bins=100);\n",
    "ax.axvline(1.4, linestyle=\"--\", c=\"red\", label=\"RUWE = 1.4\");\n",
    "ratio_large_ruwe = len(df_f[df_f[\"RUWE\"] > 1.4]) / len(df_f)\n",
    "ax.text(1.6, ax.get_ylim()[1] * 0.7, f\"RUWE > 1.4s: {100* ratio_large_ruwe:.1f}%\")\n",
    "ax.legend();\n",
    "ax.set_xlabel(\"RUWE\");\n",
    "ax.set_ylabel(\"count\");\n",
    "ax.set_xlim(0, 5);"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "lkv2_1_dev",
   "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.9.10"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "cf14734aefd952c101fa7a730b7e28a3f566541dc5a8dc5204d96dcd505be75c"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}