An approximate error function for the discretization error on a given mesh is obtained by projecting (via the energy inner product) the functional residual onto the space of continuous, piecewise quadratic functions which vanish on the vertices of the mesh. Conditions are given under which one can expect this hierarchical basis error estimator to give efficient and reliable function recovery, asymptotically exact gradient recovery and convergent Hessian recovery in the square norms. One does not find similar function recovery results in the literature. The analysis given here is based on a certain superconvergence result which has been used elsewhere in the analysis of gradient recovery methods. Numerical experiments are provided which demonstrate the effectivity of the approximate error function in practice for both the square norms and the absolute value norms.