This dissertation serves as a general introduction to Wigner functions, phase space, and quantum metrology but also strives to be useful as a how-to guide for those who wish to delve into the realm of using continuous variables, to describe quantum states of light and optical interferometry. We include many of the introductory elements one needs to appreciate the advantages of this treatment as well as show many examples in an effort to make this dissertation more friendly.

In the initial segment of this dissertation, we focus on the advantages of Wigner functions and their use to describe many quantum states of light. We focus on coherent states and squeezed vacuum with a Mach Zehnder Interferometer for many of our examples, also used by experiments such as advanced LIGO. Later, we will also analyze this setup in more detail with a full example including the effects of many noise sources such as phase drift, photon loss, inefficient detectors, and thermal noise. In this setup, we also show the optimal measurement scheme, which is currently not employed in experiment. Throughout our metrology discussions, we will also discuss various quantum limits and use quantum Fisher information to show optimal bounds. When applicable, we also discuss the use of quantum Gaussian information and how it relates to our Wigner function treatment.

The remainder of our discussion focuses on investigating the effects of photon addition and subtraction to various states of light and analyze the nondeterministic nature of this process. We use examples of photon additions to a coherent state as well as discuss the properties of an photon subtracted thermal state. We also provide an argument that this process must always be a nondeterministic one, or the ability to violate quantum limits becomes apparent. We show that using phase measurement as one’s metric is much more restrictive, which limits the usefulness of photon addition and subtraction. When we consider SNR however, we show improved SNR statistics, at the cost of increased measurement time. In this case of SNR, we also quantify the efficiency of the photon addition and subtraction process.


A Dissertation

Submitted to the Graduate Faculty of the

Louisiana State University and

Agricultural and Mechanical College

in partial fulfillment of the

requirements for the degree of

Doctor of Philosophy


The Department of Physics and Astronomy


Bryan Tomas Gard

B.S, Louisiana State University, 2012

May 2016


I must first thank Jonathan P. Dowling, for the proper balance of the “carrot and stick” principle, which likely has some quantum uncertainty relation. Jon provides the ideal level of encouragement or threats, as the situation requires and keeps progress moving, while managing not to be unreasonable. If Jon does not have the answer, he can always point you in the right direction all while not making you feel inferior, for asking possibly “dumb questions”.

I also thank Hwang Lee for always having the insight to see the clever solution, when my attempt at an explanation confuses, more than explains. I would also like to comment that Hwang’s presence at meetings always has a calming effect (perhaps to offset Jon’s).

I thank the other members of my committee, Mark Wilde, Thomas Corbitt, and dean’s representative Evgueni Nesterov, for offering their time and input to my work.

I appreciate the guidance and many many helpful discussions with Emanuel Knill at NIST-Boulder, where I spent two summers working on quantum research. The vastness of his knowledge in many sub disciplines of physics is truly humbling.

I would also like to acknowledge support from the National Physical Science Consortium (NPSC) and the National Institute of Standards and Technology (NIST) for supporting my graduate career through fellowship and summer internships. This support allowed me to focus on my studies, while taking classes, and focus on research, instead of juggling these duties along with teaching or grading duties. In large part, I believe this is why I am able to graduate in a somewhat brief time as a graduate student.

I would also like to thank the many enriching collaborators I have worked with. Researchers at the Boeing Corporation, Barbara Capron, Claudio Parazzoli, and Ben Koltenbah, have provided much assistance in the form of providing many discussions on an industry project in quantum metrology which served as invaluable tool to further my understanding of quantum optics and metrology. Also included in this project were many helpful discussions with Christopher Gerry.

I also must thank my wife, Lynn Gard, who somehow put up with my travels, long nights, never-ending typing, no money, complaining, headaches, and too many other issues to list. Always a source of comfort, I know I couldn’t have accomplished this much without her.

Last but certainly not least, my parents, James and Janet Gard, who were always encouraging and at least pretended to listen to my many “exciting” realizations in the fun world of quantum mechanics. I would also like to thank my brothers, Darin and Chris Gard for acknowledging that there can be only one true nerd in the family and accepting defeat.

Chapter 1 Introduction

In order to analyze various quantum metrology configurations, we require a quantum mechanical description of light and the effects of common optical elements. There are many mathematical models that accomplish this, through the use of wave vectors [16, 32, 59], density matrices [39, 72, 80, 66, 33], and Wigner functions [84, 48, 35, 49], to name a few. In this dissertation, we will discuss the use of continuous variables in phase space [28, 67], their advantages and potential issues. We will also use this treatment to describe common interferometer setups that involve parameter estimation [52, 77, 68], quantify their photon statistics in terms of signal to noise ratio, and discuss the effects of many quantum optics techniques. Specifically, we will investigate an interferometric setup like LIGO [63, 30, 85, 6, 79] and also describe the effects of relatively exotic operations like photon addition and subtraction [8, 9, 88, 89, 17, 22, 14, 87, 56].

With the recent, first ever, direct detection of gravitational waves [3], many large interferometers around the world continually attempt to measure further gravitational wave events [63, 85, 5, 2, 79, 6]. The Laser Interferometer Gravitational-Wave Observatory (LIGO) in Livingston, Louisiana and Hanford, Washington, are two examples of such interferometers. The initial configuration of LIGO was comprised of a coherent state and vacuum coupled in a Michelson interferometer [5] (henceforth, we refer to this as the classical setup). This scheme is a classical strategy and is limited to a classical bound on the phase variance measurement, the Shot Noise Limit (SNL) [41, 31]. The objective is to measure a relative phase shift induced in one arm of the interferometer by a passing gravitational wave. Recently, the first direct measurement of gravitational waves has been shown [3]. With this amazing accomplishment, comes the need for further measurements. Despite the remarkable precision obtained by this method, improvements are still possible. One such improvement for Advanced LIGO consists of input states of a coherent state and squeezed vacuum state, a configuration first proposed by Caves [24] and shown to achieve a superior phase variance measurement as compared to the previous classical input states [1, 43, 29, 15].

While there are many technical challenges in using a true quantum setup such as this, we show here that some of the measurement techniques previously used in the classical setup, are no longer optimal and even may exhibit problems with effects such as phase drift and thermal noise (another source of noise typically found in configurations like LIGO, radiation pressure noise, is not considered here). In order to achieve an optimal measurement scheme and investigate their resistance to phase drift, we turn to quantum measurements such as homodyne and parity measurements and compare them to a standard intensity measurement. We show that, under ideal conditions, the parity measurement achieves the smallest phase variance, but under noisy conditions, the parity measurement suffers greatly, while the homodyne measurement continues to achieve superior phase measurement. In general, we divide our results into two regimes, the low power regime (), in which different detection schemes can lead to significantly different phase variances, and the high power regime (), which applies to Advanced LIGO, and where all detection schemes approach the optimal bound.

