The procedure of quantization of a classical theory has its crucial passage in the promotion of the dynamical variables to operators in the Hilbert space of the quantum states of the system. The quantization is successful if these operators form a complete set of commuting observables (this set is called ICOC). In this case, every quantum state of the system can be expressed as a sum of eigenstates, and we can write the probability amplitudes, and take the averages of all the dynamical variables. The quantized description of the electromagnetic field in the vacuum is obtained by expressing the field hamiltonian in a finite volume $V$, and then taking the limit for $V \to \infty$. The promotion to operators is performed by the position and conjugate momentum operators, which are introduced by writing the field hamiltonian in terms of a sum, possibly infinite, of uncoupled harmonic oscillators. We start by expressing the hamiltonian of the field in a volume of space $V$. $$ H = \int_{V}dV \dfrac{1}{2}\varepsilon_0 |\mathbf{E}|^2 + \dfrac{1}{2\mu_0}|\mathbf{B}|^2. $$ Now we remember that it is possible to express the fields $\mathbf{E}$, $\mathbf{B}$ using the magnetic vector potential $\mathbf{A}$ in this way, choosing the Coulomb gauge: $$ \begin{align} \mathbf{E} &= - \dfrac{\partial}{\partial t} \mathbf{A} \\ \mathbf{B} &= \nabla \times \mathbf{A}. \end{align} $$ We recall that, from the Maxwell equations, the wave equation for the field $\mathbf{A}$ $$ \left(\nabla^2 - \dfrac{1}{c^2} \dfrac{\partial^2}{\partial t^2} \right)\mathbf{A} =0. $$ Since we are in a limited volume of space, we can impose a periodic boundary condition, and express the field $\mathbf{A}$ according to its Fourier series, by summing complex exponentials. We choose to express the component plane waves as linearly polarized waves. This does not affect the generality of the procedure, as other polarizations can be expressed as sums of linearly polarized waves, with appropriate phase shifts. Let therefore $\mathbf{k}$ be the wave vector, and $s$ the linear polarization ($s \in {\alpha, \beta}$ where $\alpha, \beta$ are two possible orthogonal polarizations given the wave vector $\mathbf{k}$). We have that
$$ \mathbf{A} = \sum_{\mathbf{k}s} A_{\mathbf{k}s}(t) e^{i \mathbf{k} \cdot \mathbf{r}} , \hat{\epsilon}{\mathbf{k}s}, $$ where $\hat{\epsilon}{\mathbf{k}s}$ is the unit vector of the field, parallel to the polarization. Furthermore the field $\mathbf{A}$ must be real, so we have the condition $$ A_{-\mathbf{k}s}(t) = A_{\mathbf{k}s}^(t). $$ By imposing this condition, unless we redefine the terms $A_{\mathbf{k}s}$ by halving them, we can rewrite the sum as $$ \mathbf{A} = \sum _{\mathbf{k}s} A _{\mathbf{k}s}(t) e^{i \mathbf{k} \cdot \mathbf{r}} + A^ {\mathbf{k}s}(t) e^{-i \mathbf{k} \cdot \mathbf{r}} \hat{\epsilon} {\mathbf{k}s}. $$ Substituting the following condition $$ \ddot{A} {\mathbf{k}s}(t) + \omega \mathbf{k}^2 A {\mathbf{k}s}(t) = 0, $$ with the dispersion relation $$ \omega\mathbf{k} = c |\mathbf{k}|. $$ Since $\omega$, in the vacuum, depends only on the modulus $|\mathbf{k}| =: k$ we can also write $\omega_k$ for simplicity. So the equation is solved as $$ A{\mathbf{k}s}(t) = A{\mathbf{k}s}(0) e^{-\omega_k t}. $$ We now computer the Hamiltonian in terms of this Fourier expansion of the fields $$ H = \int{V}dV \dfrac{1}{2}\varepsilon_0 |\dfrac{\partial}{\partial t} \mathbf{A}|^2 + \dfrac{1}{2\mu_0}|\nabla \times \mathbf{A}|^2. $$ In this way we obtain $$ H = \sum{\mathbf{k}s} \varepsilon_0 \omega_k^2 ( A_{\mathbf{k}s}(t) A_{\mathbf{k}s}^* (t) + A_{\mathbf{k}s}^{} (t) A_{\mathbf{k}s}(t)) , $$ we notice that this Hamiltonian is independent of time, since the terms $A_{\mathbf{k}s}(t)A_{\mathbf{k}s}^ (t)$ and the complex conjugate are. For brevity we can simply write $A_{\mathbf{k}s} A_{\mathbf{k}s}^* $. We introduce the adimensional variables $a_{\mathbf{k}s}$, and $a_{\mathbf{k}s}^$ in such a way that $$ \begin{align} A_{\mathbf{k}s} &= \sqrt{\dfrac{\hslash}{2\epsilon_0 \omega_k}} a_{\mathbf{k}s} \\ A_{\mathbf{k}s}^ &= \sqrt{\dfrac{\hslash}{2\epsilon_0 \omega_k}} a_{\mathbf{k}s}^*. \end{align} $$
...