Astropy cosmology - how to plot particle and event horizon?

users69S 0
Viewed 20

I've been playing with astropy.cosmology, trying to make a plot similar to the second plot in this link: http://www.mysearch.org.uk/website1/html/512.Horizons.html If the link doesn't work, this is the plot I am trying to replicate:

plot to replicate

The code:

from astropy.cosmology import Planck18, z_at_value
from astropy import units as u
import numpy as np
from matplotlib import pyplot as plt
from scipy import interpolate
import astropy.constants as cc

def Hradius(zed):
    c = cc.c
    return (c/cosmo.H(zed)).to(u.Glyr)

#define cosmology
cosmo = Planck18

# Z and time
redshift = np.arange(0,1089,0.1)
time=cosmo.age(redshift)

# horizons
hubbleRadius = [Hradius(i).value for i in redshift]
lightPath = [cosmo.scale_factor(i) * cosmo.comoving_distance(i).to(u.Glyr).value for i in redshift]
particleHorizon = [cosmo.comoving_distance(i).to(u.Glyr).value for i in redshift]

#plot
fig, ax1 = plt.subplots()



ax1.plot(time, hubbleRadius, 'g-', label = 'Hubble radius')
ax1.plot(time, lightPath, 'y-', label = 'Light path')
ax1.plot(time, particleHorizon, 'b-', label = 'Particle horizon')

ax1.invert_xaxis()
ax1.legend()


ax1.set_xlabel('Time (Glyr)')
ax1.set_ylabel('Proper distance (Glyr)')

The result: enter image description here

Don't know why the particle horizon is reversed ... any ideas? Also how do I add the event horizon?

astropy python
1 Answers

I think it works ok now. I've added this:

redshift_R = redshift[::-1] #reverse order of redshift

term1 = [cosmo.scale_factor(i) for i in redshift]
term2 = [cosmo.comoving_distance(i).to(u.Glyr).value for i in redshift_R]
particleHorizon = [term1[i]*term2[i] for i in range(len(redshift))]

The result:

enter image description here

users69S0
Related