In [3]:

cd ~/Dropbox/notebooks/ghost_halos/




/Users/lia/Dropbox/notebooks/ghost_halos




In [4]:

import numpy as np
import matplotlib.pyplot as plt
#rcParams['font.size'] = 12
font = {'size'   : 15}
matplotlib.rc('font', **font)



## Load the scattering halo intensity (normalized by point source flux)



In [5]:

from astropy.io import fits




In [6]:

try:
hdulist = fits.open(fname)
except:
print "Cannot open "+fname
return
tdata = hdulist[1].data
hdulist.close()
return tdata['Obsangle'], tdata['Intensity'] # arcsec, arcsec^-2




In [7]:



## Surface brightness from a RASS BSC quasar

Minimal flux to get into RASS Bright Sources Catalog: $5 \times 10^{-13}$ cgs

From IGM Dust Calculation (evernote): Chandra 0.04 ct/s = $3 \times 10^{-13}$ cgs @ 1 keV

Estimate of Chandra count rate: $$R_{Chandra} = F \ \frac{0.04\ {\rm ct/s}}{3 \times 10^{-13}\ {\rm cgs}}$$

The typical Chandra 0.5-2 keV background is 0.1 ct/s/chip = $4 \times 10^{-7}$ ct/s/arcsec$^2$



In [8]:

RASSBSC_flux = 5.e-13 # cgs
RChandra     = 0.04 / 3.e-13 # ct/s per flux (cgs)

SXRB         = 4.e-7 # ct/s/arcsec^2

src_rate  = RASSBSC_flux * RChandra  # Minimally bright quasars
bsrc_rate = src_rate * 10.0 # Maximally bright quasars, photon flux around 5.e-12 cgs
min_rate  = 0.04 # ct/s from Table 1 of thesis, Ch 5

print "Chandra count rate for min RASS BSC flux: ~", src_rate, "ct/s"




Chandra count rate for min RASS BSC flux: ~ 0.0666666666667 ct/s




In [9]:

def plot_halos_vs_bkg(ax, xlim=(0.1,1.e3), ylim=(3.e-11,1.e-3)):

# Plot IGM dust halos
# Brightest
ax.plot(alpha, grey_igm_halo * bsrc_rate, 'k-', lw=3, label='Grey IGM')
ax.plot(alpha, big_igm_halo  * bsrc_rate, 'k-', lw=1, label=r'$1\ \mu m$ IGM')

# Moderately bright enough
ax.plot(alpha, grey_igm_halo * min_rate, color='0.6', ls='-', lw=3, label='')
ax.plot(alpha, big_igm_halo  * min_rate, color='0.6', ls='-', lw=1, label='')

# Plot the typical background
ax.axhline(SXRB, color='k', ls=':', label='')
ax.text(40.0, 6.e-7, '0.5-2 keV background', color='0.3', size=12)

plt.loglog()
ax.set_xlim(xlim)
ax.set_ylim(ylim)

plt.legend(loc='upper right', frameon=False, fontsize=15)
plt.xlabel('Observation Angle [arcsec]', size=15)
plt.ylabel('Halo Brightness [ct s$^{-1}$ arcsec$^{-2}$]', size=15)

plt.text(1.0, 1.e-8, r'$F_{\rm src} =\ 3 \times\ 10^{-13}\ {\rm cgs}$', color='0.3', size=12)
plt.text(1.0, 1.5e-4, r'$F_{\rm src} =\ 5 \times\ 10^{-12}\ {\rm cgs}$', color='k', size=12)

return




In [10]:

ax = plt.subplot(111)
plot_halos_vs_bkg(ax)






## Time delay associated with such a halo

The basic formula and some calculations

Corrales & Paerels (2012): $$\Delta \chi_{\rm max} \approx 4.5\ {\rm kpc}\ \chi_{Gpc} E^{-2}_{keV} a_{0.1}^{-2}$$

And in general (Corrales & Paerels 2012): $$\Delta \chi = \chi_{OS} \alpha^2 (1-x)$$

Keep in mind that the minimal x value is: $$x_{min} \approx \frac{\alpha}{2 \sigma}$$ where $$\sigma = \frac{10.4'}{a(0.1\ \mu{\rm m})\ E({\rm keV})}$$

Fade time for a total halo ($\alpha = \sigma$)

The furthest you can look back along a sight line is $x_{min} \approx \alpha/2\sigma$

The largest observable halo angle is $\sim \sigma$, so the total time it takes for a halo to fade is: $$t_{fade} \approx \frac{\chi_{OS}}{c} \frac{\sigma^2}{2}$$

Fade time as a function of angle

Using the fact that $x_{min} \approx \alpha/2\sigma$, I can calculate the time delay for the longest ($\chi_{max}$) value at a particular sight line ($\alpha$).

$$\Delta t_{max} = \frac{\chi_{OS}}{c} \alpha^2 (1 - \frac{\alpha}{2 \sigma})$$


In [11]:

import cosmology as cosm




In [12]:

## Units
C      = 3.e10 # cm/s
GPC    = 1.e9 * 3.e18 # cm
YEAR   = 3.e7 # seconds

## Other params
ZS = 2.0
COSMO  = cosm.Cosmology(d=3.e-6)




In [13]:

def charsig(a, E):
""" returns the characteristic scattering angle, in arcseconds
a: [um] Grain size
E: [keV] Photon energy"""
charsig0 = 1.04 * 60 # arcsec per 1 um grain size per keV
return charsig0 / a / E

def tmax(alpha, a=0.1, E=1.0, zs=ZS):
""" Calculate the maximum time delay associated with
observation angle and typical grain population
alpha: [arcmin] Observation angle
a: [micron] Grain size
E: [keV] Photon energy
zs: Redshift of source"""
time_os   = cosm.DChi(zs, 0.0) * (GPC/C) / YEAR
xterm     = 1 - alpha/(2*charsig(a,E))
return time_os * alpha_rad**2 * xterm




In [21]:

def plot_time_delay(ax, xlim=(0.1,1.e3), ylim=(-3,5)):

tbig = tmax(alpha, a=1.0)
tsma = tmax(alpha, a=0.1)

ipos_big = np.where(tbig > 0.0)[0]
ipos_sma = np.where(tsma > 0.0)[0]

ax.plot(alpha[ipos_sma], np.log10(tsma[ipos_sma]), ls='-', color='k', lw=1, label=r'$0.1\ \mu m$ grains')
ax.plot(alpha[ipos_big], np.log10(tbig[ipos_big]), ls='--', color='k', lw=2, label=r'$1.0\ \mu m$ grains')

plt.semilogx()
ax.set_xlim(xlim)
ax.set_ylim(ylim)

plt.legend(loc='upper left', frameon=False, fontsize=12)
plt.xlabel("Observation angle [arcsec]", size=15)
plt.ylabel(r"$\log\ \delta t_{\rm max}$  [yrs]", size=15)

return




In [22]:

ax = plt.subplot(111)
plot_time_delay(ax)






## Final plot for paper



In [16]:

import matplotlib.gridspec as gridspec




In [24]:

XLIM = (1.0,1.e3)

fig = plt.figure(figsize=(6,6))
gs  = gridspec.GridSpec(2, 1, height_ratios=[2,1])

ax1 = plt.subplot(gs[0])
plot_halos_vs_bkg(ax1)

plt.xlabel('')
ax1.xaxis.set_ticklabels('')

ax2 = plt.subplot(gs[1])
plot_time_delay(ax2, ylim=(-2,5.5))

gs.update(hspace=0.0)

#plt.savefig('figure03.pdf', format='pdf')







In [17]: