Dynamical modeling of pulsed twophoton interference
Abstract
Singlephoton sources are at the heart of quantumoptical networks, with their uniquely quantum emission and phenomenon of twophoton interference allowing for the generation and transfer of nonclassical states. Although a few analytical methods have been briefly investigated for describing pulsed singlephoton sources, these methods apply only to either perfectly ideal or at least extremely idealized sources. Here, we present the first complete picture of pulsed singlephoton sources by elaborating how to numerically and fully characterize nonideal singlephoton sources operating in a pulsed regime. In order to achieve this result, we make the connection between quantum Monte–Carlo simulations, experimental characterizations, and an extended form of the quantum regression theorem. We elaborate on how an ideal pulsed singlephoton source is connected to its photocount distribution and its measured degree of second and firstorder optical coherence. By doing so, we provide a description of the relationship between instantaneous source correlations and the typical experimental interferometers (Hanbury–Brown and Twiss, Hong–Ou–Mandel, and Mach–Zehnder) used to characterize such sources. Then, we use these techniques to explore several prototypical quantum systems and their nonideal behaviors. As an example numerical result, we show that for the most popular singlephoton source—a resonantly excited twolevel system—its error probability is directly related to its excitation pulse length. We believe that the intuition gained from these representative systems and characters can be used to interpret future results with more complicated source Hamiltonians and behaviors. Finally, we have thoroughly documented our simulation methods with contributions to the Quantum Optics Toolbox in Python (QuTiP) in order to make our work easily accessible to other scientists and engineers.
Dynamical modeling of pulsed twophoton interference
 1 Introduction
 2 Prototypical singlephoton sources
 3 Modeling of source dynamics
 4 Singlephoton source photocount distribution
 5 Interferometers for observing pulsewise twophoton interference
 6 Conclusions
 7 Acknowledgements
 Appendix: Normalization of the Hong–Ou–Mandel and Mach–Zehnder correlations
1 Introduction
The development of the quantum singlephoton source has ushered in the field of optical quantum information technology [1]. Such sources, serving as generators of flying photonic qubits [2], lie at the heart of nearly every quantumoptical technology [3], including photonic logic gates [4], quantum networking [5], and highly nonclassical NOON state generation [6]. Critical to the usefulness of a singlephoton source is knowledge of its temporal profile [7, 8, 9], which has been accomplished via two means: heralding of the emission through a projective measurement on a continuous entangled pair source [10] or ondemand generation through pulsed excitation of a system that releases only one photon per pulse [11]. In this paper, we thoroughly discuss the properties of ondemand pulsed singlephoton sources.
Before we begin with our detailed analysis, we provide a brief summary of the characteristics of ideal singlephoton sources. Generally, a singlephoton source creates a pulse containing at most one photon. More specifically, ideal and ondemand singlephoton sources are characterized by their:
Given that singlephoton pulses exist as wavepackets, they may also display interference, and hence their spatiotemporal profiles are critical in this respect: From this realization came one of the great discoveries of the 20th century, the concept of an identical photon [8]. When a pair of photonic wavepackets share all three of the aforementioned characteristics, they are said to be identical or indistinguishable. Meanwhile, the identical nature of photons has provided quantum optics with a rich class of experiments to explore. Perhaps the most famous is twophoton interference (seen in the Hong–Ou–Mandel interferometer), whereby two identical photons incident on a beamsplitter always exit a chosen output port together. First demonstrated in 1987 with photons produced via parametric downconversion [13], twophoton interference has been observed with photons generated via trapped ions [14], artificial atoms [11], and even circuit QED systems [15].
While writing down the wavefunction for such photon pulses is now fairly wellestablished (at least from a quantumfield theoretic point of view [16]), simulating nonideal singlephoton sources has proven extremely challenging. Only a few papers ever attempted to tackle this type of problem, and even so they investigated highly restrictive cases [17, 18]. Here, we present a more general simulation technique that allows for the quantification of singlephoton emission from systems characterized by Hamiltonians with nearly arbitrary timedependence. This work also includes techniques for directly calculating the effective pulsewise twophoton interference in both a Hong–Ou–Mandel interferometer [13, 14] and the often misunderstood unbalanced Mach–Zehnder interferometer [11, 19]. Importantly, we have worked to make this technique easily accessible to experimentalists wishing to model their ondemand sources by contributing to the underlying code of the opensource package Quantum Toolbox in Python (QuTiP) [20]. In addition, we have provided detailed examples to the QuTiP repository demonstrating many of the simulations discussed in this paper.
2 Prototypical singlephoton sources
Although the technique we present is fairly general, we have chosen to discuss only simplistic model systems that can easily be understood and run by a broad range of quantum scientists. Similarly, while we study our systems only under drive by Gaussian pulses, our technique is trivially extensible to more complicated pulse shapes and system Hamiltonians. (For instance, see our work on singlephoton emission from Jaynes–Cummings type systems [19, 21, 22].) Instead in this paper, our three choices of model systems (figure 1) will capture various behaviors of singlephoton sources with regard to their coherences, as a function of driving pulse length. These models have been chosen as a somewhat representative set of possible ondemand singlephoton sources.

