To understand quantum optics experiments, we must perform calculations that consider the principal sources of noise, such as losses, spectral impurity and partial distinguishability. In both discrete and continuous variable systems, these can be modeled as mixed Gaussian states over multiple modes. The modes are not all resolved by photon-number measurements and so require calculations on coarse-grained photon-number distribution. Existing methods can lead to a combinatorial explosion in the time complexity, making this task unfeasible for even moderate sized experiments of interest. In this work, we prove that the computation of this type of distributions can be done in exponential time, providing a combinatorial speedup. We develop numerical techniques that allow us to determine coarse-grained photon number distributions of Gaussian states, as well as density matrix elements of heralded non-Gaussian states prepared in the presence of spectral impurity and partial distinguishability.