State-of-the-art numerical simulations of quantum electrodynamical (QED) processes in strong laser fields rely on a semiclassical combination of classical equations of motion and QED rates, which are calculated in the locally constant field approximation. However, the latter approximation is unreliable if the amplitude of the fields, 𝑎 0 , is comparable to unity. Furthermore, it cannot, by definition, capture interference effects that give rise to harmonic structure. Here we present an alternative numerical approach, which resolves these two issues by combining cycle-averaged equations of motion and QED rates calculated in the locally monochromatic approximation. We demonstrate that it significantly improves the accuracy of simulations of photon emission across the full range of photon energies and laser intensities, in plane-wave, chirped and focused background fields.