## Abstract

Probing optical excitations with nanometer resolution is important for understanding their dynamics and interactions down to the atomic scale. Electron microscopes currently offer the unparalleled ability of rendering spatially resolved electron spectra with combined meV and sub-nm resolution, while the use of ultrafast optical pulses enables fs temporal resolution and exposure of the electrons to ultraintense confined optical fields. Here, we theoretically investigate fundamental aspects of the interaction of fast electrons with localized optical modes that are made possible by these advances. We use a quantum optics description of the optical field to predict that the resulting electron spectra strongly depend on the statistics of the sample excitations (bosonic or fermionic) and their population (Fock, coherent, or thermal), whose autocorrelation functions are directly retrieved from the ratios of electron gain intensities. We further explore feasible experimental scenarios to probe the quantum characteristics of the sampled excitations and their populations. In particular, we present realistic simulations for electron beams interacting with optical cavities infiltrated with optically pumped quantum emitters, which we show to undergo a varied temporal evolution in the cavity mode statistics that causes radical modifications in the transmitted electron spectra depending on pump-electron delay.

© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## Corrections

17 December 2019: A typographical correction was made to the final paragraph of Section 2.

27 January 2021: A correction was made to the funding section.

## 1. INTRODUCTION

Electron energy-loss spectroscopy (EELS) performed in electron microscopes is a fertile source of information on the dielectric properties of materials down to the nanometer scale [1–3]. This technique is widely used to identify chemical species with atomic resolution through their characteristic high-energy core losses [4,5]. Additionally, low-loss EELS provides insight into the spatial and spectral distributions of plasmons in metallic nanostructures [6–10] and, more recently, also of phonons in polaritonic materials [3,11] thanks to remarkable advances in instrument resolution. In a parallel effort, the ultrafast dynamics of nanostructured materials and their influence on optical near-fields can be studied by synchronizing the arrivals of fs light and electron pulses at the sample [12–14]. Indeed, although photons and electrons interact extremely weakly in free space due to the lack of energy–momentum matching, the evanescent field components produced upon light scattering by material structures breaks the mismatch, giving rise to efficient light–electron interaction, and effectively producing exchanges of multiple quanta between the electron and the optical field, accompanied by complex sub-fs dynamics [12,15–27]. Based on this principle, photon-induced near-field electron microscopy (PINEM) is performed by analyzing the resulting multiple gain and loss features in the electron spectra.

PINEM experiments have so far relied on coherent light sources (i.e., lasers), for which the measured spectra are well reproduced by assuming sample bosonic excitations that are coherently populated with a large number of quanta. The probability of each electron spectral peak associated with a net exchange of $ \ell $ quanta is then simply given by the squared Bessel function $ J_\ell ^2(2|\beta |) $, where a single parameter $ \beta = (e/\hbar \omega )\int\! \text{d}z {{\cal E}_z}(z){{e}^{ - {i}\omega z/v}} $ captures the strength of the electron–light interaction, mediated by the optical electric field component $ {{\cal E}_z}(z) $ along the direction of electron propagation $ z $ for an electron velocity $ v $ and photon frequency $ \omega $ [25]. For nanometer-sized samples (e.g., $ \Delta z \sim 100 \,\, \text{nm} $) illuminated at optical frequencies ($ \hbar \omega \sim 1 \,\, \text{eV} $), a field amplitude $ {\cal E} \sim {10^7} \,\, \text{V}/\text{m} $ renders $ |\beta | \sim e\Delta z{\cal E}/ \hbar \omega \sim 1 $. Eventually, even the zero loss peak (corresponding to $ \ell = 0 $) is fully depleted for $ |\beta | \approx 1.2 $. The underlying physics is thus described in terms of a classical optical field interacting with the electron through sample mediated harmonic evanescent fields. However, we expect new physics to arise when departing from this regime by considering anharmonic states of the illuminated sample, such as those associated with fermionic excitations [28–30], or when the external light source is not in a coherent state, such as that provided by a laser. As an interesting avenue in this direction, electron–photon entanglement has been recently predicted to influence the interaction with an electron when the sample is previously excited by a trailing electron [31]. In a related context, quantum aspects of fermionic two-level excitations have been probed in the cathodoluminescence signal emitted from single atomic defects [28–30]. In combination with external illumination, efficient excitation of bosonic and fermionic systems could be achieved in order to investigate the difference in quantum behavior between them and how this affects their interaction with electron beams.

In this work, we discuss the interaction of electron beams with individual optical modes and predict nontrivial characteristics of this interaction when the modes are excited through external illumination depending on the mode nature and population statistics. Specifically, the electron spectra resulting from the interaction with bosonic and fermionic excitations exhibit a radically different dependence on the external light intensity. Additionally, the electron spectra for bosonic modes depend dramatically on the photon statistics, giving rise to a varied phenomenology of asymmetric gain and loss peaks at low fluences and distinct intensity dependences under strong pumping. Interestingly, the autocorrelation functions can be directly retrieved from ratios of measured electron gain intensities. We further propose a feasible experimental realization of these ideas based on a sample consisting of an optical cavity that is fed by optically pumped three-level quantum emitters (QEs).

## 2. INTERACTION BETWEEN A BEAM ELECTRON AND A BOSONIC EXCITATION

We consider a sample characterized by a single boson mode of frequency $ {\omega _0} $ interacting with a focused beam electron of relativistic momentum and total energy tightly peaked around $ \hbar {{\bf k}_0} $ and $ {E_0} $, respectively. Assuming the sample has an extension along the beam direction sufficiently small as to preserve the transversal beam profile in the sample region, we can write the incident electron wave function as $ {\psi _0}({\bf r},t) = {{ e}^{{ i}({{\bf k}_0} \cdot {\bf r} - {E_0}t/\hbar )}}{\phi _0}({\bf r} - {\bf v}t) $, where $ {\phi _0} $ is a slowly varying function of the moving-frame position $ {\bf r} - {\bf v}t $. Further adopting the nonrecoil approximation ($ \hbar {\omega _0} \ll {E_0} $) and neglecting inelastic boson losses ($ {\gt} \text{ps} $ lifetimes) during the interaction time (in the fs range), we can write the system Hamiltonian $ \hat {\cal H} = {\hat {\cal H}_0} + {\hat {\cal H}_1} $ as the sum of unperturbed and interaction parts (see Appendix A for a self-contained derivation from the Dirac equation)

Finally, we note that if one is interested in the amplitudes of Eq. (1), the correct form of the coefficients $f_\ell^n$ for an electron interacting with an initial separable state of light $|\psi\rangle=\sum_i\alpha_i|i\rangle$ is readily obtained by substituting $\sqrt{p_{n+\ell}}$ by $\alpha_{n+\ell}$ in the previous definitions. In contrast, if the initial light state is a pure mixture described by the density matrix $\rho=\sum_i|\alpha_i|^2|i\rangle\langle i|$ (e.g., a thermal state), the amplitudes are obtained by incoherently averaging each coefficient $f_\ell^n$ (calculated with the condition $p_{n+\ell}=\delta_{n+\ell,i}$) over the distribution $|\alpha_i|^2$. Incidentally, in calculating electron spectra, identical results are obtained if we impose $p_{n+\ell}=|\alpha_{n+\ell}|^2$ in $f_\ell^n$, and therefore, for the sake of simplicity, this correspondence is implicitly employed throughout this work.

## 3. WEAK-COUPLING LIMIT

