  We study thermal axion production around the confinement scale. At higher
temperatures, we extend current calculations to account for the masses of heavy
quarks, whereas we quantify production via hadron scattering at lower
temperatures. Matching our results between the two opposite regimes provides us
with a continuous axion production rate across the QCD phase transition. We
employ such a rate to quantify the axion contribution to the effective number
of neutrino species.
