We compare GW150914 directly to simulations of coalescing binary black holes in full general relativity,including several performed specifically to reproduce this event. Our calculations go beyond existingsemianalytic models, because for all simulations—including sources with two independent, precessing spins—we perform comparisons which account for all the spin-weighted quadrupolar modes, and separatelywhich account for all the quadrupolar and octopolar modes. Consistent with the posterior distributionsreported by Abbott et al. [Phys. Rev. Lett. 116, 241102 (2016)] (at the 90% credible level), we find the dataare compatible with a wide range of nonprecessing and precessing simulations. Follow-up simulationsperformed using previously estimated binary parameters most resemble the data, even when all quadrupolarand octopolar modes are included. Comparisons including only the quadrupolar modes constrain the totalredshifted mass Mf,z in the range 64.0 M⊙ − 73.5 M⊙ and the final black hole’sdimensionless spin parameter is consistent with af = 0.62–0.73. As our approach invokes no intermediateapproximations to general relativity and can strongly reject binaries whose radiation is inconsistent with thedata, our analysis provides a valuable complement to Abbott et al. [Phys. Rev. Lett. 116, 241102 (2016)].