The use of photon addition or subtraction is an implementation of noiseless amplification. First proposed by Agarwal and Tara [8, 9], noiseless amplification can be used to enhance a general signal with no added noise, but with the requirement that it does so nondeterministically. If one desires an amplification of signal, it must either come with additional noise (e.g. a deterministic squeezer [82, 39]) or it must be probabilistic. Either of these cases ensures consistency with fundamental conditions such as no super-luminal communication.

Here we discuss the use of photon addition and subtraction as a probabilistic amplifier and its effects on various sources, including thermal and coherent light [39, 12, 62, 81, 83]. Unlike many past discussions of this implementation, we consider the case of photon addition and subtraction at the output of a Mach-Zehnder Interferometer (MZI) [86]. Since we are using an MZI model, we are then in the realm of metrology and can therefore use many previously developed techniques from this field. The reasoning behind using the probabilistic amplification operation at the output is simply a model of the limit of control over a specified system. In the case of an externally measured source, meaning a source one has no direct control over, deterministic amplification proves useless for a phase estimation problem, as the added noise always kills any benefit of the amplification. This limit extends to the metric of signal to noise ratio (SNR), where a deterministic amplifier always amplifies signal along with its noise, leaving SNR invariant at best. More concretely, this restriction means any modification to a standard MZI must be done after the phase shifter . With this restriction in mind, the question then remains, since a deterministic amplifier doesn’t provide any benefit, is there any hope for a probabilistic amplification process?

Recent discussions by Caves [26] show the use of post selection schemes and their place in quantum metrology protocols. As we discuss later, we also show that post selection schemes alone do not allow for increased phase information and also discuss some of the pitfalls when using post selection schemes that can lead to deceivingly positive results. This result however, does not invalidate the usefulness of post selection schemes in metrology, when other metrics are considered.

Chapter 2 Wigner Functions in Phase Space

2.1 Phase Space

The use of continuous variables in phase space serves many purposes. While the choice of mathematical treatment is ultimately a choice of preference, here we will discuss the advantages of using continuous variables in phase space. A perhaps more standard approach to quantum metrology is with the use of wave vectors or density matrices. Mathematically, there are many choices available when considering which method to work in, all of which offer a full description of quantum mechanics, but a specific choice may offer computational simplicity or, as we argue here in the case of Wigner functions, offer a visual aspect as well as some connections to known measurements, such as Parity measurement.

We first begin with a visual description of phase space, in terms of the conjugate variables, . This brief introduction is used as a simple illustration of various common states of light, in phase space, and we discuss a more rigorous mathematical approach in the following section. Shown in Figure 2.1, we see many different states of light depicted in phase space. We note that the absence of photons, the vacuum state (black circle), and thermal state (checkerboard) partially overlap at the center of phase space, while a Fock state (red ring), a purely quantum mechanical state which contains exactly photons, is represented by a thin ring whose radius is determined by the chosen number state. A displaced vacuum state, or coherent state (blue circle), quantum mechanical description of laser light, is also shown along with a general uncertainty in its quadrature values due to the limits of quantum mechanics. The amount of displacement in this state, given by is related to its average photon number by . In general, all of these states can be squeezed, which trades uncertainty between its quadratures. The angle of this squeezing process determines which quadrature is enhanced, while the other suffers. A squeezed coherent state (orange) is shown and has been squeezed along the axis which enhances the quadrature while increasing the uncertainty in the quadrature, in compliance with the uncertainty principle relating these quadratures, . Unless otherwise stated, for the remainder of this document, we use the natural units convention of . A review of these states of light can be found in [39].

Figure 2.1: Various states of light shown in phase space in terms of position, momentum space, (x,p). Fock state (red) shown as a ring. The radius of this ring depends on the photon number chosen. Vacuum state (black circle) shown at the center. Thermal state (checkerboard) shares partial overlap with the vacuum but is always strictly larger. Coherent state (blue circle) shown in top right quadrant with phase angle and uncertainties in each quadrature also shown. This state is displaced from the center by an amount given by . A squeezed coherent state (orange) is shown in the bottom right quadrant and is squeezed to reduce the uncertainty in the quadrature.

These displayed states of light compose a typical set of the most commonly described forms of light. More exotic forms of quantum light, such as photon added or subtracted states, two mode squeezed vacuum, Schrödinger cat states, etc. can be visually depicted with various combinations of the states shown in Figure 2.1. For example, Schrödinger cat states can be shown by a superposition of two coherent states, while photon subtracted thermal states have the vacuum portion of a thermal state removed. These pictures in phase space of various forms of light can be very instructive when we consider things such as the various statistics of these states. We can see that the Fock state reduces to the vacuum state for , but is not allowed to reduce to a single point, as this would violate the uncertainty principle. A similar comparison can be made for the thermal state. A thermal state of zero average photon number, also reduces to the vacuum state, as does a coherent state with no displacement. In this way, one may say that the vacuum state is the principal state of light, which other states are modifications of, through various optical processes.

2.2 Wigner Functions

Wigner functions, first introduced by Eugene Wigner in 1932 [84], is a quasi-probability distribution for a given state of light, in phase space. The term “quasi” is used since these distributions may take on negative values, which also means they are not typical (classical) probabilities. One can show that any state whose Wigner function obtains a negative value, is a quantum state, but this statement is not an if and only if statement, meaning all negative Wigner functions correspond to quantum states, but not all quantum states have Wigner functions which attain negative values. This treatment of light is a full mathematical description and is connected to that of density matrices by [7],


where y are eigenvectors of the quadrature operators satisfying and can also be connected to the so called characteristic function with the following ordered relations,


with standing for the Euclidean norm on and for Wigner functions, , Husimi Q-functions, and P-functions, . For our purposes, we will solely focus on Wigner functions for the remainder of this document. As mentioned earlier, the coherent state is also known as displaced vacuum by virtue of, . Note that the Wigner functions can be defined for any two conjugate variables, not always position, momentum space . Another typical representation is in complex phase space . These two bases are connected by the relations,


The position and momentum operators obey the bosonic commutation relations , while the creation and annihilation operators obey


where the commutator is defined by . Now that the mathematical construction of Wigner functions had been covered, we turn to the practicality of using Wigner functions in quantum optical metrology.

