Probing the Purcell effect without radiative decay: Lessons in the frequency and time domains

The microscopic processes underlying electro-optic sampling of quantum-vacuum fluctuations are discussed, leading to the interpretation of these experiments in terms of an exchange of virtual photons. With this in mind it is shown how one can directly study the Purcell effect, i.e. the changes induced by cavities upon the quantum vacuum, in the frequency and time domains. This forges a link between electro-optic sampling of the quantum vacuum and geometry-induced vacuum effects.

Ground-state fluctuations of the electromagnetic field can be seen as responsible for observable effects such as the Lamb shift [1], Casimir forces [2] or spontaneous emission [3]. Just as the ground-state fluctuations can alter states of matter, matter can in turn influence the quantum vacuum: ground-state fluctuations in close proximity to macroscopic objectsthe so-called polaritonic or medium-assisted quantum vacuum-are inherently different from their free-space counterpart. This is the essence of the Purcell effect, see Fig. 1 (b). Shaping the quantum vacuum using optical environments, such as cavities or surface plasmon-polaritons, can thus be exploited to e.g. enhance or suppress spontaneous emission rates [4] or resonant energy transfer [5]. In the strong coupling regime it even allows the alteration of chemical properties of molecules, for example their reaction rates [6,7].
A novel route to studying ground-state fluctuations of the electromagnetic field has been introduced by means of electro-optic sampling (EOS) experiments [8][9][10], see Fig. 1 (a). In the setup described in Ref. [9] two linearly-polarised, ultra-short laser pulses propagate through a nonlinear crystal, in which they couple to vacuum fluctuations of the electromagnetic field via the nonlinear susceptibility inside the crystal. This coupling leads to a change of the polarisation direction of the laser pulses, which can be measured to infer information about the quantum vacuum inside the crystal on a sub-cycle time scale. This has been used to measure bare [8,11] as well as squeezed [12,13] quantum-vacuum noise. By tuning the temporal and spatial shifts between the two laser pulses, one can additionally detect vacuum correlations between distinct spatio-temporal regions and this way access the spectrum of the electromagnetic ground state [9], making EOS a promising tool for an in-depth study of medium-assisted vacuum fluctuations [14,15].
In this work it is shown how EOS can be used to directly access environment-induced changes of the quantum vacuum. By interpreting ground-state correlations as the exchange of virtual photons, it is shown how one can observe the polaritonic quantum vacuum in the time domain. This allows one to study the dynamical formation of cavity modes in the quantum vacuum via multiple reflections, leading to a time-frequency uncertainty relation between the frequency resolution of the observed vacuum correlations and the time during which the quantum vacuum is observed. To do so, we build upon previous theoretical results which have been introduced and compared in Ref. [14,15] to experimental data without consideration of the Purcell effect. The EOS signal accounting for absorption and dispersion inside the nonlinear crystal, as well as allowing for general optical environments and pump-pulse profiles, was found to be given by: × F (r, r , Ω, Ω , δr , δt) Ê vac,x (r, Ω)Ê † vac,x (r , Ω ) .
Here, F (r, r , Ω, Ω , δr , δt) is a filter function depending on the spatio-temporal shape of the two laser pulses as well as on their lateral and temporal shifts δr and δt, respectively, see Ref. [16], and V C is the volume of the crystal. Ê vac (r, Ω)Ê † vac (r , Ω ) is the two-point ground-state correlation function of the electromagnetic field in an absorptive and dispersive   Figure 1. (a) Electro-optic sampling. Basic setup: Two laser pulses with tuneable spatial offset δr and delay δt mix with vacuum fluctuations inside a nonlinear crystal through its nonlinear susceptibility χ (2) such that the outgoing pulses contain information about the quantum vacuum inside the crystal. (b) Purcell effect. Spectrum of the coincidence limit of the vacuum-correlation function in the bulk of a ZnTe-crystal with infinite extension (solid line), at a distance d = 17 µm to a plate with a surface-plasmon polariton at Ω = 1.85 2πTHz (dotted-dashed) and inside a cavity of length L = 0.5 mm (dashed). arXiv:2107.02647v1 [quant-ph] 6 Jul 2021 Im(rp) a.u. d→∞ found when the plate is removed, and we use L = 1 mm. optical environment, provided by macroscopic QED as Ê vac (r, Ω)Ê † vac (r , Ω ) =h Ω 2 c 2 ε 0 π δ(Ω − Ω )ImG(r, r , Ω) [17]. Here, ε 0 is the free-space permittivity, c the speed of light in free space and G(r, r , Ω) is the dyadic Green's function describing the propagation of a photon at frequency Ω from r to r [16]. G and thus Ê 2 vac depend on the geometry of the optical environment [ Fig. 1 (b)]. Equation (1) suggests that by tuning the laser pulses and thus F , one can use EOS to access various characteristics of the two-point correlation functions of the electromagnetic ground state.
Throughout this paper we assume that all involved laser pulses are linearly polarised and Gaussian, with beam waist w = 80 µm, central frequency ω c = 375 × 2π THz and pulse duration ∆t = 80 fs [16]. Unless stated otherwise, we consider the nonlinear crystal of length L = 0.1 mm to be a ZnTe crystal whose optical characteristics are summarized in the supplementary material [16].
a. Theory In order to show how medium-induced changes of the quantum vacuum can be found in EOS experiments, we consider a plate and a cavity with different orientations, attached to the nonlinear crystal (see bottom right of Fig. 1a). In the presence of additional surfaces the Green's tensor splits into its bulk part G (0) and scattering part G (1) with G = G (0) + G (1) . The term G (1) accounts for all reflection effects, such that restricting G to G (0) is equivalent to neglecting all influences of any macroscopic object near to the nonlinear crystal as well as reflections at the surfaces of the crystal itself. We focus our attention on the changes these objects induce in the quantum vacuum correlations, and therefore neglect their influence on the laser pulses (aside from obscuring part of the beams). This can be justified by assuming that the reflection coefficients of the plate and the cavity are close to zero in the frequency range of the laser, but different from zero for the resolved frequencies of the vacuum field. Inserting the full Green's tensor into Eq. (1) we get two contributions, one stemming from the bulk G (0) and one from the scattering part G (1) . The latter describes the change of the correlation function due to the presence of the macroscopic plate(s). Neglecting absorption effects inside the nonlinear crystal (e.g. ε(Ω) real-valued) and applying the 'laser paraxial' approximation suitable in the parameter range considered here [14][15][16] the corresponding bulk (j = 0) and scattering (j = 1) contributions take the form [16] The full signal is given by g = g (0) + g (1) . In Eq. (2), √ C ∝ χ (2) LI (L: length of the crystal, χ (2) : nonlinear susceptibility, I: intensity of the laser pulses) determines the sampling efficiency [9] and E 2 vac (Ω) gives the strength of the vacuum fluctuations in a bulk crystal at frequency Ω. The integral over the parallel wave vector of the vacuum field q (q = n(Ω)Ω/c = q 2 ⊥ + q 2 ) describes the propagation of the virtual photon from one laser pulse to the other. Here, O (j) (q, δr ) accounts for the obscuring of the beam due to the presence of the plate. Most importantly, this integral contains the propagation factor p (j) (q, δr ) which for the bulk contribution is simply given by p (0) = q 2 /(q ⊥ q)e iq·δr , whereas p (1) depends on the chosen geometry of the attached plate(s). For a plate at x = −d for example one finds p (1) = p (0) e 2iq x d R p , where R p is the p-polarised Fresnel reflection coefficient. The additional factor in p (1) compared to p (0) accounts for the additional propagation of an exchanged virtual photon to the plate and back, see Fig. 3 (a), and further discussion below. Similar expressions can be found for other geometries and a list of all propagation factors p (1) considered here is found in the Supplementary Material [16]. Furthermore, the integral in Eq. (2) contains the known response function R(q) = e −(q 2 [14,15] accounting for the averaging over vacuum modes inside the finite lateral pulse profile, phase-matching and which contains the spectral autocorrelation function f (Ω) = e −πΩ 2 ∆t 2 /4 . b. Observing the Purcell effect We start by considering a reflecting plate attached to the nonlinear crystal in either the x = −d or y = −d planes, which is thus parallel to the propagation direction of the laser pulses offset by a distance d, see inset in the bottom right of Fig. 1 (a). The contribution to the EOS signal with 'coincident' pulses (δt = δr = 0) is shown in Fig. 2 as a function of d. First, in Fig. 2 (a), we assume a perfectly reflecting plate described by reflection coefficients for p and s polarised waves R p = 1 and R s = −1, respectively. The electro-optic sampling signal changes with the beamplate distance as a result of competition between two effects. On the one hand, the signal decreases when the beam is closer to the plate, since a larger fraction of the beam becomes obscured by it. On the other hand, the effects of reflection upon the vacuum field start contributing significantly at a distance of approximately d = 4w. The opposite signs of these additional, plate-induced contributions ('scattering contributions') for the cases where the plate is in the x = −d or y = −d plane can be understood by realising that they arise from image fluctuations: in the former (latter) geometry the image-fluctuations are parallel (antiparallel) with respect to the x-polarised fluctuations the laser pulse singles out. This leads to same (opposite) signs of the image fluctuating field compared to the bare fluctuations and thus to an enhancement (reduction) of the total fluctuating field. In both cases, the influence of the scattering contributions and thus of the plate-induced changes upon the quantum-vacuum correlations is clearly visible in the predicted full EOS signal.
A second model for the optical response of the plate is a Drude-Lorentz model permittivity defined by 2π THz, ω c = 0.04 × 2π THz and Γ = 0.056 × 2π THz. These parameters are chosen such that the material's surface-plasmon polariton resonance coincides with the frequencies that the filter function picks out from the vacuum. Consequently there is a peak in the imaginary part of the Fresnel reflection coefficient R p for p-polarized light as can be seen in Fig. 2c, which corresponds to modes evanescent at the interface between the plate and the crystal. These evanescent modes dominate the spectrum of the vacuum's contribution to the variance and lead to an increase of the EOS signal by up to a factor of 5.5 when the beam gets close to the surface, cf. Fig. 2 We thus have revealed in which way one can identify the Purcell effect, i.e. the changes of the medium-assisted quantum vacuum due to the coupling of the electromagnetic field to media, in the statistics of EOS experiments. This makes the latter a suitable tool for an in-depth study of the sculpted quantum vacuum in different optical environments.
c. Time domain perspective Thus far we have restricted the discussion to a frequency domain picture: the modes of the quantum-vacuum fluctuating at different frequencies are accessed by averaging them over the finite space-time volume of the laser pulses, cf. Eq. (1). Revealing the microscopic processes involved leads to a complementary time-domain picture. As we discuss in more detail in the supplementary material [16], the basic mechanism leading to the EOS signal in Eq. (1) are two successive, nonlinear processes which are correlated via the quantum vacuum: in the first process a (virtual) photon is generated into the vacuum mode at position r [∝Ê † vac,x (r)] which is annihilated in a second process at position r [∝Ê vac,x (r )]. One example of two such correlated processes is displayed in Fig. 3 (a) and is given by spontaneous parametric down-conversion and subsequent sum-frequency generation. The two photons which arise from this process are hence correlated and this correlation, which is proportional to the vacuum correlation function between the points r and r , is measured in EOS experiments of the quantum vacuum. On a microscopic level the observed process is hence an exchange of a virtual photon between two points within the nonlinear medium. The emerging time-domain picture of these experiments goes beyond their previous interpretation as a means to measure static, pre-existing vacuum-noise which can be squeezed or shaped by e.g. a cavity. Rather, the signal is interpreted as arising from the propagation of a virtual photon, which thus should experience retardation effects. To illustrate this, we study the signal as a function of the time delay between the pulses δt using different values for the beam separations δy, see Fig. 3 (b). In Fig. 2 (b) one sees that when δy becomes greater than the beam waist w we find maximal correlations for space-time regions which are shifted also in time by δyc n as dictated by the finite velocity of the virtual photons. This is a clear signature of retardation effects in the quantumvacuum correlations which also underlie all vacuum induced phenomena such as the Purcell effect and e.g. leads to the causal behaviour of Van der Waals forces between atoms mediated by such virtual photons [19,20].

d. Observing the Purcell effect in time domain
In the virtual-photon picture of the quantum vacuum outlined in the last section, the changes of the quantum vacuum due to the presence of additional surfaces can be understood as follows: the virtual photon can not only propagate directly from r to r , but also propagate along a path which includes reflections from boundaries, see right hand side of Fig. 3 (a). The former (latter) is described by the bulk (scattering) Green's tensor G (0) (G (1) ). Again, the dynamics of this process are not arbitrarily fast and, thus, in order to resolve it, the laser pulses which accesses the quantum vacuum must scan the quantum correlations for a time interval long enough such that (multiple) reflections can occur. This leads to a type of timefrequency uncertainty relation. To show this dynamical aspect of the medium-induced quantum vacuum, we consider the effect of reflections from the front and back surfaces of the electro-optic crystal: these already constitute a planar cavity structure, see inset in the bottom right of Fig. 1 (a). The electro-optic signal can be computed via Eq. (2) as before with O (0) = O (1) = 1 and p (j) given in the supplementary material [16]. By definition, the bulk contribution remains unaffected by the surfaces, whereas for the scattering contribution the phase-matching condition in the response function selects only certain propagation paths. In order to be phase-matched, the virtual photon before and after reflection must propagate in the same direction as one of the two pulses such that it can be 'picked up' by that laser pulse. This means, for laser pulses propagating in the z-direction, only virtual photons which have been reflected at least twice can be phase-matched, see illustrations in Fig. 3 (c).
The resulting signal as a function of δt is shown in the first row of Fig. 3 (c). One finds that the scattering contribution is very small for δt = 0, since the

Time domain Frequency domain
Fourier transfrom pulses have already left the crystal (cavity) before the virtual photons have been reflected at least twice (this takes 2L/c n = 2.1 ps, pulse duration ∆t = 80 fs). This is again a clear signature of retardation effects in the vacuum correlation function. However, if the delay between the pulses δt is a multiple of 2L/c n the second pulse arrives at the crystal precisely when the virtual photon generated by the first pulse has been reflected an even number of times. Thus the peaks for even values of δtc n /L in Fig. 3 (c) display its propagation in time domain. Fourier transforming the signal with respect to δt one can obtain the signal's spectrum g(Ω) [g(δt = 0) = dΩg(Ω)] [15], which shows the expected mode structure [see Fig. 1 (b)] with resonances at multiples of c n L/π, compare right hand side of Fig. 3 (c). However, in order to obtain g(Ω) experimentally one has to perform a measurement of g(δt) for a wide range of δt, i.e. one needs to resolve the correlations arising from different numbers of reflections individually in time. In the case of a single measurement at δt = 0, the positive and negative contributions in g (1) (Ω) mostly cancel each other such that the cavity field remains unseen when it is averaged over a single space-time region. This can be seen as a time-frequency uncertainty relation.

Scattering Scattering coated Bulk
To further improve the visibility of the cavity contribution we consider two beams propagating in opposite directions, i.e. the first into z, the second into −z. As a result, the bulk contribution is phase-mismatched and reduced considerably [16]. The scattering contribution is dominated by that stemming from virtual photons which are reflected an odd number of times: they are generated by the first laser pulse propagating into z direction and are 'picked up' while propagating into −z direction by the second laser pulse. Since these virtual photons only have to propagate over a shorter distance compared to the configuration in which both pulses propagate into the same direction, this improves the visibility of the cavity modes, cf. second row of Fig. 3 (c).
Lastly, we consider the case where the durations of the laser pulses are much longer then the time a photon needs to propagate back and forth between the front and back side of the crystal, i.e. ∆t L/c n . In this case, multiple reflection can in principle occur but the spectral autocorrelation function restricts the accessed quantum vacuum to frequencies much shorter than the lowest resonant mode, i.e. f (Ω) = 0 for all Ω > c n π/L. Hence, all the accessed modes interfere destructively so that in this case the vacuum field is completely suppressed such that g(0, 0) ≈ 0, compare lowest row of Fig. 3 (c).
e. Conclusion In this work we have shown how electro-optic sampling can be exploited to measure environment-induced changes of the electromagnetic vacuum fluctuations (i.e. the Purcell effect) in the frequency and time domains. Interpreting vacuum correlations as arising from the exchange of virtual photons leads to a time domain picture of how vacuum correlations evolve which reveals retardation effects and a time-frequency uncertainty relation. This dynamical picture might allow one to reveal other space-time properties of the quantum vacuum in the future [21][22][23][24]. It also forges a missing link between EOS and well established quantum-vacuum effects which can also be seen as arising from the exchange of virtual photons such as e.g. van der Waals forces, Casimir forces or resonant energy transfer and the influence of different optical environments upon them. Future extensions might include studying the vacuum field in more complex geometries such as plasmonic cavities [25].

A. Bulk Green's tensor
The bulk Green's tensor solves Eq. (3) for an isotropic permittivity, i.e. (r, Ω) = (Ω). In a (2+1)-dimensional Weyl decomposition relative to a plane whose normal direction is denoted by r ⊥ it reads: [1] Here, we have defined the wave vector k = (Ω)Ω/c, which can be split into perpendicular (k ⊥ ) and parallel

B. Scattering Green's tensor
The scattering part of the Green's tensor depends on the geometry under consideration. Here, we are firstly interested in the geometry of a plate attached to the nonlinear crystal. Note, that we here neglect reflections from the surfaces of the crystal. We assume that the plate is thick enough such that it can be approximated by a semi-infinite half space with refractive index n (Ω) whose interface is placed at r ⊥ < −d, with r ⊥ being the plate's normal vector. For the refractive index used to characterize the attached plate, see Sec. II of this Supplementary Matrial. For r ⊥ > −d we thus find the nonlinear crystal with refractive index n(Ω). In part b of the main text, we consider precisely this configuration with r ⊥ = x or r ⊥ = y. The scattering part of the Green's tensor inside the nonlinear crystal (r ⊥ , r ⊥ > −d) reads [1] As in the main text, R σ are the reflection coefficients at the nonlinear crystal/plate interface. In case of perfectly reflecting plates they are given by R p = 1, R s = −1 and in case of a plate with finite permittivity they are the usual Fresnel reflection coefficients which are given by Figure 4. Real and imaginary parts of the refractive index in the THz as measured in Ref. [3]. The dashed red lines show an interpolation line of the measured data whereas the black solid line indicates the refractive index used in the simulations. The latter differs from the measured data in frequency ranges where the measured data is not reliable, i.e. for small frequencies and in region where Im[n(Ω)] < 0.
The other geometry considered in this work is that where reflections at the surfaces of the crystal are taken into account. This can be done by including the scattering part of the Green's tensor for a configuration where the refractive index for −L/2 < z < L/2 is given by the one of the nonlinear crystal and otherwise it is defined to be the vacuum refractive index, i.e. n (ω) = 1. For this geometry one finds [1] G (1) for r ⊥ , r ⊥ ∈ [−L/2, L/2]. Note that here we have chosen r ⊥ = z and k ⊥ = k z = k − k 2 x − k 2 y and the terms D σ = 1 − R σ R σ e 2ik ⊥ L in the denominators account for multiple reflections. In the case that we assume that the crystal's surfaces are coated in order to increase the reflectivity in the THz range, we assume R p = 0.95 and R s = −0.95 and otherwise the reflection coefficients are given by Eqs. (8) and (9) with n (ω) = 1.

II. PARAMETERS
In this section we give all optical parameters of the nonlinear crystal and its optical surroundings used in the main text to simulate the signal of electro-optic sampling experiments.
The ordinary and group refractive indices of the nonlinear crystal at ω c are n = 2.85 and n g = 3.2 as measured in Ref. [2]. For the nonlinear refractive index we neglect its dispersion and use [2] with r 41 = 1.17 × 10 −21 CV −2 [2].
In the THz range we use the data for n which was measured using time domain spectroscopy in Ref. [3]. Its resulting real and imaginary part are shown in Fig. 4.

III. THEORY OF ELECTRO-OPTIC SAMPLING IN GENERAL ENVIRONMENTS
In this section we derive analytic expressions for the electro-optic signal in the different geometries considered in the main text. The general approach is always the same and build upon Eqs. (1) and (2) of the main text: One has to find the filter function as well as the Green's tensor for the different configurations, insert these expressions into Eqs. (1) and (2) of the main text and solve as many of the resulting integrals as possible.

A. General equations
As shown in Refs. [4][5][6], the electro-opitc signal is given by : . (13) Here, η(ω) is the detector efficiency which is assumed to be one over the frequency range of the laser pulses, S is the symmetrization operator which is defined by its action onto a product of two operatorsÔ 1,2 : SÔ 1Ô2 = (Ô 1Ô2 +Ô 2Ô1 )/2, : · · · : denotes normal ordering, . . . means that the ground state expectation value is taken, andÊ 1,2 is the electric field at detector 1 and 2, respectively.Ê 1,2 can be perturbatively expanded in terms of the nonlinear susceptibility χ (2) and consists of the vacuum fieldÊ vac , the classical probe field 1 (E p,1 ) or 2 (E p,2 ) and the field consisting of the nonlinear mixing of the vacuum field with eitherÊ p,1 orÊ p,2 , respectively. This means, the electric field emerging from the crystal consists of the free fields (vacuum and coherent laser pulses) as well as of their mixing via the nonlinear coupling inside the crystal. See Refs. [5,6] for details. Allowing for absorption inside the crystal as well as general optical environments it was shown [5,6] that g is given by Eqs. (1) and (2) of the main text which are repeated here for clarity Ê vac (r, Ω)Ê † vac (r , Ω ) =h Ω 2 c 2 ε 0 π δ(Ω − Ω )ImG(r, r , Ω), =⇒ g(δt, δr ) =h dΩF (r, r , Ω, Ω , δr , δt)ImG(r, r , Ω).
The filter function is given by with To further simplify this expression, we insert the Gaussian laser pulse defined in the main text for E p,i into Eq. (18), i.e. we use E p,i (r, t) = dωE p,i (r, ω)e iωt with Note, that here r = (x, y, 0) T , and k = n(ω)ω/c. Also we insert G xx (ω) = G (0) xx (ω) with G (0) given in Eq. (4) into Eqs. (17) and (18). This is justified since throughout this work we assume that the optical environments do not affect the near-infrared laser pulses apart from obscuring them as discussed in the main text. Neglecting absorption in the frequency range of the laser pulses and applying the laser-paraxial approximation introduced in Refs. [5,6] one finds F (r, r , Ω) = 2 2χ (2) µ 0 cN ω p f (Ω) w 2 n(ω c ) 2 e −in g Ω(z−z )/c × e −2 r 2 +(r +δr ) 2 /w 2 e iΩδt + e −2 r 2 +(r +δr ) 2 /w 2 e −iΩδt . (22) Here, n g is the group refractive index at the central frequency ω c of the laser pulse and ω p and f (Ω) are the averaged detected frequency and the spectral autocorrelation function, respectively, given by After having obtained the filter function, as well as the Green's tensor for the different geometries (see Sec. I) we have all the ingredients needed in order to calculate the electro-optic sampling signals in the different configurations using Eqs. (14) and (15) which is done in the next two sections.

B. Plate with different orientations parallel to the propagation direction of the laser pulses
Here, we consider a plate which is parallel to the propagation direction of the laser pulses with two different orientations, i.e. in the x < −d and y < −d half spaces. We neglect absorption inside the nonlinear crystal by assuming (Ω) to be real. One can split the signal g in Eq. (16) into contributions stemming from the bulk and the scattering Green's tensors. Note, that although the bulk Green's tensor inside the nonlinear crystal is unaffected by the additional optical plates considered here, its contribution to the EOS signal still changes due to the fact that the laser pulses might get obscured by the plate whenever the pulse/plate distance becomes comparable to the beam waist. This is included in Eq. (16) by the restriction of the spatial integrals to the volume of the crystal −∞ dx for the two different orientations of the plate considered here, respectively. Inserting the bulk or scattering Green's tensors in Eqs. (4) and (7) with r ⊥ chosen perpendicular to the applied surface and the filter function in Eq. (22) into Eq. (16) one finds after a tedious but straight forward calculation that the signal can always be brought into the form This expression is identical to Eq. (3) of the main text. Here, √ C = 2χ (2) Lω p N (d, δr )/n 0 c with N (d, δr ) being the total number of detected photons given by N is the total number of detected photons without the plate obscuring parts of the laser pulses, and x 0 dt e −t 2 is the error function. E 2 vac (Ω) is the coincidence limit of the bulk two-point correlation function of the electric field operator neglecting absorption effects and it reads Also, note that the wave vector q has been split into a component which is parallel to the surface of the attached plate (q ) and one which is perpendicular to it (q ⊥ = q 2 − q 2 ), with q = n(Ω)Ω/c. Finally, as also stated in the main text, the response function R(q) is given by ≈ e −(q 2 x +q 2 Here, we restricted the response in the second row to the phase-matched contribution only. Note, that this approximation is only a good approximation for propagating modes with q < Re[q]. It thus is used to calculate the EOS signal in Fig. 2 (a) but not in Fig. 2 (b), since in the latter case the main contribution to the scattering part of the signal stems from evanescent modes with q > Re[q]. The obscuring factors O (j) (q, δr , d) in case of a plate in x < −d are given by The expressions for O (j) (q, δr ) in case of a plate in y < −d can be obtained from Eqs. (31) and (32) by replacing q x ↔ q y and δx ↔ δy. Finally, the propagation factor p (j) differs for the bulk and scattering contribution as well as for the different geometries (plates in x < −d or y < −d) and reads for the different configurations:

C. Cavity setup: Including reflections from the crystal's front and back surfaces
Similar to the last section, we can include the effect from reflections at the front and back surfaces of the crystal by using the appropriate scattering Green's tensor in addition to the bulk Green's tensor in Eq. (16). We first again neglect absorption by assuming that n(Ω) is real. The bulk contribution is not further restricted in this configuration and thus agrees with the signal considered in previous studies in which reflection effects have been neglected [5,6]. It can be derived by using the bulk Green's tensor (Eq. (4)) in Eq. (16) and the resulting expression is again given by Eq. (25) with q ⊥ = q z , O (j) (q, δr ) = 1, N (d, δr ) = N and p (0) (q, δr ) = q The scattering contribution is obtained by inserting the Green's tensor of the cavity in Eq. (10) into Eq. (16). After some calculation very similar to the ones in the last section one obtains In the case that the crystal length is very short, e.g. set to L = 1 µm (third row of Fig. 3c in the main plot), one needs to take all these terms into account. However, in case L = 0.1 mm (as used throughout the rest of the paper) it is enough to only include the phase-matched contribution, i.e. only the term proportional to Note, that this propagation factor only includes contributions which are at least quadratic in the reflection coefficients, meaning that the virtual photon is reflected at least twice.
In the case where the two pulses are assumed to propagate in opposite direction one has to replace e ikz by e −ikz in Eq. (20). The calculation is very similar to the one before. For the bulk contribution, the result is again given by Eq. (25) with q ⊥ = q z , O(q, δr , d) = 1, N (d, δr ) = N and p (0) = q 1−q 2 x /q 2 q ⊥ e iq·δr as for the case were the two pulses propagate into the same direction but the response function now reads We see that there is no phase-matched contribution as expected since the two pulses propagate into opposite directions. The scattering contribution is also given by Eq. (25) with q ⊥ = q z , O(q, δr , d) = 1, N (d, δr ) = N , the response function in Eq. (30) and Note, that this propagation factor now includes contribution which are proportional to R σ , meaning that the virtual photon can this time also be reflected only once and still be phase-matched as expected.
Lastly, we calculate the signal accounting not only for reflections at the front and back surfaces of the crystal but also for absorption in the THz region. To do so we allow n(Ω) to have a nonvanishing imaginary part (compare right hand side of Fig. 4). For the bulk contribution with both pulses propagating into the same direction one finds as previously also obtained in Ref. [6] Here, only the phase-matched contribution is included.
For the scattering contribution we repeat the same calculation leading to Eq. (33) but allowing for a imaginary part of n(Ω). This leads to In the case where the two pulses propagate into opposite direction and absorption effects are considered, the bulk contribution to the signal is given by

IV. MICROSCOPIC INTERPRETATION
We discuss the microscopic processes leading to the contribution to the signal's variance found in Eq. (14). We first realize that, on the one hand, terms of the structureÊ x (r , Ω)Ê y (r , ω − Ω)G xx (r, r , ω) correspond to a processes generally referred to as sum-frequency generation [7]. Here an atom absorbs two photons with polarization x, y and frequency Ω, ω − Ω respectively and the excited atom subsequently emits an x polarized photon with frequency ω. On the other hand, terms of the structureÊ † x (r , Ω)Ê y (r , ω + Ω)G xx (r, r , ω) describe the process of parametric down-conversion [7], where an atom is excited by one photon of frequency ω + Ω which subsequently deexcites by a two-photon emission process where both photons are x-polarised and have frequencies Ω and ω, respectively. Now let us turn our attention to Eq. (13), i.e.
Note, that in the lowest order in χ (2) ,Ê x (r, ω) is given by [5] E x (r, ω) = χ (2) such that it is apparent from the previous discussion that this field stems from the nonlinear process of either parametric down-conversion or sum-frequency generation depending on the sign of the frequencies. When taking the ground state expectation value only terms proportional toÊ vac,x (Ω)Ê † vac,x (Ω) contribute where Ω is now positive. This means that on a microscopic level only those nonlinear processes contribute where a (virtual) photon at frequency Ω is generated in a first nonlinear process and subsequently absorbed by a second one. To illustrate this, we insert Eq. (42) into Eq. (13) and consider one of the resulting terms given bŷ E x (r , Ω)Ê y (r , Ω − ω )G xx (r, r , ω )Ê † x (r , Ω)Ê y (r , Ω + ω)G xx (r, r , ω).
It is apparent that this describes the processes of spontaneous parametric down-conversion at position r and sum-frequency generation of one of the generated photon together with a photon of the laser pulse at position r .