In the context of ultra-relativistic nuclear collisions, we present a fast method for calculating the final particle spectra after the direct decay of resonances from a Cooper–Frye integral over the freeze-out surface. The method is based on identifying components of the final particle spectrum that transform in an irreducible way under rotations in the fluid-restframe. Corresponding distribution functions can be pre-computed including all resonance decays. Just a few of easily tabulated scalar functions then determine the Lorentz invariant decay spectrum from each space-time point, and simple integrals of these scalar functions over the freeze-out surface determine the final decay products. This by-passes numerically costly event-by-event calculations of the intermediate resonances. The method is of considerable practical use for making realistic data to model comparisons of the identified particle yields and flow harmonics, and for studying the viscous corrections to the freeze-out distribution function.