import numpy as np
import scipy.integrate as integrate
import matplotlib.pyplot as plt
plt.rcParams['font.family'] = "serif"
plt.rcParams['font.serif'] = ["Source Han Serif JP", "DejaVu Serif", ]
plt.rcParams["mathtext.fontset"] = "cm"
plt.rcParams["font.size"] = 12
c_per_H0 = 2.998e3 #[Mpc/h]
Omega_mat = 0.315
Omega_rad = 9.14E-5
Omega_lam = 1. - Omega_mat - Omega_rad
z_max = 7
def E(z):
'''Hubble parameter per Hubble constant H(z)/H0'''
return np.sqrt( Omega_lam + Omega_mat*(1.+z)**3 + Omega_rad*(1.+z)**4 )
z = np.linspace(0, z_max, 64)
def comoving(z):
result, err = integrate.quad(lambda x: 1./E(x), 0, z)
return c_per_H0 * result
comoving = np.vectorize(comoving)
x = comoving(z)
dl = (1. + z) * x
da = x / (1. + z)
fig = plt.figure(figsize=(6,4), )
ax = fig.add_subplot(111)
ax.plot(z, x, c="#005aff", ls="-", label="共動距離 $x$")
ax.plot(z, dl, c="#ff4b00", ls="-.", label="光度距離 $d_L$")
ax.plot(z, da, c="#03af7a", ls=":", label="角径距離 $d_A$")
ax.set_xlim([0, z_max])
ax.set_ylim([0, 7e3])
fontsize = 13
ax.set_xlabel(r"赤方偏移 $z$", fontsize=fontsize)
ax.set_ylabel(r"距離 [$\mathrm{Mpc}/h$]", fontsize=fontsize)
ax.grid()
ax.legend(loc=7)
left=0.14
plt.subplots_adjust(bottom=0.15, left=left, right=1-left)
plt.plot()
#plt.show()
plt.savefig("redshift_vs_cosmological_distance_jpn.svg")
plt.close()