We begin with the equation for classical rotational kinetic energy obtained in section Types of motion in molecules and rewrite it as follows

$$ T_\text{rot} = \frac{1}{2}\sum_i m_i (\boldsymbol{\omega}\times\mathbf{r}_i)^2 = \frac{1}{2}\sum_i m_i \left(\begin{array}{c} \omega_yz_i - \omega_zy_i \\ \omega_zx_i - \omega_xz_i \\ \omega_xy_i - \omega_yx_i \\ \end{array}\right)^2 = \\ = \frac{1}{2}\sum_i m_i\left[ (\omega_yz_i - \omega_zy_i)^2 + (\omega_zx_i - \omega_xz_i)^2 + (\omega_xy_i - \omega_yx_i)^2 \right] = \\ = \sum_i m_i \left[ \omega_x^2(y_i^2+z_i^2) + \omega_y^2(x_i^2+z_i^2) + \omega_z^2(y_i^2+x_i^2) - \right. \\ -\left. 2\omega_x\omega_yx_iy_i - 2\omega_x\omega_zx_iz_i - 2\omega_y\omega_zy_iz_i \right] = \\ = \boldsymbol{\omega}^T\cdot \mathbf{I}\cdot \boldsymbol{\omega}. $$

Here, $\mathbf{I}$ is the moments of inertia tensor with elements

$$ \mathbf{I} = \left(\begin{array}{ccc} \sum_i m_i(y_i^2+z_i^2) & -\sum_i m_ix_iy_i & -\sum_i m_ix_iz_i \\ -\sum_i m_ix_iy_i & \sum_i m_i(x_i^2+z_i^2) & -\sum_i m_iy_iz_i \\ -\sum_i m_ix_iz_i & -\sum_i m_iy_iz_i & \sum_i m_i(x_i^2+y_i^2) \end{array}\right). $$

The classical angular momentum operator $\mathbf{J} =\mathbf{I}\cdot\boldsymbol{\omega}$, from which the angular frequency vector can be written as $\boldsymbol{\omega}=\mathbf{I}^{-1}\cdot\mathbf{J}$, and the expression for rotational kinetic energy becomes $T_\text{rot} = \frac{1}{2}\mathbf{J}^T\cdot\mathbf{I}^{-1}\cdot\mathbf{J}$. We get the quantum mechanical formulation by substituting the classical angular momentum operator with its quantum mechanical equvalent,

$$ \hat{T}\text{rot} = \frac{1}{2}\sum{\alpha,\beta=x,y,z}\hat{J}\alpha I{\alpha\beta}\hat{J}_\beta. $$

As mentioned in the section Types of motion in molecules, one of the reasons for choosing the molecular frame embedding is to simplify the rotational kinetic energy operator. For example, we can compute Cartesian coordinates of atoms $r_{i\alpha}$ in some arbitrary center-of-mass frame and then 'rotate' the coordinate system such that the moment of inertia tensor $\mathbf{I}$ becomes diagonal. In this case $\hat{T}\text{rot}=\frac{1}{2}\left(I{aa}^{-1}\hat{J}a^2 +I{bb}^{-1}\hat{J}b^2 +I{cc}^{-1}\hat{J}_c^2 \right)$. Such coordinate system is called the principal axes system (PAS). The rotation of Cartesian coordinates of atoms from the reference (arbitrary) to the PAS frame can be expressed as follows:

$$ r_{i\alpha} = \sum_{\beta=x,y,z}d_{\alpha\beta}\bar{r}_{i\beta}, $$

where $\bar{r}{i\alpha}$ and $r{i\alpha}$ are coordinates of atoms in the reference and PAS frames, respectively, and $\mathbf{d}$ is a $3\times 3$ frame rotation matrix ($\mathbf{d}\mathbf{d}^T=\mathbf{I}$). The molecular frame embedding conditions can now be expressed in terms of three equations for the elements of rotation matrix $\mathbf{d}$. For example, the PAS frame constraints on the diagonal form of the moment of inertia tensor provide the following equation for the rotation matrix:

$$ [\mathbf{d}\mathbf{u}\mathbf{d}^T]{\alpha\beta} = 0~~~\text{for}~\alpha\neq \beta, \\ u{\alpha\beta} = \sum_i m_i \bar{r}{i\alpha}\bar{r}{i\beta}. $$

The Eckart condition results in the following equation:

$$ \mathbf{u}\mathbf{d}^T-\mathbf{d}\mathbf{u}^T=0, \\ u_{\alpha\beta} = \sum_i m_i r_{i\alpha}^\text{(eq)}\bar{r}_{i\beta}. $$

In the principal axes system the inertia tensor is diagonal, resulting in simplified expression for the kinetic energy operator:

$$ \hat{T}\text{rot}=\frac{1}{2}\left(I{aa}^{-1}\hat{J}a^2 +I{bb}^{-1}\hat{J}b^2 +I{cc}^{-1}\hat{J}_c^2 \right). $$

By introducing the rotational constants (in units of Hz),

$$ A=\frac{\hbar}{4\pi I_{aa}},~~~B=\frac{\hbar}{4\pi I_{bb}},~~~C=\frac{\hbar}{4\pi I_{cc}}, $$

we can rewrite the kinetic energy operator as

$$ \hat{T}_\text{rot}=A\hat{J}_a^2 +B\hat{J}_b^2 +C\hat{J}_c^2. $$