  We derive the light-front wave function (LFWF) representation of the
$\gamma^* \gamma^* \to \eta_c(1S)\,,\eta_c(2S)$ transition form factor
$F(Q_1^2, Q_2^2)$ for two virtual photons in the initial state. For the LFWF,
we use different models obtained from the solution of the Schr\"odinger
equation for a variety of $c \bar c$ potentials. We compare our results to the
BaBar experimental data for the $\eta_c(1S)$ transition form factor, for one
real and one virtual photon. We observe that the onset of the asymptotic
behaviour is strongly delayed and discuss applicability of the collinear and/or
massless limit. We present some examples of two-dimensional distributions for
$F (Q_1^2,Q_2^2)$. A factorization breaking measure is proposed and
factorization breaking effects are quantified and shown to be almost model
independent. Factorization is shown to be strongly broken, and a scaling of the
form factor as a function of $\bar Q^2 = (Q_1^2 + Q_2^2)/2$ is obtained.
