from yprof_halomodel.py import * 

ell, Dl, ybar, ybartab = powerspectrum()
ybar=np.array([ybar],dtype=np.float64)
np.savez('halomodel.npz',ell=ell,Dl=Dl,ybar=ybar,ybartab=ybartab)
print("ybar, ybartab = ",ybar,ybartab)
import matplotlib.pyplot as plt
plt.loglog(ell,Dl/ybar**2,label='power spectrum using Tinker08 + websky y(theta|M,z')
plt.gca().set_ylabel(r'$10^{12}c_{\ell}^{yy}\ell(\ell+1)/(2\pi)/<y>^2$')
plt.gca().set_xlabel(r'$\ell$')
plt.show()