Many familiar properties in terms of density matrices convert to Wigner functions, but the advantageous aspects are apparent that we transition from an infinite dimensional discrete sum to a continuous variable integral of size , where is the number of spatial modes. Specifically some of these properties are,

where this property is seen as the normalization requirement of any quantum state (which enforces that probabilities sum to one). From this simple constraint, we see that, instead of performing the trace of an infinite sum, we instead integrate our Wigner function over a finite set of variables. While both of these techniques are typically straightforward, computationally we find that integrals are typically much more easily managed, without the need of typical “truncation” tricks as used with sums. In general this comment applies to all measurements with Wigner functions, meaning we are typically able to obtain analyctical results, while working with density matrices frequently (though, not always) results in numerical answers. The purity of a quantum state is also commonly used to classify states. In terms of Wigner functions this is simply,

where the state is pure if and is mixed if . This condition also has a pleasing visualization in terms of the Bloch sphere, where pure states lie on the surface of the unit sphere, while mixed states lie inside the volume of the unit sphere. In general, we can see that a trace over a density matrix corresponds to an integral of our Wigner function. This idea extends to that of partial traces. For example, consider a two mode density matrix, . This state has the property, , where we have traced over the “B” mode. Similarly, our Wigner function has the property and we have integrated out the “B” mode.

2.3 Gaussian States

We can now discuss another particular strength of working in phase space and the use of Wigner functions, that of Gaussian form, which classify many typical states of light. Any Wigner function that is Gaussian in form, has many simplifications that can be made. In this section we will review many of the properties of such states. A general Gaussian function can be written as,

where, and is an positive definite matrix [7] and ensures normalization, such that . In terms of Wigner functions then, the simplest example of a Gaussian form is the Wigner function of the vacuum state, given by,


which one can notice has the promised Gaussian form. Shown in Figure 2.2, we see that we also may visualize these various states of light easily from the use of Wigner functions. Compared to our phase space picture shown in Figure 2.1, we can notice that the phase space view is simply a projection (or slice) of the full Wigner function.

Figure 2.2: Wigner function for the vacuum state as a function of the phase space quadratures, . We can easily observe this state’s Gaussian form.

Another example of a typical state of light, in terms of Wigner functions is the coherent state,


where is amplitude of the coherent state and is the phase angle. In this form, it is instructive to notice that the form is similar to that of Eq. (2.0) but is displaced in both the and directions. The amount of displacement is controlled by the size of and the direction of displacement is controlled by . We can also show the form of a thermal state is given by,


where is the average photon number in the thermal state. Again, one can connect this to the vacuum state for .

All the previous Wigner functions adhere to the Gaussian form, but as an example of a non-Gaussian form, we turn to the simplest example of a Fock state, the single-photon state. All Wigner function Fock states can be described by


where is the Laguerre polynomial. Note that for this reduces to Eq. (2.0). For , we then have the single-photon Fock state, which is necessarily non-Gaussian and shown in Figure 2.3.

Figure 2.3: Wigner function for the single-photon Fock state as a function of the phase space quadratures, . We can observe that not only is this a non-Gaussian form but also attains negative values near the origin.

It is important to note that this Gaussian discussion is not specific to Wigner functions and is, in no way, a limitation on these usefulness of this treatment, but merely a choice of simplification. While, in general there is no problem with representing non-Gaussian forms with Wigner functions, the remainder of this section will discuss strategies that one can utilize when restricting to Gaussian only states. The advantage in using this restriction is in properties of the Gaussian form itself. Any Gaussian function can be fully described by its first and second moments. With our choice of basis, , this amounts to only needing to specify and . In practice, the quantities we are really interested in are the mean and covariance. The mean and covariances given by


where , corresponding to the two conjugate phase space variables . One can notice that for , twice the variance of the quadrature. This factor of two is merely a convention and due to the definitions discussed earlier, but should be noted when comparing to references with other definitions. Now that we have established definitions of our mean and covariance, we can connect them back to Wigner functions with the relation,


where is the number of spatial modes, X,d are vectors of phase space variables and means, respectively, and is the full covariance matrix for the desired spatial modes. One may also notice that this construction lends itself to one measurement choice in particular, homodyne measurement, as we work specifically in the first and second moments of the phase space quadratures, which is exactly what a homodyne process attempts to measure. We will discuss this aspect more thoroughly, in later sections. This, along with a treatment of Quantum Fisher Information (QFI) are the main benefits from using a Gaussian-only treatment, but this discussion is also left for later sections.

Chapter 3 Interferometer Model

3.1 Mach Zehnder Interferometer Model

In order to fully model the interferometric process, we must have a description of the effect of various optical elements on the various states of light presented earlier. There are many choices of how to model states of light, but there are also choices in how to describe the propagation of this light through optical elements. In general, one can describe the movement of the the entire state through various the optical elements, a Schrödinger picture, or describe the effect of these elements on the mode operators, a Heisenberg picture. While both are mathematically complete descriptions, they are not necessarily computationally equivalent [38] and therefore we choose to describe this propagation in the Heisenberg picture.

We have seen some examples of common types of light used in theoretical quantum optics, in terms of Wigner functions; here we will discuss how we model the evolution of these states through various optical elements, in terms of general Wigner functions, as well as the Gaussian restriction discussed earlier. First we will show how the mode operators, specifically our basis choice of , evolve through various optical elements. In general, each optical element is represented by a symplectic () matrix, of dimension and acts on a vector of phase space variables of length , where is the number of spatial modes, by the following input-output relation.


Perhaps the simplest transform to begin with is that of the displacement operator, briefly mentioned in the previous section. This operator can take a vacuum state into a coherent state and can be described by the transformation,


where determines the magnitude of displacement and , the direction. For another example, a typical optical element, a beam splitter, splits the amplitude of two incident electric fields, according to a transmission parameter intrinsic to the beam splitter. The action of this device can be described by the matrix,


where is the transmissivity of the beam splitter. Another common optical element is a phase shifter. This single mode element can be described by,


where is the phase imparted to a single mode, usually treated as an unknown parameter to be estimated. While this is a typical way to model a phase shifter, it is also instructive to construct a so-called symmetric phase shifter, which is used to form a balanced phase between two modes. This symmetric phase shifter is a two mode transform of the form,


Squeezing is a quantum operation that is commonly used to outperform classical only treatments. A single mode squeezing operator can be described by,


where is the squeezing parameter and is the squeezing angle. The squeezing strength can be related to the gain of the squeezer by the relation , for . In some circumstances, it may be possible to assume and change variables using to obtain,


which, if applicable, can greatly improve computation time. A two mode version of this squeezer is described similarly as,


