I originally asked this on the physics Stack Exchange site, but perhaps it could be more easily answered here.
Given the definition of the correlation function for CMB temperature fluctuations as
$$ C\left(\theta\right) = \left\langle \frac{\delta T}{T}\left(\hat{n}_1\right) \frac{\delta T}{T}\left(\hat{n}_2\right) \right\rangle_{\hat{n}_1\cdot \hat{n}_2 = \cos\theta} ,$$
I should be able to simplify it to
$$ C\left(\theta\right) = \frac{1}{4\pi} \sum_{l=0}^\infty (2l + 1) \, C_l \, P_l\left(\cos\theta\right) $$
(where $P_l \left(x\right)$ are the Legendre polynomials) by decomposing the temperature fluctuations into spherical harmonics like this
$$ \frac{\delta T}{T} = \sum_{l=0}^\infty \sum_{m=-l}^l a_{lm} Y_{lm}. $$
I think the first step of this procedure should look like this
$$ C\left(\theta\right) = \left\langle \sum_{l_1=0}^\infty \sum_{m_1=-l_1}^{l_1} a_{l_1 m_1} Y_{l_1 m_1}\left(\hat{n}_1\right) \sum_{l_2=0}^\infty \sum_{m_2=-l_2}^{l_2} a_{l_2 m_2} Y_{l_2 m_2}\left(\hat{n}_2\right) \right\rangle_{\hat{n}_1\cdot \hat{n}_2 = \cos\theta} .$$
I understand that the spherical harmonics can be written in the form
$$ Y_{lm}(\theta,\phi) \propto P_{lm} \left(\cos\theta\right) e^{i m \phi} $$
(where $P_{lm}(x)$ are the associated Legendre polynomials) and that $C_l$ should come out as
$$ C_l = \frac{1}{2l + 1} \sum_{m=-l}^l a_{lm} a_{l-m} $$
(though I could be off on this last piece). However, I am unsure of the mathematical steps involved in simplifying the four sums down to one. What identities, properties, or other insights will allow me to make this simplification?
Thanks!