The first-order perturbation solution of Eq. (2) readily leads to loss and gain probabilities (see Appendix C) $ {P_{ - 1}} = ({\bar n} + 1)|{\beta _0}{|^2} $ and $ {P_1} = {\bar n} {|{\beta _0}|^2} $, respectively, where $ {\bar n} = \sum\nolimits_{n = 1}^\infty n{p_n} $ is the average mode population and $ {\beta _0} $ is given by Eq. (5). We then find that both loss and gain peaks increase in strength with ${|{\beta _0}|^2} $, but their difference is independent of $ \bar n $. In fact, the electron–boson interaction strength is determined by $ {\bar n} |{\beta _0}{|^2} $, which allows us to separate the regimes of weak and strong coupling depending on $ {|{\beta _0}|^2} $ and $ \bar n $, as shown in Fig. 1. We remark that this separation between weak- and strong-coupling regimes is consistent with the commonly used criterion based on the existence of a perturbative solution to the dynamics of the system, which in the present case depends on the combined value of the mode population $ \bar n $ and coupling parameter $ |{\beta _0}| $. Additionally, we observe that the ratio of gains to losses approaches 1 in the $ \bar n \gg 1 $ limit, which is consistent with the behavior of the weak-coupling ratio $ \bar n/(\bar n + 1) $.

Continuing the perturbation series, we find that the leading order contribution to the $ \ell \gt 0 $ gain peak is directly proportional to the $ \ell $th-order autocorrelation function at zero delay $ {{ g}^{(\ell )}} = \sum\nolimits_{n = \ell }^\infty n(n - 1) \ldots (n - \ell + 1){p_n}/{{\bar n}^2} $ [36]; we find the powerful result $ {P_\ell }/P_1^\ell = {{ g}^{(\ell )}}/(\ell {!)^2} $ (see Appendix 10), which holds for $ \sqrt{\bar{n}}|\beta_0| \lesssim 0.1 $ (low-excitation limit; see Fig. S1 in Supplement 1) and allows us to extract the autocorrelation functions from ratios of measured gain intensities.

#### A. Fermionic Excitations

For a fermionic sample mode (e.g., a two-level system, for which the operators $ a $ and $ {a^\dagger } $ satisfy anticommutation relations instead of boson commutation relations) under external cw illumination, we have $ f_\ell ^{n \gt 1} = 0 $, $ {p_1} = \bar n $, and $ {p_0} = 1 - \bar n $. Now, to first-order perturbation in the electron–fermion interaction, these populations readily lead to electron energy loss and gain probabilities $ {P_{ - 1}} = (1 - \bar n)|{\beta _0}{|^2} $ and $ {P_1} = {\bar n} |{\beta _0}{|^2} $, respectively. Simulating the fermionic mode as a two-level system with losses (see details in Appendix E), we find $ \bar n = 1/2 $ in the high light-intensity limit. More details on the intensity dependence of these peaks and their comparison with those of a bosonic system are provided in Fig. S2 (see Supplement 1). Nevertheless, free-electron interactions with fermionic systems are generally weak. Additionally, because only two states are involved, the dependence on incident photon statistics is trivial. In contrast, bosonic systems can realistically reach higher coupling to electrons and further exhibit a nontrivial dependence on photon statistics (see below).

## 4. INTERACTION WITH A DIPOLAR EXCITATION

For completeness, we discuss a mode described by a transition dipole $ {\bf p} $ at the origin so that the single-mode electric field reduces to $ {\vec {\cal E}_0} = ( {{\bf p}k_0^2 + {\bf p} \cdot \nabla \nabla } ){{ e}^{{ i}{k_0}r}}/r $, where $ {k_0} = {\omega _0}/c $, $ r = \sqrt {{b^2} + {z^2}} $, and $ b $ is the electron impact parameter. This is an excellent approximation for plasmonic particles and Mie resonators in a spectral region dominated by a dipolar excitation [6]. Inserting this field into Eq. (5), we find the coupling parameter $ {\beta _0} = ( - 2e{\omega _0}/\hbar {v^2}\gamma ) [ {{ i}{p_x}{K_1}(\zeta ) + ({p_z}/\gamma ){K_0}(\zeta )} ] $, where the modified Bessel functions $ {K_m} $ are evaluated at $ \zeta = {\omega _0}b/v\gamma $, we introduce the relativistic Lorentz factor $ \gamma = 1/\sqrt {1 - {v^2}/{c^2}} $, and the electron beam is taken to move along $ z $ and cross the $ x $ axis. Reassuringly, the resulting lowest-order loss probability for the ground state sample $ {|{\beta _0}|^2} $ coincides with the integrated EELS probability for a dipolar scatterer (see Appendix F), thus corroborating that the expression assumed for $ {\vec {\cal E}_0} $ has the correct single-dipole-mode normalization. The values of $ {\beta _0} $ considered in this study are in the range estimated via this model for the dipolar response of metallic (plasmonic) and dielectric (Mie resonances) optical cavities interacting with few to 100 s keV electrons (see Appendices G and H, respectively).

## 5. DEPENDENCE ON BOSON POPULATION STATISTICS

The boson mode can be populated in different ways, and in particular when pumped by external illumination, the photon statistics are directly imprinted on it. Here, we study three representative distributions with average occupation $ \bar n $: (1) a Fock state described by $ {p_n} = {\delta _{n,{\bar n}}} $, as one would obtain by coupling to quantum emitters (see below); (2) a coherent state characterized by a Poissonian distribution $ {p_n} = {{ e}^{ - \bar n}}{\bar n^n}/n! $ [37], as provided by external laser illumination; and (3) a thermal distribution, as produced by illuminating with chaotic light or by heating the sample mode at a temperature $ T $ commensurate with $ \hbar {\omega _0}/{k_{B}} $. The latter is characterized by $ {p_n} = ( {1 - {{ e}^{ - \theta }}} ) {{ e}^{ - n\theta }} $, corresponding to a Bose–Einstein average occupation $ \bar n = 1/( {{{ e}^\theta } - 1} ) $, where $ \theta = \hbar {\omega _0}/{k_{ B}}T $. Examples of these distributions are plotted in Figs. 2(a)–2(c) for $ {\bar n} = 2 $, 10, and 50. Incidentally, a thermal population for low-energy modes (e.g., $ \hbar {\omega _0} \sim 10 \,\, \text{meV} $) requires temperatures $ \sim {\bar n} \hbar {\omega _0}/{k_{ B}} \sim {\bar n} \times 100\,\, \text{K} $ that are for example attainable in graphene $ \pi $ electrons under ultrafast optical pumping [38]; alternatively, this level of thermal population could be reached under external intense chaotic light irradiation.

Figure 1 already shows a clear dependence on boson population statistics by examining the gain-to-loss ratio $ \big({\sum\nolimits_{\ell \gt 0} {P_\ell }} \big)/\big( {\sum\nolimits_{\ell \lt 0} {P_\ell }} \big)$. In particular, this quantity displays oscillations for a Fock state when either $ \bar n $ or $ |{\beta _0}| $ is varied, presumably as a remnant of similar oscillations observed in the associated electron spectra (see below). In contrast, averaging over the different Fock states involved in coherent and thermal distributions produces more monotonic ratios. Interestingly, Fock states lead to ratios closer to 1 within the range of parameters explored in the figure, indicating that the probabilities for climbing or descending the $ n $ ladder (i.e., $ \propto n + 1 $ and $ n $) quickly approach similar values as $ n $ increases compared with coherent and thermal distributions, which involve a substantial contribution from low $ n $’s up to relatively large $ \bar n $.

We further present in Figs. 2(d)–2(o) the evolution of the transmitted electron spectra for each of the distributions as a function of the coupling parameter $ \sqrt {\bar n} |{\beta _0}| $ (vertical scales). The spectra become more asymmetric for smaller $ \bar n $ because the number of gains cannot substantially exceed $ \bar n $ (i.e., the electrons cannot absorb more bosons than present already in the sample; see also Fig. 1), as observed in recent experiments [39], while the number of losses increases indefinitely with $ |{\beta _0}| $. Incidentally, the similarity between spectra obtained for coherent and thermal states for low $ \bar n $ stems from the fact that their respective mode populations $ {p_n} $, although different, also bear some similarity in that limit, and the resulting spectra are fully determined by $ {p_n} $, as shown by Eqs. (3) and (4).

