This has to be in the lit somewhere. Can someone point me to this in any accessible book or lit? It's just a big trig, sinusoidal, Fourier series thing:
$$\begin{align}
 y(t) &= \sum_{k=0}^{K} a_k \big( A \cos(\omega t) \big)^k \\
\\
     &= \sum_{n=0}^{K} b_n \cos(n \omega t) \\
\end{align}$$
The result I got (I think I fixed my error with the limits in swapping summations) is:
apparently, the DC component ($n=0$) of the Fourier cosine series is:
$$ b_0 =  a_0 + \sum_{k=1}^{\left\lfloor \tfrac{K}{2} \right\rfloor} a_{2k} \frac{(2k)!}{(k!)^2} \left(\frac{A}{2}\right)^{2k}  $$
The coefficient for the even harmonic ($n$ even) terms:
$$ b_{n} = 2 \sum_{k=\frac{n}{2}}^{\left\lfloor \tfrac{K}{2} \right\rfloor} a_{2k} \frac{(2k)!}{\left(k+\frac{n}{2}\right)! \left(k-\frac{n}{2}\right)!} \left(\frac{A}{2}\right)^{2k}  $$
And the coefficient for the odd harmonic ($n$ odd) terms:
$$ b_{n} = 2 \sum_{k=\frac{n-1}{2}}^{\left\lfloor \tfrac{K-1}{2} \right\rfloor} a_{2k+1} \frac{(2k+1)!}{(k+\frac{n+1}{2})!(k-\frac{n-1}{2})!} \left(\frac{A}{2}\right)^{2k+1}  $$
where $\lfloor u \rfloor$ is the floor() function, which is the largest integer no larger than the argument $u$.
$$ \lfloor u \rfloor \le u  < \lfloor u \rfloor + 1 \qquad \qquad \lfloor u \rfloor \in \mathbb{Z}, u \in \mathbb{R}$$
or
$$  u-1 < \lfloor u \rfloor \le u  \qquad \qquad \lfloor u \rfloor \in \mathbb{Z}, u \in \mathbb{R}$$
 A: Since the OP wants a reference, I would say the use the sums found in Prudnikov, Integrals and Series, vol. 1. 4.4.1.16 and 4.4.1.17. Ultimately they are not any more profound than the binomial theorem. In this particular situation they are, with a little algebra to clean them up
$$ \cos^{2n}(x)=2^{1-2n}\sum_{k=1}^{n}\binom{2n}{n-k}\cos{(2\,k\,x)} + 2^{-2n}\binom{2n}{n}$$
$$\cos^{2n+1}(x)=2^{-2n}\sum_{k=0}^{n}\binom{2n+1}{n-k}\cos{((2\,k+1)x)} $$
The OP's procedure basically amounts to inserting these expressions in the original sum and interchanging the summations.
A: Note that the Chebyshev polynomials $T_n$ are such that $T_n(\cos\theta) = \cos n\theta$; see the wikipedia article for details.  So you are given a polynomial $A(x)$ with known coefficients and want to find $b_k$ so that $A(x)=\sum_k b_k T_k(x)$.  A recipe is given here where it says "with inverse $x^n = 2^{1-n}\sum\dots$".  Multiply both sides of that by your $a_n$ and sum;
the coefficient of $T_k(x)$ will be your desired $k$.
Added: A reference for the needed identity is Corollary 12 in this paper; I'm sure there are others.
