{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analyzing your Runs\n",
    "\n",
    "In this tutorial you will learn: \n",
    "\n",
    "1. How to visualize and analyze your runs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from bokeh.io import output_notebook \n",
    "from bokeh.plotting import show, figure\n",
    "from bokeh.palettes import Colorblind\n",
    "output_notebook()\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import astropy.units as u\n",
    "\n",
    "#here is pyeddy\n",
    "import virga.justdoit as jdi\n",
    "import virga.justplotit as jpi\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This time let's load in a Hot Jupiter PT profile and go through the motions of the analysis. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mieff_dir = '/data/virga'\n",
    "\n",
    "metallicity = 1 #atmospheric metallicity relative to Solar\n",
    "mean_molecular_weight = 2.2 # atmospheric mean molecular weight\n",
    "\n",
    "#set the run \n",
    "a = jdi.Atmosphere(['MnS','Cr','MgSiO3'],\n",
    "                  fsed=1,mh=metallicity,\n",
    "                 mmw = mean_molecular_weight)\n",
    "\n",
    "#set the planet gravity\n",
    "a.gravity(gravity=7.460, gravity_unit=u.Unit('m/(s**2)'))\n",
    "\n",
    "#Get preset pt profile for testing\n",
    "a.ptk(df = jdi.hot_jupiter())\n",
    "\n",
    "#get full dictionary output \n",
    "all_out = jdi.compute(a, as_dict=True,\n",
    "                      directory=mieff_dir)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PT Profile with Condensation Curves\n",
    "\n",
    "Thick lines in this figure signify the gases that were supplied in the run, NOT those that were recommended by `pyeddy.recommend`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "show(jpi.pt(all_out,plot_height=450))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cumulative Optical Depth By Gas \n",
    "\n",
    "This is the cumulative optical depth for the top of the atmosphere down. Remember that this is the optical depth for conservative geometric scatterers. They won't identically match up with the `opd_per_layer` output shown below. This quantity, therefore, is just to help you visualize about where and how optically thick the cloud is. The actual quantities in `opd_per_layer` should be used for any real calculations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "show(jpi.opd_by_gas(all_out))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualize Optical Depth, Asymmetry, and Single Scattering\n",
    "\n",
    "The key to analyzing these plots: \n",
    "\n",
    "1. Start with the middle plot, locate the cloud decks\n",
    "2. Scan over to the single scattering and asymmetry plots and get a sense for what the optical properties are doing where the cloud decks are"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "show(jpi.all_optics(all_out))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also see the optics but in a single dimension (averaged over wavelength)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "show(jpi.all_optics_1d(all_out, wave_range=[1,2]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Particle Radii \n",
    "\n",
    "Currently, all particle distributions in `pyeddy` are log-normal. There is an input (`sig`) in the run call if you'd like to alter the width.\n",
    "\n",
    "This figure is important for two things: \n",
    "1. Determining the full profile of particle radii in your atmosphere \n",
    "2. Assessing whether or not the grid that was used to compute Mie properties was sufficient. **If your distribution is cut off at the large or small end, consider changing `rmin` and `nradii` in `pyeddy.calc_mie_db`.** "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "fig, dndr = jpi.radii(all_out,at_pressure=1e-3)\n",
    "show(fig)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mean Mass Mixing Ratio of the Condensates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "show(jpi.condensate_mmr(all_out))"
   ]
  },
  {
   "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.5"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": true,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}