Spatially correlated photonic qutrit pairs using pump beam modulation technique

Higher dimensional quantum systems have a very important role to play in quantum information, computation as well as communication. While the polarization degree of freedom of the photon is a common choice for many studies, it is restricted to only two orthogonal states, hence qubits for manipulation. In this paper, we theoretically model as well as experimentally verify a novel scheme of approximating photonic qutrits by modulating the pump beam in a spontaneous parametric down conversion process using a three-slit aperture. The emerging bi-photon fields behave like qutrits and are found to be highly correlated in the spatial degree of freedom and effectively represent spatially correlated qutrits with a Pearson coefficient as high as 0.9. In principle, this system provides us a scalable architecture for generating and experimenting with higher dimensional correlated qudits.


INTRODUCTION
The ability to generate and detect correlated photon pairs by means of spontaneous parametric down conversion (SPDC) offers lots of experimental possibilities ranging from fundamental tests of quantum theory [1] to practical applications in quantum information [2] and quantum communication [3]. The influence of the pump beam spatial mode on the characteristics of the resulting photon pairs has been analyzed in [4] and in terms of optimization of the collecting of the photon pairs in single mode fibres in [5]. Recently it was analysed how the pump beam spatial mode influences the spatial characteristics of the resulting photon pairs [3]. In Ref. [6] the authors analyze the case where the pump beam is prepared in Bessel-Gauss mode. In turn, Pugh et al [3] show the way to control the pump beam for the application for long distance quantum communication with satellites.
Along these lines, there is an extensive research towards encoding of quantum information in a single photon's degree of freedom. Spatial degree of freedom offers a variety of possibilities such as orbital angular momentum [7][8][9][10] or spatial qudits [11][12][13][14]. The latter framework proved to be useful when testing fundamentals of quantum theory [15][16][17]. In earlier work [16], it has been shown that a triple slit aperture placed in the path of a photon leads to the generation of a spatial qutrit. Note that in practice the transmission of the system, which is fundamentally very low, limits the robustness in terms of counting rates. The use of calcite beam displacer can reduce this problem for a qubit implementation [18,19]. However the scaling towards higher dimensional states is challenging within this framework.
Here we study the case where the pump beam is prepared such that its spatial mode resembles three slits, and the spatial structure is carried through to the resultant signal and idler photons. We investigate the correlations between the signal and idler photons in the spatial domain and effectively verify our spatially correlated qutrits have a Pearson coefficient as high as 0.9.
In this work, we show that the structure of the pump field is preserved in the process of SPDC and we find very convincing match between theoretical predictions and experimental results for the correlation between the signal and the idler photons. We also lay down the recipe for the optimal choice of experimental parameters which can be exploited in future architecture to maximize the correlations.

