diff --git a/docs/notebooks/Develop_Rates_Simulations.ipynb b/docs/notebooks/Develop_Rates_Simulations.ipynb index e2b95d8..aa5bf2f 100644 --- a/docs/notebooks/Develop_Rates_Simulations.ipynb +++ b/docs/notebooks/Develop_Rates_Simulations.ipynb @@ -1,42 +1,11 @@ { "cells": [ - { - "cell_type": "code", - "execution_count": null, - "id": "378a73b4", - "metadata": {}, - "outputs": [], - "source": [ - "#FIXME: remove after done developing\n", - "import highz_funcs\n", - "import importlib\n", - "importlib.reload(highz_funcs)\n", - "\n", - "# imports \n", - "import numpy as np1\n", - "import json\n", - "from astropy import units as u\n", - "from scipy.stats import norm\n", - "from astropy.cosmology import Planck18 as cosmo\n", - "\n", - "# my functions\n", - "from highz_funcs import get_volumetric_BTS_rates\n", - "from highz_funcs import calculate_rates_vs_redshift\n", - "from highz_funcs import plot_rates_vs_redshift\n", - "from highz_funcs import plot_rates_grid\n", - "from highz_funcs import get_BTS_rates_from_filtered\n", - "from highz_funcs import apply_z_dependent_rates\n", - "from highz_funcs import get_z_limit\n", - "from highz_funcs import calculate_rubin_detectable_rates\n", - "from highz_funcs import plot_observed_vs_rubin" - ] - }, { "cell_type": "markdown", "id": "0353c712", "metadata": {}, "source": [ - "## This notebook simulates expectations for Rubin observations of bright transients at high redshifts.\n", + "# Simulate expectations for Rubin observations of bright transients at high redshifts.\n", "\n", "The main sections are:\n", "\n", @@ -55,7 +24,7 @@ "id": "eb5a3e46", "metadata": {}, "source": [ - "Assumptions: \n", + "Notes and caveats: \n", "- We account for modeled transient rates as a function of redshift. These may have high uncertainty.\n", "- We currently use luminosities distributions of local transient, which may differ from the high redshift population\n", "- Rubin's unique cadence and colors could change detection rates in unmodeled ways\n", @@ -66,28 +35,43 @@ { "cell_type": "code", "execution_count": null, - "id": "7c780453", + "id": "8d7788a4", "metadata": {}, "outputs": [], "source": [ - "# FIXME NOW\n", - "# by hand simulate TDE, SN lightcurves. AT DIFFERENT Z. understand completely each step\n", - "# use redback to simulate the above\n", - "# use snana to simulate the above\n", - "# produce plots of the above for comparison" + "#FIXME: remove after done developing\n", + "\n", + "import importlib\n", + "import thor.utils.rates_sims\n", + "importlib.reload(thor.utils.rates_sims)" ] }, { "cell_type": "code", "execution_count": null, - "id": "6cab4ad8", + "id": "e22d0546", "metadata": {}, "outputs": [], "source": [ - "# FIXME\n", + "# imports \n", + "import numpy as np1\n", + "import json\n", + "from astropy import units as u\n", + "from scipy.stats import norm\n", + "from astropy.cosmology import Planck18 as cosmo\n", "\n", - "# plasicc rates compare to BTS rates\n", - "# apply cuts based on recovery due to cadence etc." + "# my functions\n", + "from thor.utils.rates_sims import (\n", + " get_volumetric_BTS_rates,\n", + " calculate_rates_vs_redshift,\n", + " plot_rates_vs_redshift,\n", + " plot_rates_grid,\n", + " get_BTS_rates_from_filtered,\n", + " apply_z_dependent_rates,\n", + " get_z_limit,\n", + " calculate_rubin_detectable_rates,\n", + " plot_observed_vs_rubin\n", + ")" ] }, { @@ -622,7 +606,7 @@ "outputs": [], "source": [ "# stacked histogram\n", - "from highz_funcs import plot_rates_layered\n", + "from thor.utils.rates_sims import plot_rates_layered\n", "\n", "plot_rates_layered(\n", " BTS_calculated_rates_vs_z,\n", @@ -732,7 +716,7 @@ "source": [ "# plot current and future tdes\n", "\n", - "with open(\"data/tde_redshifts.json\") as f:\n", + "with open(\"../../data/tde_redshifts.json\") as f:\n", " raw = json.load(f)\n", "tde_z_obs = np.array([z for z in raw if z is not None], dtype=float)\n", "\n", @@ -782,30 +766,8 @@ "# lyman alpha forest at 1216 Angstroms in rest frame, redshifted to 1216*(1+z) in observed frame\n", "# transmission model: np.exp(-0.0036 * (lam[laf_mask] / lya) ** 3.46)\n", "\n", - "from highz_funcs import plot_tde_sed_with_filters\n", - "plot_tde_sed_with_filters(temperature_k=30000, redshifts=[0, 0.5, 1.0, 1.5, 2.0, 2.5])" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "lh2679qxnx", - "metadata": {}, - "outputs": [], - "source": [ - "# SNIa SED — blackbody at 10,000 K\n", - "hf.plot_snia_sed_with_filters(temperature_k=10000, redshifts=[0, 0.5, 1.0, 1.5, 2.0, 2.5])" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "yhi6hh83ix", - "metadata": {}, - "outputs": [], - "source": [ - "# CCSN SED — blackbody at 6,000 K\n", - "hf.plot_ccsn_sed_with_filters(temperature_k=6000, redshifts=[0, 0.5, 1.0, 1.5, 2.0, 2.5])" + "from thor.utils.lightcurve_sims import plot_blackbody_sed_with_filters\n", + "plot_blackbody_sed_with_filters(temperature_k=30000, redshifts=[0, 0.5, 1.0, 1.5, 2.0, 2.5])" ] }, { @@ -816,18 +778,6 @@ "# Simulated Lightcurves" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "aa7b70ca", - "metadata": {}, - "outputs": [], - "source": [ - "import highz_funcs as hf \n", - "importlib.reload(hf)\n", - "import matplotlib.pyplot as plt " - ] - }, { "cell_type": "code", "execution_count": null, @@ -837,246 +787,21 @@ "source": [ "# simulate LSST lightcurves\n", "# FIXME: clarify object evolution (with correct time dilation) \n", - "# FIXME: improve rubin cadence modeling \n", + "# FIXME: improve rubin cadence modeling \n", + "\n", + "from thor.utils.lightcurve_sims import (simulate_lsst_lightcurve, plot_lsst_lightcurve)\n", " \n", "# Single TDE at z=1\n", - "lc = hf.simulate_lsst_tde_lightcurve( \n", + "lc = simulate_lsst_lightcurve( \n", " redshift=1.0, \n", " t_peak_mjd=60000.0, \n", " L_peak_erg_s=1e44, \n", " T_peak_k=30000, \n", - " t_rise_days=30, \n", + " t_rise_days=30,\n", + " transient_class='TDE' \n", ") \n", - "hf.plot_lsst_tde_lightcurve(lc, t_peak_mjd=60000.0, redshift=1.0, T_k=30000, L_peak=1e44, ymax=20, xmin=-70, xmax=50) " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "0059286b", - "metadata": {}, - "outputs": [], - "source": [ - "# Grid of redshifts — same TDE moved to different distances\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=True)\n", - "fig.suptitle('Custom TDE Simulations', fontsize=20, fontweight='bold', y=1.01)\n", - "\n", - "params = [\n", - " (0.1, 1e43), (0.5, 1e44), (1.0, 1e44),\n", - " (1.5, 1e45), (2.0, 1e45), (2.5, 1e45),\n", - "]\n", - "for ax, (z, Lpeak) in zip(axes.flat, params):\n", - " lc = hf.simulate_lsst_tde_lightcurve(\n", - " redshift=z, t_peak_mjd=60000.0,\n", - " L_peak_erg_s=Lpeak, T_peak_k=30000,\n", - " t_rise_days=30,\n", - " )\n", - " hf.plot_lsst_tde_lightcurve(lc, t_peak_mjd=60000.0,\n", - " redshift=z, T_k=30000, L_peak=Lpeak, ax=ax, ymax=20, xmax=300)\n", - "plt.tight_layout()\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "71849b5a", - "metadata": {}, - "outputs": [], - "source": [ - "# simulate snia lightcurves\n", - "\n", - "lc = hf.simulate_lsst_snia_lightcurve( \n", - " redshift=0.5,\n", - " t_peak_mjd=60000.0, \n", - ") \n", - "hf.plot_lsst_snia_lightcurve(lc, t_peak_mjd=60000.0, redshift=0.5, ymax=20, xmax=50)\n", - " \n", - "# Grid across redshifts\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=True) \n", - "for ax, z in zip(axes.flat, [0.1, 0.3, 0.5, 0.75, 1.0, 1.5]): \n", - " lc = hf.simulate_lsst_snia_lightcurve(redshift=z, t_peak_mjd=60000.0) \n", - " hf.plot_lsst_snia_lightcurve(lc, t_peak_mjd=60000.0, redshift=z, ax=ax, ymax=15, xmax=50) \n", - "plt.tight_layout(); plt.show() " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "aqdz1gio46c", - "metadata": {}, - "outputs": [], - "source": [ - "# Grid of redshifts — Custom CCSN (IIP) moved to different distances\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=True)\n", - "fig.suptitle('Custom CCSN Simulations', fontsize=20, fontweight='bold', y=1.01)\n", - "\n", - "for ax, z in zip(axes.flat, [0.1, 0.3, 0.5, 0.75, 1.0, 1.5]):\n", - " lc = hf.simulate_lsst_ccsn_lightcurve(redshift=z, t_peak_mjd=60000.0)\n", - " hf.plot_lsst_ccsn_lightcurve(lc, t_peak_mjd=60000.0, redshift=z,\n", - " T_k=6000, L_peak=1e43, ax=ax, ymax=18, xmax=250)\n", - "plt.tight_layout()\n", - "plt.show()" - ] - }, - { - "cell_type": "markdown", - "id": "713db385", - "metadata": {}, - "source": [ - "## package simulators" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "9487c195", - "metadata": {}, - "outputs": [], - "source": [ - "import highz_funcs as hf \n", - "importlib.reload(hf)\n", - "import matplotlib.pyplot as plt " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "69f27e71", - "metadata": {}, - "outputs": [], - "source": [ - "# redback TDE lightcurve simulation using tde_analytical model\n", - "\n", - "lc_rb = hf.simulate_redback_tde_lightcurve(\n", - " redshift=0.1,\n", - " t_peak_mjd=60000.0,\n", - " l0=1e55,\n", - " t_0_turn=50.0,\n", - " mej=1.0,\n", - " vej=1e4,\n", - " kappa=0.1,\n", - " kappa_gamma=10.0,\n", - " temperature_floor=5000.0,\n", - " ebv=0.1,\n", - " apply_igm=True,\n", - ")" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "44c10f80", - "metadata": {}, - "outputs": [], - "source": [ - "hf.plot_redback_tde_lightcurve(lc_rb, t_peak_mjd=60000.0, redshift=0.1, l0=1e55, t_0_turn=50.0)" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "xe47rcw0pb", - "metadata": {}, - "outputs": [], - "source": [ - "# Grid of redback TDE lightcurves at different redshifts\n", - "redshifts = [0.1, 0.5, 1.0, 1.5, 2.0, 2.5]\n", - "\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=False)\n", - "fig.suptitle('Redback TDE Simulations', fontsize=20, fontweight='bold', y=1.01)\n", - "\n", - "for ax, z in zip(axes.flat, redshifts):\n", - " lc = hf.simulate_redback_tde_lightcurve(\n", - " redshift=z,\n", - " t_peak_mjd=60000.0,\n", - " l0=1e55,\n", - " t_0_turn=50.0,\n", - " mej=1.0,\n", - " vej=1e4,\n", - " kappa=0.1,\n", - " kappa_gamma=10.0,\n", - " temperature_floor=5000.0,\n", - " ebv=0.1,\n", - " apply_igm=True,\n", - " )\n", - " hf.plot_redback_tde_lightcurve(\n", - " lc, t_peak_mjd=60000.0, redshift=z, l0=1e55, t_0_turn=50.0,\n", - " ax=ax, ymin=32, ymax=14,\n", - " )\n", - " for band, m5 in hf._LSST_M5.items():\n", - " ax.axhline(m5, color=hf._BAND_COLORS[band], lw=0.8, ls='--', alpha=0.5)\n", - "\n", - "plt.tight_layout()\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "9p7srbp3nlc", - "metadata": {}, - "outputs": [], - "source": [ - "# Grid of redback SNIa lightcurves at different redshifts\n", - "redshifts = [0.1, 0.5, 1.0, 1.5, 2.0, 2.5]\n", - "\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=False)\n", - "fig.suptitle('Redback SNIa Simulations', fontsize=20, fontweight='bold', y=1.01)\n", - "\n", - "for ax, z in zip(axes.flat, redshifts):\n", - " lc = hf.simulate_redback_snia_lightcurve(\n", - " redshift=z,\n", - " t_peak_mjd=60000.0,\n", - " f_nickel=0.6,\n", - " mej=1.4,\n", - " kappa=0.1,\n", - " kappa_gamma=10.0,\n", - " vej=1e4,\n", - " temperature_floor=3000.0,\n", - " ebv=0.1,\n", - " apply_igm=True,\n", - " )\n", - " hf.plot_redback_snia_lightcurve(lc, t_peak_mjd=60000.0, redshift=z, ax=ax, ymin=32, ymax=14)\n", - " for band, m5 in hf._LSST_M5.items():\n", - " ax.axhline(m5, color=hf._BAND_COLORS[band], lw=0.8, ls='--', alpha=0.5)\n", - "\n", - "plt.tight_layout()\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "bvzcn1whszc", - "metadata": {}, - "outputs": [], - "source": [ - "# Grid of redback CCSN lightcurves at different redshifts\n", - "redshifts = [0.1, 0.5, 1.0, 1.5, 2.0, 2.5]\n", - "\n", - "fig, axes = plt.subplots(2, 3, figsize=(18, 10), sharey=False)\n", - "fig.suptitle('Redback CCSN Simulations', fontsize=20, fontweight='bold', y=1.01)\n", - "\n", - "for ax, z in zip(axes.flat, redshifts):\n", - " lc = hf.simulate_redback_ccsn_lightcurve(\n", - " redshift=z,\n", - " t_peak_mjd=60000.0,\n", - " f_nickel=0.07,\n", - " mej=8.0,\n", - " kappa=0.1,\n", - " kappa_gamma=10.0,\n", - " vej=8e3,\n", - " temperature_floor=3000.0,\n", - " ebv=0.1,\n", - " apply_igm=True,\n", - " )\n", - " hf.plot_redback_ccsn_lightcurve(lc, t_peak_mjd=60000.0, redshift=z, ax=ax, ymin=32, ymax=14)\n", - " for band, m5 in hf._LSST_M5.items():\n", - " ax.axhline(m5, color=hf._BAND_COLORS[band], lw=0.8, ls='--', alpha=0.5)\n", - "\n", - "plt.tight_layout()\n", - "plt.show()" + "\n", + "plot_lsst_lightcurve(lc, t_peak_mjd=60000.0, redshift=1.0, T_k=30000, L_peak=1e44, ymax=20, xmin=-70, xmax=50) " ] }, { @@ -1109,8 +834,10 @@ "#FIXME intrinsic z=0 values look off. use abs M_r -19.1 for SN1a, -19 for CCSN, and -19 for TDE. use sn1a rise 10 days, ccsn rise 15 days, tde rise 25 days\n", "#FIXME plot distributions not median values\n", "\n", - "# Single panel, r-band, IGM on \n", - "hf.plot_transient_magnitude_timescale( \n", + "# Single panel, r-band, IGM on \n", + "from thor.utils.lightcurve_sims import plot_transient_magnitude_timescale\n", + " \n", + "plot_transient_magnitude_timescale( \n", " redshifts=[0.05, 0.1, 0.2, 0.5, 0.75, 1.0, 1.5, 2.0, 2.5, 3.0], \n", " band='r', \n", " apply_igm=True, \n",