  Extracting information from cosmic surveys is often done in a two-step
process, construction of maps and then summary statistics such as two-point
functions. We use simulations to demonstrate the advantages of a general
Bayesian framework that consistently combines different cosmological
experiments on the field level, and reconstructs both the maps and cosmological
parameters. We apply our method to jointly reconstruct the primordial CMB, the
integrated Sachs-Wolfe effect, and six tomographic galaxy density maps on the
full sky on large scales along with several cosmological parameters. While the
traditional maximum a posterior estimator has both two-point level and
field-level bias, the new approach yields unbiased cosmological constraints and
improves the signal-to-noise ratio of the maps.