THREE-SLIT-BASED QUTRITS
Let us consider a setup depicted in Fig. 1 consisting of a Type-1 non-linear crystal. We are using the collinear geometry for our down conversion process. The pump beam spatial mode is prepared by transferring a Gaussian beam through a set of three slits (centered at the middle slit) and imaging the result on the crystal using a lens in the 2f -2f configuration. The coordinate system orientation is chosen such that the propagation is along z axis, and x and y axes are parallel to the shorter and longer length of the slits, respectively.
We model the pump beam as three box functions weighted by a Gaussian as the pump profile in the crystal to generate the bi-photon qutrits and their resultant correlations so that the simulations ended in finite time.
We followed [20] to solve for phase matching in the Type I SPDC process. The crystal that we have used is BBO with a cut Fig. 1. Schematic of the experimental set-up. Horizontal pump beam is made incident on a triple slit aperture. Lens L 1 is used to transfer the image of the pump beam on the Type-1 BBO crystal. After appropriate filtering of the blue pump beam, another lens L 2 is used to transfer the signal and idler spatial profiles to actuated detectors placed on either side of a beam splitter. The spatial profiles of the signal and idler photons are measured using detectors D 1 and D 2 and the spatial correlation is measured using an appropriate co-incidence logic unit.
angle of 29.3°for non collinear SPDC process. For type I SPDC, we define the phase mismatch as Here, ω p , ω s and ω i refer to the frequency of the pump, signal and idler photons respectively, and n e and n o refer to extraordinary and ordinary refractive indices. The pump wave vector inside the crystal depends on the angle α it makes with the optic axis of the crystal.
The SPDC process considered here is e → oo. For collinear degenerate SPDC, we expect ∆ k = 0. Using parameters of refractive index from BBO Sellmeier equations, we have the collinear phase matching at 28.81°.
Given a phase mismatch ∆ k, we assign the intensity weight to it as where L x is the transverse length of the crystal and L z is the thickness of the crystal along beam propagation direction. We have ignored the transverse length along the height i.e. y-axis for ease of computation. Following [20,21], we can find the Hamiltonian for the SPDC process where the pump is treated classically and the signal and idler are treated as perturbations. For degenerate downconversion ω p = 2ω s = 2ω i , and the probability amplitude to get a coincidence for degenerate photons generated from collinear SPDC will be: Here, A p (x, y, z) is the amplitude of the pump as a function of space. z is assumed to be the direction of propagation of pump, k p , k s and k i are the pump, signal and idler wave vectors respectively. The integral over dk s depends on the angle that the detector subtends to the crystal. For point detectors we have a single value of k s,i . However, the limits of integration for finite sized detectors for the k s will be of the form k µ s ± k ∆ s . Here k µ s is determined by the position of the detectors and k ∆ s is determined by the size of the detectors.
We now attempt to compute A numerically. However, integration of the complex exponential is numerically cumbersome and does not converge. Hence we assume that k p does not depend on position inside the crystal, so that we can integrate the exponential independently and replace it with a sinc function. We use the approximation Here Π i=x,y,z sinc(∆k i (r)L i ) should be interpreted as the weight associated with a point classical pump and here the phase matching is weighted by the length of the crystal. Then we consider that the pump is composed of several such point pumps and we assume that different points in the crystal are independent as far as wave vectors are considered. Note that phase mismatch depends on ∆r, the difference between crystal and detector positions (and not just on the crystal position).
In the simulation code, we integrate over two coordinates (x and z). We are interested in the transverse correlation (along x). The height of the slit is assumed to be infinity, pump wave vector is assumed to be along z-direction only, assuming a thin crystal and each point in the crystal is assumed to have SPDC independently with the weight sinc(∆k i (∆r)L i ).
We follow the following steps to compute |A(k s , k i )| 2 .
• Compute ∆k for a pair of signal and idler points and for a given point on the crystal.
• Compute the point to point correlation as the product of intensity of the pump at the given crystal point |A p | 2 and the weight function |Π i=x,y,z sinc(∆k i (r)L i )| 2 • Integrate over the signal/idler position about the mean detector position and with range equal to the size of detector.
• Repeat the above for different idler points to get the spatial correlation between a signal position and the entire idler profile.
• Sum the idler profile computed for different signal positions to get the signal profile. We moved the signal detector with a step size of 3µm to attain convergence in this code.
The above simulations were carried out for different choice of simulation parameters to obtain an optimal set of parameters that gives us a high value of Pearson correlation coefficient ρ [22]. In these simulations, we have assumed point detectors to enable faster simulations. In the result of the numerical simulation we get an approximation of the probability distribution of the coincidence detection of signal photon at position x 1 and idler photon at position x 2 . Based on the numerical data one can easily compute correlation coefficient as a ratio of co-variance over a product of respective variances for signal and idler photon detector positions: ρ =< x 1 x 2 > / < x 2 1 >< x 2 2 >. We perform simulations for both Type I and Type II SPDC using the same crystal, lens and slit parameters and found that the Pearson coefficient for Type II process is 0.985(2) whereas for Type I, it is 0.966 (2). We selected an experiment using a Type I crystal as we were concerned that transverse walk off which could happen to one of the two orthogonal polarization exhibited by signal and idler photons could mask the otherwise high correlations that have been predicted from theory.
We have chosen the three slit system to have slit width of 30µm, inter-slit distance of 100µm, slit height of 300µm. For this, we have chosen an incident pump beam with a Gaussian RMS width of 300µm. This is an optimal choice as we wish to strike a balance between throughput and uniform distribution. For a wider beam width, the number of singles and coincidence counts is expected to drop significantly whereas a smaller beam width would result in the three slits not all being centered close to the peak of the Gaussian. We have done the simulations for each of the slit parameters keeping other parameters constant and ascertained that the Pearson coefficient is high for our final choice. Table 1 shows the comparison between different slit and crystal parameters in terms of the Pearson coefficient. Further details of these parameter optimization simulations are discussed in the Appendix.

