  We present a method to numerically estimate the densities of a discretely
sampled data based on binary space partitioning tree. We start with a root node
containing all the particles and then recursively divide each node into two
nodes each containing roughly equal number of particles,until each of the nodes
contains only one particle. The volume of such a leaf node provides an estimate
of the local density. We implement an entropy-based node splitting criterion
that results in a significant improvement in the estimation of densities
compared to earlier work. The method is completely metric free and can be
applied to arbitrary number of dimensions. We apply this method to determine
the phase space densities of dark matter halos obtained from cosmological
N-body simulations. We find that contrary to earlier studies, the volume
distribution function $v(f)$ of phase space density $f$ does not have a
constant slope but rather a small hump at high phase space densities. We
demonstrate that a model in which a halo is made up by a superposition of
Hernquist spheres is not capable in explaining the shape of $v(f)$ vs $f$
relation, whereas a model which takes into account the contribution of the main
halo separately roughly reproduces the behavior as seen in simulations. The use
of the presented method is not limited to calculation of phase space densities,
but can be used as a general-purpose data-mining tool and due to its speed and
accuracy it is ideally suited for analysis of large multidimensional data sets.
