光子（photon）屬於Bose-Einstein粒子，各狀態的光子個數不受包例不相容原理限制（關於受包立不相容原理限制的Fermi-Dirac粒子的統計物理性質，可參考 自由電子海的能量分布與Fermi能階）。在黑體中假如有M個狀態、N個相同無法辨明的光子，則系統總共有$W=H^M_N = C^{N+M-1}_{N}$種可能的組態，此時系統的亂度為
$$ S = k_B \ln(W) = k_B \ln \left( \frac{(N+M-1)!}{N!(M-1)!}\right) $$
利用Stirling近似，當$M$和$N$皆趨近於無限時，系統亂度有底下的漸近關係
$$ S \rightarrow k_B \left( (N+M) \ln(N+M) - N \ln(N) - M \ln(M) \right) = k_B M \left( (1+n) \ln(N+M) - n \ln(N) - \ln(M) \right) = k_B M \left( (1+n) \ln(1+n) - n \ln(n) \right) $$
其中$n=\frac{N}{M}$為每個狀態平均光子數。
定溫、定容，但粒子數可改變的系統中，我們定義總位能（grand potential）$\Phi_G$為
$$ \Phi_G = U - TS - \mu N $$
其中$U$為系統內能、$T$為系統溫度、$\mu$為系統化學能。根據熱力學第二定律，在熱力學平衡的情況下，光子在各個狀態的分布$N_i$（$i$表不同狀態）會使得總位能有最小值。定義$n_i = \frac{N_i}{M}$，則系統總位能可以寫成
$$ \Phi_G = M \sum_i (\varepsilon_i - \mu) n_i - k_B T \left( (1+n_i) \ln(1+n_i) - n_i \ln(n_i) \right) $$
其中$\varepsilon_i$為光子在狀態$i$時的內能。
當系統總位能有最小值時，其對各個$n_i$的偏微分應當為零：
$$ \left. \frac{\partial \Phi_G}{\partial n_i} \right|_{T,V, \mu} = M \left( \varepsilon_i - \mu - k_B T \left( \ln(1+n_i) - \ln(n_i) \right) \right) = 0 $$
$$ \ln \left( \frac{1 + n_i}{n_i} \right) = \frac{\varepsilon_i - \mu}{k_B T} $$
$$ \frac{1 + n_i}{n_i} = \exp \left(\frac{\varepsilon_i - \mu}{k_B T} \right) $$
$$ n_i = \frac{1}{\exp \left(\frac{\varepsilon_i - \mu}{k_B T} \right) - 1} $$
此即Bose-Einstein分布。而黑體中光子可以自由地被表面吸收獲釋放，因此該系統內並不需要考慮化學能，即$\mu_i=0$。代入單一光子所具有的內能$\varepsilon_i = h|\vec{\nu}_i| = c |\vec{p}_i|$（這裡$h$為普朗克常數、$\vec{\nu}_i$和$\vec{p}_i$各為狀態$i$之光波頻率與光子動量），則光子的在不同狀態的分布為