Satisfying the Einstein-Podolsky-Rosen criterion with massive particles

In 1935, Einstein, Podolsky and Rosen (EPR) questioned the completeness of quantum mechanics by devising a quantum state of two massive particles with maximally correlated space and momentum coordinates. The EPR criterion qualifies such continuous-variable entangled states, where a measurement of one subsystem seemingly allows for a prediction of the second subsystem beyond the Heisenberg uncertainty relation. Up to now, continuous-variable EPR correlations have only been created with photons, while the demonstration of such strongly correlated states with massive particles is still outstanding. Here, we report on the creation of an EPR-correlated two-mode squeezed state in an ultracold atomic ensemble. The state shows an EPR entanglement parameter of 0.18(3), which is 2.4 standard deviations below the threshold 1/4 of the EPR criterion. We also present a full tomographic reconstruction of the underlying many-particle quantum state. The state presents a resource for tests of quantum nonlocality and a wide variety of applications in the field of continuous-variable quantum information and metrology.

I n their original publication 1 , Einstein, Podolsky and Rosen describe two particles A and B with correlated position x B ¼ x A þ x 0 and anti-correlated momentum p B ¼ À p A (Fig. 1a). When coordinates x A and p A are measured in independent realizations of the same state, the correlations allow for an exact prediction of x B and p B . EPR assumed that such exact predictions necessitate an 'element of reality', which predetermines the outcome of the measurement. Quantum mechanics, however, prohibits the exact knowledge of two noncommuting variables like x B and p B , since their measurement uncertainties are subject to the Heisenberg relation Dx B Dp B ! 1 2 . EPR thus concluded that quantum mechanics is incomplete-under their assumptions that are today known as 'local realism'. Later, the notion of EPR correlations was generalized to a more realistic scenario, yielding a criterion 2,3 for the uncertainties Dx inf B , Dp inf B of the inferred predictions for x B and p B . The EPR criterion is met if these uncertainties violate the Heisenberg inequality for the inferred uncertainties Dx inf B Dp inf B ! 1 2 . The EPR criterion also certifies steering, a concept termed by Schrödinger 4,5 in response to EPR, which has attracted a lot of interest in the past years 6 . An experimental realization of states satisfying the EPR criterion is not only desirable in the context of the fundamental questions raised by EPR, but also provides a valuable resource for many quantum information tasks, including dense coding, quantum teleportation 7 and quantum metrology 8 . Some quantum information tasks specifically require the strong type of entanglement that is tested by the EPR criterion, as for example one-sided device independent entanglement verification 9 .
Up to now, the creation of continuous-variable entangled states satisfying the EPR criterion was only achieved in optical systems. In a seminal publication 10 , the EPR criterion was met by a twomode squeezed vacuum state generated by optical parametric down-conversion. In this experiment, and in more recent investigations 11,12 , continuous variables are represented by amplitude x A/B and phase p A/B quadratures, satisfying the commutation relation [x A/B , p A/B ] ¼ i. These quadratures can be measured accurately by optical homodyning. The correlations are captured by the four two-mode variances V AE Þ . These variances were proven to fulfil a symmetric form of Reid's inequality 3 V À x ÂV þ p o 1 4 , which is a sufficient EPR criterion since V À x ¼ ðDx inf B Þ 2 and V þ p ¼ ðDp inf B Þ 2 . In recent years, continuous-variable entangled optical states have been applied for proof-of-principle quantum computation and communication tasks 7 . Despite these advances with optical systems, an experimental realization of EPR correlations with massive particles is desirable, because of the similarity to the original EPR proposal and since massive particles may be more tightly bound to the concept of local realism 2,3 .
Entangled states of massive particles have been generated in neutral atomic ensembles, promising fruitful applications in precision metrology due to the large achievable number of entangled atoms [13][14][15][16] . They have been created by atom-light interaction at room temperature 14,17 , in cold samples [18][19][20][21][22] and by collisional interactions in Bose-Einstein condensates 13,16,[23][24][25] . For Gaussian states of two collective atomic modes, the inseparability criterion 26,27 has been used to demonstrate entanglement 14,17,28 , but the strong correlations necessary to meet the more demanding EPR criterion V À x ÂV þ p o 1 4 have not been achieved so far. Here we report on the creation of an entangled state from a spinor Bose-Einstein condensate (BEC), which meets the EPR criterion. We exploit spin-changing collisions to generate a twomode squeezed vacuum state in close analogy to optical parametric down-conversion. The phase and amplitude quadratures are accessed by atomic homodyning. Their correlations yield an EPR entanglement parameter of 0.18 (3), which is 2.4 s.d. below the threshold 1/4 of the EPR criterion. Finally, we deduce the density matrix of the underlying many-particle state from a maximum likelihood reconstruction.

