  We apply for the first time a two-dimensional fitting statistic, $\tau^2$, to
rotational evolution models (REMs) of stars (0.1 to 1.3 $M_{\odot}$) on the
period-mass plane. The $\tau^2$ statistic simultaneously considers all cluster
rotation data to return a goodness of fit, allowing for data-driven improvement
of REMs. We construct data sets for Upper Sco, the Pleiades and Praesepe, to
which we tune our REMs. We use consistently determined stellar masses
(calculated by matching $K_\textrm{s}$ magnitudes to isochrones) and literature
rotation periods. As a first demonstration of the $\tau^2$ statistic, we find
the best-fitting gyrochronology age for Praesepe, which is in good agreement
with the literature. We then systematically vary three parameters which
determine the dependence of our stellar wind torque law on Rossby number in the
saturated and unsaturated regimes, and the location of the transition between
the two. By minimising $\tau^2$, we find best-fit values for each parameter.
These values vary slightly between clusters, mass determinations and initial
conditions, highlighting the precision of $\tau^2$ and its potential for
constraining REMs, gyrochronology, and understanding of stellar physics. Our
resulting REMs, which implement the best-possible fitting form of a broken
power-law torque, are statistically improved on previous REMs using similar
formulations, but still do not simultaneously describe the observed rotation
distributions of the lowest masses, which have both slow and fast rotators by
the Praesepe age, and the shape of the converged sequence for higher masses.
Further complexity in the REMs is thus required to accurately describe the
data.