In contrast, in the $ {\bar n} \gg 1 $ limit, the electron spectra become symmetric with respect to $ \ell = 0 $ [Figs. 2(g), 2(k), and 2(o)] when $ |\ell | \ll {\bar n} $. We can then approximate the square roots of Eq. (2) by $ \sqrt {n + \ell } $, leading to the same equation as for PINEM [16,25,32,40], whose solution reads (see Appendix D)

## 6. INTERACTION WITH AN OPTICAL CAVITY POPULATED THROUGH PUMPED QES

As a feasible system to explore the above ideas, we consider an optical cavity (e.g., a Mie resonator, similar to those already probed by EELS [42]) hosting a spectrally isolated mode (frequency $ {\omega _0} $, inelastic damping rate $ \kappa $) fed by a number $ N $ of three-level QEs, as illustrated in Fig. 3(a). A similar scheme applies to four-level systems, which are also extensively used in experimental realizations of QEs coupled to optical cavities [43]. We note that optical cavities hosting quantum emitters are customarily used as gain media, for example, in infiltrated silica nanoparticles [44] and photonic crystals [45]. The emitters are initialized in their excited state at $ t = 0 $ (e.g., by optical pumping with a $ \pi $ pulse, typically requiring moderate laser fluences that depend on the kind of QE), from which they decay to an intermediate state by coupling to the cavity at a rate $ g $ (the same for all QEs for simplicity) with a resonant transition frequency $ {\omega _0} $. We further assume fast internal decay from the intermediate to the ground state so that each QE interacts with the cavity only once. Incidentally, if inhomogeneous irradiation of the emitters takes place so that only a fraction of them are excited, the system can then be understood as consisting of a smaller effective number of QEs. The combined probability $ p_n^m $ for a cavity Fock state $ |n\rangle $ with $ m $ remaining excited emitters follows the equation of motion $ dp_n^m/dt = g [ {n(m\, +\, 1)p_{n - 1}^{m + 1}\, -\, (n + 1)mp_n^m} ] \,+\, \kappa [ {(n + 1)p_{n + 1}^m - np_n^m} ] $, which we solve numerically to obtain the time-dependent distribution $ {p_n} = \sum\nolimits_{m = 0}^N p_n^m $.

Examples of the evolution of the resulting average mode population $ \bar n $ and second-order autocorrelation $ {{ g}^{(2)}} $ are shown in Figs. 3(b) and 3(c), respectively. The latter starts at $ {{ g}^{(2)}} = 2(1 - 1/N) $ at $ t = 0 $ and in the absence of damping evolves toward $ 1 - 1/N $ at long times (see Appendix I), as expected for an assembly of $ N $ single-photon emitters. For finite cavity damping, $ \bar n $ reaches a maximum $ \lt N $, from which it exhibits an exponential decay, while $ {{ g}^{(2)}} $ eventually jumps to large values when $ \bar n $ becomes very small. For a sufficiently large number of QEs, we find a substantial average population while $ {{g}^{(2)}} $ varies from nearly 2 down to a quantum regime characterized by $ \lt 1 $ values. This evolution strongly affects the resulting electron spectra as a function of the time $ t $ at which the electron–boson interaction occurs after pumping the QEs, as shown in Figs. 3(g)–3(i) under the assumption that the interaction time is small compared with both $ 1/g $ and $ 1/\kappa $. The spectra initially resemble those of the thermal distributions of Fig. 2, as expected from the $ {{g}^{(2)}} \approx 2 $ values, and gradually become similar to those of Fock states; for finite cavity damping, they undergo an attenuation similar to Fig. 2(g) as $ \sqrt {\bar n} |{\beta _0}| $ decreases [the right part of Figs. 3(h) and 3(i)]. In particular, as the average number of photons in the cavity goes down with time, gain and loss peaks concentrate at increasingly lower orders $ \ell $. Incidentally, under the conditions of Fig. 3(i), the oscillations in the spectra as a function of $ \ell $ are attenuated compared with those of less lossy cavities [Figs. 3(g) and 3(h)], and additionally, an oscillation-free regime is reached at long times when the population is severely reduced and $ {{ g}^{(2)}} $ takes high values [see also the dashed green curves in Figs. 3(b) and 3(c)].

In order to illustrate the feasibility of the assumed parameters, we consider a Mie resonance in a silicon sphere ($ \epsilon = 12 $, radius $ a $) for which we estimate (see Appendix H) a sharp resonance of quality factor $ {\omega _0}/\kappa \sim {10^4} $ for a size parameter $ {\omega _0}a/c = 2.6775 $, causing an increase in the decay rate of QEs embedded in it (Purcell effect [46]) by an enhancement factor $ g/{g_0} \sim 100 $, so for natural QE decay rates $ {g_0} \sim \text{GHz} $ and optical frequencies $ {\omega _0} \sim 100 \,\,\text{THz} $, we have $ \kappa /g \sim 0.1 $, while the coupling parameter for this mode is $ |{\beta _0}| \sim 0.03 $ with 100–200 keV electrons.

## 7. CONCLUSIONS

In summary, we have shown that the interaction of electron beams with near optical fields depends on both the quantum nature of the sample excitations (fermionic versus bosonic) and the statistics of their populations. For bosonic modes, the spectral distribution of losses and gains varies dramatically when comparing Fock, coherent, and thermal distributions. Our simulations reveal that these regimes can be explored by populating an optical cavity through optically pumped quantum emitters, for which we have elaborated a model based on realistic optical cavities (e.g., Mie resonators, such as those used in a recent PINEM experiment [47]) infiltrated with optically pumped quantum emitters (e.g., gain atoms or molecules like rhodamine). We predict that the autocorrelation functions of the population distributions are directly retrievable from the peak intensities in the electron spectra. An implementation of this idea to probe the statistics of light with electron beams might be facilitated by schemes in which light–electron coupling is enhanced, for example, by using mirrors [25], photonic cavities [47,48], or aloof interaction with total internally reflected light in a prism [49]. Our results hold strong potential for resolving the nature of many-body excitations supported by complex materials, including strongly correlated systems, and open an unexplored avenue for the study of ultrafast plasmon, hot-electron, and phonon dynamics in optically pumped nanostructures using time-resolved electron microscope spectroscopy.

## APPENDIX A. INTERACTION BETWEEN FAST ELECTRONS AND OPTICAL MODES: A RELATIVISTIC QUANTUM OPTICS APPROACH

The interaction of a beam electron with the evanescent field surrounding an illuminated nanostructure has been extensively studied under the assumption of a classical light field [15,16,25,32,40]. Here, we present a self-contained derivation of the theory needed to describe the interaction with nonclassical evanescent fields. But first, we start by revisiting the interaction with a classical field in a fully relativistic approach [40]. We represent the optical field by the scalar and vector potentials $ \varphi $ and $ {\bf A} $ and describe the interaction with the electron through the Dirac equation [50,51]

We focus in this work on incident electrons prepared in a state involving wave vectors $ {\bf k} $ tightly packed around a central value $ {{\bf k}_0} $ and further adopt the nonrecoil approximation by assuming that the interaction with the optical field effectively produces $ {\bf k} $ components also satisfying $ |{\bf k} - {{\bf k}_0}| \ll {k_0} $. Under these conditions, we approximate the energy in Eq. (A4) by the first-order Taylor expansion $ {E_k} \approx {E_0} + ({\hbar ^2}{c^2}/{E_0}){{\bf k}_0} \cdot ({\bf k} - {{\bf k}_0}) $, where $ {E_0} = {E_{{k_0}}} $. We further notice that $ {\bf k} $ inside the $ {\bf k} $ sum can be substituted by $ - { i}\nabla $ outside it. Putting these elements together, we can recast Eq. (A1) as

We now extend Eq. (A7) to account for nonclassical fields by working in a gauge in which $ \varphi = 0 $ and substituting the vector potential by the operator [52]

## APPENDIX B. ANALYTICAL SOLUTION OF EQ. (3) IN THE MAIN TEXT

Equation (2) is obtained from Eqs. (A8)–(A10) by assuming a single mode $ j = 0 $ and expressing the wave function in terms of amplitudes $ f_\ell ^n $ through Eq. (1). Importantly, Eq. (2) is formally equivalent to the Schrödinger equation for the coupling of a quantum harmonic oscillator to a classical perturbation $ g(t) $. In order to show this in a tutorial way, we review some textbook physics and consider simpler forms of the free and interaction Hamiltonians $ {\hat {\cal H}_0} = \hbar {\omega _0}{\hat a^\dagger }\hat a $ and $ {\hat {\cal H}_1} = {g^*}(t){\hat a^\dagger } + g(t)\hat a $, and write the Schrödinger equation $ ({ i}\hbar \partial /\partial t)|{\psi _{\text{S}}}\rangle = ( {{{\hat {\cal H}}_0} + {{\hat {\cal H}}_1}} )|{\psi _{S}}\rangle $ with the wave function $ |{\psi _{S}}\rangle = \sum\nolimits_n {\alpha _n}{{ e}^{ - { i}n{\omega _0}t}}|n\rangle $ expanded in terms of number states $ |n\rangle $. The time-dependent amplitudes $ {\alpha _n} $ satisfy the equation

In order to obtain an analytical solution, it is convenient to move to the interaction picture, where the wave function $ |{\psi _{ i}}\rangle = {{ e}^{{ i}{{\hat {\cal H}}_0}t/\hbar }}|{\psi _{S}}\rangle\! $ satisfies $\! ({ i}\hbar \partial /\partial t)|{\psi _{ i}}\rangle = {\hat {\cal H}_{1{ i}}}|{\psi _{ i}}\rangle $ and $\! {\hat {\cal H}_{1{ i}}} = {{ e}^{{ i}{{\hat {\cal H}}_0}t/\hbar }}\!\!$ ${\hat {\cal H}_1}{{ e}^{ - { i}{{\hat {\cal H}}_0}t/\hbar }} = {g^*}(t){{ e}^{{ i}{\omega _0}t}}{\hat a^\dagger } + g(t){{ e}^{ - { i}{\omega _0}t}}\hat a $. In terms of number states, we have

The sought-after analytical solution is then obtained as [33,34] where is the time-evolution operator,## APPENDIX C. ELECTRON–BOSON INTERACTION IN THE WEAK-COUPLING LIMIT

A perturbative solution can be produced for Eq. (2) of the main text in the weak-coupling limit, provided the variations of all amplitudes $ f_\ell ^n $ are small during electron–boson interaction. When preparing the incident electron in $ \ell = 0 $ and the boson in state $ |n\rangle $, the nonvanishing elements of the perturbation series $ f_\ell ^n = \sum\nolimits_{s = 0}^\infty f_\ell ^{n,s} $ satisfy the equations

In the above series expansion, the lowest-order contribution to the $ \ell \gt 0 $ gain peak corresponds to the amplitude $ f_\ell ^{n - \ell ,\ell } $, which satisfies the equation

For coherent states (i.e., a Poissonian distribution), we have $ {{g}^{(\ell )}} = 1 $ (with $ \ell \gt 0 $) [36], leading to gain peak intensity ratios $ {P_\ell }/P_1^\ell = (1/\ell {!)^2} $. In contrast, for a thermal distribution one has $ {{ g}^{(\ell )}} = \ell ! $ [36], which produces more intense gain peaks with $ {P_\ell }/P_1^\ell = 1/\ell ! $. We stress that these results are valid only for weak interactions, as we are assuming that $ |\beta | = \sqrt {\bar n} |{\beta _0}| \ll 1 $.

Incidentally, if we also assume $ \bar n \gg 1 $ then $ {P_\ell } $ is given by analytical expressions for coherent and thermal mode populations [see Eq. (6) in the main text and Section Appendix D below], for which, using the small argument approximation $ \approx {(z/2)^\ell }/\ell ! $ of both $ {J_\ell }(z) $ and $ {I_\ell }(z) $ with $ \ell \ge 0 $ [54], we find $ {P_\ell } \approx |\beta {|^{2\ell }}/(\ell {!)^2} $ and $ {P_\ell } \approx |\beta {|^{2\ell }}/\ell ! $, respectively, therefore directly recovering the above results for the $ {P_\ell }/P_1^\ell $ ratio. Additionally, we note that Fock states lead to the same ratio as coherent states in the $ \bar n \gg 1 $ limit because they are characterized by $ {{ g}^{(\ell )}} = {\bar n}({\bar n} - 1) \ldots ({\bar n} - \ell + 1)/{{\bar n}^\ell } \approx 1 $ (with $ \ell \gt 0 $).

## APPENDIX D: ELECTRON–BOSON INTERACTION IN THE $ {\bar N} \gg 1 $ LIMIT

In the $ {\bar n} \gg 1 $ limit, the bulk of the electron–boson interaction involves high $ n $’s, which we consider to be much larger than the net number of quanta exchanges $ \ell $. This condition is satisfied if the interaction-strength parameter is small ($ |{\beta _0}| \ll 1 $), which is still compatible with a high total effective interaction $ {\bar n} |{\beta _0}{|^2} \sim 1 $ for sufficiently large $ \bar n $. We can then approximate both $ \sqrt n $ and $ \sqrt {n + 1} $ by $ \sqrt {n + \ell } $ in Eq. (2) of the main text; for each value of $ n + \ell $, which is conserved during propagation of the electron amplitudes $ f_\ell ^n $ along $ z $, the resulting equation coincides with Eq. (3) of Ref. [32] for the PINEM interaction with an optical field $ {{\cal E}_z} = {{\cal E}_{0z}}\sqrt {n + \ell } $, and therefore, we take from that reference the solution $ f_\ell ^n(\infty ) = \sqrt {{p_{n + \ell }}} {{ e}^{{ i}\ell \arg \{ - {\beta _0}\} }}{J_\ell }(2\sqrt {n + \ell } |{\beta _0}|) $, where $ {\beta _0} $ is the electron-mode interaction parameter defined by Eq. (5) in the main text and $ {p_{n + \ell }} $ is the population distribution of the mode Fock state $ |n + \ell \rangle $ before interaction with the electron. Because $ n \gg |\ell | $, we can further approximate $ n + \ell \approx n $ and write the probability associated with a net number $ \ell $ of quanta exchanges [Eq. (3) in the main text] as

For a boson prepared in the Fock state $ |n\rangle $, the interaction probabilities reduce to $ {P_\ell } = J_\ell ^2(2|\beta |) $, where $ \beta = \sqrt {\bar n} {\beta _0} $ and obviously the average mode population is $ \bar n = n $. Likewise, for a coherent state the population distribution approaches a Gaussian $ {p_n} \approx {{ e}^{ - {{(n - \bar n)}^2}/2\bar n}}/\sqrt {2\pi \bar n} $ in the $ {\bar n} \gg 1 $ limit [55], the width of which ($ \sim \sqrt {\bar n} $) becomes increasingly small compared with the average population $ \bar n $ as this one increases; we can thus approximate $ n \approx \bar n $ inside the Bessel function of Eq. (D1), which leads to $ {P_\ell } \approx J_\ell ^2(2|\beta |) $, that is, the same result as for the Fock state. We thus conclude that in the large average population limit both Fock and coherent states of the boson mode produce the same types of electron spectra as observed in PINEM experiments.The situation is however different for a chaotic thermal distribution $ {p_n} = (1 - {{ e}^{ - \theta }}) {{ e}^{ - n\theta }} $ with average population $ {\bar n} = 1/({{ e}^\theta } - 1) $, where $ \theta = \hbar {\omega _0}/{k_{B}}T $ and $ T $ is the mode temperature. We approach the $ {\bar n} \gg 1 $ limit at high temperatures, for which $ \theta \ll 1 $, and consequently, $ \theta \approx 1/{\bar n} $ and $ {p_n} \approx {{ e}^{ - n/{\bar n}}}/{\bar n} $. Inserting these expressions into Eq. (D1) and approximating the sum as an integral with the change of variable $ n/{\bar n} = {x^2} $, we find

## APPENDIX E. STEADY-STATE POPULATION OF OPTICALLY PUMPED BOSONIC AND FERMIONIC MODES

## 1. FERMIONIC MODE

A two-level system (states $ j = 0 $, 1 of energies $ \hbar {\varepsilon _j} $) coupled to a monochromatic light field $ {\bf E}(t) = {{\bf E}_0}{{ e}^{ - { i}\omega t}} + {\bf E}_0^*{{ e}^{{ i}\omega t}} $ constitutes a textbook example of light–matter interactions, commonly described through the optical Bloch equations [57]. The density matrix of the system satisfies the equation of motion

The interaction with a fast electron can be described following exactly the same formalism discussed in the main text for a bosonic mode but replacing the commutating bosonic operator $ \hat a $ by the anticommutating fermionic operator $ \hat \sigma $. This prescription leads to Eq. (2) with $ f_\ell ^n $ vanishing unless $ n = 0 $ or 1. In the weak electron–fermion coupling regime, this results in loss and gain probabilities

## 2. BOSONIC MODE

A bosonic mode is described by an expression identical to Eq. (E1) with $ \hat \sigma $ and $ {\hat \sigma ^\dagger } $ replaced by the annihilation and creation operators $ a $ and $ {a^\dagger } $, respectively, which now satisfy the commutation relation $ [\hat a,{\hat a^\dagger }] = 1 $. Additionally, $ j $ must be replaced by a boson occupation index $ n = 0,1, \ldots $ that labels Fock states $ |n\rangle $, while $ {\varepsilon _j} $ must be substituted by the ladder frequencies $ (n + 1/2){\omega _0} $. The boson density matrix admits a rigorous analytical solution [33,58]: $ \hat \rho = |\xi (t)\rangle \langle \xi (t)| $, where $ |\xi \rangle = {{ e}^{ - |\xi {|^2}/2}}\sum\nolimits_{n = 0}^\infty ( {{\xi ^n}/\sqrt {n!} } ){{ e}^{ - { i}n{\omega _0}t}}|n\rangle $ is a coherent state of amplitude $ \xi $ (in the interaction picture) satisfying $ a|\xi \rangle = {{ e}^{ - { i}{\omega _0}t}}\xi |\xi \rangle $ [37] and describing a Poisonian distribution with occupation probabilities $ {p_n} = \langle n|\hat \rho |n\rangle = {{ e}^{ - |\xi {|^2}}}|\xi {|^{2n}}/n! $ and average population $ {\bar n} = |\xi {|^2} $. By inserting this expression of $ \hat \rho $ into the equation of motion, one finds the solution $ \xi (t) = ( - { i}/\hbar )\int_{ - \infty }^t dt^\prime g(t^\prime ){{ e}^{{ i}{\omega _0}t^\prime - (t - t^\prime )\kappa /2}} $ for the mode amplitude driven by and external classical perturbation $ g(t) $. In particular, for the monochromatic light field considered above, this integral leads to

## APPENDIX F. INTERACTION WITH A DIPOLAR MODE AND ENSUING EELS PROBABILITY

We consider a dipolar mode of frequency $ {\omega _0} $ characterized by a transition electric dipole moment $ {\bf p} $ placed at the origin. As an ansatz, we write the single-mode electric field as the one produced by this dipole

In the weak-coupling regime (see Appendix E), the integral of the EELS probability over the mode spectral peak when the mode is initially depleted ($ {\bar n} = 0 $) reduces to $ {P_{ - 1}} = |{\beta _0}{|^2} = (2e{\omega _0}/\hbar {v^2}\gamma {)^2} [ {|{p_x}{|^2}K_1^2(\zeta ) + {{(|{p_z}|/\gamma )}^2}K_0^2(\zeta )} ] $. For an isotropic particle characterized by a triply degenerate mode of electric dipoles along the Cartesian directions, the probability is given by the sum over the three polarization directions, which amounts to setting $ {p_x} = {p_y} = {p_z} = p $; this leads to

In order to corroborate the correctness of the normalization of $ {\vec {\cal E}_0} $ in the above ansatz, we compare $ P_{ - 1}^{\text{isotropic}} $ with the result derived from the classical EELS probability for an isotropic dipolar particle [6], $ {\Gamma _{\text{EELS,dip}}}(\omega ) = (1/\hbar \pi )(2e\omega /{v^2}\gamma {)^2}f(\omega b/v\gamma ) \text{Im}\{ \alpha (\omega )\} $, where $ \alpha (\omega ) $ is the polarizability. Linear response theory allows us to write the latter as [59] $ \alpha (\omega ) = (|p{|^2}/\hbar ) [ {1/({\omega _0} - \omega - { i}{0^ + }) + 1/({\omega _0} + \omega + { i}{0^ + })} ] $ in terms of the mode dipole $ p $ and frequency $ {\omega _0} $, which upon insertion into the spectral integral $ {P_{ - 1}} = \int_0^\infty \text{d} \omega \;{\Gamma _{\text{EELS,dip}}}(\omega ) $ reproduces Eq. (F1), therefore confirming the ansatz.

## APPENDIX G. ELECTRON–BOSON COUPLING STRENGTH FOR PLASMONIC CAVITIES

Plasmons in metallic nanoparticles constitute excellent candidates to explore the interaction between electrons and optical cavities. Here, we estimate the coupling parameter $ |{\beta _0}| $ [see the leading factor in Eq. (F1)] for two types of metallic nanoparticles in which the aspect ratio allows one to tune their frequency, also affecting the values of $ |{\beta _0}| $; in particular, we consider prolate ellipsoids and spherical shells made of silver (permittivity $ \epsilon (\omega ) \approx { \epsilon _{ b}} - \omega _{ p}^2/\omega (\omega + { i}/\tau ) $ with $ { \epsilon _{b}} \approx 4 $, $ \hbar {\omega _{ p}} = 9.17 \,\,\text{eV} $, and $ \hbar /\tau = 21 \,\,\text{meV} $ [60]). Following similar methods as those of Ref. [61], a prolate ellipsoid of volume $ V $ is found to exhibit a normal-to-the-symmetry-axis resonance frequency $ {\omega _0} = {\omega _{ p}}/\sqrt {{ \epsilon _{ b}} - { \epsilon _0}} $ and an effective transition dipole $ p \approx (1 - { \epsilon _0})\sqrt {\hbar {\omega _0}V/8\pi ({ \epsilon _{b}} - { \epsilon _0})} $, where $ { \epsilon _0} = 1 - 1/L $, $ L = ({r^2}/2){\Delta ^{ - 3}} [ {\pi /2 - \text{arc} \tan (1/\Delta ) - \Delta /{r^2}} ] $ is the depolarization factor for a height-to-diameter aspect ratio $ r \gt 1 $, and $ \Delta = \sqrt {{r^2} - 1} $. For a spherical metal shell of small thickness $ t $ compared with the radius $ a $, filled with a core dielectric $ { \epsilon _{ c}} $, and treated in the $ t \ll a $ limit, we find $ {\omega _0} \approx ( {{\omega _{\text{p}}}/\sqrt {{ \epsilon _{\text{c}}} + 2} } )\sqrt {2t/a} $ and $ p \approx \sqrt {3\hbar {\omega _0}{a^3}/2({ \epsilon _{ c}} + 2)} $. Numerical inspection of these two types of particles yields optimum coupling values that can reach $ |{\beta _0}| \sim 0.1 $ with plasmon energies in the 1 eV region, particle diameters of $ \sim 20 \,\text{nm} $, and electron energies $ \sim 10 \,\text{keV} $ when considering either disk-like prolate ellipsoids (aspect ratio $ r \sim 5 $) or thin shells ($ t/a \sim 0.2 $) filled with silica ($ { \epsilon _{ c}} = 2 $). The values adopted in the main paper are commensurate with these parameters.

## APPENDIX H. INTERACTION STRENGTH OF QUANTUM EMITTERS AND BEAM ELECTRONS WITH DIELECTRIC OPTICAL CAVITIES

High-index dielectric nanostructures can trap light with small radiative leakage. For example, a Si sphere of radius $ a $ modeled with a permittivity $ \epsilon = 12 $ exhibits a narrow resonance at a size parameter $ {\rho _0} = {\omega _0}a/c \approx 2.6775 $ with a quality factor (frequency divided by width) $ Q = {\omega _0}/\kappa \sim {10^4} $. Analytical EELS calculations based on a previously published formula [6] predict a peak-integrated excitation probability $ \sim {10^{ - 3}} $ (i.e., $ |{\beta _0}| \sim 0.03 $) for grazingly passing electrons of 100–200 keV kinetic energy (see the illustrative calculation in Fig. S3 of Supplement 1). Besides this relatively high probability, the large $ Q $ value of the Mie resonance under consideration produces a high Purcell enhancement in QEs when they are embedded inside the structure, implying nearly perfect QE–cavity coupling and negligible radiative losses. Indeed, following a quantum optics formalism for dispersionless and lossless dielectrics [62], we can express the electromagnetic Green tensor as $ G({\bf r},{\bf r}^\prime ,\omega ) = \sum\nolimits_j {{\bf f}_j}({\bf r}) \otimes {\bf f}_j^*({\bf r}^\prime )/({\omega ^2} - \omega _j^2 + { i}{0^ + }) $, where the sum extends over photon modes $ j $ in the presence of the cavity, and the mode functions are orthonormalized according to $ \int \text{d}^{3}{\bf r}\;{{\bf f}_j}({\bf r}) \cdot {\bf f}_{j^\prime }^*({\bf r}) \epsilon ({\bf r}) = {\delta _{jj^\prime }} $; we now argue that the modes contributing to the Mie resonance $ {\omega _0} $ should have similar spatial profiles inside the cavity, so we approximate $ G({\bf r},{\bf r}^\prime ,\omega ) \approx {{\bf f}_0}({\bf r}) \otimes {\bf f}_0^*({\bf r}^\prime )/(\omega (\omega + { i}\kappa ) - \omega _0^2) $ by phenomenologically introducing the resonance width $ \kappa $; the Purcell enhancement factor (EF) is then proportional to the local density of optical states normalized to the vacuum value [63], which can be calculated from $ G $ as $ \text{EF} \approx ( - 6\pi {c^3}/\omega )\text{Im}\{ \hat{\bf n} \cdot G({\bf r},{\bf r},\omega ) \cdot \hat{\bf n}\} $, with $ {\bf r} = {\bf r}^\prime $ corresponding to the position of the QE; for resonant coupling $ \omega = {\omega _0} $, arguing from the normalization condition that $ {|{{\bf f}_0}|^2} \sim 1/ \epsilon V $, where $ V = 4\pi {a^3}/3 $ is the cavity volume, we find $ \text{EF} \sim 9Q/2 \epsilon \rho _0^3 \sim 250 $ for the cavity under consideration. We now envision QEs with a natural decay rate $ {g_0} \sim \text{GHz} $, whose coupling rate increases to $ g = \text{EF} \times {g_0} \sim {10^2}\,\, \text{GHz} $; for an optical frequency in the $ {\omega _0} \sim 100 \,\, \text{THz} $ range, the cavity damping rate is $ \kappa \sim {\omega _0}/Q \sim 10\,\, \text{GHz} $, thus leading to small values of $ \kappa /g \sim 0.1 $ similar to those assumed in the main text.

## APPENDIX I. POPULATION OF A CAVITY COUPLED TO AN ENSEMBLE OF QUANTUM EMITTERS

We consider an optical cavity hosting a bosonic excitation mode of frequency $ {\omega _0} $ coupled to $ N $ three-level QEs that are prepared in their excited state at time $ t = 0 $ [see Fig. 3(a) in the main text]. The QEs can decay to their intermediate states by coupling to the cavity at a rate $ g $, which we assume to be the same for all QEs. The transition frequency in this decay is taken to coincide with the cavity mode $ {\omega _0} $. Once in their intermediate states, the QEs experience internal decay at a fast rate compared with $ g $ so that they no longer interact with the cavity. The present analysis also applies to four-level atoms, such as those used in lasers with their intermediate-states transition matching $ {\omega _0} $. We further incorporate an inelastic decay of the cavity mode at a rate $ \kappa $, associated with either intrinsic losses in the materials or radiative emission.

We characterize the QE–cavity system in terms of the time-dependent probabilities $ p_n^m(t) $ for the cavity to be in the Fock state $ |n\rangle $ while $ m $ emitters are still in their excited states. The equation of motion for these probabilities can be readily written as

We are interested in the cavity occupation distribution $ {p_n} = \sum\nolimits_m p_n^m $, the average population $ \bar n = \sum\nolimits_n n {p_n} $, and how they influence the coupling to a passing electron. As a rough way of characterizing the population distribution, we consider the instantaneous second-order correlation $ {{ g}^{(2)}} = ({\overline{n^2}} - \bar n)/{\bar n^2} $, where $ {\overline{ n^2}} = \sum\nolimits_n {n^2}{p_n} $. At small times right after pumping the QEs to their excited states, $ {{ g}^{(2)}} $ can be obtained analytically by Taylor expanding $ p_n^m $; we find $ {\bar n} = Ngt + O({t^2}) $ and $ {\overline{n^2}} - {\bar n} = 2N{(N - 1)(gt)^2} + O({t^3}) $, from which $ {{g}^{(2)}} = 2(1 - 1/N) $ at $ t{ = 0^ + } $. By numerically solving Eq. (I1), we observe that $ {{g}^{(2)}} $ is generally a decreasing function of $ gt $ (see Fig. 3 in the main text). In particular, for $ \kappa = 0 $ the number of excitations in the system must be conserved (i.e., $ n + m = N $), and eventually all of them are transferred to the cavity, leading to the asymptotic solution $ p_n^m(\infty ) = {\delta _{nN}}{\delta _{m0}} $, which in turn produces $ {{ g}^{(2)}} = 1 - 1/N $. Incidentally, a large increase in $ {{ g}^{(2)}} $ is eventually observed when $ {\bar n} $ drops to very small values in lossy cavities.

We present additional numerical results for cavities with $ N = 100 $ QEs in Fig. S4 and $ N = 1,2 $, and 5 in Fig. S5 (see Supplement 1), showing stronger side bands as both $ N $ and $ {\beta _0} $ are increasing, and with the gain side limited to $ \ell = N $.

## Funding

Ministerio de Economía y Competitividad (MAT2017-88492-R, SEV2015-0522); European Research Council (Advanced Grant 789104-eNANO); Catalan CERCA (CERCA); Fundación Cellex (Fundació Privada Cellex); European Commission (Marie Skłodowska-Curie Grant 713729).

## Acknowledgment

We thank Vahagn Mkhitaryan, Matteo Lostaglio, and Sophie Meuret for helpful and enjoyable discussions.

See Supplement 1 for supporting content.

## REFERENCES

**1. **A. Howie, “Valence excitations in
electron microscopy: resolved and unresolved issues,”
Micron **34**,
121–125 (2003). [CrossRef]

**2. **K. A. Mkhoyan, T. Babinec, S. E. Maccagnano, E. J. Kirkland, and J. Silcox, “Separation of bulk and
surface-losses in low-loss eels measurements in stem,”
Ultramicroscopy **107**,
345–355 (2007). [CrossRef]

**3. **O. L. Krivanek, T. C. Lovejoy, N. Dellby, T. Aoki, R. W. Carpenter, P. Rez, E. Soignard, J. Zhu, P. E. Batson, M. J. Lagos, R. F. Egerton, and P. A. Crozier, “Vibrational spectroscopy in
the electron microscope,” Nature **514**, 209–214
(2014). [CrossRef]

**4. **D. A. Muller, L. Fitting Kourkoutis, M. Murfitt, J. H. Song, H. Y. Hwang, J. Silcox, N. Dellby, and O. L. Krivanek, “Atomic-scale chemical imaging
of composition and bonding by aberration-corrected
microscopy,” Science **319**, 1073–1076
(2008). [CrossRef]

**5. **G. Zhu, G. Radtke, and G. A. Botton, “Bonding and structure of a
reconstructed (001) surface of SrTiO_{3} from
TEM,” Nature **490**,
384–387 (2012). [CrossRef]

**6. **F. J. G. de Abajo, “Optical excitations in
electron microscopy,” Rev. Mod. Phys. **82**, 209–275
(2010). [CrossRef]

**7. **D. Rossouw and G. A. Botton, “Plasmonic response of bent
silver nanowires for nanophotonic subwavelength
waveguiding,” Phys. Rev. Lett. **110**, 066801
(2013). [CrossRef]

**8. **M. Kociak and O. Stephan, “Mapping plasmons at the
nanometer scale in an electron microscope,”
Chem. Soc. Rev. **43**,
3865–3883 (2014). [CrossRef]

**9. **A. T. A. Hörl and U. Hohenester, “Full three-dimensional
reconstruction of the dyadic green tensor from electron energy loss
spectroscopy of plasmonic nanoparticles,” ACS
Photon. **2**,
1429–1435 (2015). [CrossRef]

**10. **G. Guzzinati, A. Beche, H. Lourenco-Martins, J. Martin, M. Kociak, and J. Verbeeck, “Probing the symmetry of the
potential of localized surface plasmon resonances with phase-shaped
electron beams,” Nat. Commun. **8**, 14999 (2017). [CrossRef]

**11. **M. J. Lagos, A. Trügler, U. Hohenester, and P. E. Batson, “Mapping vibrational surface
and bulk modes in a single nanocube,”
Nature **543**,
529–532 (2017). [CrossRef]

**12. **B. Barwick, D. J. Flannigan, and A. H. Zewail, “Photon-induced near-field
electron microscopy,” Nature **462**, 902–906
(2009). [CrossRef]

**13. **A. Ryabov and P. Baum, “Electron microscopy of
electromagnetic waveforms,” Science **353**, 374–377
(2016). [CrossRef]

**14. **M. Kozák, J. McNeur, K. J. Leedle, H. Deng, N. Schönenberger, A. Ruehl, I. Hartl, J. S. Harris, R. L. Byer, and P. Hommelhoff, “Optical gating and streaking
of free electrons with sub-optical cycle precision,”
Nat. Commun. **8**, 14342
(2017). [CrossRef]

**15. **F. J. G. de Abajo, A. A. Garcia, and M. Kociak, “Multiphoton absorption and
emission by interaction of swift electrons with evanescent light
fields,” Nano Lett. **10**, 1859–1863
(2010). [CrossRef]

**16. **S. T. Park, M. Lin, and A. H. Zewail, “Photon-induced near-field
electron microscopy (PINEM): theoretical and
experimental,” New J. Phys. **12**, 123028 (2010). [CrossRef]

**17. **F. O. Kirchner, A. Gliserin, F. Krausz, and P. Baum, “Laser streaking of free
electrons at 25 keV,” Nat. Photonics **8**, 52–57
(2014). [CrossRef]

**18. **L. Piazza, T. T. A. Lummen, E. Quiñonez, Y. Murooka, B. Reed, B. Barwick, and F. Carbone, “Simultaneous observation of
the quantization and the interference pattern of a plasmonic
near-field,” Nat. Commun. **6**, 6407 (2015). [CrossRef]

**19. **A. Feist, K. E. Echternkamp, J. Schauss, S. V. Yalunin, S. Schäfer, and C. Ropers, “Quantum coherent optical phase
modulation in an ultrafast transmission electron
microscope,” Nature **521**, 200–203
(2015). [CrossRef]

**20. **T. T. A. Lummen, R. J. Lamb, G. Berruto, T. LaGrange, L. D. Negro, F. J. G. de Abajo, D. McGrouther, B. Barwick, and F. Carbone, “Imaging and controlling
plasmonic interference fields at buried interfaces,”
Nat. Commun. **7**, 13156
(2016). [CrossRef]

**21. **K. E. Echternkamp, A. Feist, S. Schäfer, and C. Ropers, “Ramsey-type phase control of
free-electron beams,” Nat. Phys. **12**, 1000–1004
(2016). [CrossRef]

**22. **A. Feist, N. Bach, T. D. N. R. da Silva, M. Mäller, K. E. Priebe, T. Domräse, J. G. Gatzmann, S. Rost, J. Schauss, S. Strauch, R. Bormann, M. Sivis, S. Schäfer, and C. Ropers, “Ultrafast transmission
electron microscopy using a laser-driven field emitter: femtosecond
resolution with a high coherence electron beam,”
Ultramicroscopy **176**,
63–73
(2017). [CrossRef]

**23. **K. E. Priebe, C. Rathje, S. V. Yalunin, T. Hohage, A. Feist, S. Schäfer, and C. Ropers, “Attosecond electron pulse
trains and quantum state reconstruction in ultrafast transmission
electron microscopy,” Nat. Photonics **11**, 793–797
(2017). [CrossRef]

**24. **E. Pomarico, I. Madan, G. Berruto, G. M. Vanacore, K. Wang, I. Kaminer, F. J. G. de Abajo, and F. Carbone, “meV resolution in
laser-assisted energy-filtered transmission electron
microscopy,” ACS Photon. **5**, 759–764
(2018). [CrossRef]

**25. **G. M. Vanacore, I. Madan, G. Berruto, K. Wang, E. Pomarico, R. J. Lamb, D. McGrouther, I. Kaminer, B. Barwick, F. J. G. de Abajo, and F. Carbone, “Attosecond coherent control of
free-electron wave functions using semi-infinite light
fields,” Nat. Commun. **9**, 2694 (2018). [CrossRef]

**26. **W. Cai, O. Reinhardt, I. Kaminer, and F. J. G. de Abajo, “Efficient orbital angular
momentum transfer between plasmons and free
electrons,” Phys. Rev. B **98**, 045424
(2018). [CrossRef]

**27. **Y. Morimoto and P. Baum, “Attosecond control of electron
beams at dielectric and absorbing membranes,”
Phys. Rev. A **97**,
033815 (2018). [CrossRef]

**28. **L. H. G. Tizei and M. Kociak, “Spatially resolved quantum
nano-optics of single photons using an electron
microscope,” Phys. Rev. Lett. **110**, 153604
(2013). [CrossRef]

**29. **S. Meuret, L. H. G. Tizei, T. Cazimajou, R. Bourrellier, H. C. Chang, F. Treussart, and M. Kociak, “Photon bunching in
cathodoluminescence,” Phys. Rev.
Lett. **114**, 197401
(2015). [CrossRef]

**30. **R. Bourrellier, S. Meuret, A. Tararan, O. Stéphan, M. Kociak, L. H. G. Tizei, and A. Zobelli, “Bright UV single photon
emission at point defects in h-BN,” Nano
Lett. **16**,
4317–4321 (2016). [CrossRef]

**31. **O. Kfir, “Entanglements of electrons and
cavity photons in the strong-coupling regime,”
Phys. Rev. Lett. **123**,
103602 (2019). [CrossRef]

**32. **F. J. G. de Abajo, B. Barwick, and F. Carbone, “Electron diffraction by
plasmon waves,” Phys. Rev. B **94**, 041404
(2016). [CrossRef]

**33. **P. Carruthers and M. M. Nieto, “Coherent states and the forced
quantum oscillator,” Am. J. Phys. **33**, 537–544
(1965). [CrossRef]

**34. **A. A. Lucas, E. Kartheuser, and R. G. Badro, “Electron-phonon interaction in
dielectric films. application to electron energy loss and gain
spectra,” Phys. Rev. B **2**, 2488–2499
(1970). [CrossRef]

**35. **S. Franke, S. Hughes, M. K. Dezfouli, P. T. Kristensen, K. Busch, A. Knorr, and M. Richter, “Quantization of quasinormal
modes for open cavities and plasmonic cavity quantum
electrodynamics,” Phys. Rev. Lett. **122**, 213901
(2019). [CrossRef]

**36. **R. Loudon, *The Quantum Theory of Light*
(Oxford University,
2000).

**37. **R. J. Glauber, “Coherent and incoherent states
of the radiation field,” Phys. Rev. **131**, 2766–2788
(1963). [CrossRef]

**38. **I. Gierz, J. C. Petersen, M. Mitrano, C. Cacho, I. C. E. Turcu, E. Springate, A. Stöhr, A. Köhler, U. Starke, and A. Cavalleri, “Snapshots of non-equilibrium
Dirac carrier distributions in graphene,” Nat.
Mater. **12**,
1119–1124 (2013). [CrossRef]

**39. **P. Das, J. D. Blazit, M. Tencé, L. F. Zagonel, Y. Auad, Y. H. Lee, X. Y. Ling, A. Losquin, O. S. C. Colliex, F. J. G. de Abajo, and M. Kociak, “Stimulated electron energy
loss and gain in an electron microscope without a pulsed electron
gun,” Ultramicroscopy **203**, 44–51
(2019). [CrossRef]

**40. **S. T. Park and A. H. Zewail, “Relativistic effects in
photon-induced near field electron microscopy,”
J. Phys. Chem. A **116**,
11128–11133 (2012). [CrossRef]

**41. **J. H. Eberly, B. Shore, Z. Białynicka-Birula, and I. Białynicki-Birula, “Coherent dynamics of n-level
atoms and molecules. I. Numerical experiments,”
Phys. Rev. A **16**,
2038–2047 (1977). [CrossRef]

**42. **J. K. Hyun, M. Couillard, P. Rajendran, C. M. Liddell, and D. A. Muller, “Measuring far-ultraviolet
whispering gallery modes with high energy electrons,”
Appl. Phys. Lett. **93**,
243106 (2008). [CrossRef]

**43. **A. Reiserer and G. Rempe, “Cavity-based quantum networks
with single atoms and optical photons,” Rev.
Mod. Phys. **87**,
1379–1418 (2015). [CrossRef]

**44. **R. Barbosa-Silva, A. F. Silva, A. M. Brito-Silva, and C. B. de Araújo, “Bichromatic random laser from
a powder of rhodamine-doped sub-micrometer silica
particles,” J. Appl. Phys. **115**, 043515
(2014). [CrossRef]

**45. **Y. Nishijima, K. Ueno, S. Juodkazis, V. Mizeikis, H. Fujiwara, K. Sasaki, and H. Misawa, “Lasing with well-defined
cavity modes in dye-infiltrated silica inverse opals,”
Opt. Express **17**,
2976–2983 (2009). [CrossRef]

**46. **E. M. Purcell, “Spontaneous emission
probabilities at radio frequencies,” Phys.
Rev. **69**, 37–38
(1946). [CrossRef]

**47. **O. Kfir, H. Lourenço-Martins, G. Storeck, M. Sivis, T. R. Harvey, T. J. Kippenberg, A. Feist, and C. Ropers, “Controlling free electrons
with optical whispering-gallery modes,”
arXiv:1910.09540.

**48. **K. Wang, R. Dahan, M. Shentcis, Y. Kauffmann, S. Tsesses, and I. Kaminer, “Coherent interaction between
free electrons and cavity photons,”
arXiv:1908.06206.

**49. **R. Dahan, S. Nehemia, M. Shentcis, O. Reinhardt, Y. Adiv, K. Wang, O. Beer, Y. Kurman, X. Shi, M. H. Lynch, and I. Kaminer, “Observation of the stimulated
quantum Cherenkov effect,”
arXiv:1909.00757.

**50. **A. Messiah, *Quantum Mechanics*
(North-Holland,
1966).

**51. **J. J. Sakurai, *Modern Quantum Mechanics*
(Addison-Wesley,
1994).

**52. **P. W. Milonni, *The Quantum Vacuum: An Introduction to
Quantum Electrodynamics* (Academic
Press, 1994).

**53. **A. A. Abrikosov, L. P. Gor’kov, and I. Y. Dzyaloshinskii, *Quantum Field Theoretical Methods in
Statistical Physics*
(Pergamon,
1965).

**54. **M. Abramowitz and I. A. Stegun, *Handbook of Mathematical
Functions* (Dover,
1972).

**55. **W. Feller, *An Introduction to Probability Theory
and Its Applications* (Wiley,
1968).

**56. **I. S. Gradshteyn and I. M. Ryzhik, *Table of Integrals, Series, and
Products* (Academic,
2007).

**57. **M. O. Scully and M. S. Zubairy, *Quantum Optics*
(Cambridge University,
1997).

**58. **F. J. G. de Abajo, “Multiple excitation of
confined graphene plasmons by single free electrons,”
ACS Nano **7**,
11409–11419 (2013). [CrossRef]

**59. **D. Pines and P. Nozières, *The Theory of Quantum Liquids*
(W. A. Benjamin,
1966).

**60. **P. B. Johnson and R. W. Christy, “Optical constants of the noble
metals,” Phys. Rev. B **6**, 4370–4379
(1972). [CrossRef]

**61. **E. J. C. Dias and F. J. G. de Abajo, “Fundamental limits to the
coupling between light and 2D polaritons,” ACS
Nano **13**,
5184–5197 (2019). [CrossRef]

**62. **R. J. Glauber and M. Lewenstein, “Quantum optics of dielectric
media,” Phys. Rev. A **43**, 467–491
(1991). [CrossRef]

**63. **F. J. G. de Abajo and M. Kociak, “Probing the photonic local
density of states with electron energy loss
spectroscopy,” Phys. Rev. Lett. **100**, 106804
(2008). [CrossRef]