Emission Measure Calculations



In [2]:

    
#%matplotlib notebook



In [3]:

    
# imports
from importlib import reload
import numpy as np
from matplotlib import pyplot as plt

from astropy import units
from astropy import constants

from frb import em
from frb import halos

EM from H$\alpha$

Follow Reynolds 1977



In [4]:

    
obs_Ha = 1e-17 * units.erg / units.cm**2 / units.s / units.arcsec**2



In [5]:

    
Ha = 6564. * units.Angstrom
E_Ha_photon = constants.c * constants.h / Ha
E_Ha_photon.to('erg')









    Out[5]:





$3.0262733 \times 10^{-12} \; \mathrm{erg}$



In [6]:

    
(obs_Ha*units.ph/E_Ha_photon).decompose().to('rayleigh') * 68









    Out[6]:





$120.1326 \; \mathrm{R}$



In [7]:

    
(7.96e8 * E_Ha_photon / units.s / units.m**2 / units.sr).to('erg/s/cm**2/sr')









    



WARNING: UnitsWarning: 'erg/s/cm**2/sr' contains multiple slashes, which is discouraged by the FITS standard [astropy.units.format.generic]






    Out[7]:





$2.4089136 \times 10^{-7} \; \mathrm{\frac{erg}{s\,sr\,cm^{2}}}$

FRB 121102



In [8]:

    
reload(em)
em_121102 = em.em_from_halpha(6.8e-16*units.erg/units.cm**2/units.s/units.arcsec**2, 0.1927)
em_121102









    Out[8]:





$668.58868 \; \mathrm{\frac{pc}{cm^{6}}}$

FRB 180924



In [9]:

    
Ha_total = 28.1 * 1e-17 * units.erg/units.s/units.cm**2

Assume 5% of the average surface brightness at the FRB; this should be conservative



In [10]:

    
Ha_180924 = 0.05 * Ha_total / units.arcsec**2
Ha_180924









    Out[10]:





$1.405 \times 10^{-17} \; \mathrm{\frac{erg}{s\,{}^{\prime\prime}^{2}\,cm^{2}}}$



In [11]:

    
EM_180924 = em.em_from_halpha(Ha_180924, 0.3214)
EM_180924









    Out[11]:





$20.813195 \; \mathrm{\frac{pc}{cm^{6}}}$

DM from EM -- Reynolds and Cordes

FRB 121102



In [12]:

    
reload(em)
DM_s = em.dm_from_em(em_121102, 1*units.kpc)
DM_s









    Out[12]:





$408.52143 \; \mathrm{\frac{pc}{cm^{2}}}$

FRB 180924



In [13]:

    
DM_s_180924 = em.dm_from_em(EM_180924, 0.1*units.kpc)
DM_s_180924/(1+0.32)









    Out[13]:





$17.267553 \; \mathrm{\frac{pc}{cm^{2}}}$

Fooling around in the Halo

Emissivity -- $j_\nu \sim n_e^2 \, \exp[-(h\nu - Z^2 h \nu_0/n^2)/kT] / T^{3/2}$



In [14]:

    
mw = halos.MilkyWay()

Plot $n_e^2$



In [15]:

    
xyz = np.zeros((3,100))
Z = np.linspace(10,300,100)
xyz[2,:] = Z



In [16]:

    
nesq = mw.ne(xyz)**2



In [17]:

    
# 
plt.clf()
ax = plt.gca()
ax.plot(Z, nesq)
# Labels
ax.set_xlabel('Z Distance (kpc)')
ax.set_ylabel(r'$n_e^2 \; (\rm cm^{-3})$')
plt.show()

Plot Emission per radial shell



In [18]:

    
plt.clf()
ax = plt.gca()
ax.plot(Z, nesq* Z**2)
# Contrast with ISM (ignoring Temperature which reduces the halo further)
#ax.plot([0., 300], [0.1**2 * 10.**2]*2, 'b--')
# Labels
ax.set_xlabel('Z Distance (kpc)')
ax.set_ylabel(r'$n_e^2 \, Z^2 \; (\rm cm^{-1})$')
#
plt.show()



In [ ]: