  This paper reviews and analyses methods used to identify neighbours in 6D
space and estimate the corresponding phase-space density. It compares SPH
methods to 6D Delaunay tessellation on statical and dynamical realisation of
single halo profiles, paying attention to the unknown scaling, S_G, used to
relate the spatial dimensions to the velocity dimensions. The methods with
local adaptive metric provide the best phase-space estimators. They make use of
a Shannon entropy criterion combined with a binary tree partitioning and with
SPH interpolation using 10-40 neighbours. Local scaling implemented by such
methods, which enforces local isotropy of the distribution function, can vary
by about one order of magnitude in different regions within the system. It
presents a bimodal distribution, in which one component is dominated by the
main part of the halo and the other one is dominated by the substructures.
While potentially better than SPH techniques, since it yields an optimal
estimate of the local softening volume (and the local number of neighbours
required to perform the interpolation), the Delaunay tessellation in fact
poorly estimates the phase-space distribution function. Indeed, it requires,
the choice of a global scaling S_G. We propose two methods to estimate S_G that
yield a good global compromise. However, the Delaunay interpolation still
remains quite sensitive to local anisotropies in the distribution. We also
compare 6D phase-space density estimation with the proxy, Q=rho/sigma^3, where
rho is the local density and 3 sigma^2 is the local 3D velocity dispersion. We
show that Q only corresponds to a rough approximation of the true phase-space
density, and is not able to capture all the details of the distribution in
phase-space, ignoring, in particular, filamentation and tidal streams.
