The accuracy of gravitational-wave models of compact binaries has traditionally been addressed by the mismatch between the model and numerical-relativity simulations. This is a measure of the overall agreement between the two waveforms. However, the largest modelling error typically appears in the strongfield merger regime and may affect subdominant signal harmonics more strongly. These inaccuracies are often not well characterised by the mismatch. We explore the use of a complementary, physically motivated tool to investigate the accuracy of gravitational-wave harmonics in waveform models: the remnant’s recoil, or kick velocity. Asymmetric binary mergers produce remnants with significant recoil, encoded by subtle imprints in the gravitational-wave signal. The kick estimate is highly sensitive to the intrinsic inaccuracies of the modelled gravitational-wave harmonics during the strongly relativistic merger regime. Here we investigate the accuracy of the higher harmonics in four state-of-the-art waveform models of binary black holes. We find that the SEOBNRv4HM_ROM, IMRPhenomHM, IMRPhenomXHM and NRHybSur3dq8 models are not consistent in their kick predictions. Our results enable us to identify regions in the parameter space where the models require further improvement and support the use of the kick estimate to investigate waveform systematics. We discuss how numerical-relativity kick estimates could be used to calibrate waveform models further, proposing the first steps towards kick-based gravitational-wave tuning.