Coherently excited twolevel system: This source was chosen to illustrate performance degradation primarily in its secondorder coherence with increasing pulse length relative to its spontaneous emission rate. This system is described by the following Hamiltonian [23]
(1) where is the transition frequency, the laser frequency, the pulse driving strength, and the atomic lowering operator. We will consider the dynamics of this system with only one collapse operator, corresponding to spontaneous emission at a rate , i.e. [23]. Example experimental systems of this type can be found in coherently excited quantum dots [24, 25], trapped ions [14], and circuit QED platforms [15].

Incoherently excited threelevel ladder system: This source was chosen to model performance degradation only in its firstorder coherence with increasing pulse length relative to its spontaneous emission rate. First, we assume that the system is initialized into the state each emission cycle. From there, we will simulate only dissipative evolution through a cascade driven by two collapse operators of rates and , i.e. . Example experimental systems of this type can be found in electrically injected or incoherently excited quantum dots [26] and the polaritonphonon cascades of solidstate Jaynes–Cummingslike systems [22].

Coherently excited threelevel lambda system: This source was chosen to illustrate lack of performance degradation in its ideal limit. Specifically, it always has the perfect values of its first and secondorder coherences for a singlephoton source. Therefore, such a source can generate singlephotons with nearly arbitrary wavefunctions. This system is described by the following Hamiltonian
(2) where is the transition frequency, the laser frequency, the pulse driving strength, and an atomic lowering operator. We will consider the dynamics of this system after having been initialized each cycle to the state and under the influence of two collapse operators, corresponding to spontaneous emission from the excited state towards the system’s two ground states. These give rise to dynamics governed by two collapse operators with rates and , i.e. . To study the ideal case of arbitrary singlephoton generation, we will set . Example experimental systems of this type can be found in trapped ions [5, 27].
For each of these systems, its firstorder coherence can be decreased by the effect of any pure dephasing terms as well. These terms are modeled by the inclusion of collapse operators with the form , which may potentially be timedependent.
Finally, we note that it is trivial to model coupling to different spatial modes by decomposing a given collapse operator into independent collapse operators. For example, consider the replacement where : Then the emission occurs at the same total rate but into two separate spatial modes (with potentially different polarizations). Because this phenomenon is well understood [7], we will implicitly assume that all of our sources emit into the same spatial mode for this paper.
3 Modeling of source dynamics
The dynamics of the above systems are governed by the quantumoptical master equation [23, 28]
(3) 
where is the Liouvillian superoperator that characterizes the timedependent system evolution. More specifically, in the Markovian limit of systemreservoir interactions
(4) 
where are the collapse operators referenced above ( are the operators through which the system couples to the environmental modes). As discussed in the introduction, singlephoton sources are characterized by potentially nontrivial first and secondorder optical coherences. Therefore, we will be interested to use the system dynamics to calculate correlations of the form
(5) 
where , , and are each some combination of , , and the identity matrix. Although this calculation is often performed with the quantum regression theorem [29], its formal statement excludes timedependent Liouvillians in all references to the authors’ knowledge. Yet, when studying a system under pulsed excitation, a timedependent Liouvillian invariably arises due to the dynamical driving. Therefore, we have extended the quantum regression theorem to timedependent Liouvillians below.
When the quantum regression theorem is applied for timedependent Liouvillians, it inherits the approximations from the quantumoptical master equation and yields the following result
(6) 
where is governed by the evolution equation
(7) 
and is subject to the initial condition
(8) 
Although this equation only can evolve forward in time, the correlators discussed in this paper either give rise to physical measurements and inherently require [30] or possess a conjugate timereversal symmetry about [23]. We note that the authors have added this algorithm as a routine to QuTiP, so a simple call to the correlators automatically calculates twotime correlations with timedependent Liouvillians.
This evolution equation has a particularly nice interpretation when the probability for three successive photodetection events is negligible and
(9) 
which then represents the probability density of detecting an excitation at time followed by another excitation at time . The quantum regression theorem clearly captures this through: first a reduction of the density of matrix by one excitation
(10) 
and by then computing the probability density of a second reduction
(11) 
Finally, we comment on the connection between the system correlators and the free radiationfield correlators. Because singlephoton detectors make measurements on the free radiation field, we physically measure correlations having to do with the fieldflux operators, e.g. . Fortunately, Gardiner and Zoller’s input–output theory [29] provides a direct connection between the internal system operators and external radiation mode operators, such that measuring flux correlations is equivalent to measuring the internal system correlators. More specifically, if a system operator is coupled to the external radiation modes represented by and we wish to calculate a physical correlator involving , then we may simply make the replacement .
4 Singlephoton source photocount distribution
As discussed in the introduction, an ideal pulsed singlephoton source would contain only a singlephoton per pulse and hence only a single quanta of energy in its wavepacket—in this section we will fully explore this concept. Consider photon flux incident on a photon counter with finite timing resolution. Mathematically, the probability that the flux results in a detection event between and is
(12) 
where is the total detection and collection efficiency and is the instantaneous photon flux operator (analogous to classical intensity) [31]. Notably, the instantaneous detection probability must be negligible such that multiple photoionizations in the detector cannot occur. Integrating over the photon pulse duration, , there exists a classical probability distribution that governs the number of expected photocounts. A singlephoton source requires
(13) 
where is a classical random variable that represents the number of photocounts.
Thus, it is sufficient to characterize a singlephoton source by a measure which roughly corresponds to the probability of two events occurring within : the secondorder factorial moment of its photocount distribution [28]
(14) 
(Because the photocount distribution is over the classical random variable , we use the notation to denote the classical expectation value.) While the first detection still depends on , the second detection depends on because the first detection subtracts a quantum of energy from the distribution [31]. Therefore, for an ideal pulsed singlephoton source .
Unfortunately, the detection and collection efficiency can be quite difficult to experimentally measure, which makes this correlation challenging to directly estimate. Instead we will focus on determining its normalized version
(15) 
This normalized secondorder factorial moment of is referred to as the measured degree of secondorder coherence at zero time delay and doesn’t depend on the efficiency [31]. This quantity is useful because any value is disallowed by classical physics [29] and because approximates the error probability of the source to produce two detection events relative to the number of single detection events.
Because the photocount distribution can only be estimated from the outcome of many pulsewise experiments, a typical experimental cycle will involve periodic generation of the photon wavepacket and its subsequent detection every seconds. The number of photodetections from a given pulse at the time bin can then be represented by the classical random variable . We can then extend the pulsewise definition of to nonzero time delays, i.e.
(16) 
which we refer to as the unnormalized secondorder intensity correlation at time delay . We note that its time origin is irrelevant due to the quasistationary nature of the pulses. If the probability for three photodetections over a given pulse is negligible, then this correlation represents the probability that two photons are detected between pulses separated by the time difference . From this definition, we can arrive at another definition of which will turn out to be experimentally most useful:
(17) 
provided that we have correctly chosen longer than the correlation time of any system operators (where ). Experimental realities may sometimes inhibit the realization of this criterion and one can then simply take
(18) 
While this definition may seem a trivial extension, this form is most useful (due to the limitations of legacy measurement instrumentation) in discussing experimental setups that estimate . Thus, we have discussed how experimental singlephoton sources can be readily characterized even in the face of unknown detection and collection efficiencies.
4.1 Hanbury–Brown and Twiss interferometer
While one can easily estimate from the detection record of ideal photon counters, this is not possible for most experimental detectors and singlephoton sources. Experimental photodetectors have a socalled “dead time”, for which the detector cannot register a second count following the first, that usually is much longer than the temporal length of the photon pulse. With a maximum of one registered count per pulse, such a detector could never distinguish between single and multiphoton sources. Fortunately, it turns out that under the right approximations, an experimental setup known as the Hanbury–Brown and Twiss interferometer (figure 2) is still capable of precisely estimating by using two photon counters, and hence is still capable of characterizing realistic singlephoton sources [32].
Because our focus is on pulsed singlephoton sources, we will only discuss the Hanbury–Brown and Twiss (HBT) interferometer operated in a pulsed mode of operation, where we repeatedly excite our source every seconds [32]. In the HBT setup, every pulse is pathentangled between two channels by a beamsplitter and each channel is fed into a singlephoton detector. The periodic photon absorption events registered by the two detectors can be represented as the classical random variables and where ; They may take on unity values at times to represent photon detections. Here, the square brackets will indicate the discrete nature of the random variables and their intensity correlations, due to both their periodic pulsed nature and any timing uncertainty in their detection (either explicit by detector jitter or implicit through purposeful erasure of timing information). Although an estimate for could be built up through many periodic trials and used to directly compute , traditionally has been extracted from a histogram of the timecorrelated detection records, i.e.
(19) 
where , is the integration time, and is the time difference between detections [33]. Each histogrammed timebin is an independent and binomiallydistributed random variable, whose standard deviation estimator is given by for most experiments ( is again the combined detection and collection efficiency).
After careful consideration of the detector nonidealities such as dead time and dark counts, as well as of the expected statistics of a singlephoton source [33], one can arrive at the approximation
(20) 
Here, the factor of accounts for the action of the beamsplitter to halve the signal at each detector. In making the connection back to the secondorder factorial moment at zero delay that was discussed in the previous section, consider
(21) 
While it may be surprising that this expression is equivalent to , i.e. , the key insight is that the random variables and are not independent because a detection event by either detector pulls one quantum of energy from the total pathentangled field.
As previously mentioned, however, it is difficult to experimentally estimate due to unknown setup efficiencies. Fortunately, since all correlations in are lost at long times such that
(22) 
then can serve as an intensity reference. Thus, we can obtain an estimate for the measured degree of secondorder coherence from the following ratio
(23) 
Importantly, several approximations were required to arrive at this result [33]:

The net detection probability per pulse must be very small relative to the dead time.

The net detection probability per pulse must be very small relative to the repetition time .

The probability of manyphoton detection must be moderately low; Analogously, the higherorder factorial moments of the photocount distribution must be of order unity, i.e. for .

The detectors’ dark count rates are minimal, such that
(24)
As an example of such a histogram, we reproduce data from reference [21] in figure 3. In this experiment, a singlephoton source was run in pulsed mode and timecorrelated using an HBT setup. Because all pulsewise correlations have decayed already after one repetition time in this experiment,
(25) 
and therefore
(26) 
which yields .
Finally, we note that many experimental HBT interferometers do not exactly histogram . Rather, they only approximate by electronically timecorrelating detections onthefly. This is done by taking the first detector as the signal to start timing and the second detector as the signal to stop timing: Each startstop sequence generates a count in the timebin of the timer value. This way, the histogram is built up in realtime as new measurement correlations are recorded. The downside to this method is that each successively longer timebin requires more failed detections to register a count, where the actual histogram constructed is
(27) 
Here the product term means that at long timecorrelations, approximately and hence it decays to zero. Therefore, this electronic method of timecorrelating to estimate may not always work well if the shows correlated behavior for large . Now that we have a good understanding of experimentally how to characterize the photocount distribution of a singlephoton source, we will discuss the theoretical connection to the instantaneous correlations of the photon wavepacket’s fields.
4.2 Connection to correlations of instantaneous fields
From a numerical modeling perspective, estimating for a given system Liouvillian is fairly straightforward using Monte–Carlo wavefunction techniques [34]. Here, evolution of the system wavefunction is modeled, conditioned on the detection events of an ideal singlephoton detector with infinite bandwidth. Such a detector may absorb one quantum of energy from the emitted field at a time, in proportion to the photon flux at a given instant, and can distinguish detection events with infinite timing resolution. Through this process, detection records are generated that build up an estimate for and hence . However, the Monte–Carlo method often requires an extremely large number of simulated evolutions (trajectories) to obtain an acceptable approximation of the measured degree of secondorder coherence.
Fortunately, for systems with reasonablydimensioned Hamiltonians, there exists a faster and more intellectually satisfying way of computing —as used for pulsed nonclassical light sources, this algorithm was first implemented numerically in our previous work [21]. This method directly relies on using the quantum regression theorem, as outlined in section 3, to compute correlations between the instantaneous system operators previously discussed. These correlations in turn are related to the instantaneous correlations of the continuousmode freefield operator that describes the emitted photonfield flux. In this section, we now fully elaborate on this explicit connection of the measured degree of coherence with its associated instantaneous field correlations.
We begin with a more directly quantum mechanical model of our ideal singlephoton detector that formalizes and combines the stories from references [16, 31] and [32]. In analogy to the classical integrated mean intensity, the quantum mechanical operator
(28) 
represents the total photon number arriving at an ideal detector over the time interval (where is again over the duration of the photon pulse). As such, its quantum mechanical expectation value yields . Comparing our semiclassical definition of to the quantum mechanical operator we have
(29) 
Here, denotes the quantum mechanical expectation value of the normallyordered second moment of the photon number operator . Writing out this moment explicitly
(30) 
where the operators indicate the timeordering required of a physical measurement [35] (operators with higher time indices towards the center of the expression).
We can also compute the quantum mechanical version of the normalized secondorder moment , with
(31) 
and its explicit form
(32) 
Written another way by defining
(33) 
then
(34) 
Therefore, the measured actually represents the sum of all fieldflux correlations over the detection time. This expression agrees with intuition when higherorder correlations are negligible, as it then represents the probability to detect two photons at every possible pair of times, normalized to the total photon number squared. As previously discussed, input–output theory provides a direct connection between the internal system operators and external mode operators such that measuring zero delay flux correlations is equivalent to calculating the mode correlators. Thus, we may simply calculate
(35) 
where could be an atomic lowering operator in the case of radiation from a fewlevel system or a cavity mode operator in the case of radiation from a cavity. Hence, we have finished outlining our novel method of numerical simulation that will be used for the rest of this paper in modeling the dynamics of pulsed singlephoton sources.
Finally, we note that the above expressions will later also be expanded to directly consider correlations between two field operators labeled, for instance, and with
(36) 
(37) 
and
(38) 
4.3 Simulated secondorder optical coherences of singlephoton sources
Now that we have fully elaborated on the theory behind the HBT setup, we can put all of the aforementioned pieces together and begin to simulate the three different singlephoton sources discussed in section 2. These systems were chosen because even singlephoton sources with more complicated Hamiltonians and characters may potentially be mapped onto the source behaviors we will present in the following simulation sections.
The first feature we will consider is the energetic shape of the wavepackets emitted from the three systems. In fact for a given pulse, the three systems generate nearly identically shaped wavepackets when excited to emit, on average, one photon. Here, the shape of the wavepacket refers to the profile of the average energy density at a given point in space and time. Importantly, for long pulses, we opted against the standard definition [36] of using constant pulse area to define our pulse lengths. Instead, we chose to determine our pulse lengths such that each system emits, on average, one photon per pulse. For short pulses the two definitions agree, but for long pulses this difference allows us to focus exclusively on photon statistics for comparable average photon flux.
Using these definitions, we describe the disparate situations when the three sources are excited by temporally short or long pulses. When excited by a short Gaussian pulse, where short is relative to the systems’ characteristic emission times, then the resulting wavepackets have exponentially decaying energy densities with time. An example of this behavior is shown in figure 4a (red), with a pulse length of in energy. On the other hand, consider the systems when excited by a long Gaussian pulse. Then the resulting wavepackets have almost Gaussian shapes, as shown in figure 4b (blue), with a pulse length of in energy. Wavepackets of these two shapes will be considered for the rest of the paper and their shapes will always be denoted by the colors of red (exponential) or blue (Gaussian).
Given that the energy density profiles match so closely between our prototypical sources excited with the same pulses, the differences between their wavepackets must lie elsewhere. As we will see, such differences will be seen in their quantum statistics (optical coherences). In this section, we make comparisons between their secondorder optical coherences, while in the next simulation section 5.3 we will compare values based on firstorder optical coherences. The simulated values for wavepackets emitted by both threelevel sources, systems (ii) and (iii), are quite trivial: They are manifestly zero for all possible excitation pulses. This result is quite simple to understand—the act of initializing the system to a third state [state for system (ii) and state for system (iii)] means only a single quanta of energy can leave the system before it is reinitialized. The initial and final states for both sources are unconnected so that the systems cannot be reexcited to emit multiple photons [27]. These results are summarized in the legends of figures 4a and 4b.
On the other hand, because the initial and final states are the same for the twolevel singlephoton source [system (i)], reexcitation may occur. Therefore, its values exhibit a more complicated character and may be nonzero. For example, consider the effects of pulse length on for the twolevel system (shown in figure 4c). If a short excitation pulse drives photon emission, then the excitation occurs over a very short timescale as compared to the actual emission and reexcitation is very unlikely to occur. Therefore, the system possesses low secondorder coherence and acts as a singlephoton source (over the red side of the curve). In fact, this reexcitation probability is linear with pulse length at low powers and therefore we can fit
(39) 
(Surprisingly, some papers discussing experimental results on stateoftheart quantum dot sources quote values for that are lower than this limit, e.g. those in reference [37].) Because the reexcitation probability is always finite, also has a limiting value even for arbitrarily low powers. This effect can be seen in figure 4d, where as a function of pulse area is shown for a short pulse (pulse length of in energy). Therefore, even in the absence of all other nonidealities, pulsed twolevel systems will never exhibit perfect singlephoton emission. We note that our modeling agrees very well with recent experimental results regarding the emission from a neutral quantum dot [25].
Meanwhile, if a long excitation pulse drives photon emission from system (i), then reexcitation is highly likely and the system inherits the photon statistics of the laser pulse (=1). This regime can be seen in the blue side of figure 4c and is tabulated in figure 4b. From these simulations, we can conclude that there are regimes where a pulsed twolevel system does not operate as singlephoton emitter at all.
5 Interferometers for observing pulsewise twophoton interference
In the previous section, we discussed how the HBT interferometer characterizes a single photon by measuring the secondorder intensity interference of its pathentangled wavepacket. There the spatiotemporal profiles were irrelevant to the pulsewise results—now we will consider interferometers whose outputs depend on these profiles. However, we will implicitly assume that the wavepackets occupy comparable spatial modes so we can focus on their temporal characteristics. Next, we will introduce several important findings and concepts that we will elaborate on and prove in the subsequent sections.
Specifically, we will discuss two interferometers that compare two independent photon wavepackets and depend both on secondorder (intensity) and firstorder (field) interference: the Hong–Ou–Mandel interferometer [7] and the unbalanced and doublyexcited Mach–Zehnder interferometer [11]. For these interferometers in the limit of identical singlephoton inputs, the intensity correlations equal zero just like in the HBT setup. Interestingly, this means that the remaining information is encoded in the mutual firstorder optical interference between wavepackets. Therefore, this correlation could in principle be measured through a traditional Michelson interferometer and is more recognizable as interference in classical optics. However, a Michelson interferometer can only measure the ensemble average of firstorder coherence over extremely long timescales, while the Hong–Ou–Mandel and Mach–Zehnder interferometers can directly integrate the pulsewise firstorder coherence of singlephoton wavepackets.
First, we will discuss the Hong–Ou–Mandel interferometer, which measures the indistinguishability of two singlephoton inputs by comparing the likeness of their wavefunctions directly [7, 8]. In the case of perfect identical photons at the inputs, detection of a photon by the first detector projects the second photon into the same channel as the first photon—this phenomenon is known as twophoton interference. Unfortunately, producing a single source of indistinguishable photons is often quite challenging (let alone two), but one would still like to characterize the source’s instantaneous degree of firstorder coherence. Or in an alternative source architecture, an individual singlephoton source has recently been utilized to generate several coincident indistinguishable photons through timemultiplexing, for use in a complicated quantumoptical experiment known as Boson sampling [38, 39]. In this case, the indistinguishability of photons generated at different times, but by the same system, needs to be compared.
For these purposes, we will discuss the unbalanced and doublyexcited Mach–Zehnder interferometer. Using this interferometer, one can quantify how a source would behave in a Hong–Ou–Mandel experiment [19] by interfering pulses emitted at different times. Notably, unlike in the continuous excitation case, the two interferometers (Hong–Ou–Mandel and Mach–Zehnder) do not necessarily measure the same interference visibilities for comparable sources, nor do they share the same threshold for nonclassicality.
5.1 Hong–Ou–Mandel interferometer
Here, we detail interference in a Hong–Ou–Mandel (HOM) interferometer. In such an interferometer, two identically independent sources are interfered on two detectors by the action of a single beamsplitter. This setup is shown schematically in figure 5: Consider two systems periodically driven by the pulses and , respectively. Their outputs, the Heisenberg freefield operators and , are subsequently fed into the beamsplitter which mixes the two according to the unitary transformation
(40) 
The detection events are then correlated electronically to arrive at a temporal coincidence histogram
(41) 
In terms of the instantaneous correlations discussed in section 4.2, we wish to calculate which is . To perform this calculation, we first decomposed the underlying instantaneous correlations of , i.e. , into more manageable components based on and that allowed for the averaging over all quicklyvarying phase terms. These phase terms are difficult to observe experimentally since they depend on femtosecond phase locking and are, regardless, not required for the observation of twophoton interference.
It has been shown that after performing this procedure [7] one arrives at
(42) 
where represents the sources’ firstorder optical coherences and . Then we again remark that in its pulsewise form
(43) 
This formulation can easily represent any form of intersource distinguishability, possibly arising due to a frequency difference, temporal or amplitude mode mismatch, multiphoton pulse, or an excitation timing jitter. However, we need a valid intensity reference to compare this correlation against in experiment due to unknown detection and collection efficiencies (see the Appendix for a detailed discussion on this normalization choice). Similar to the longtime correlation in the HBT setup, here the only true reference is in proportion to the statistically independent crosscorrelation
(44) 
Since we’re only interested in determining the pulsewise correlations, we consider
(45) 
whose estimate from the HOM setup is
(46) 
While the analytical expression for is somewhat complicated once expanded, one can derive insight by investigating limiting cases. For instance, consider two (potentially distinguishable) singlephoton sources as the inputs to the HOM interferometer. Then,
(47) 
which could represent the interference between two singlephoton sources with different temporal delays, center frequencies, or pulse shapes. This expression can be understood as measuring the overlap of the sources’ singlephoton wavefunctions: When the sources perfectly share the aforementioned characteristics, then they have identical firstorder coherences. When the firstorder coherences match, so do the wavefunctions and thus is measured [8]. Importantly, we now have the second requirement for a good singlephoton source—a good source must be able to interfere with other sources, which requires both the wavepackets to be identical and hence in pure states [7, 8]. This is met only if both states have complete firstorder coherence. Unlike for the HBT interferometer, it is also this dependence on the firstorder coherence that makes Monte–Carlo simulation of the HOM histogram much more difficult.
Finally, in preparation for our comparison with the Mach–Zehnder interferometer, we look to make a few more simplifications. We note that if the sources are both identically independent and temporally matched, but we allow for nonideal secondorder coherence, then equation (45) simplifies significantly to
(48) 
Here, the source correlations between and have terms only coming from a single source (the indices again could trivially be switched as and are identically independent). Now the proportionality of the square of the experimental intensity reference reduces to
(49) 
Therefore, we are interested in the computation of the integrated and normalized correlation
(50) 
where
(51) 
One may immediately notice the nonclassical threshold of due to the classical limit of . Further, for a singlephoton source with no error rate, is precisely the trace purity of the singlephoton emission. But given a finite error rate, this equivalence no longer holds. To see this, consider the resonantly driven twolevel system with no dephasing: the emitted state is a pure state (with unity trace purity) even though .
Considering Eq. 50, the observable pulsewise HOM interference for identically independent sources simply comprises terms involving the total first and secondorder coherences of the source. As with the continuous form of HOM interference [7], whether this pulsewise interference is truly just twophoton interference therefore depends on the nature of the source outputs. Additionally, we note that in many real systems the firstorder coherence extracted here is wildly different than the firstorder coherence extracted from a Michelson interferometer due to the longtime averaging action of the Michelson.
5.2 Unbalanced and doublyexcited Mach–Zehnder interferometer
While the HOM interferometer compares two single photons between a pair of different sources, the goal of the Mach–Zehnder (MZ) interferometer (shown schematically in figure 6) is to ascertain how a single source would behave were it and its identical twin fed into an HOM interferometer. To this end, it was previously realized that some aspect of twophoton interference can be observed for a single source at the input to a MZ interferometer [11]. However, we note that further formalism was required beyond previous analyses, which we will provide here.
In order to observe this type of interference, the system must be doubly excited at a time interval that matches the temporal delay between the two paths of the MZ interferometer. Additionally, this interval must be long enough that the emission resulting from the two excitations is identically independent. For perfect singlephoton input the MZ output correlation is in fact twophoton interference and thus equivalent to HOM interference. However, in many cases the literature has incorrectly compared the two when the source has some finite probability of emitting nonsingle photons. Thus, we now outline the derivation of the correct measured correlation. Consider the schematic in figure 6: A single system is periodically driven with the coherent driving field and has fieldmode operator output , where is long enough so that the state of the system is reset between excitations. The first beamsplitter simply splits in two—we have dropped the vacuum operator in anticipation that our detectors only measure a normallyordered moment of their impinging fields. Next, one path is delayed by the excitation interval and the second beamsplitter then mixes with the time delayed version of itself. The time delay is absolutely critical: For time delays much less than , can interfere with itself, and therefore reproduce a correlation similar to the HOM crosscorrelation. Finally, the detection apparatus generates the histogram
(52) 
Following the methods outlined in the previous section, we similarly calculate in terms of the instantaneous correlations. Although computation of in terms of source correlations certainly appears daunting due to the doubly excited source, the independence of for times of order or larger dramatically simplifies calculation of . This expansion is further simplified by consideration of the correlations centered around time delays that are integer multiples of . It is fairly trivial to show that the majority of terms for delays larger than are zero, with only a few additional phasedependent terms arising. We note that while one might expect the correlations to be phase locked by the MZ interferometer, typical integration times quickly destroy the phase interference; Besides, this phase dependence is unwanted for twophoton interference [7].
Consider the correlations about zero delay: and its time delayed version are statistically independent so equation (48) holds but with the source correlations of instead of . Because the source is doubly excited the ratio of the correlations in terms of is altered, which can easily be seen by applying the above rules. Now,
(53) 
Again, we need a longtime reference for the square of the experimental intensity, which is proportional to
(54) 
and so we consider the normalized pulsewise correlation
(55) 
Thus, we arrive at
(56) 
whose estimate from the MZ setup is given by
(57) 
Importantly, we note that this different interference visibility as compared with equation (50) for the HOM setup holds only for the pulsed case [22]. Additionally, the threshold for a nonclassical correlation has now changed to . On the other hand, if the MZ interferometer is used for twophoton interference under continuous excitation, its interference visibility and nonclassical threshold actually match that of the HOM interferometer.
Finally as an experimental example, we discuss the observed twophoton interference in a MZ setup when interfering singlephoton generation from a cavityquantum electrodynamical source (figure 7); The data is reproduced from reference [22]. The MZ histogram shows groups of five peaks with nontrivial correlations. By applying the above analysis to each delay, i.e. , we can compute in terms of source correlations. However, the experimental beamsplitters often have a transmittivity () to reflectivity () ratio deviating from . Therefore, to fit the observed fivepeak pattern in correlation measurements, we need to include the individual transmittivities and reflectivities of the first beamsplitter () and second beamsplitter (). Recently, a very intuitive set of pictorial rules of this process were given in the supplementary information of reference [40]. Combining these rules with our formalism, the amplitudes of the five peaks in correlation measurements are given by:
(58) 
5.3 Simulated firstorder optical coherences of singlephoton sources and hence twophoton interference
In the previous simulation section 4.3, we examined the secondorder optical coherences of our three prototypical systems. We saw that under the right conditions, each system could act as a singlephoton source that emitted precisely only one photon per pulse. However in this simulation section, we will discuss the other requirement of an ideal singlephoton source, that it possess complete firstorder optical coherence. We will calculate the expected firstorder coherence for each of our model systems, and will directly show the expected twophoton interference between photons of the prototypical systems.