We introduce a method to quantify the quality-of-fit between data and observables depending on the large scale Galactic magnetic field. We combine WMAP5 polarized synchrotron data and rotation measures of extragalactic sources in a joint analysis to obtain best fit parameters and confidence levels for GMF models common in the literature. None of the existing models provide a good fit in both the disk and halo regions, and in many instances best-fit parameters are quite different than the original values. We note that probing a very large parameter space is necessary to avoid false likelihood maxima. The thermal and relativistic electron densities are critical for determining the GMF from the observables but they are not well constrained. We show that some characteristics of the electron densities can already be constrained using our method and with future data it may be possible to carry out a self-consistent analysis in which models of the GMF and electron densities are simultaneously optimized.