Rational Krylov model-order reduction techniques are used to accelerate the characterization of frequency selective surfaces (FSSs) over broad spectral ranges, including frequencies at which a higher order Floquet mode stops evanescing and begins to propagate. The procedure comprises three main stages: (1) construction of the spectral Galerkin system at a small set of frequencies, (2) linearization of that system, and (3) reduction of the linearized system using the rational Krylov technique. The inclusion of blazing frequencies in the band of interest complicates the second and third of these steps because of the branch point singularity in the periodic Green's function. This difficulty is avoided by removing the blazing modes contributions to the spectral Galerkin matrix using the Woodbury formula for low-rank updates. The algorithm results in a set of small linear systems producing outputs that are combined to approximate the reflection and transmission coefficients of all propagating modes. The technique is applied to three different frequency selective surfaces and is shown to be accurate and efficient in all cases.