

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]:

    
def load_halo(fname):
    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]:

    
alpha, grey_igm_halo = load_halo('halo_uIGM_grey_gra.fits')
alpha, big_igm_halo  = load_halo('halo_uIGM_gra.fits')

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
ARCSEC = 5.e-6 # radians
ARCMIN = 3.e-4 # radians
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
    alpha_rad = alpha * ARCSEC
    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]: