Angular momentum and spin precession are expected to be generic features of a significant fraction of binary black hole systems. As such, it is essential to have waveform models that faithfully incorporate the effects of precession. Here, we assess how well the current state of the art models achieve this for waveform strains constructed only from the ell = 2 multipoles. Specifically, we conduct a survey on the faithfulness of the waveform models SEOBNRv5PHM, TEOBResumS, IMRPhenomTPHM, IMRPhenomXPHM to the numerical relativity (NR) surrogate NRSur7dq4 and to NR waveforms from the SXS catalog. The former assessment involves systems with mass ratios up to six and dimensionless spins up to 0.8. The latter employs 317 short and 23 long SXS waveforms. For all cases, we use reference inclinations of zero and 90 degrees. We find that all four models become more faithful as the mass ratio approaches unity and when the merger-ringdown portion of the waveforms are excluded. We also uncover a correlation between the co-precessing (2, +/-2) multipole mismatches and the overall strain mismatch. We additionally find that for high inclinations, precessing (2, +/-1) multipoles that are more faithful than their (2, +/-2) counterparts, and comparable in magnitude, improve waveform faithfulness. As a side note, we show that use of uniformly-filled parameter spaces may lead to an overestimation of precessing model faithfulness. We conclude our survey with a parameter estimation study in which we inject two precessing SXS waveforms (at low and high masses) and recover the signal with SEOBNRv5PHM, IMRPhenomTPHM and IMRPhenomXPHM. As a bonus, we present preliminary multidimensional fits to model unfaithfulness for Bayesian model selection in parameter estimation studies.