Self-calibrating Quantum State Tomography

We introduce and experimentally demonstrate a technique for performing quantum state tomography on multiple-qubit states despite incomplete knowledge about the unitary operations used to change the measurement basis. Given unitary operations with unknown rotation angles, our method can be used to reconstruct the density matrix of the state up to local sigma-z rotations as well as recover the magnitude of the unknown rotation angle. We demonstrate high-fidelity self-calibrating tomography on polarization-encoded one- and two-photon states. The unknown unitary operations are realized in two ways: using a birefringent polymer sheet---an inexpensive smartphone screen protector---or alternatively a liquid crystal wave plate with a tuneable retardance. We explore how our technique may be adapted for quantum state tomography of systems such as biological molecules where the magnitude and orientation of the transition dipole moment is not known with high accuracy.


Introduction
Quantum state characterization is essential to the development of quantum technologies, such as quantum computing [1,2], quantum information [3] and quantum cryptography [4]. The measurement of multiple copies of the quantum state and subsequent reconstruction of the state's density matrix is known as quantum state tomography (QST) [5,6].
QST relies on performing measurements in a number of different bases to gain access to all the information required to reconstruct the density matrix. Different measurement bases are typically obtained by applying unitary operations to the state before measurement. For systems that couple to the electromagnetic field via their transition dipole, precise knowledge of the magnitude and orientation of the dipole moment is vital to performing such measurement-basis changes. For example, QST of trapped ion qubits requires precise control of the orientation of the ions via magnetic fields [7].
These techniques can not be immediately applied to systems that are either not well characterized, such as biological molecules, e.g. those found in photosynthetic systems [8,9], where the magnitude and orientation of the transition dipole moment is not known with high accuracy; or systems that experience variability during the fabrication process, such as colloidal quantum dots [10].
With such systems in mind we introduce a method for self-calibrating tomography (SCT), which successfully reconstructs the state of the system despite incomplete knowledge about the unitary operations used to change the measurement basis. Other efforts to characterise photosynthetic systems are also underway, e.g. quantum process tomography using ultrafast spectroscopy was proposed by Yuen-Zhou et al. [11,12] and Hamiltonian tomography was proposed by Maruyama et al. [13].
In Section 2, we show that given local unitary operations with unknown rotation angles, but known and adjustable rotation axes, it is possible to reconstruct the density matrix of a state up to localσ z rotations, as well as recover the magnitude of the unknown rotation angles.
We demonstrate SCT in a linear-optical system using polarized photons as qubits in Section 3. An inexpensive smartphone screen protector, i.e. an uncharacterized birefringent polymer sheet, is used to change the measurement basis. We go on to investigate the technique's robustness to measurement noise and retardance magnitude, and demonstrate SCT of a two-qubit state using liquid crystal wave plates with tuneable retardances.
In section 4, we investigate the application of this technique to quantum state tomography of photosynthetic systems.

Theory
In this section, we give a brief review of quantum state tomography and introduce our technique for self-calibrating tomography.

Quantum State Tomography
The state of a qubit, given by the density matrixρ, can be decomposed into a sum of orthogonal operators σ i : whereσ 0 is the identity operator andσ 1−3 are the Pauli matrices. The coefficients λ j are given by the expectation values of the orthogonal operators, Given an unknown state,ρ can be reconstructed by performing a number of projective measurementsμ j on subsequent copies of the state. The measurement statistics are given by where N j is a constant that depends on the duration of data collection, detector efficiency, loss etc. Inserting (1) into (2) gives The equations given by (3) can be solved for λ i , which can then be inserted into (1) to reconstruct the unknown density matrixρ. The operatorsμ j must be chosen such that the equations in (3) form a linearly independent set. If N j are known, three different measurements suffice, given the normalization condition λ 0 = 1. Typically, N j is unknown, but assumed to be equivalent for all measurements, in which case, four measurements are required to reconstruct the state. In practice, due to fluctuations in measurement statistics, it is advantageous to employ a maximum likelihood estimation method [5] to reconstructρ rather than directly solving (2), which can lead to unphysical density matrices. We can define the measurement operatorsμ j in terms of a fixed known measurement operatorμ 0 and different unitary operationsÛ j as follows.
In terms of the measurement statistics, this is equivalent to applying the unitary operator U j to the state prior to making measurements given byμ 0 : where we used the cyclic properties of the trace.

