import yprof_halomodel as yh
import yprof_battaglia as yb
import numpy as np
import matplotlib.pyplot as plt

theta = 0.5
M200c = 3e15
z     = 1.0

yhf = np.vectorize(yh.y)
ybf = np.vectorize(yb.y)

y_table = yhf(theta,M200c,z)[0]
y_batta = ybf(theta,M200c,z)

theta1d = np.logspace(-4,1,100)

y1d_table = yhf(theta1d,M200c,z)[0]
y1d_batta = ybf(theta1d,M200c,z)

print('y_integration / y_table at center: ',y1d_batta[0]/y1d_table[0])

plt.loglog(theta1d,y1d_table,label='websky table')
plt.loglog(theta1d,y1d_batta,label='integration')

plt.title(r'$M_{\rm 200c}=$3e15, $z=1$')

plt.gca().set_xlabel(r'$\theta$ [arcmin]')
plt.gca().set_ylabel(r'$y$')

plt.legend()
plt.savefig('comparison_3e15_1.png',dpi=300)