which can also be written in terms of the gain of this squeezing process, if so desired. These squeezing processes have many different implementations, including the use of non-linear crystals and atomic clouds. The specific differences and challenges with each of these implementations is not the focus of our discussion however [78, 74]. We have described several typical optical elements and how they transform phase space variables; therefore we can now model various states of light and the action of many optical elements on these states. Since we have given various single and dual mode transforms, a natural question is how to deal with the case of mixing various combinations of single and dual mode transforms, as these are of different dimension. In treatments of density matrices, this is given by the tensor product, , but in phase space it is given by the direct sum, . For example, if we wish to pass a two mode Wigner function through a 50/50 beam splitter and a single mode displacement on the “upper mode”, we would use a full transform of,


where it is clear that we have simply inserted the displacement operator into the top-left block and the identity matrix into the bottom right block, as the bottom mode does not incur a transformation at this point. We also use the convention that the first entry listed in a direct sum labels the first mode, etc. Some caution should be used at this point as one can note that the direct sum of two single mode squeezers, .

Now as an illustrative example, let us consider a full model of a typical Mach-Zehnder Interferometer (MZI), which is mathematically equivalent to a Michelson Interferometer (MI), shown for a specific choice of input states, in Figure 3.1.

Figure 3.1: A general Mach-Zehnder interferometer with coherent () and squeezed vacuum () states as input. A phase shift, , represents the phase difference between the two arms of the MZI, due to a path length difference between the two arms. Our goal is to estimate , which in the case of LIGO, would be caused by a passing gravitational wave.

For reasons discussed later, let us consider input states of a coherent state and squeezed vacuum. Since, in general, these two states are assumed to be initially independent, the two mode Wigner function for this choice is simply the product of the individual Wigner functions. The full Wigner function can then be written as,


where X is a list of the mode variables , corresponding to the position and momentum components of each spatial mode. The average photon number in the coherent state is and in the squeezed vacuum state . Both states have equal phases, as this gives rise to the optimal phase sensitivity (discussed later) and are taken to be .

The propagation of this Wigner function is accomplished by a simple transformation of the phase space variables through the MZI, dictated by its optical elements, discussed earlier. These transformations are described by


where both beam splitters are fixed to be 50/50 and we have chosen to use a symmetric phase model, shown in Eq. (3.0), in order to simplify calculations as well as agree with many other references [80, 51]. Using these transforms, the total transform for phase space variables is given by,


From here, the final variables (denoted by etc.) are inserted to the initial Wigner function to obtain the Wigner function at the output. This is a straightforward task but we do not show the result here as the state is fairly cumbersome. This process can be followed for the various choices of input states and showcases the basic method of Wigner function evolution, in the Heisenberg picture.

We can again discuss the use of a Gaussian restriction and show an alternative method to evolve a given Gaussian state of light through the MZI. This method utilizes Quantum Gaussian Information (QGI) [7] and lends itself to simplifying some measurement schemes, as well as allowing a particularly useful calculation of Quantum Fisher Information (QFI), discussed later. As we showed earlier, once a Gaussian form is assumed, one need only be concerned with the mean and covariance. It is these two quantities that we will evolve through our MZI. Initially, we require the mean and covariance of the various states of light, and luckily most have a very simple form, which showcases the usefulness of this method. The mean of a given state of light is tied to its displacement value, which for any state centered at the origin, in phase space, is clearly zero. Therefore, states such as the thermal state and vacuum state (along with the squeezed versions of these) have a zero mean. For completeness, we show all the means in vector form as,


where it is understood that the first entry corresponds to and the second entry, . The covariance matrix for each of these states is given by,


Note that the coherent state carries its statistics only in its mean, while the thermal state carries its statistics only in its covariance. Much like the process described earlier, we now need a way to evolve these two parameters (mean and covariance) through various optical elements. The transforms that act on these parameters are the same as we used earlier, with the only difference being in how the covariance evolves. Both mean and covariance evolve according to,


where is given by the same optical element transforms shown in Eq. (3.0)-(3.0). At this point, our mean and covariance have evolved through various optical elements and we now have the statistics at the output. If so desired, we could now construct the Wigner function at the output with Eq. (2.0) which would give the same Wigner function with the previous evolution method, described above. For completeness, we will again use the example of our MZI from Figure 3.1 and list the mean and covariance, at the output, as our example input states, a coherent state and squeezed vacuum both maintain Gaussian form. This choice of input states and the fixed topology of our MZI evolve the mean and covariance to,


where and , where again we have fixed .

Chapter 4 Quantum Measurement

4.1 Measurements using Wigner functions

Now that we have shown ways in which one can evolve various states of light through an MZI, we can discuss the role of the detectors. While there are many different types of detectors, instead of restricting to specific physical processes of detection, we describe the detection process in terms of quantum operators, shown in Eq. (4.0). In general, measurement operators are computed, in terms of Wigner functions, by utilizing [7],


where is a symmetrically ordered function of phase space operators, . A significant difference between the density matrix approach and Wigner functions is the requirement of symmetric ordering. Using Wigner functions, all measurements are assumed to be symmetrically ordered in their field operators and therefore specific measurements must take this symmetric ordering into account. For example, a typical measurement operator, the number operator (also defined as an intensity measurement), given by is not symmetrically ordered. In order to symmetrize this operator, we write it in symmetric form as,

and utilize the commutation relations shown in Eq. (2.2) to obtain,

and therefore have that a measurement of intensity on a single mode Wigner function, in terms of symmetric ordering is given by


and its second moment, which requires similar, but more complicated symmetric ordering,


While this symmetrization is always required when using Wigner functions, it is a straightforward process. Note that this requirement is needed for any products of the field operators, but a typical homodyne measurement, i.e. , is already symmetrically ordered.

Also of typical use as a state characterization tool is constructing the photon number distribution for a given state. This is physically done by running many trials of an experiment and performing number counting at either detector, each trial giving a certain number of photons. Over many trials then, we can construct the probabilities of the state having any number of photons. When using density matrices, if one works in the number basis, then this distribution is essentially calculated for free but luckily using Wigner functions, we can also obtain this distribution with a relatively simple calculation given by constructing a generating function and differentiating it according to,


where is our generating function, our photon number distribution and it’s clear that , since is necessarily discrete. This construction allows us to calculate the photon number distribution for any state of light and use it to characterize the properties of this state, which we will utilize when we analyze our full examples later in this dissertation. This construction also demonstrates that we can recover the discrete nature of the quantized nature of photons, even when describing them in a continuous variable space.

Investigating the so-called balanced homodyne measurements further, we see that it is typically discussed in terms of its implementation, that is, the target beam, is incident on a 50/50 beam splitter along with a coherent state (or commonly referred to as a local oscillator) of the same frequency as the target beam (typically this coherent state is derived from the same source as the input state). After the beam splitter, the two outputs are collected on detectors and an intensity difference measurement is performed. This process is shown in Figure 4.1.

Figure 4.1: A typical balanced homodyne measurement is performed with a target beam incident on a 50-50 beam splitter along with a coherent state. The two detectors then perform an intensity difference measurement.

In practice, the phase of the coherent state is adjusted, allowing a homodyne measurement along any arbitrary direction in phase space, . Therefore this process amounts to nothing but an attempt to measure the mean value of any arbitrary superposition of the state in phase space. This allows us to notice that we do not need to fully model this measurement scheme as in Figure 4.1, but can instead simply calculate any homodyne measurement as , on the state represented by . For simplicity, if we consider one particular homodyne measurement, along the direction then a specific homodyne measurement, along with its variance, can be calculated, using Wigner functions with,


This technique also gives us some insight into the usefulness of a homodyne measurement, in terms of our Gaussian information techniques described earlier. Since this technique works directly in terms of mean and covariance, it’s clear that homodyne measurements are essentially calculated for free with this treatment, as the output state’s mean vector is exactly what a homodyne measurement attempts to capture.

We can also notice that with a bit of work, we can use the Gaussian techniques to easily calculate an intensity measurement. This can be done by utilizing


where we have used the relation . In the proper, symmetrical form, this can be written in terms of mean and covariance as,


Note that are quantities already calculated in this Gaussian treatment, and so we need only square them. Also of significance is that is only a 22 matrix (not infinite like for a density matrix) and its trace is trivial. Therefore, this detection scheme is also calculated very simply when using this Gaussian treatment. Calculating the second moment of this operator using Gaussian techniques, while possible, is not as simple as one would perhaps desire. If we assume for the chosen input state of light, this calculation becomes tractable, but its advantages over using the Wigner becomes less apparent. With this observation, we then can say, that Gaussian information techniques can greatly simply some calculations, when particular input states are assumed (Gaussian) and particular measurement schemes are considered (mainly homodyne). Contrasted with this Gaussian only treatment, we can then see that Wigner functions provide a completely general approach to handle any types of input state and any types of measurements.

Another typical measurement choice is one that we briefly mentioned earlier, an intensity difference measurement defined by , where it is clear this is now a measurement on two spatial modes. However, one can easily see that this detection scheme is a direct adaptation of the single mode intensity measurement. Specifically this measurement can be calculated, in terms of Wigner functions, as,


As a final example of another measurement choice, we consider the single mode parity measurement, defined as . Under this choice, we are able to utilize another benefit of describing our system in terms of Wigner functions, since it has been shown [70] that the parity measurement satisfies , or in words, the expectation of the parity measurement is given by the value of the Wigner function, at the origin of phase space. From this property along with the property, , so , the parity measurement, is perhaps the simplest to calculate of the choices presented here, as an integral is not required.

These cover some of the common types of measurements and while, in principle, there are many more, the methods discussed here showcase some of the benefits and properties of utilizing Wigner functions and Gaussian information. It seems only fair then to discuss some of the difficulties with using Wigner functions. Likely, the same difficulties exist for other, analytical methods, but these are computational issues and not problems with the construction itself. We have shown how a general quantum measurement can be calculated using Wigner functions, which involves integrating the quantum operator against the output Wigner function. In general we have seen that these Wigner functions are typically exponentials and can, depending on the input states, have very complicated exponents. We assume one would calculate these integrals using software such as Mathematica and the following uses formatting and language to follow this software, but many other software likely have a similar form. In the case of these complicated integrals, even with sufficient assumptions within the software, the computational time required to complete these integrals can be lengthy (greater than a day on a relatively powerful desktop PC). Instead of moving to a higher end PC (or access to a super computer), we note a technique that allows us to calculate such intractable integrals. This is only due to inefficiencies in how the software handles complicated integrals but is useful for our purposes and therefore relevant to include in our discussion. Instead of integrating something such as, , a generic quantum measurement, where is given, for example, by , we instead form a list of replacement rules. First we group our Wigner functions exponent into terms dependent on the phase space operators , with Mathematica’s “Collect” command. For example, if our Wigner function is named “Wig01” then we would perform

where “[[2,2]]” simply selects only the exponent of Wig01 (assuming the exponent is in the list position 2,2). With the form of this collection, we then see the form that our Wigner function takes and can integrate a general exponential of this form and create a replacement rule for this integration. In principle, integrating a general exponential of the same form as our Wigner function shouldn’t provide any benefit, as the only difference is in what we consider constants, but the software seems to not understand this difference well. Once we have the form of the exponential, we integrate a general form and create a replacement rule for this form. Then, instead of integrating our specific Wigner function, we instead use the replacement rule with the formatted Wigner function. The collection of this process with a specific choice of a coherent state and squeezed vacuum, as input, and homodyne measurement is shown in Figure 4.2. Here we show an example of this code as well as the times (in seconds) shown for the various calculations and while this calculation is easily manageable in either treatment, one can see the benefit of using our replacement method when the calculation times begin to become cumbersome. Note that the total time to calculate the general integral, perform the replacement, and simplify is , while a direct integration and simplification is , approximately 80 times longer. The main drawback of this method is the organization required for composing the specified Wigner function to a general form, but one can imagine forming a large database of these various forms which could alleviate some of the overhead of this process. In even more complicated cases, this replacement method still maintains its benefit over the direct integration method but typically requires more prep-work, depending on the chosen measurement scheme.

Figure 4.2: Example of replacement method instead of directly integrating our Wigner function. Outputs show the time taken to calculate (in seconds) as well as the output. The fourth line from the bottom shows the direct integration time, along with its result, while the last line shows our replacement methods time and its result.

4.2 Phase Measurement

In the previous section we discussed how to utilize Wigner functions to perform quantum measurements. In this section we will show how these measurements apply to a typical goal in metrology, performing a phase measurement [13, 76, 47]. As shown in Figure 3.1, we can use a MZI, along with an unknown phase , to represent the unknown effect of a material (perhaps a material with an unknown index of refraction) interrogated inside the MZI. While the exact effect of this material can vary, depending on the physical situation, we can model the effects of this material as an unknown variable, . If this is the model one uses, then there are some typical benchmarks we can discuss.

The so-called Shot Noise Limit (SNL) is typically defined for the MZI in terms of the total average photon number which enters it [1, 73]. First developed by Caves [24] in a fixed setup of sending a coherent state and vacuum into an MZI, the SNL is defined by . Using the example from previous sections, if we consider sending a coherent state and squeezed vacuum into an MZI, we then have, . This sets a boundary on the variance of our unknown phase measurement distinguishing classical techniques from quantum techniques. A similar limit, the Heisenberg limit, while not a hard limit, is typically also used as a reference and defined by . This Heisenberg limit should be considered with caution as it is obtained from the uncertainty relation and assuming that , so that, . However, this assumption is violated for certain states of light, such as the two-mode squeezed vacuum, and it has been shown that this state indeed can achieve a phase estimate below the HL [11]. Spefically, the HL is technically only a hard limit for states of definite photon number (such as Fock states [39, 36]) but not necessarily a hard limit for states of average photon number (all Gaussian states [7] as well as some non-Gaussian). So while we may use the HL as a reference, it should be considered carefully.

There are two principal ways in which we can calculate the variance of our phase measurement. The first way one can show the phase variance of a chosen measurement is to use error propagation [42, 57, 25] in the form of,


where is a chosen quantum operator. The use of this formula connects the variance of our unknown phase, to any quantum operator, with the utilization of the Taylor expansion according to,


where is a non-linear function, and we have assumed that is the only parameter to be estimated and is uncorrelated with any other variables. This gives us a variance of , a kind of quality of our measurement. In general, one particular quantum operator can outperform, giving a smaller phase variance, another quantum operator, and therefore the process of searching for the optimal measurement scheme can be very challenging. Note that in order to utilize Eq. (4.0), we require the first and second moments of the chosen operator, since . This requirement may be trivial, as in the case of the parity operator or may be fairly complicated, as in the case of the intensity difference operator.

An alternative to using this error propagation treatment is with the use of Classical Fisher Information (CFI) [71, 23]. Instead of considering a specific quantum measurement, we instead consider probabilities of events (though, certain quantum measurements lead to probabilistic events). For example, the intensity measurement gives a measurement of the average number of photons entering a detector. A particular question we could ask in this case is,“What is the probability of the detector receiving exactly one photon?” The probability distribution for this case can then be used with the relation,


with the condition that , that is the ’s represent complete probabilistic events. Returning to our brief example, if we let be the probability of our detector detecting a single photon, then,

where we have used that . In order to obtain the probability distributions in question, we again use our Wigner functions and form projective measurements. In this case, the probability of a detector receiving exactly one photon is given by the projection of the Wigner function onto the single photon subspace. This is achieved by performing the measurement of [64],


where is the Wigner function for the single photon Fock state, shown in Eq. (2.0). This implementation shows that we have the option of constructing probability distributions through projective measurements, which may be simpler or more physical than a particular quantum operator. Another important note is that of so called click detection (on/off detection), which a typical APD performs, meaning it responds to the presence of any number of photons, but cannot discriminate the number of photons that is present. In this case, we can model this APD by the projection onto the subspace of all photon numbers states, other than vacuum, by using the projection . This would give us the probability distribution of a detector receiving any number of photons, other than zero, just as the APD measures. While projection treatment is perhaps advantageous in that it does not require the second moment calculations, it instead requires probability distributions and, typically, projective measurements. It is worth noting that some quantum operators themselves are also probability distributions, such as in the case of the parity operator and therefore either method can be used to obtain the phase variance. This treatment also suffers the same problem as in the error propagation; however, it has many possible choices of probabilities and thus a reasonable question to ask is, which one gives the smallest phase variance?

In order to show that we have obtained the smallest phase variance possible, we construct the Quantum Cramér Rao Bound (QCRB), through the Quantum Fisher Information (QFI) [18, 80, 65, 53, 54]. This treatment differs from the previous in that it will not depend on a measurement choice, and instead gives the best phase variance bound possible with any possible measurement choice allowed by quantum mechanics (meaning it is described by a positive operator valued measure (POVM)). In general this quantity is difficult to calculate and while it gives us the ultimate lower bound on a phase variance measurement, as mentioned, it does not depend on the measurement choice and therefore does not tell us the “optimal” measurement directly, merely the best possible with any measurement. In order to show the optimality of a measurement then, we separately calculate the phase variance obtained from a specific measurement and compare it to the QCRB. While it is not an exclusive bound (meaning multiple measurement may achieve the QCRB), if a chosen measurement achieves the QCRB, then no other measurement can outperform this. A particular measurement which reaches the QCRB is then said to be optimal, but it may have many technical challenges in actually implementing and a full noise model is useful at this point to characterize the effect of various noise sources on various measurements, which we discuss in later sections.

While there are a variety of ways to calculate the QFI and QCRB, most apply easily to some specific class of states but not others, or may be widely applicable but very calculation intensive when exotic states of light are considered. Here we will focus on the calculation of the QFI through Gaussian information as well as directly from the Wigner function. The Gaussian information method applies to any Gaussian state, pure or mixed, while using the Wigner function applies only to pure states. In terms of Wigner functions, we may calculate the QFI of any pure state according to [69],


where is the number of spatial modes. It should be noted, as we mentioned earlier with our rule replacement method, integrals of Wigner functions can be computationally intensive at times, but the rule replacement method works equally well for this calculation. While this method seems deceivingly simple, in practice, even with the discussed methods, computing this integral with a relatively complicated Wigner function can be quiet challenging. If we further assume that our state is pure and Gaussian then we can find that,


where are the mean and covariance, respectively. Note that this definition only differs slightly from [69], as we use different definitions of quadratures, . If we relax the condition that the state must be pure, but retain the requirement of Gaussian form, we can still find that the QFI is manageable using Gaussian information. Specifically, the QCRB in this case takes the form [37],


where is the symplectic matrix defined by


We should also note that this construction is valid after a change of basis, as it proves to be to a computational advantage, with the change of basis according to,


This calculation of the QCRB is particularly useful since the only requirement is the Gaussian form, but can accommodate mixed states, which include the classical thermal state. This thermal state is crucial in some noise models [44, 90] and necessarily takes any pure state to a mixed state when thermal noise is considered in the noise model and therefore this calculation is particularly useful when considering realistic noise models.

4.3 Noise Modeling

4.3.1 Photon Loss

For any accurate model, one must always consider the effect of various sources of noise. Here, we will discuss the modeling of various sources of loss and noise, typically found in interferometers. This includes photon loss to the environment [60], inefficient detectors [46], phase drift, and thermal noise [44, 90]. In general, each of these can be mitigated through various techniques, but not completely eliminated and therefore, in our attempt to model a realistic interferometer, we must have a way to model these unavoidable effects.

First we consider photon loss in the model by way of two mechanisms, photon loss to the environment inside the interferometer and photon loss at the detectors, due to inefficient detectors [39]. Both of these can be modeled by placing a fictitious beam splitter, of variable transmissivity, inside the interferometer with vacuum and a interferometer arm as input and tracing over one of the output modes, to mimic loss of photons to the environment. This process is shown in Figure 4.3.

Figure 4.3: A model of linear photon loss by passing a general state of light, , through a variable transmissivity beam splitter, with vacuum as the second input. One of the output modes is then traced over to represent the loss of photons which then results in a lossy state, .

These linear photon loss mechanisms manifest themselves as a simple change of variables with respect to the average photons. Using our example input states from before, this amounts to a change of variables of and , where is the detector efficiency (in decimal), is assumed to be equal for all detectors and is the photon loss inside the MZI, assumed to be equal in both arms. This variable replacement greatly simplifies calculations over the full model of inserting many fictitious beam splitters and our assumption of equal losses and equal detector efficiencies is fairly reasonable if one reasons that both arms of the interferometer are in similar media and both detectors are identical. It is clear that the value of the transmissivity of the variable beam splitter, , controls the amount of loss, in decimal, to the environment. This process of variable replacement can be generalized to other input states with a similar variable replacement condition of simply, , where is the average photon number in each spatial mode. If our assumption of equal losses in both arms and identical detector inefficiencies is relaxed, then a full calculation becomes necessary.

4.3.2 Phase Drift

Another common effect in MZI’s is the random drift of phase. The way in which we treat phase drift, to the best of our knowledge, has not been presented elsewhere. Normally, we assume the unknown phase to be a fixed value, but in practice it may vary slightly over many experimental trials; we call this effect, phase drift. Generally, when we calibrate our MZI, we would place a control phase in one arm of the interferometer, which allows us to tune the interference between the two arms. We could assume we have infinite precision with our control phase, but in practice, whatever mechanism controls the value of our control phase, it can drift slightly over many experimental trials. In principle, we attempt to set this control phase to an optimal value; in order to give rise to a measurement with the best statistics, however, the control phase value will vary around this optimal phase setting. For this reason we aim to show this phase draft in a more mathematical way and therefore we can use the analytical forms of the various measurement phase variances, as a function of unknown phase, , and simulate phase drift by computing a running average of the phase variance, with a pseudo-randomly chosen phase, near the optimal phase, for each measurement. This is accomplished in the following way: we find the true optimal phase (typically a multiple of ) and allow a pseudo-random number generator to choose a phase near this optimal phase (within 20% above or below the optimal value), which we use in place of the optimal phase. We perform this pseudo-random process over many trials, each trial forming the average of the choices from previous trials. In this way, after many trials, our random choosing approaches the true optimal phase. This mimics the idea of the experiment in that, a single measurement provides very little information on the unknown phase, , and it is only with many measurements that we can say we have obtained a good estimate of the unknown parameter. The behavior of various detection schemes under the effects of phase drift are not necessarily identical, as we show in later sections, some measurement schemes, specifically those that utilize multiple spatial modes, perform better than single mode measurement schemes, under phase drift.

4.3.3 Thermal Noise

In addition to photon loss, detector efficiency, and phase drift, we also model the inevitable interaction with thermal noise from the environment [44, 90]. This is accomplished much in the same way as a photon loss model, but here we consider a thermal state incident on a fictitious beam splitter, on both arms of the interferometer, inside the interferometer and trace out one of its output modes, shown in Figure 4.4. This allows a tunable amount of thermal noise (by changing the average photon number in the thermal state), into the interferometer. Again, as we will see in a later section, the effects of this injection of thermal photons can vary, depending on the measurement scheme that is considered.

Figure 4.4: Model of the interaction of a general state , inside the MZI, with the environment, represented by a thermal state, .

Chapter 5 Photon Addition and Subtraction

5.1 Photon Addition and Subtraction

5.1.1 Photon Addition and Subtraction with Wigner functions

We have discussed much of the general way one can use Wigner functions to model an MZI and the various benchmarks one uses to qualify many of the typical measurement choices. We will now discuss the modeling of an interesting, but less typical consideration in MZIs, the use of addition or subtraction of photons from a beam of light. First proposed by Agarwal and Tara [8, 9], an implementation of noiseless amplification, photon addition and subtraction has received much attention over the past decade [88, 89, 17, 22, 14, 87, 56]. While there is some choice on how to exactly model the operation of photon addition and subtraction, we first suggest a word of caution. As this is the realm of probabilistic noiseless amplification, we should not overlook the fact that we must consider the probabilistic nature of this process.

One way that some have chosen to model the operator of photon addition and subtraction is with the use of the photon creation and annihilation operators. Using these operators, we can write the addition of a photon, in terms of a density matrix, as,


where and is the average photon number, prior to photon addition, and shows that this process must also be properly re-normalized. Similarly, we can write photon subtraction as,


where . However, we would of course like this representation, in terms of Wigner functions, therefore with some work we obtain, [17]


with representing the Wigner function to be photon added or subtracted in the spatial mode. All the previously presented methods are valid for any state of light, but we shall refer to this construction as the mathematical treatment of photon addition or subtraction, as it is mathematically correct, but does not necessarily faithfully reproduce all the aspects of photon addition and subtraction, as we will show in the following discussion.

The use of the mathematical description of photon addition and subtraction may be sufficient, depending on the specific model one is considering, but if a full model of the effects of this process are desired, then an alternative method should be used. This alternative method consists of using a physical process, along with projective measurements, to implement photon addition and subtraction.

Figure 5.1: SPDC model of photon addition. A general state, and vacuum are incident on a pumped non-linear crystal. Occasionally, depending on the pump strength, this crystal emits photon pairs into the two output modes. One output mode is used to herald that the other mode has been photon added.
Figure 5.2: Beam splitter model of photon addition. A general state, and a single photon Fock state are incident on a variable transmissivity beam splitter. On the condition that one output receives no photons, the remaining mode is known to be photon added.
Figure 5.1: SPDC model of photon addition. A general state, and vacuum are incident on a pumped non-linear crystal. Occasionally, depending on the pump strength, this crystal emits photon pairs into the two output modes. One output mode is used to herald that the other mode has been photon added.

For photon addition we can consider a few different implementations, but here briefly highlight two implementations. One uses Spontaneous Parametric Down Conversion (SPDC) [21, 40, 58] to facilitate photon addition while the second model uses a variable transmissivity beam splitter, along with a single photon source. While these two processes involve physically different optical elements, both of these models accurately reproduce the mathematical treatment described above, in the limit of vanishing interaction. Both implementations of these models are shown in Figure 5.2 and Figure 5.2.

Shown in Figure 5.2, we see that the first of these models uses an active optical element, the SPDC process, to probabilistically generate photon pairs. One individual from the pair strikes a detector, which confirms that the other spatial mode has the remaining of the pair. This outgoing beam is then known to have been photon added. We have presented the framework needed to model this process, in terms of Wigner functions, but for clarity the specific transforms one requires are given by, the preperation of the initial two mode state, given by,


where is the Wigner function of the vacuum state in spatial mode one. The action of the two mode squeezer is then,


where , and we use the convention that the upper spatial mode is labeled by the subscript one. The action of the condition that a single photon must be emitted into one of the output modes and detected is then modeled as a projection of that mode onto the single photon subspace, given by,


where is the Wigner function for a single photon Fock state, given by Eq. (2.0), in the lower mode, as we also use the convention that mode labels transfer through optical elements in the transmitted labels are kept; i.e. , across an optical element, the lower mode becomes the upper and the upper becomes the lower. It is also important to note that since we have performed a projective measurement, it now requires proper normalization, performed simply by,


We use this normalization constant to enforce the condition that our Wigner function is normalized at all times,

It is also important to note that this normalization constant is also the probability that the detector in the herald mode, detects a single photon. This is particularly useful later when we account for the probabilistic nature of this scheme.

The transforms for modeling the photon addition process with the beam splitter model follows much the same construction, shown in Figure 5.2. Specifically, the initial state,


the action of the beam splitter of arbitrary transmissivity,


and the projection,


while the normalization constraint is much the same as previously described.

Clearly these two physical processes appear different in the way in which they function, as they have different initial states, different parameters and projections. However, one will find that in the limit of vanishing interaction, both of these models exactly match that of the mathematical treatment of photon addition. Specifically this means that, for the SPDC model, as the final state will exactly match the action of the mathematical treatment. Similarly, for the beam splitter model, this condition is satisfied for . This leads us to why the mathematical model is somewhat limiting. Also in these limits, as we will see later, the probability of these events (given by ) approach zero. This means that while these treatments reduce to the mathematical model, physically the process is extremely unlikely to ever occur, something not very useful in experiments!

Figure 5.3: Beam splitter model of photon subtraction. A general state, and a vacuum state are incident on a variable transmissivity beam splitter. On the condition that one output receives one photon, the remaining mode is known to be photon subtracted.

In the case of photon subtraction, we use a similar model to the beam splitter model of addition, though with different conditions. Shown in Figure 5.3, this process utilizes the same beam splitter of variable transmissivity, but different initial states and projective measurements. Following the same process as before we would perform this photon subtraction following the regimen of preparing the initial state,


utilizing the beam splitter transform shown in Eq. (5.0), a similar projection as in Eq. (5.0) and proper normalization as in Eq. (5.0). As with the models of photon addition, this model also agrees with the mathematical model in the limit that , but again , in this limit.

These processes can both be generalized to the case of photon addition and subtraction by generalizing the initial state and post selected projective measurement in a straight forward way. Specifically, for photon additions we modify the SPDC model to the case where photons are detected in the herald mode, confirming that photons have been added to the remaining mode. In the case of the beam splitter photon addition model, we consider the photon Fock state as the initial state but still condition on receiving no photons in the output mode, confirming photons have been added to the remaining output mode. For the case of photon subtraction we consider a projection onto the photon Fock state, which projects the remaining state into the photon subtracted subspace. As one may expect, while these generalize fairly easily, the probabilities of these photon additions and subtractions also decrease, with increasing , as we will show later.

5.1.2 Photon Addition and Subtraction Statistics

Now that we have described the ways in which to model photon addition and subtraction, we can investigate some of the interesting properties of these processes. While their models are perhaps fairly intuitive to follow, we will show that their effects are not exactly as one may expect. As a first step, shown in Figure 5.4 and Figure 5.5 are the Wigner functions for a photon added coherent state and a photon subtracted thermal state, respectively. Note that a coherent state is invariant under photon subtraction since it obeys the relation , meaning a coherent state under photon subtraction, returns the same coherent state. These states are generated following the prescription described above. We can notice that in the case of the Single Photon Added Coherent State (SPACS), even though we have only added a single photon, the Wigner function has changed drastically as a result, attaining negativity near the origin, confirming that this is a quantum state.

In the case of a Single Photon Subtracted Thermal State (SPSTS), shown in Figure 5.5, we see a drastically different distribution, which is positive everywhere, but has a small dip at the origin, caused by the photon subtraction process.

Figure 5.4: Wigner function for a SPACS. Notice that a typical coherent state has a particularly simple Gaussian Wigner function, but our SPACS is clearly non-Gaussian and even attains negativity in its Wigner function. The relevant parameters have been fixed at values of,
Figure 5.5: Wigner function for a SPSTS. Notice that a typical thermal state has a particularly simple Gaussian Wigner function, but our SPSTS is clearly non-Gaussian and contains a “dip” at the origin. The relevant parameters have been fixed at values of, and
Figure 5.6: Photon number distributions for normal coherent and thermal states and their photon added and subtracted transformations. Parameters have been fixed at, .

While these Wigner functions serve as a useful way to observe these states, we can analyze them further by constructing their photon number distributions as discussed in Eq. (4.0). In Figure 5.6, we show the photon number distributions for a, normal coherent state and thermal state, along with the SPACS and SPSTS. We can notice that in both cases, the photon number distribution has shifted to the higher photon numbers and the probabilities of lower photon number has decreased. In the case of photon addition, this is fairly expected, as we are adding photons to the state. However, in the case of photon subtraction, the photon number has also increased, a somewhat surprising claim. This can be explained by the fact that, under the assumption of photon subtraction, we gain partial knowledge of the photon number distribution. Specifically, consider the thermal state of low average photon number (); if the photon subtraction process has succeeded, then we know that the state must have contained at least a single photon, which partially collapses the photon number distribution near the origin. This has the effect of shifting the distribution to the higher photon numbers. A similar effect happens in the case of photon addition. Under the assumption of successful photon addition, we are sure that the distribution cannot contain zero photons; thus our distribution shifts to the higher photon states. This is a visual description of the effects of the simple cases of photon addition and subtraction, but in general it is advantageous to have a qualitative way to describe the effect of photon addition and subtraction. Using the process described earlier to model the process of photon addition and subtraction, we can see that a pattern emerges in the average photon numbers for the addition of photons to a coherent state as well as subtraction of photons from a thermal state. These expressions can be found to be,