  We devise fast and provably accurate algorithms to transform between an
$N\times N \times N$ Cartesian voxel representation of a three-dimensional
function and its expansion into the ball harmonics, that is, the eigenbasis of
the Dirichlet Laplacian on the unit ball in $\mathbb{R}^3$. Given $\varepsilon
> 0$, our algorithms achieve relative $\ell^1$ - $\ell^\infty$ accuracy
$\varepsilon$ in time $O(N^3 (\log N)^2 + N^3 |\log \varepsilon|^2)$, while
their dense counterparts have time complexity $O(N^6)$. We illustrate our
methods on numerical examples.