EXPERIMENT
The Type I BBO crystal is cut for non-collinear phase matching at 405 nm to 2 × 810 nm at 29.3°which translates to collinear phase matching at 28.8°. The parameter of the slits and pump beam are as mentioned above and the focal length of the lens performing image transfer from slit plane to center of non linear crystal is 146mm for blue incident beam. The focal length of the lens performing the image transfer of the signal and idler photons to the detector plane is 150mm for IR wavelength.
The detectors D1 (D2) are mounted on motorized stages (actuated with stepper motorized actuators ZST225B and ZST213 from Thorlabs respectively) allowing control of their position in the plane orthogonal to the propagation of signal (idler) photon. The accuracy of the motors is sub micron level. We measured single counts and coincidences of detection at both detectors with a coincidence time window of 1024 ps using FPGA electronics (UQD LOGIC-16). We scanned the characteristic range in the direction orthogonal to the slits' longer dimension with both the detectors. While one of the detectors (D1) scanned the signal photon spatial profile, D2 scanned the idler photon profile. By keeping D1 fixed at different positions of the signal profile, we scanned the detector D2 to yield correlations between the signal position and the entire idler profile. We decided on 13 fixed detector positions for sufficient statistics. These 13 positions correspond to peaks, dips and asymmetrically chosen slope positions to give us maximum information as per a Nyquist sampling criterion. For each D1 fixed position, D2 was moved with a step size of 10µm. The data acquisition time was 180 sec at each point, resulting in one complete D2 run taking close to three hours. We repeat the measurement 5 times for better averaging. Thus, for each fixed signal position, we take close to 15 hours to generate the idler profile and resultant correlations.

RESULTS
The result of coincidence counting is plotted in Fig. 2(a). The experimental correlations between signal and idler positions are appropriately captured by the measurement of coincidences. For the measure of correlation we take the Pearson's coefficient.
Here the random variables are detectors' positions, x 1 and x 2 and the probability distribution is estimated directly by our measurement. The probability of getting a coincidence detection at positions x 1 and x 2 is proportional to the coincidence counts measured and depicted in Fig. 2(a). The estimate of the Pearson's correlation coefficient is 0.9(2). We estimated the uncertainty of the coefficient by simulating 10 5 probability distributions based on the measured statistics and assuming Poisson statistics of the counts.
Red dots in Fig. 2(b), (c) and (d) correspond to coincidences, R c , measured as the idler detector is moved while the signal detector is kept fixed at first peak (slit A), second peak (slit B) and third peak (slit C) respectively. Error bars for both position and number uncertainty have been included. The blue lines represent the theoretically simulated correlation profiles.
There are three ways in which one can generate the commensurate theory graphs with respect to the modeling of the pump profile at the crystal. Two methods, discussed in the Appendix , involve using a detailed image transfer formalism or using three box functions weighted by a Gaussian. In order to capture the correlations due to the actual pump beam profile that has been transferred to the crystal, we measured the pump beam profile that is transferred to the position of the center of the crystal using a lens in 2f-2f configuration and used this profile itself to generate the signal and idler profiles as well as calculating the spatial correlations between the two profiles. Slight difference in magnification between the experimental and theoretically simulated profiles has been accounted for in the simulations. An example of a pump beam profile at the crystal position is given in the supplementary material. The use of the experimentally measured pump profile to generate the correlation has the advantage that any non-idealness that may exist in the experiment in terms of alignment or otherwise would then be captured in the theory and the comparison would not suffer from comparing experimental data with ideal theory conditions. While Fig. 2(a) shows R c for all measured detector positions, in Fig. 2(b), (c) and (d), we have chosen to highlight the correlations at the peak positions in a 2-D format to enable representation of the error bars in terms of position and number uncertainty and also indicate the extent of overlap between theoretical predictions and experimental results. Fig. 2(e) shows the comparison between the single photon counts R s measured as a function of detector position for the signal photon and the theoretically generated profile. Experimental and theoretically generated R c and R s have been appropriately normalized by their respective maxima.  Actual measurements indicate that the slit widths are not exactly 30µm and 100µm respectively as taken in simulations. Slit C is the widest at about 36µm, slit B the thinnest at around 29µm while slit A is in between at around 33µm. This is reflected in the singles profile where one can see that slit C has maximum counts while B has minimum. The data has been normalized with the maximum (here C) counts. As the coincidence counts are in principle proportional to the singles, when fixed detector is at position of slit C, we do measure maximum coincidences there. However, in measurement of coincidences, a couple more factors also come in. While the R S is generated with a single detector (D1) motion controlled by one actuator, R C involves both the detectors controlled by two different actuators. Slight nonangular misalignment in the lens positioning and/or difference in behavior between the two actuators as well as slightly offset detector positioning affects the R C data and this is reflected in figures b) and c) where slightly less coincidences are measured when fixed detector is at A than at B in spite of B being thinner than A by 4µm. Point to be noted is that theoretical simulations assume all slits to have equal widths and the slit B being centered on a perfectly Gaussian incident laser beam. Thus, in theory, slit B always has maximum R S as well as R C . Thus difference in slit dimensions as well as small lens and detector misalignment cause some difference between the normalized theory and experimental graphs. However, the spatial correlation exhibited by our scheme is emphatic and as high as 0.9 of a possible maximum of 1.0 and is almost independent of the small inherent and unavoidable experimental non-idealness.

DISCUSSION
We show a novel approach of generating photon pairs in SPDC which have intrinsic spatial correlation. We demonstrate a configuration which approximates the behavior of two threedimensional quantum systems. Our approach can be used for implementation quantum information protocols which require higher dimensional quantum systems. While it has previously been shown how to encode a qutrit using a system of three slits similarly in the Ref. [16], our approach provides direct access to two correlated qutrits by modulating the laser profile of the SPDC setup.
We demonstrated a novel and simple approach to generate higher dimensional entanglement paving a possible route towards quantum communication and information processing and fundamental entanglement studies, using higher dimensional entangled photon states.

A. Image transfer formalism
There are certain assumptions which go into the image transfer formulation. All calculations are done in two dimensions which are the beam propagation direction "z" and the transverse "x" direction. We are assuming that the slit height can be approximated to be infinitely long compared to the slit width. Scalar field paraxial approximations for a thin lens are used. The center of the crystal is assumed to be perfectly at twice the focal length from the lens.
The system of lenses and mirrors transfers the image of the slits to the crystal. The center of the crystal is at 2f from the thin lens. Its transfer function is given by: Here (x 0 , z 0 ) are the coordinates of the object whose image is to be transferred to a location (x i , z i ) and λ is the wavelength of the incident beam. Thus, in order to transfer the scalar field U(x 0 ; z 0 ) to U(x i ; z i ), we would use the following transformation equation and the scalar field after the lens would then be given as: It yields the image transferred at crystal position. The integration in Eq. (5) is done analytically whereas the final integration in Eq. (6) is evaluated by numerical means using Mathematica 11. Thus instead of using uniform top hat functions as representations of the slit profiles, we have used the lens transfer formulation to transfer the image of the slits to the center of the crystal when the slits are illuminated by a Gaussian.

B. Simulations to determine the optimal set of parameters for high resultant spatial correlations
First, we varied the slit width keeping other parameters constant. For slit widths ranging from 5µm to 40µm, the Pearson coefficient remained around 0.96, see Fig. 3 red curve, which indicates that for these conditions, choosing a slit width in the above range should be sufficient. We decided to choose 30µm. Next, we varied the inter-slit distance from 50µm to 200µm. The Pearson coefficient is found to increase with increasing distance between the slits kept at a constant slit width as can be seen in Fig. 3 blue curve. The conclusion is that when the slits are more separated, the overlap between them goes down, as a result of which the point to point correlation increases. However, a 200µm interslit distance would entail a much bigger incident pump beam which could again lead to less throughput so we decided to choose the 100µm interslit distance as a compromise between throughput and correlation coefficient. When slit width w is varied, inter-slit distance is kept fixed at 100µm whereas when interslit distance is varied, slit width w is kept constant at 30µm. The crystal length L Z has been kept fixed at 10mm for both these simulations.
Before choosing the crystal length, we also simulated for different crystal lengths keeping slit and beam parameters constant. If we consider the intensity weight function associated with phase matching, we have along the longitudinal axis sinc(∆k z × L z /2). As L z increases the momentum uncertainty of the photon pairs inside the crystal also increases. As a result the correlation between the two down-converted photons decreases. This was substantiated by the simulations which showed a steady decrease in Pearson coefficient as crystal length was increased from 5mm to 100mm. On the other hand the thinner crystals give lower pair production rate. Thus, we selected a crystal length which is not too short but also yields an expected high correlation coefficient i.e. 10mm. The transverse length of the crystal needs to be larger than the extent of the transverse pump profile and was chosen to be 5mm.  When the lens is used for image transfer experimentally, a magnification is introduced in the system, which has also been incorporated in theory. While the y-axis denotes the crystal length along beam propagation direction, the x-axis denotes the image along the transverse crystal direction. The figure on the right is the experimentally measured image of the modulated pump at the position corresponding to center of the crystal using a lens.