Self-calibrating Tomography
Now consider the case where the unitary operator, used to change the measurement basis, is characterized by an unknown parameter α. Recall that a unitary operation is simply a rotation on the Bloch sphere. In this paper, we restrict ourselves to rotations about axes that lie on the equatorial plane of the Bloch sphere, i.e. the x−y plane, given byÛ where the index j refers to a particular ϕ and ν, where ϕ is the known angle between the axis of rotation and the x−axis and ν is real number that, along with α, determines the rotation angle. As we see later, the x−y plane restriction arrises naturally in the systems we consider. The measurement operators are defined aŝ The relationship between the measurement statistics and the parameters which characterize the density matrix will be given by By making an additional measurement, it is possible to solve the series of equations given by (9) for λ i and |α|. As before, measurements are chosen such thatμ j (α) generate a linearly independent set of equations. If this condition is met, the equations yield two sets of solutions, corresponding to ±α, i.e. either the actual state or the phase-flipped version of the state. It is straightforward, but non-trivial, to extend this formalism to multiple-qubit systems. However this is not required as in practice one may use a maximum likelihood estimation method adapted for SCT to reconstruct the density matrix, which we will introduce in the next section.

Maximum Likelihood
The solutions for λ i in (3) and (9) are typically quite sensitive to fluctuations in measurement statistics, which can lead to unphysical density matrices. It is therefore favourable to use a maximum likelihood algorithm to reconstruct the density matrix.
Assuming the noise on the measurements has a Gaussian probability distribution, the probability of obtaining a set of N measurement results {n j } is [5] where A is a normalization constant and σ j is the standard deviation of the jth measurement (given approximately by √ n j ). The expectation value of each measurmeent is given bȳ whereρ is the state we are trying to solve for. To ensure a normalized, positive and hermitian density matrix, we restrictρ to the form where for a single qubit and for a two qubit system. Maximising the likelihood that the physical density matrix ρ gave rise to the data {n j }-i.e. maximising (10)-can be reduced to finding the minimum of the function For large uncertainty in the measurement statistics, the likelihood function can have multiple local minima and the maximum likelihood algorithm may converge to a state and retardance that do not correspond to the actual state. This can be overcome given a priori knowledge of the purity of the state, or by taking more measurements.

Linear-optical implementation
For a linear-optical test of this method, we performed SCT on polarization-encoded oneand two-qubit states.

Single qubits
To prepare single-qubit states, we first pumped a 1mm long type-I down-conversion beta-barium borate (BBO) crystal with a 405nm diode laser, as shown in Fig. 1  a), which leads to spontaneous generation of pairs of photons. The detection of a horizontally polarized photon in one spatial mode of the down-converted photon-pair heralds the presence of a horizontally polarized photon in the other mode, i.e. the single qubit. Using the method described in Bocquillon et al. [31], and based on an estimated average rate of 2-photon events in the signal conditional on a detected photon in the heralding arm of 0.96s −1 , we calculated the g (2) (0) of our heralded single-photon source Figure 1. Experimental scheme: a) A single-qubit state is prepared, unknown unitary operations are implemented using BPSs or an LCWP and the state is measured using a trusted, fixed measurement; b) An entangled two-qubit state is prepared, unknown local unitary operations are implemented using LCWPs and the two-qubit state is measured using a trusted, fixed measurement.
to be g (2) (0) = 0.067, indicating a high quality single-photon source. Quarter-and half-wave plates can then prepare the qubit in arbitrary polarization states.
For the trusted, fixed projector, we choseμ 0 = |R R|, implemented with a quarterwave plate and a polarizing beamsplitter, followed by coincidence-counting, as shown in Fig. 1 a). To change the measurement basis, a unitary operation was applied to the state prior to measurement, implemented using either one or two identical pieces of a birefringent polymer sheet (BPS). The BPS is an optical element that alters the polarization state of light in the same way as a wave-plate, with the crucial difference that the retardance of the BPS is not known. The effect of the BPS on the state of the photon can be described byR φ (α) given in (6), where α is the retardance of the BPS, ϕ is the alignment of the BPS's optical axis with respect to the polarization of the light, andσ i are in the |R /|L basis. Using different alignments of the BPS, we constructed the following operatorsμ j (α) =Û j (α) †μ 0Ûj (α) using different implementations ofÛ j (α): Note that rotations of 2α were required to satisfy the linear-independence condition on (3). These rotations were realized with a double-pass through the BPS-the sheet was cut in half and the pieces placed in succession. The beam passed through sections of the pieces that were within close proximity before partition, to ensure consistency in the birefringence. One piece was then removed to perform rotations of α, before removal of the remaining piece to measureμ 0 . This measurement sequence eliminated the need for realignment between measurements.
In this particular linear-optical implementation of the scheme, both output ports of the PBS are easily accessible. We therefore measured the complimentary projectors, µ j (α) =Û j (α) †μ 0Û j (α) whereμ 0 = |L L|, to determine N j . Using this method, SCT was performed on 14 states: six that lie at the ends of the three axes and eight other states defined by all combinations of θ = ± π 4 and φ = {± π 4 , ± 3π 4 } where |ψ = cos θ 2 |R + e iφ sin θ 2 |L . Each state was reconstructed using the maximum-likelihood estimation method described in Section 2.3 in conjunction with the YALMIP solver [32].
For comparison, standard tomography (ST) was also performed on each of the states. ST required fewer projectorsμ j =Û † jμ 0Ûj , given by the unitary operators: and their complimentary projectors,μ j =Û † jμ 0Û j whereμ 0 = |L L|. We used the fidelity F = Tr[ √ρ stρsct √ρ st ] 2 to characterize the success of SCT, whereρ st andρ sct are the density operators reconstructed using ST and SCT respectively. We achieved an average fidelityF = 0.99 ± 0.01 and determined the retardance of the BPS (at 405nm) to be α = 0.58 ± 0.03 (close to π/6). Individual fidelities for each state are shown in Fig. 2 a).

Dependence on retardance
To determine the effectiveness of SCT as a function of the retardance, we repeated the experiment using BolderVision Optik liquid crystal wave plates (LCWP), as shown in Fig. 1 a). In an LCWP, the birefringence is a function of the applied voltage, allowing for the retardance to be set arbitrarily. A single LCWP was used to implement rotations of α and 2α using different voltage settings. Fidelities for different setting of α are shown in Fig. 2 b). The average fidelity decreases with the size of the retardance. The fidelity for states |ψ = |R and |ψ = |L always remains close to unity. This is a consequence of our choice of measurement projectors. Using the LCWP, we find a slightly lower average fidelity with a greater spread over the different states, even for values of α larger than that of the BPS. We attribute this to errors associated with the non-linear dependence of the retardance on the voltage passed through the LCWP. This influences how accurately 2α can be set with respect to α.

Dependence on noise
To investigate the dependence of the fidelity on the amount of noise, we performed SCT on the state |ψ = |H using a retardance α = π/32. The noise scales inversely with the square root of the number of counts and by collecting from 150 to 3000 photons we varied the noise between roughly 1% and 10%. For each noise setting, SCT was repeated 100 times. The resulting fidelities and predicted retardances are presented in Fig. 2 c). Notice that the majority of predicted fidelities and retardances correspond to the actual values. However, recall that the likelihood function can have multiple local minima, making it possible for the maximum likelihood algorithm to converge to a state and retardance that do not correspond to the actual state. This shows up in the high-noise, or low photon-number, region and can be seen as a small bump in the fidelity histogram (and to a lesser extent in the retardance histogram) in Fig. 2 c). As the noise decreases, only high-fidelity results corresponding to the actual state and retardance remain. Fig.  3 shows the distribution of states on the Bloch sphere as predicted by SCT for high and low noise amounts. The darker red cluster of states in Fig. 3 a) also correspond to a local minimum in the likelihood function.

Entangled qubits
For multi-qubit tomography, entangled two-qubit states, |ψ = a|HH + b|V V , were prepared using two 1mm long type-I down-conversion BBO crystals with optical axes aligned in perpendicular planes [33] pumped by a 405nm diode laser, as shown in Fig.  1 b). The parameters a and b which characterise the state were tuned by changing the polarization of the pump beam with a HWP. LCWPs in each arm, implemented the unknown unitary operations, followed by detection ofμ 0 = |R R| in both modes. For this two-qubit experiment, LCWPs were chosen over the BPSs due to realignment issues that accompany the insertion (as opposed to removal) of optical elements in the  Fig. a), which corresponds to a local minimum in the likelihood function. There are 100 states shown in each sphere.
beam-path. The concurrence [34] for states reconstructed using SCT and ST, as well as the fidelity between the two, are shown for a variety of entangled states in Table 1.

Potential applications to photosynthetic systems
In this section, we give some detail about how this technique may be applied to performing quantum state tomography on biological molecules, such as chromophores in photosynthetic systems, where the transition dipole moment and its orientation are not known with high accuracy.

Unitary operation
Consider a molecule with two states of interest-the ground state |g and the excited state |e -that undergoes an interaction with an EM field, e.g. a laser pulse, given by  Table 1.

the interaction Hamiltonian
where E(t) is the electric field,d is the transition dipole moment operator and d eg = e|d|g = d * ge . An ultra-fast pulse with a carrier frequency ω c can be described by E(t) = E(t)e iωct e + h.c., where E(t) is the temporal profile of the pulse and e is the polarization vector. Under the rotating wave approximation, we can write the Hamiltonian in the interaction picture aŝ H = ( d eg · e)E(t)e i(ωc−ωe)t |e g| + ( d eg · e) * E(t) * e −i(ωc−ωe)t |g e| , where ω e is the |g → |e transition energy. The unitary governing the interaction is given byÛ [Ω]σy)/2 = e −iνα(cos(ϕ)σx+sin(ϕ)σy)/2 , (26) where Ω = 2 ( d eg · e) dtE(t)e i(ωc−ωe)t and Comparing with (6), we can see that (26) constitutes a rotation of angle να, where ϕ-determined by the relative phase of E(t)-defines the rotation axis. Notice from (27) that ν is linearly proportional to the square root of the intensity of the pulse at the transition frequency. We can therefore execute rotations of 2α by changing the intensity of the pulse, without any knowledge of d eg · e.

Projective measurement
We can perform a measurement corresponding to the operator by preparing multiple copies of the state and performing statistics on the number of photons emitted. A molecule in the excited state will-on a long enough time-scaleemit a photon, while one in the ground state will not; a molecule in a superposition of ground and excited states a|g + b|e will emit a photon with probability |b|. The measurement statistics will be given by whereÛ j (α) is given in (26) and N j depends on the solid angle of the detector, the position of the detector with respect to the dipole moment of the molecule and losses. Measurement statistics can be collected either by repeated preparation and measurement of a single molecule, recently demonstrated by Hildner et al. [35,36], which is analogous to the linear-optical implementation discussed above; or by preparing an ensemble of molecules in the same state and measuring the fluorescence of the sample, where the time-integrated intensity of the emitted light is directly proportional to the number of photons. An open question is whether it will be possible for the selfcalibrating technique to compensate for an ensemble of randomly oriented molecules. This is currently being investigated.
In principle, the extension to multiple entangled molecules should be straightforward: even if individual molecules can not be spatially isolated, independent rotations on each qubit could be executed with appropriately shaped pulses, as long as there are different transition energies for each molecule. However, in the case of coupled molecules, there are technical difficulties associated with addressing the transition energies of the individual molecules-the system rather absorbs light at the eigen-energies of the coupled Hamiltonian. In this case, a system of two coupled two-level systems can be approximated as a three-level system (where the fourth level is typically disregarded). Self-calibrating quantum state tomography of coupled two-level system will require an extension of our formalism to a three-dimensional Hilbert space. This is also currently being investigated.

Conclusion
In conclusion, we have introduced and demonstrated a technique for performing quantum state tomography of multi-qubit systems which does not rely on complete knowledge of the unitary operation used to change the measurement basis. Our technique characterizes the state of interest with high fidelity as well as recovering the unknown unitary operation. We find that in the context of QST of polarization-encoded qubits, it is possible to do away with a well-calibrated half-wave plate and replace it with an inexpensive uncharacterized piece of birefringent material.
A future extension of this work will attempt to incorporate unknown rotation axes and non-unitary transformations, such as those that introduce noise to the system, thereby reducing the purity of the state by an unknown amount.
We anticipate our method to be applicable to quantum state tomography of systems that do not have well-characterized transition dipole moments, such as biological molecules like those found in photosynthetic systems, as well as systems that experience variability during the fabrication process, such as colloidal quantum dots.
AMB, AD and DFVJ thank DARPA (QuBE) for support. DHM, LAR and AMS thank NSERC, CIFAR and QuantumWorks for support. We also thank Xingxing Xing for the serendipitous discovery of the birefringence of his smartphone screen protector.