  We provide a complete theoretical description of the two-neutrino electron
capture in $^{124}$Xe, improving both the nuclear and the atomic structure
calculations. We improve the general formalism through the use of the Taylor
expansion method, leading to higher order terms in the decay rate of the
process. The nuclear part is treated with pn-QRPA and interacting shell model
(ISM) methods. The nuclear matrix elements (NMEs) are calculated with the
pn-QRPA method with spin restoration by fixing the input parameters so that the
experimental decay rate is reproduced, resulting in values significantly lower
than in previous calculations. The validity of the pn-QRPA NMEs is tested by
showing their values to be comparable with the ones for double-beta decay with
emission of two electrons of $^{128,130}$Te, which have similar pairing
features. Within the ISM, we reproduce the total experimental half-life within
a factor of two and predict the capture fraction to the KK channel of about
74\%. We also predict the capture fractions to other decay channels and show
that for the cumulative decay to the $\rm{KL_{1}}$-$\rm{KO_{1}}$ channels, a
capture fraction of about 24\% could be observed experimentally. On the atomic
side, calculations are improved by accounting for the Pauli blocking of the
decay of innermost nucleon states and by considering all $s$-wave electrons
available for capture, expanding beyond the K and L$_1$ orbitals considered in
previous studies. We also provide improved atomic relaxation energies of the
final atomic states of $^{124}$Te, which may be used as input for background
modeling in liquid Xenon experiments.
