Abstract: We present a comprehensive analysis of thermal QCD axion production across various mass thresholds, focusing on the KSVZ and DFSZ frameworks. We extend current calculations to account for heavy quark masses at higher temperatures and quantify axion production via hadron scattering at lower temperatures. A continuous axion production rate is derived across the QCD phase transition, which is used to estimate the axion contribution to the effective number of neutrino species, ΔNeff. By revisiting cosmological constraints, we provide updated bounds on axion mass. Our results offer robust limits on ΔNeff and axion mass, with significant implications for future CMB surveys. Additionally, we perform a rigorous momentum-space analysis to track dark radiation production via two-body decays and scatterings, highlighting the importance of non-instantaneous decoupling and quantum statistical effects in determining ΔNeff, with errors within the sensitivity of upcoming CMB experiments. Finally, thermal axion production in low-reheating scenarios will be briefly discussed.