  Binary population synthesis calculations and associated predictions,
especially event rates, are known to depend on a significant number of input
model parameters with different degrees of sensitivity. At the same time, for
systems with relatively low formation rates, such simulations are heavily
computationally demanding and therefore the needed explorations of the
high-dimensional parameter space require major -- often prohibitive --
computational resources. In the present study, to better understand several key
event rates involving binary evolution and binaries with two compact objects in
Milky Way-like galaxies and to provide ways of reducing the computational costs
of complete parameter space explorations: (i) we perform a methodical parameter
study of the \emph{StarTrack} population synthesis code ; and (ii) we develop a
formalism and methodology for the derivation of {\em multi-dimensional fits}
for event rates. We significantly generalize our earlier study, and we focus on
ways of thoroughly assessing the accuracy of the fits. We anticipate that the
efficient tools developed here can be applied in lieu of large-scale population
calculations and will facilitate the exploration of the dependence of rate
predictions on a wide range binary evolution parameters. Such explorations can
then allow the derivation of constraints on these parameters, given empirical
rate constraints and accounting for fitting errors. Here we describe in detail
the principles and practice behind constructing these fits, estimating their
accuracy, and comparing them with observations in a manner that accounts for
their errors.