Results
Two-mode squeezed vacuum. In our experiments, a BEC with 2 Â 10 4 87 Rb atoms in the Zeeman level (F, m F ) ¼ (1, 0) generates atom pairs in the levels (1, ±1) due to spin-changing collisions (Fig. 1b), ideally yielding the two-mode squeezed state where x ¼ Ot is the squeezing parameter, which depends on the spin dynamics rate O ¼ 2p Â 5.1 Hz and the spin dynamics duration t ¼ 26 ms. The notation n; m j i represents a two-mode Fock state in the two Zeeman levels (1, ±1). The generated two-mode squeezed state can be characterized by the quadratures A=B þ a A=B Þ and p A=B ¼ i ffiffi 2 p ða w A=B À a A=B Þ for the two levels (1, ±1). These exhibit EPR correlations, since the variances The state fulfills Reid's EPR criterion for x4 1 2 ln 2 ð Þ % 0:35 which corresponds to a spin dynamics duration of 411 ms. In the limit of large squeezing, our set-up presents an exact realization of the perfect correlations with massive particles envisioned by EPR.
Quadratures and the EPR criterion. The quadratures in the two modes are simultaneously detected in our experiments by unbalanced homodyne detection (see Methods). Atomic homodyne detection was first demonstrated in ref. 28, and later applied to discriminate between vacuum and few-atom states in a quantum Zeno scenario 29 . A small radiofrequency pulse couples 15% of the BEC in the level (1, 0) (the local oscillator) symmetrically to the two modes (1, ±1). The local oscillator phase is represented by the BEC phase relative to the phase sum of the two ensembles in (1, ±1). It can be varied in our experiments by shifting the energy of the level (1, À 1) for an adjustable time. From the measured number of atoms in both levels, we obtain a linear combination of the quadratures according to X A=B y ð Þ ¼ x A=B cos y À p 4 À Á þ p A=B sin y À p 4 À Á . Figure 2a shows two-dimensional histograms of these measurements for y ¼ 3 4 p and y ¼ 5 4 p, corresponding to the x-and p-quadratures. The histograms demonstrate the strong correlation and anticorrelation of these two quadratures, as expected for the EPR case. The variances along the two diagonals, represented by Þ , are shown in Fig. 2b and reveal the expected two-mode squeezing behaviour. From these measurements, we quantify the EPR entanglement by Reid's criterion, yielding V þ x ÂV À p ¼ 0:18 3 ð Þ, which is 2.4 s.d. below the limit of 1 4 . In addition, the data also fulfil the inseparability criterion as V þ x þ V À p ¼ 0:85 8 ð Þ, which is 15 s.d. below the classical limit of 2 (Fig. 2d), and meets the criterion for a symmetric ('two-way') steering between the systems 6 . We estimate that the product value could be reduced to V þ x ÂV À p ¼ 0: 13 3 ð Þ if the radiofrequency intensity noise was eliminated by stabilization or postcorrection. The experimental creation of entangled massive particles that satisfy the continuous-variable EPR criterion presents the main result of this publication.
Squeezing dynamics. Figure 3 shows the squeezing dynamics due to the spin-changing collisions. For these measurements, we fix the local oscillator phase to the values yE3p/4 and yE5p/4 to record only the x-and p-variances. As a function of the evolution time, the variances V À x ; V þ p are squeezed below the vacuum reference of 1, while the variances V þ x ; V À p exhibit an antisqueezing behaviour (Fig. 3a,b). From these data, we extract the EPR parameter V À x ÂV þ p , as a function of evolution time (Fig. 3c). The EPR parameter is quickly pushed below 1 and follows the prediction for an ideal-squeezed state. It eventually reaches a minimum at the optimal squeezing time of 26 ms, as used for the data in Fig. 2. The data were well reproduced by a simple noise model, which includes a radiofrequency intensity noise of 0.4% and a local oscillator phase noise of 0.044p (see Methods). ARTICLE Full state reconstruction. The total of 2,864 homodyne measurements obtained for different local oscillator phases at the optimal evolution time allow for a full reconstruction of the underlying many-particle state. Previously, tomography of an atomic state was demonstrated either by reconstruction of the Wigner function 30 or the Husimi Q-distribution 21,25 . However, both methods yield a characterization of the state's projection on the fully symmetric subspace only. The well-developed methods in quantum optics 31 allowed for a full reconstruction of an optical two-mode squeezed state by homodyne tomography 11,32 . Despite the beautiful tomography data, the optical state reconstruction assumed either Gaussian states or averaged over all phase relations, such that the coherence properties could not be resolved.
In contrast, we obtain an unbiased, positive semidefinite density matrix by maximum likelihood reconstruction 31,33 of the experimental data, free of any a priori hypothesis. This represents the second major result of this publication. The recorded data for each local oscillator phase are binned in two-dimensional histograms (Fig. 2a) presenting the marginal distributions for the x A/B and p A/B variables. The reconstructed state is the one that optimally reproduces the measured histograms by a superposition of harmonic oscillator wave functions 31 . The coefficients of this superposition are estimates of the density matrix elements of the underlying quantum state (see Methods). Figure 4 shows the result of the reconstruction. The diagonal matrix elements (Fig. 4a) witness the predominant creation of atom pairs. The two-particle twin Fock state 1; 1 j i shows the strongest contribution besides the vacuum state. Likewise, the twin Fock states 2; 2 j i to 4; 4 j i have the strongest contribution for a given total number of particles. The strong nonclassicality of the reconstructed state becomes also apparent in the distributions of the difference and the sum of the particles (Fig. 4b,c). The distribution of the number difference is strongly peaked at zero and is much narrower than a Poissonian distribution with the same number of particles. The distribution of the total number of atoms shows an indication of the characteristic even/odd oscillations, which is caused by the pair production in the underlying spin dynamics.

Discussion
For an evaluation of the created state, we have extracted a logarithmic negativity of 1.43±0.06 from the reconstructed density matrix. This value is above the threshold of zero for separable states and signals non-separability of the reconstructed state. The quantum Fisher information 34 F Q for the state projected on fixed-N subspaces reveals that F Q = n ¼ 1: n is the average number of particles. Since F Q = n41 the state is a resource for quantum enhanced metrology 34 . Furthermore, we fit an ideal two-mode squeezed state with variable squeezing parameter x to the reconstructed two-mode density matrix with maximal fidelity. With a fidelity of 78.4%, the experimentally created state matches a two-mode squeezed state with a squeezing parameter of x fit ¼ 0.63. The fidelity increases to 90% if we include local oscillator phase noise and statistical noise. The unwanted contributions in the density matrix, including the off-diagonal terms in Fig. 4a, can be well explained by four effects. First, the purity of the reconstructed state is limited by the finite number of homodyne measurements. Second, small drifts in the microwave intensity of the dressing field (on the order of 0.1%), which shifts the level (1, À 1), result in a small drift of the local oscillator phase. Third, a small drift of the radiofrequency coupling strength during homodyning virtually increases the variance in the (x A þ x B ) and the (p A þ p B ) directions. Finally, we did not correct for the detection noise of our absorption imaging.
Our experimental realization of the EPR criterion demonstrates a strong form of entanglement intrinsically connected to the notion of local realism. In the future, the presented atomic twomode squeezed state allows to demonstrate the continuousvariable EPR paradox with massive particles. Since the two modes A and B are Zeeman levels with an opposite magnetic moment, the modes can be easily separated with an inhomogenous magnetic field to ensure a spatial separation. The nonlocal EPR measurement could then be realized by homodyning with two spatially separated local oscillators. These can be provided by splitting the remaining BEC into the levels (2, ±1) which have the same magnetic moment as the two EPR modes. Furthermore, this set-up can be complemented by a precise atom number detection to demonstrate a violation of a Clauser-Horne-Shimony-Holt-type inequality. Such a measurement presents a test of local realism with continuous-variable entangled states. In this context, neutral atoms provide the exciting possibility to investigate the influence of increasingly large particle numbers and possible effects of gravity. (1, þ 1) to (2, þ 2) and from (1, À 1) to (2, À 2) followed by another light pulse. This cleaning procedure ensures that no thermal or other residual excitations are present in the two output modes, which may destroy the EPR signal 35 . Figure 5 shows a schematic overview of the following experimental sequence. A microwave frequency which is red-detuned to the transition between the levels (1, À 1) and (2, À 2) by about 208 kHz is adiabatically ramped on within 675 ms. The microwave shifts the level (1, À 1) by about 500 Hz, depending on the chosen detuning, to compensate for the quadratic Zeeman effect of q ¼ 491 Hz, such that multiple spin dynamics resonances can be addressed 16,36 . Each resonance condition belongs to a specific spatial mode of the states (1, ±1) to which the atoms are transferred. If the energy of the level (1, À 1) is reduced, then more internal energy is released, and higher excited spatial modes are populated (for details, see ref. 36). Here we choose the first resonance, where spin dynamics leads to a population of the levels (1, ±1) in the ground state of the effective potential. This ensures an optimal spatial overlap between the atoms in the three contributing levels. This resonance condition is reached, when the input state (two atoms in the BEC in the level (1, 0) at the energy of the chemical potential) is exactly degenerate with the output state (two atoms in the levels (1, ±1) including dressing, trap energy and mean-field shift). Due to this degeneracy, the phase relation between the initial condensate and the output state stays fixed during the spin dynamics evolution time. For this configuration, we have checked that spin dynamics is the only relevant process, which produces atoms in the state (1, ±1) (see ref. 29, Fig. 3). Subsequently, the microwave dressing field is ramped down within 675 ms, stays off for a variable duration between 25 and 1,150 ms and is quickly switched on again. The variable hold time allows for an adjustment of the local oscillator phase relative to the output state.

Methods
For the atomic homodyning, a radiofrequency pulse with a frequency of 1.834 MHz and a duration of t ¼ 30 ms couples the level (1, 0) with the levels (1, ±1). The microwave dressing field is chosen such that both radiofrequency transitions are resonant, but the resonance condition for spin dynamics is not fulfilled. Afterwards, the dipole trap is switched off to allow for a ballistic expansion. After an initial expansion of 1.5 ms to reduce the density, a strong magnetic-field gradient is applied to spatially separate the atoms in the three Zeeman levels. Finally, the number of atoms in the three clouds is detected by absorption imaging on a charge-coupled device camera with a large quantum efficiency. The statistical uncertainty of a number measurement is dominated by the shot noise of the photoelectrons on the camera pixels and amounts to 16 atoms. We estimate the uncertainty of the total number calibration to be o1%.
To calculate the mode transformation, we use H; where c ¼ cos(Ot/2), s ¼ sin(Ot/2), and are rescaled Rabi frequencies. Below, we illustrate how the measurement of the number of particles in the m F ¼ ±1 mode after the radiofrequency coupling, gives access to the number conserving quadratures where n 0 ¼ ha w 0 a 0 i is the average number of particles in the condensate before homodyne (similarly, n A;B ¼ hâ w A;BâA;B i). In our experiment, n A þ n B ð Þ =n 0 $ 10 À 4 . We thus neglect fluctuations of the number of particles in the m F ¼ 0 mode, replacing a y 0 a 0 with its mean value n 0 En A þ n 0 þ n B ¼ N tot .
Number difference. The quadrature difference can be experimentally obtained by measuring the difference of the number of particles in the ±1 modes. From equation (2) we can directly calculate N A À N B . To the leading order in n 0 , we have Since O þ 1 and O À 1 only differ by 1.7% in our experiments, and cÕ 2 þ 1 ÀÕ 2 À 1 À Á ( 2Õ 2 AE 1 , we can simplify this equation and obtain TakingÕ À 1 %Õ þ 1 % 1, we have Finally, the mean transfer of particles from m F ¼ 0 to m F ¼ ±1 and the mean number difference is used to calculate Observing a transfer of 15% of the atoms from m F ¼ 0 to m F ¼ ±1 we deduce c 2 E0.85. To summarize, equations (6) and (8) are used to experimentally obtain p A ±p B from the measurement of the number of particles in the output modes. The quadratures x A ±x B are obtained with the same method, by applying a relative p/2 phase between the pump and side modes before homodyne detection.
Entanglement criteria for continuous variables. Criteria for identifying continuous-variable entanglement between the systems A and B, with no assumption on the quantum state of the local oscillator, have been discussed in ref. 37.
Þare the variances of quadrature sum and difference, respectively. A violation of equation (10) signals non-separability, that is, rar sep . Equation (10) generalizes the criterion of refs 26,27 that was derived for standard quadrature operators (that is, when the m F ¼ 0 mode is treated parametrically, the operator a 0 being replaced by ffiffiffiffi ffi n 0 p ).
EPR criterion. Reid's EPR criterion corresponds to a violation of the Heisenberg uncertainty relation on system B, when measurements are performed on system A. This requires the two-mode state to be non-separable and to have strong correlations between the sum and difference of position and momentum quadratures, x A ±x B and p A Ç p B . We point out that not all non-separable states fulfil Reid's criterion. The position-momentum quadratures for the B mode satisfy the commutation relation x B ; p B ½ ¼ Àiða w B a B À a w 0 a 0 Þ=n 0 . The corresponding Heisenberg uncertainty relation is D 2 x B D 2 p B ! 1 4 ð1 À nB n0 Þ 2 . Let us introduce the quantities x ext (x A ) and p ext (p A ), which are the estimate of x B and p B on system B, respectively, given the results of quadrature measurements on the system A. We then indicate as D 2 x inf B the squared deviation of the estimate from the actual value, averaged over all possible results x A , and similarly for D 2 p inf B , where P(x A , x B ) is the joint probability. Reid's criterion thus reads 37 where the bar indicates statistical average, Reid's criterion translates into a condition for the product of quadrature variances: In our case, n A =n 0 $ n B =n 0 % 10 À 4 . Therefore, corrections in equations (10) and (12) due to finite number of particles in the m F ¼ 0 are negligible. We are thus in a continuous-variable limit. We point out that the above EPR criterion-consistent with the analysis of the experimental data presented in Figs 2b and 3-uses quadrature variances with symmetric contributions from A and B. In this case the EPR threshold is 1/4. The above inequalities and entanglement criteria can be generalized (and optimized) for asymmetric contributions, see refs 3,12.
Quantum-state tomography. Here we discuss the protocol used for quantum state tomography and, very briefly, its theoretical basis. A more detailed discussion can be found in refs 31,33. We point out that our state reconstruction is performed without any assumtpions neither on the state nor on the experimental quadrature distribution, in particular we do not assume our state to be a Gaussian state.
We have collected a total of N ¼ 2,864 measurements of the quadratures x A and x B at different values of the local oscillator phase y relative to the side modes. The measurement results are binned in 2D histograms (see Fig. 2a, where the typical bin width is dx ¼ 0.25) such that we take x A,B to have a discrete spectrum.
To simplify the notation, let us indicate as x the square bin [x A , x A þ dx], [x B , x B þ dx]. Given a quantum state r (unknown here), the probability to observe a certain sequence of results (n x, y measurement in the bin x, when the phase value is set to y, with P x;y indicated as likelihood function. In equation (13), P r x; y ð Þ¼P r x y j ð ÞP y ð Þ is the joint probability, P r x y j ð Þ ¼ x h jU w y rU y x j i is the conditional probability, with ð Þ , and P(y) is the fraction of measurements done when the phase is equal to y. The maximum likelihood (ML) is the state that maximizes L n x;y È É r j À Á on the manifold of density matrices. To find the ML we use the chain of inequalities 31,33 L n x;y where a x,y are arbitrary positive numbers (a ¼ {a x, y } is the corresponding vector), f x, y ¼ n x, y /N are relative frequencies (f ¼ {f x, y } is the corresponding vector), and is a non-negative operator with largest eigenvalue l(f, a). The second inequality is saturated by taking r ¼ r ML with support on the subspace corresponding to the maximum eigenvalue of R: R(f, a)r ML ¼ l(f, a)r ML . The first inequality is a Jensen's inequality between the geometric and the arithmetic average (which follows from the concavity of the logarithm). It is saturated if and only if a x,y ¼ P(x, y), which also implies Tr[R(f, a)r ML ] ¼ 1 and thus l(f, a) ¼ 1. In conclusion, the search for the ML can be recast in the operator equation Rr ML ¼ r ML or, equivalently (since R and r ML are Hermitian operators), Numerically, this equation is solved iteratively: we start the protocol from a unit matrix r 0 ð Þ ¼ N 1 ½ and apply repetitive iterations according to equation (17), being the kth step of the algorithm, where N denotes normalization to unit trace. The convergence (which is not guaranteed in general) is checked. The method guarantees that r ML is a non-negative definite operator. In practical implementations, it is most convenient to work in the atomnumber representation and write r ¼ P ncut nA;nB ;mA;mB¼0 r nA;nB;mA;mB n A ; n B j im A ; m B h j , where n cut is a cutoff (in our case n cut E10). We use n x A;B ¼ e À x A;B =2 p 1=4 ffiffiffiffiffiffiffi ffi 2 n n ! p H n x A;B À Á , where H n is the Hermite polynomial of order n.
Simulation of ideal-state reconstruction. To check the consistency of the used tomography method, we have simulated the reconstruction of an ideal two-mode squeezed vacuum state x j i, equation (1). The simulation follows three steps: (i) we generate distributions for the quadratures x A, B at different values of the phase shift, according to the probability P x A;B y j À Á ¼ x h jU y x A;B 2 ; (ii) we generate p random quadrature data for each y value (for a total of N ¼ p Â n y , where n y is the number of y values considered). This simulates, via Monte Carlo sampling, the acquisition of experimental data. (iii) We perform a ML reconstruction, using the same numerical code and method used for the analysis of the experimental data. In Fig. 6, we plot the quantum fidelity between the reconstructed state, r ML , and the two-mode squeezed vacuum state, . When the number of measurements p per y value is increased, the fidelity converges to an asymptotic value F 1 t1. The asymptotic fidelity F 1 tends to 1 when decreasing the bin size dx.
Furthermore, to characterize the entanglement of the reconstructed state, we have evaluated the logarithmic negativity and the quantum Fisher information (QFI). The logarithmic negativity is defined as E r ½ ¼ log 2 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r ppt r T ppt q , where r ppt is the partial transpose of r ML . Mode-entanglement is obtained for 39 E[r]40. The QFI for the state projected over subspaces of a fixed number of particles n, r ¼ P n Q n r n , is given by 40 where r n ¼ P kn p kn k n j i k n h j is in diagonal form and J r is the collective pseudo-spin operator (pointing along an arbitrary direction r in the three-dimensional space). The QFI is then maximized over r, for further details see ref. 8. Particle entanglement, useful for sub-shot-noise metrology, is obtained for 40 n corresponds to the average number of particles in the two-mode state. Similarly to the results of simulations shown in Fig. 6 we obtain that, in the limit p-N and dx-0, the logarithmic negativity and QFI converge to E x j i ½ ¼ 2x=log2 and F Q x j i ½ ¼ 4sinh 2 xcosh 2 x, respectively, which are analytical values calculated for the two-mode squeezed vacuum state.
Noise model and simulation of noisy state reconstruction. The main sources of noise in our apparatus are phase fluctuations and noise of the radiofrequency coupling strength. The phase noise is assumed to have a Gaussian distribution P s y ð Þ ¼ e À y 2 =2s 2 = ffiffiffiffiffiffiffiffiffiffi 2ps 2 p and we estimate a width sE0.36. Correlated fluctuations of O þ 1 and O À 1 affect (to first order) only measurements of the quadrature sum. We have evaluated that this effect systematically increases the variance by 0.12. Both these effects are included in the solid line of Fig. 3.
We have simulated the state reconstruction in presence of these noise effects. We model the state in presence of phase noise as r pn s ð Þ ¼ where x; y j i¼ P þ 1 n¼0 e À iny tanhx ð Þ n coshx n; n j i. The systematic shift of the quadrature sum is included in the calculation of the quadrature distributions used to generate random data. Results are shown in Fig. 7. We see that statistical noise (that is, the limited sample size) and phase noise are responsible for the appearance of off-diagonal terms, very similar to the ones observed in Fig. 4. Note that phase noise alone is not responsible for the appearance of off-diagonal terms in the density matrix. This can be seen by rewriting equation (19) as r pn s ð Þ ¼ P þ 1 n;m¼0P s n À m ð Þ tanhx ð Þ n þ m cosh 2 x n; n j i m; m h j, wherẽ P s n À m ð Þ¼ R p À p dyP s y ð Þe i n À m ð Þ y . Figure 7 shows a slight asymmetry of the reconstructed state due to the systematic shift of the variance sum: this effect is also observed in Fig. 4. The quantitative agreement between the simulated density matrix r sim and the experimental density matrix r exp is excellent, with a quantum fidelity Tr ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi ffiffiffiffiffiffiffiffi r sim p r exp ffiffiffiffiffiffiffiffi r sim p p h i % 0:9. Here we show the quantum fidelity between the reconstructed state r ML and the ideal two-mode squeezed vacuum state, as a function of the number of measurements p per y value (here we consider n y ¼ 29 phase values). The dark blue circles are obtained for a binning of dx ¼ 0.25, the light blue diamonds for dx ¼ 0.1. The dashed lines are a guide to the eye, the error bars (one s.d.) are obtained with a bootstrap method.