Given a basis $B = \{g_0, g_1,\ldots , g_{nk}\} $, I wish to construct a set of $nk$ functions $S$ such that $\langle S_i,S_j \rangle = \delta_{i,j}$ (so that $S$ is an orthonormal set) where $\langle f,g\rangle = \int_{0}^{\infty}f g$.
Assuming the following definitions: $S_k = c_{k,0}\cdot g_0 + c_{k,1}\cdot g_1 +\cdots + c_{k,nk}\cdot g_{nk}$
$C$ is a matrix where $C_{i,j} = c_{i,j}$
$G$ is a matrix where $G_{i,j} = \langle g_i,g_j \rangle $
Then it is easy to see that $CGC^T$ is a matrix with elements equal to $\langle S_i,S_j \rangle $. Thus, if $CGC^T = I$, then $C$ is the matrix containing the correct coefficients to produce the set of functions $S$. Now, upon diagonalizing the matrix $G$, we can generate a matrix $P$ where $PGP^T$ is diagonal with eigenvalues on the main diagonal.
What changes do I need to make to $P$ in order to produce $C$? Is it as simple as making small changes to $P$ or will I need to pursue another method altogether?