In [1]:
import numpy as np
import matplotlib.pyplot as plt
font = {'size' : 15}
matplotlib.rc('font', **font)
In [2]:
cp /Users/lia/Academic/CygX-3/6601_paper1/best_fits.py /Users/lia/Dropbox/notebooks/cygx3_redo/best_fits.py
In [3]:
import sigma_scat as ss
import dust
import constants as c
In [4]:
NH, d2g = 1.e21, 0.009
dust_mass = NH * c.mp() * d2g
ERANGE = np.power( 10.0, np.arange(-0.6,1.0,0.005) )
In [5]:
RHO_SIL, RHO_GRA, RHO_AVG = 3.8, 2.2, 3.0 # g cm^-3; see Draine book
In [6]:
MRN_RANGE = np.arange(0.005,0.25001,0.05)
MRN_sil = dust.Dustdist( rad=MRN_RANGE, p=3.5, rho=RHO_SIL )
MRN_gra = dust.Dustdist( rad=MRN_RANGE, p=3.5, rho=RHO_GRA )
MRN_avg = dust.Dustdist( rad=MRN_RANGE, p=3.5, rho=RHO_AVG )
In [7]:
# choices motivated by crazy fit
BIG_RANGE = np.arange(0.005, 1.5, 0.05)
BIG_sil = dust.Dustdist( rad=BIG_RANGE, p=3.5, rho=RHO_SIL )
BIG_gra = dust.Dustdist( rad=BIG_RANGE, p=3.5, rho=RHO_GRA )
BIG_avg = dust.Dustdist( rad=BIG_RANGE, p=3.5, rho=RHO_AVG )
In [8]:
## Rayleigh-Gans plus Drude approximation
RGD_mrn = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=MRN_avg, md=dust_mass), scatm=ss.makeScatmodel('RG','Drude') )
RGD_big = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=BIG_avg, md=dust_mass), scatm=ss.makeScatmodel('RG','Drude') )
In [9]:
## Mie scattering for the small grain MRN distribution
Mie_mrn_sil = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=MRN_sil, md=dust_mass), scatm=ss.makeScatmodel('Mie','Silicate') )
Mie_mrn_gra = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=MRN_gra, md=dust_mass), scatm=ss.makeScatmodel('Mie','Graphite') )
In [10]:
## Mie scattering for the slightly bigger grain distribution
Mie_big_sil = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=BIG_sil, md=dust_mass), scatm=ss.makeScatmodel('Mie','Silicate') )
Mie_big_gra = ss.Kappascat( E=ERANGE, dist=dust.Dustspectrum(rad=BIG_gra, md=dust_mass), scatm=ss.makeScatmodel('Mie','Graphite') )
In [15]:
MD = dust_mass * 10.0 # for N_H = 1.e22
plt.plot( RGD_mrn.E, RGD_mrn.kappa * MD, '0.4', lw=2, label='RG-Drude' )
plt.plot( Mie_mrn_sil.E, Mie_mrn_sil.kappa * MD, 'g', lw=2, label='Mie-Silicate' )
plt.plot( Mie_mrn_gra.E, Mie_mrn_gra.kappa * MD, 'b', lw=2, label='Mie-Graphite' )
plt.legend( loc='upper right', fontsize=12 )
plt.plot( RGD_big.E, RGD_big.kappa * MD, '0.4', lw=1, ls='-' )
plt.plot( Mie_big_sil.E, Mie_big_sil.kappa * MD, 'g', lw=1, ls='-' )
plt.plot( Mie_big_gra.E, Mie_big_gra.kappa * MD, 'b', lw=1, ls='-' )
plt.loglog()
plt.xlim(0.3,10)
plt.ylim(1.e-2,10.0)
plt.xlabel( "Energy [keV]", size=15 )
plt.ylabel( r"Scattering Opacity [$\tau$ per $N_{\rm H}/10^{22}$]", size=15 )
plt.text( 0.5, 0.1, '$0.25\ \mu m$\ncut-off', size=12)
plt.text( 1.2, 3.0, '$1.5\ \mu m$\ncut-off', size=12)
#plt.savefig('/Users/lia/Dropbox/Writing/CygX3-paper1/figure5.pdf', format='pdf')
In [12]:
from scipy.interpolate import interp1d
In [13]:
kappa_interp = interp1d( RGD_mrn.E, RGD_mrn.kappa )
kappa_1keV = kappa_interp(1.0)
print "Tau (1 keV) for RG-Drude with MRN:", kappa_1keV
In [14]:
plt.plot( RGD_mrn.E, RGD_mrn.kappa )
plt.plot( 1.0, kappa_1keV, 'ko', markersize=10 )
plt.loglog()
Out[14]:
In [14]: