  Context. In previous work, we developed a quasi-Gaussian approximation for
the likelihood of correlation functions, which, in contrast to the usual
Gaussian approach, incorporates fundamental mathematical constraints on
correlation functions. The analytical computation of these constraints is only
feasible in the case of correlation functions of one-dimensional random fields.
  Aims. In this work, we aim to obtain corresponding constraints in the case of
higher-dimensional random fields and test them in a more realistic context.
  Methods. We develop numerical methods to compute the constraints on
correlation functions which are also applicable for two- and three-dimensional
fields. In order to test the accuracy of the numerically obtained constraints,
we compare them to the analytical results for the one-dimensional case.
Finally, we compute correlation functions from the halo catalog of the
Millennium Simulation, check whether they obey the constraints, and examine the
performance of the transformation used in the construction of the
quasi-Gaussian likelihood.
  Results. We find that our numerical methods of computing the constraints are
robust and that the correlation functions measured from the Millennium
Simulation obey them. Despite the fact that the measured correlation functions
lie well inside the allowed region of parameter space, i.e. far away from the
boundaries of the allowed volume defined by the constraints, we find strong
indications that the quasi-Gaussian likelihood yields a substantially more
accurate description than the Gaussian one.
