  Accurate modelling of redshift-space distortions (RSD) is challenging in the
non-linear regime for two-point statistics e.g. the two-point correlation
function (2PCF). We take a different perspective to split the galaxy density
field according to the local density, and cross-correlate those densities with
the entire galaxy field. Using mock galaxies, we demonstrate that combining a
series of cross-correlation functions (CCFs) offers improvements over the 2PCF
as follows: 1. The distribution of peculiar velocities in each split density is
nearly Gaussian. This allows the Gaussian streaming model for RSD to perform
accurately within the statistical errors of a ($1.5\,h^{-1}$Gpc)$^3$ volume for
almost all scales and all split densities. 2. The PDF of the density field at
small scales is non-Gaussian, but the CCFs of split densities capture the
non-Gaussianity, leading to improved cosmological constraints over the 2PCF. We
can obtain unbiased constraints on the growth parameter $f\sigma_{12}$ at the
per-cent level, and Alcock-Paczynski (AP) parameters at the sub-per-cent level
with the minimal scale of $15\,h^{-1}{\rm Mpc}$. This is a $\sim$30 per cent
and $\sim$6 times improvement over the 2PCF, respectively. The diverse and
steep slopes of the CCFs at small scales are likely to be responsible for the
improved constraints of AP parameters. 3. Baryon acoustic oscillations (BAO)
are contained in all CCFs of split densities. Including BAO scales helps to
break the degeneracy between the line-of-sight and transverse AP parameters,
allowing independent constraints on them. We discuss and compare models for RSD
around spherical densities.
