def METHOD_NAME(self,z=[],l=np.arange(3000)+1,pk_params=None,cosmo_h=None,cosmo=None,pk_func=None):
    if not cosmo_h:
        cosmo_h=self.cosmo_h
    nz=len(z)
    nl=len(l)
    pk,kh=pk_func(z=z,pk_params=pk_params)
    cls=np.zeros((nz,nl),dtype='float32')*u.Mpc**2
    for i in np.arange(nz):
        DC_i=cosmo_h.comoving_transverse_distance(z[i]).value#because camb k in h/mpc
        lz=kh*DC_i-0.5
        DC_i=cosmo_h.comoving_transverse_distance(z[i]).value
        pk_int=interp1d(lz,pk[i]/DC_i**2,bounds_error=False,fill_value=0)
        cls[i][:]+=pk_int(l)*u.Mpc*(c/(cosmo_h.efunc(z[i])*cosmo_h.H0))
    return cls