  Among dynamical modeling techniques, the made-to-measure (M2M) method for
modeling steady-state systems is among the most flexible, allowing
non-parametric distribution functions in complex gravitational potentials to be
modeled efficiently using N-body particles. Here we propose and test various
improvements to the standard M2M method for modeling observed data, illustrated
using the simple setup of a one-dimensional harmonic oscillator. We demonstrate
that nuisance parameters describing the modeled system's orientation with
respect to the observer---e.g., an external galaxy's inclination or the Sun's
position in the Milky Way---as well as the parameters of an external
gravitational field can be optimized simultaneously with the particle weights.
We develop a method for sampling from the high-dimensional uncertainty
distribution of the particle weights. We combine this in a Gibbs sampler with
samplers for the nuisance and potential parameters to explore the uncertainty
distribution of the full set of parameters. We illustrate our M2M improvements
by modeling the vertical density and kinematics of F-type stars in Gaia DR1.
The novel M2M method proposed here allows full probabilistic modeling of
steady-state dynamical systems, allowing uncertainties on the non-parametric
distribution function and on nuisance parameters to be taken into account when
constraining the dark and baryonic masses of stellar systems.
