Skip to content

Commit

Permalink
fix the notebooks
Browse files Browse the repository at this point in the history
  • Loading branch information
Joran Angevaare committed Aug 17, 2022
1 parent 4c00219 commit d2d5544
Show file tree
Hide file tree
Showing 2 changed files with 43 additions and 23 deletions.
56 changes: 37 additions & 19 deletions notebooks/Checks, plots.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,10 @@
"source": [
"try:\n",
" plt.style.use('latex_thesis')\n",
" _latex_thesis=True\n",
"except (FileNotFoundError, OSError):\n",
" print(\"Can't use my favourite plotting style... oh well\")"
" print(\"Can't use my favourite plotting style... oh well\")\n",
" _latex_thesis=False"
]
},
{
Expand All @@ -48,7 +50,8 @@
"source": [
"kms = nu.km/nu.s\n",
"\n",
"from wimprates import v_0, v_earth, v_esc \n",
"from wimprates import StandardHaloModel, v_earth\n",
"v_0, v_esc = StandardHaloModel().v_0, StandardHaloModel().v_esc\n",
"\n",
"def galactic_v_dist(v):\n",
" # Maxwell Boltzmann distribution with\n",
Expand Down Expand Up @@ -393,10 +396,14 @@
],
"source": [
"# Get mean energies and energy spectra\n",
"from laidbax import wimps\n",
"from s2only import wimprates as wr\n",
"try:\n",
" from laidbax import wimps\n",
" from s2only import wimprates as wr\n",
"except ImportError:\n",
" import wimprates as wr\n",
" wimps = None\n",
"\n",
"es = np.logspace(-2, 3.5, 500)\n",
"es = np.logspace(-2, 3.5, 100)\n",
"ms = mchi / gevcsq\n",
"rates = []\n",
"rates_bs = []\n",
Expand All @@ -407,7 +414,10 @@
"means_bs = np.zeros(len(ms), dtype=np.float)\n",
"\n",
"for i, mw in enumerate(tqdm(ms)):\n",
" r = wimps.wimp_recoil_spectrum(es, mass=mw, sigma=1e-45)\n",
" if wimps is not None:\n",
" r = wimps.wimp_recoil_spectrum(es, mass=mw, sigma=1e-45)\n",
" else:\n",
" r = wr.rate_wimp_std(es, mw, 1e-45)\n",
" rates.append(r)\n",
" rate_per_bin = (r[1:] + r[:-1])/2 * widths\n",
" if rate_per_bin.sum() != 0:\n",
Expand Down Expand Up @@ -435,9 +445,13 @@
"metadata": {},
"outputs": [],
"source": [
"from s2only.utils import load_pickle, dump_pickle\n",
"dump_pickle(dict(means=means, means_bs=means_bs, rates=rates, rates_bs=rates_bs, ms=ms),\n",
" 'dm_mean_energies.pickle')"
"try:\n",
" from s2only.utils import load_pickle, dump_pickle\n",
"except ModuleNotFoundError:\n",
" pass\n",
"else:\n",
" dump_pickle(dict(means=means, means_bs=means_bs, rates=rates, rates_bs=rates_bs, ms=ms),\n",
" 'dm_mean_energies.pickle')"
]
},
{
Expand Down Expand Up @@ -574,13 +588,13 @@
"tyear = 1e3 * 365\n",
"\n",
"plt.plot(ms, [tyear * total_rate(m, threshold=5) for m in ms], \n",
" label='NR, $E_T =\\SI{5}{keV}$', \n",
" label='NR, $E_T =\\SI{5}{keV}$' if _latex_thesis else 'NR, $E_T =5$ keV', \n",
" c='blue', linewidth=1)\n",
"plt.plot(ms, [tyear * total_rate(m, threshold=1) for m in ms], \n",
" label='NR, $E_T = \\SI{1}{keV}$', \n",
" label='NR, $E_T = \\SI{1}{keV}$' if _latex_thesis else 'NR, $E_T = 1$ keV', \n",
" c='g', linewidth=1)\n",
"plt.plot(ms, [tyear * total_rate(m, bs=True, threshold=0.18) * 1e-10 for m in ms], \n",
" label='BS, $E_T =\\SI{180}{eV}$',\n",
" label='BS, $E_T =\\SI{180}{eV}$' if _latex_thesis else 'BS, $E_T =180$ eV',\n",
" c='r', linewidth=1)\n",
"plt.yscale('log')\n",
"plt.ylim(None, 1e4)\n",
Expand All @@ -589,8 +603,12 @@
"leg = plt.legend(loc='upper right')\n",
"leg.get_frame().set_linewidth(0.0)\n",
"leg.get_frame().set_alpha(0.85)\n",
"plt.ylabel('Rate (\\si{events/(ton\\;year)})')\n",
"plt.xlabel('$m_\\chi$ (\\si{GeV/c^2})')\n",
"if _latex_thesis:\n",
" plt.ylabel('Rate (\\si{events/(ton\\;year)})')\n",
" plt.xlabel('$m_\\chi$ (\\si{GeV/c^2})')\n",
"else:\n",
" plt.ylabel('Rate')\n",
" plt.xlabel('$m_\\chi$')\n",
"plt.gca().invert_yaxis()\n",
"#plt.plot(es, total_rate(, c='b', label='$\\SI{1}{GeV/c^2}$')\n",
"plt.xlim(1e-1, 1e4)\n",
Expand All @@ -615,10 +633,10 @@
}
],
"source": [
"plt.plot(es, nearest_rate(1), c='b', label='$\\SI{1}{GeV/c^2}$')\n",
"plt.plot(es, nearest_rate(100), label='$\\SI{100}{GeV/c^2}$', c='g')\n",
"plt.plot(es, nearest_rate(1), c='b', label='$\\SI{1}{GeV/c^2}$' if _latex_thesis else None)\n",
"plt.plot(es, nearest_rate(100), label='$\\SI{100}{GeV/c^2}$' if _latex_thesis else None, c='g')\n",
"plt.plot(es, nearest_rate(1, bs=True) * 1e-10, c='b', linestyle=':',\n",
" label='$\\SI{1}{GeV/c^2}$, Bremsstrahlung' # $\\\\times 10^{10}$'\n",
" label='$\\SI{1}{GeV/c^2}$, Bremsstrahlung' if _latex_thesis else None # $\\\\times 10^{10}$'\n",
" )\n",
"leg = plt.legend(loc='upper right', frameon=False)\n",
"#leg.get_frame().set_linewidth(0.0)\n",
Expand Down Expand Up @@ -695,7 +713,7 @@
"metadata": {
"anaconda-cloud": {},
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -709,7 +727,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.7"
"version": "3.10.4"
}
},
"nbformat": 4,
Expand Down
10 changes: 6 additions & 4 deletions notebooks/DM-electron scattering.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -325,11 +325,13 @@
" \n",
" for shell, rn in drsn.items():\n",
" rn *= (1000 * nu.kg * nu.year)\n",
" plt.plot(n_el, rn, linestyle='steps-mid', label=shell)\n",
" plt.plot(n_el, rn, drawstyle='steps-mid', label=shell)\n",
"\n",
" plt.plot(n_el, np.sum(list(drsn.values()), axis=0),\n",
" label='Total',\n",
" drawstyle='steps-mid', linestyle='--', c='k')\n",
" drawstyle='steps-mid', \n",
" linestyle='--', \n",
" c='k')\n",
"\n",
" plt.title(title + (' -- SWAP 4s<->4p' if do_swap else ''))\n",
" plt.legend(loc='upper right', ncol=2)\n",
Expand All @@ -355,7 +357,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -369,7 +371,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.7"
"version": "3.10.4"
}
},
"nbformat": 4,
Expand Down

0 comments on commit d2d5544

Please sign in to comment.