I'm attempting to perform the integration that will yield the sound horizon at recombination: $$ c_s^2 = \frac{c^2}{3}\left[\frac{3}{4}\frac{\rho_{b,0}(1+z)^3}{\rho_{\gamma,0}(1+z)^4} + 1\right]^{-1} $$The present day density of baryons is reasonably easy to calculate as $\rho_b=\Omega_b\times \rho_{crit}$. Using the Planck Collaboration results, I get a present day value of $\rho_b=4.14\times 10^{-25}\space g\space m^{-3}$.
I've seen one post calculate the photon density as $$ \rho_\gamma = \frac{a_B\, T_0^4}{c^2} = 4.64511\times 10^{-31}\;\text{kg}\,\text{m}^{-3}. $$ where $$ a_B = \frac{8\pi^5 k_B^4}{15h^3c^3} = 7.56577\times 10^{-16}\;\text{J}\,\text{m}^{-3}\,\text{K}^{-4} $$ But then the author goes on to calculate the radiation density and suggests plugging the total radiation density value into the formula for sound velocity which doesn't seem right (as the electrons are bound to photons, not all radiation).
What is the right value for $\rho_{\gamma,0}$ for the purpose of calculating the velocity of sound in the pre-recombination fluid?