In this work, we address an everlasting issue in model order reduction (MOR) in electromagnetics that has remained unnoticed until now. Contrary to what has been previously done, we identify for the very first time spurious modes in MOR for time-harmonic Maxwell’s equations and propose a methodology to remove their negative influence on the reduced order model (ROM) response. These spurious modes have nonzero resonance frequency and may have shown up in the past giving rise to spikes in the frequency response, in effect, deteriorating the accuracy and efficiency of the MOR process. However, they were never characterized as spurious mode contributions, rather they were most likely considered as poor localized approximation issues in the MOR process. When we try to get further physical insights from the ROM, rather than simple frequency domain data, we cannot afford any poorly localized approximation issue, that is, any spurious mode in the band of analysis. Otherwise, these mathematical, but nonphysical, modes will mislead the physical behavior of the device under analysis. A computationally inexpensive variational divergence condition is established to identify spurious modes in the band of analysis, since any physical in-band mode must be divergence-free. In addition, once a spurious mode is identified in the band of analysis, its influence is removed from the ROM by a physics-based coupling strategy. As a result, a robust spurious mode contribution-free MOR in electromagnetics is proposed. Finally, several actual microwave circuits, such as a quad-mode filter and a triple-mode triple-band filter, will illustrate the capabilities and efficiency of the proposed approach.