{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "5a99ae8e",
   "metadata": {},
   "source": [
    "# Quick start"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8f04291",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import astropy.units as u\n",
    "\n",
    "import xdust"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7c9f1258",
   "metadata": {},
   "source": [
    "## Extinction curves\n",
    "\n",
    "To compute an extinction curve, start with a Grain Population (specifying a grian size distribution, composition, and scattering model). `SingleGrainPop` accepts short string names for the grain size distribution, composition, and scattering model. Here we use `'Powerlaw'` for the size distribution, `'Silicate'` for the composition (Draine 2003), and `'Mie'` for the scattering physics."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "81af396e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A grid of wavelengths spanning the optical/UV\n",
    "wavelength_grid = np.logspace(-1, 1, 100) * u.micron\n",
    "\n",
    "silicate_pop = xdust.SingleGrainPop('Powerlaw', 'Silicate', 'Mie')\n",
    "silicate_pop.calculate_ext(wavelength_grid)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "239958b4",
   "metadata": {},
   "source": [
    "The `plot_ext` method plots the extinction, scattering, and absorption efficiencies, integrated over the grain size distribution.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0c554355",
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = plt.subplot(111)\n",
    "silicate_pop.plot_ext(ax, 'all', frameon=False)\n",
    "plt.loglog()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8f6953c2",
   "metadata": {},
   "source": [
    "## An X-ray dust scattering halo for the same grain population\n",
    "\n",
    "X-ray scattering halos use the `xdust.halos` module, which needs an energy grid, an angular grid, and a `GrainPop` (here, our `silicate_pop`) to model. We compare two geometries:\n",
    "\n",
    "- `UniformGalHalo` -- dust distributed uniformly along the line of sight\n",
    "- `ScreenGalHalo` -- dust confined to a thin screen at a fraction `x` of the distance to the source"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68f067ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "# X-ray energy grid and angular grid for the halo calculation\n",
    "n_energy, n_theta = 50, 200\n",
    "energy_grid = np.logspace(-1, 1, n_energy) * u.keV\n",
    "theta_grid = np.logspace(-1, 4, n_theta) * u.arcsec\n",
    "\n",
    "# Initialize the Halo objects for uniform and screen geometries\n",
    "uniform_halo = xdust.halos.UniformGalHalo(energy_grid, theta_grid)\n",
    "screen_halo = xdust.halos.ScreenGalHalo(energy_grid, theta_grid)\n",
    "\n",
    "dust_pop = xdust.make_MRN_RGDrude(md=1.e-6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "de93574b",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute the scattering halo profiles\n",
    "uniform_halo.calculate(dust_pop)\n",
    "screen_halo.calculate(dust_pop, x=0.33)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2d6964ad",
   "metadata": {},
   "source": [
    "Both halos are calculated from the same grain population, so we can compare their surface brightness profiles directly. Below, we plot the normalized intensity at a single photon energy.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f3509db2",
   "metadata": {},
   "outputs": [],
   "source": [
    "i=10\n",
    "plt.plot(theta_grid, screen_halo.norm_int[i, :], label='Screen halo')\n",
    "plt.plot(theta_grid, uniform_halo.norm_int[i, :], label='Uniform halo')\n",
    "plt.loglog()\n",
    "plt.legend()\n",
    "\n",
    "plt.title(\"Photon Energy = {:.2f}\".format(energy_grid[i]))\n",
    "plt.xlabel(r\"Angular Distance from Point Source ({})\".format(theta_grid.unit))\n",
    "plt.ylabel(r\"Normalized Intensity ({})\".format(screen_halo.norm_int.unit))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "66774e9a",
   "metadata": {},
   "source": [
    "To compute the surface brightness, multiply each row [i,:] by the absorbed flux[i] in units of your choice (e.g., photon/s/cm$^2$). See the next pages in {ref}`tutorials` for a closer look at customizing grain populations, choosing scattering models, and working with halo intensities and point-source spectra."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}