Tuning between photon-number and quadrature measurements with weak-field homodyne detection

TL;DR

This paper demonstrates an experimental weak-field homodyne detector capable of continuously tuning between photon-number and quadrature measurements, enhancing quantum state analysis and engineering.

Contribution

The authors experimentally implement a tunable weak-field homodyne detector that bridges photon-number and quadrature measurements, enabling optimized quantum property testing and state preparation.

Findings

Discrete difference statistics converge to quadrature distribution with increased coherent amplitude

Controlled photon-number distribution of heralded states demonstrated

Detector enables continuous measurement operator tuning

Abstract

Variable measurement operators enable the optimization of strategies for testing quantum properties and the preparation of a range of quantum states. Here, we experimentally implement a weak-field homodyne detector that can continuously tune between measuring photon numbers and field quadratures. We combine a quantum signal with a coherent state on a balanced beam splitter and detect light at both output ports using photon-number-resolving transition edge sensors. We observe that the discrete difference statistics converge to the quadrature distribution of the signal as we increase the coherent state amplitude. Moreover, in a proof-of-principle demonstration of state engineering, we show the ability to control the photon-number distribution of a state that is heralded using our weak-field homodyne detector.

Tables

Table 1. Table 1: Measured parameters.

Parameter	Value
${\eta }_h$	0.395 $\pm$ 0.002
${\eta }_c$	0.274 $\pm$ 0.001
${\eta }_d$	0.352 $\pm$ 0.002
$\left|\lambda \right|$	0.797 $\pm$ 0.001
$ℳ$	0.800 $\pm$ 0.060

Equations

$$\Delta \hat{n}=i({\hat{a}}^{†}\hat{b}-\hat{a}{\hat{b}}^{†}),$$

$$\Delta {\hat{n}}_{\mathrm{classical}}=i|\alpha |\hat{X}(\theta ),$$

$$\begin{array}{rl}{P}^{(j,\alpha )}(\Delta n)& =\sum _{m=\max (0,\Delta n)}^{\infty }\frac{e^{-|\alpha {|}^2}j!|\alpha {|}^{2(2m-\Delta n-j)}}{2^{2m-\Delta n}m!(m-\Delta n)!}\\ & \times {\left|\sum _{k=0}^j\left(\genfrac{}{}{0ex}{}{m}{m+k-j}\right)\left(\genfrac{}{}{0ex}{}{m-\Delta n}{k}\right)(-1{)}^k\right|}^2.\end{array}$$

$$P_{\mathrm{classical}}^{(j,\alpha )}(\Delta n)=\frac{e^{-\Delta n^2/2|\alpha {|}^2}}{\sqrt{2\pi }{2}^{j}j!|\alpha |}{\left|H_j\left(\Delta n/\sqrt{2}\alpha \right)\right|}^2,$$

$$|{\Psi }_{ha}⟩=\sqrt{1-|\lambda {|}^2}\sum _{f=0}^{\infty }{\lambda }^f|f_h,f_a⟩,$$

$${\mathcal{S}}_{\mathrm{classical}}=\frac{1}{\nu }\sum _{\Delta n}|P_{\exp }^{(j,\alpha )}(\Delta n)-P_{\mathrm{classical}}^{(j,\alpha )}(\Delta n){|}^2,$$

$$\begin{array}{rl}{P}^{(j,\alpha )}(m,n)& ={\left|⟨m_c,n_d⟩j_a,{\alpha }_b\right|}^2\\ & =\frac{1}{m!n!}{\left|⟨0_c,0_d⟩{\hat{c}}^m{\hat{d}}^n{j}_a,{\alpha }_b\right|}^2,\end{array}$$

$$\begin{array}{r}\hat{c}=(\hat{a}+i\hat{b})/\sqrt{2},\\ \hat{d}=(\hat{b}+i\hat{a})/\sqrt{2},\end{array}$$

$$\begin{array}{rl}{P}^{(j,\alpha )}(m,n)& =\frac{1}{2^{m+n}m!n!}{\left|⟨0_c,0_d⟩(\hat{a}+i\hat{b}{)}^m(\hat{b}+i\hat{a}{)}^n{j}_a,{\alpha }_b\right|}^2\\ & =\frac{1}{2^{m+n}m!n!}{\left|⟨0_c,0_d⟩\sum _{l=0}^m\sum _{k=0}^n\left(\genfrac{}{}{0ex}{}{m}{l}\right)\left(\genfrac{}{}{0ex}{}{n}{k}\right)i^{l+k}{\hat{a}}^{m-l+k}{\hat{b}}^{l+n-k}{j}_a,{\alpha }_b\right|}^2.\end{array}$$

$$P^{(j,\alpha )}(m,n)=\frac{e^{-|\alpha {|}^2}j!}{2^{m+n}m!n!}|\alpha {|}^{2(m+n-j)}{\left|\sum _{k=0}^j\left(\genfrac{}{}{0ex}{}{m}{m+k-j}\right)\left(\genfrac{}{}{0ex}{}{n}{k}\right)(-1{)}^k\right|}^2,$$

$$P_{{\eta }_c,{\eta }_d}^{(j,\alpha )}(m,n)=\sum _{x=m}^{\infty }\sum _{y=n}^{\infty }\left(\genfrac{}{}{0ex}{}{x}{m}\right)\left(\genfrac{}{}{0ex}{}{y}{n}\right){\eta }_c^m{\eta }_d^n(1-{\eta }_c{)}^{x-m}(1-{\eta }_d{)}^{y-n}{P}^{(j,\alpha )}(x,y).$$

$$|{\alpha }_b⟩\to |\sqrt{\mathcal{M}}{\alpha }_{||}⟩\otimes |\sqrt{1-\mathcal{M}}{\alpha }_{\perp }⟩,$$

$$P_{{\eta }_c,{\eta }_d,\mathcal{M}}^{(j,\alpha )}(m,n)=P_{{\eta }_c,{\eta }_d}^{(j,\sqrt{\mathcal{M}}\alpha )}(m,n)\ast P_{{\eta }_c,{\eta }_d}^{(0,\sqrt{1-\mathcal{M}}\alpha )}(m,n),$$

$${\hat{\rho }}_a^{(j)}=\mathcal{X}{\mathrm{Tr}}_h\left({\hat{\Pi }}_h|{\Psi }_{ha}⟩⟨{\Psi }_{ha}|\right)=\mathcal{X}\sum _{f=j}^{\infty }\left(\genfrac{}{}{0ex}{}{f}{j}\right){\eta }_h^j(1-{\eta }_h{)}^{f-j}|\lambda {|}^{2f}|f⟩⟨f|.$$

$$P_{{\eta }_c,{\eta }_d,{\eta }_h,\mathcal{M}}^{(j,\alpha )}(m,n)=\mathcal{X}\sum _{f=j}^{\infty }\left(\genfrac{}{}{0ex}{}{f}{j}\right){\eta }_h^j(1-{\eta }_h{)}^{f-j}|\lambda {|}^{2f}{P}_{{\eta }_c,{\eta }_d,\mathcal{M}}^{(f,\alpha )}(m,n).$$

$$\boxed{P_{{\eta }_c,{\eta }_d,{\eta }_h,\mathcal{M}}^{(j,\alpha )}(\Delta n)=\sum _{m=\max (0,\Delta n)}^{\infty }{P}_{{\eta }_c,{\eta }_d,{\eta }_h,\mathcal{M}}^{(j,\alpha )}(m,m-\Delta n).}$$

$$p^{(j,\alpha )}(\Delta n)=\frac{1}{\sqrt{2\pi }|\alpha |}\frac{1}{2^{j}j!}{\left|H_j\left(\frac{\Delta n}{\sqrt{2}\alpha }\right)e^{-\Delta n^2/4{\alpha }^2}\right|}^2,$$

$$p_{{\eta }_c,{\eta }_d}^{(j,\alpha )}(\Delta n)=p^{(j,g)}(\Delta \tilde{n})\ast \frac{1}{\sqrt{2\pi }\sigma }{e}^{-\Delta {\tilde{n}}^2/2{\sigma }^2},$$

$$p_{{\eta }_c,{\eta }_d,\mathcal{M}}^{(j,\alpha )}(\Delta n)=p_{{\eta }_c,{\eta }_d}^{(j,\sqrt{\mathcal{M}}\alpha )}(\Delta n)\ast p_{{\eta }_c,{\eta }_d}^{(0,\sqrt{1-\mathcal{M}}\alpha )}(\Delta n).$$

$$\boxed{p_{{\eta }_c,{\eta }_d,{\eta }_h,\mathcal{M}}^{(j,\alpha )}(\Delta n)=\mathcal{X}\sum _{f=j}^{\infty }\left(\genfrac{}{}{0ex}{}{f}{j}\right){\eta }_h^j(1-{\eta }_h{)}^{f-j}|\lambda {|}^{2f}{p}_{{\eta }_c,{\eta }_d,\mathcal{M}}^{(f,\alpha )}(\Delta n).}$$

$$\begin{array}{rl}{M}_{x,y}^{(j)}& =\frac{2}{E_{tot}^{(j)}}\sum _{k,l}\left({\delta }_{k,x}{\delta }_{l,y}+{\delta }_{k,y}{\delta }_{l,x}\right)E^{(j)}(k,l)\\ & +\frac{1}{[E_{tot}^{(j)}{]}^2}\left(\sum _{k,l}({\delta }_{k,x}+{\delta }_{l,x})E^{(j)}(k,l)\right)\\ & \times \left(\sum _{k,l}({\delta }_{k,y}+{\delta }_{l,y})E^{(j)}(k,l)\right),\end{array}$$

$$P^{(j)}(n)=\frac{1}{E_{tot}^{(j)}}\sum _{k+l=n}{E}^{(j)}(k,l).$$

$$g^{(2)(j)}=\frac{{\sum }_n(n^2-n)P^{(j)}(n)}{{\left({\sum }_{n}n{P}^{(j)}(n)\right)}^2}.$$

Full text

Introduction. Counting particles is at the core of all measurements in quantum physics. Many detection strategies used in, e.g., atomic Gerlach and Stern (1922), nuclear Rutherford and Geiger (1908), or particle Glaser (1952) physics rely on this fundamental concept. In optics, recent technological developments have enabled counting the number of quanta in light using photon-number-resolving detectors Hadfield (2009). Since these detectors determine the exact number of photons in an optical signal, information about the complementary property, i.e., the wave nature of the signal, is lost. In order to measure its wavelike properties, it is necessary to mix the signal with a phase-reference field. This is usually achieved using homodyne detection, which uses a classical phase reference to measure the amplitude and phase of the signal Yuen and Shapiro (1980). However, because the photodetectors used in arrangements like homodyne often have significant electronic noise, the phase reference must be quite strong, typically containing billions of photons Hansen et al. (2001).

Weak-field homodyne detection (WFHD) encompasses both the strategies mentioned above by first mixing a quantum signal with a weak phase reference and then detecting the resulting fields using photon-number-resolving detectors. When the phase reference is comparable in strength to the signal, it cannot be treated as a classical field Braunstein (1990); Vogel and Grabow (1993); Raymer et al. (1995); Banaszek and Wódkiewicz (1997); Tyc and Sanders (2004). In this regime, WFHD spans the region between the complementary detection strategies of direct photon counting and homodyne. This enables WFHD to reveal fundamental properties of light such as the presence of nonclassical correlations that would be hidden from conventional detection strategies Grangier et al. (1988); Banaszek and Wódkiewicz (1998); Kuzmich et al. (2000).

The versatility of WFHD also makes it a promising tool for optical quantum information processing. On the one hand, homodyne detection performs phase-sensitive Gaussian measurements by projecting light onto continuous quadrature states. This is the standard measurement used in quantum state tomography Lvovsky and Raymer (2009) and protocols that encode quantum information in continuous-variable degrees of freedom of light Braunstein and van Loock (2005). On the other hand, photon counting performs phase-insensitive non-Gaussian measurements by projecting light onto discrete photon-number states. This measurement is key when information is encoded in discrete variables (e.g., polarization, spatial mode) and also provides the non-Gaussian resource required for universal continuous-variable quantum computing Menicucci et al. (2006). In principle, WFHD can be reconfigured to realize both of these measurements simply by controlling the strength of the phase reference. Moreover, WFHD has the unique ability to perform phase-sensitive non-Gaussian measurements since it can access the region between homodyne and photon counting Puentes et al. (2009); Zhang et al. (2012); Donati et al. (2014). Such measurements provide a powerful state preparation and characterization tool, especially in hybrid discrete- and continuous-variable protocols Andersen et al. (2015). Although the tunability of WFHD has been investigated theoretically Braunstein (1990); Vogel and Grabow (1993); Raymer et al. (1995); Banaszek and Wódkiewicz (1997); Tyc and Sanders (2004); Puentes et al. (2009), it has not yet been demonstrated experimentally due to limitations in the efficiency, noise, and dynamic range of the photon-number-resolving detectors used in previous implementations of WFHD Puentes et al. (2009); Zhang et al. (2012); Donati et al. (2014); Allevi et al. (2017); Olivares et al. (2019).

In this paper, we overcome these limitations by using state-of-the-art transition edge sensors Lita et al. (2008). We experimentally demonstrate the ability to tune between photon-number and quadrature measurements of quantum signals. We combine a heralded photon-number state (our signal) with a weak coherent state (our phase-reference) on a beam splitter. We then measure the photon-number difference between both output ports using transition edge sensors. We observe that the discrete difference statistics converge to the quadrature distribution of the signal as we increase the coherent state strength. Moreover, in a proof-of-principle demonstration of state engineering, we show the ability to control the photon-number distribution of a state that is heralded using WFHD.

Theory. In homodyne detection, a signal $\hat{\rho}_{a}$ and coherent state of amplitude $\alpha$, $\ket{\alpha_{b}}$, are mixed on a balanced beam splitter (see Fig. 1; subscripts label modes). The difference in photon number $\Delta\hat{n}=\hat{n}_{c}-\hat{n}_{d}$ between both outputs of the beam splitter is given by

[TABLE]

where $\hat{a}$ and $\hat{b}$ are the input annihilation operators in modes $a$ and $b$, respectively. When $\alpha=0$, measuring the observable $\Delta\hat{n}$ with ideal number-resolving detectors projects $\hat{\rho}_{a}$ onto non-Gaussian photon-number states. In the limit of large $|\alpha|$, one can invoke the classical field approximation $\hat{b}\rightarrow|\alpha|e^{i\theta}$, where $|\alpha|$ and $\theta$ give the strength and phase of the coherent state, respectively Tyc and Sanders (2004). With this approximation, Eq. (1) becomes:

[TABLE]

where $\hat{X}(\theta)=\hat{a}^{\dagger}e^{i\theta}-\hat{a}e^{-i\theta}$ is the quadrature operator. This is the usual treatment of homodyne detection found in textbooks Leonhardt (1997). In contrast to $\Delta\hat{n}$, measuring $\Delta\hat{n}_{\mathrm{classical}}$ projects $\hat{\rho}_{a}$ onto Gaussian quadrature states that are eigenstates of $\hat{X}(\theta)$ Yuen and Shapiro (1980).

The transition from $\Delta\hat{n}$ to $\Delta\hat{n}_{\mathrm{classical}}$ implies the ability to tune between photon-number and quadrature measurements by varying $|\alpha|$. To exploit this tunability in an experiment, one must be able to make $|\alpha|$ sufficiently large to ensure that the classical field approximation is valid. References Tyc and Sanders (2004); Skotiniotis et al. (2017) estimate the required $|\alpha|$ to be $|\alpha|\gg N$, where $N$ is the average number of photons in $\hat{\rho}_{a}$. This condition ensures that quantum fluctuations in the coherent state are larger than $N$. If these fluctuations were smaller than $N$, i.e. $|\alpha|<N$, a measurement of $\Delta\hat{n}$ would reveal some information about $\hat{a}^{\dagger}\hat{a}$, the photon number of $\hat{\rho}_{a}$. This information gained about $\hat{a}^{\dagger}\hat{a}$ comes at the cost of disturbing the measurement of $\hat{X}(\theta)$ since these two operators do not commute.

We investigate this transition in the particular case of a signal in a photon-number state, i.e. $\hat{\rho}_{a}=\ket{j_{a}}\bra{j_{a}}$ Braunstein (1990); Vogel and Grabow (1993); Raymer et al. (1995). Given some $|\alpha|$, the probability of measuring a difference of $\Delta n$ photons at the output of the beam splitter is given by (see the Supplemental Material SM for a derivation)

[TABLE]

Based on the arguments made above we expect that for $|\alpha|\gg j$, $P^{(j,\alpha)}(\Delta n)$ should converge to the measurement statistics of the observable $\Delta\hat{n}_{\mathrm{classical}}$, that is

[TABLE]

where $H_{j}$ is the Hermite polynomial of degree $j$. Notably, $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ is simply the quadrature distribution of a photon-number state $\ket{j}$ scaled by $|\alpha|$ Leonhardt (1997). In the Supplemental Material SM , we generalize Eqs. (3) and (4) to include imperfections that arise in our experiment such as optical losses, detection inefficiency, and mode mismatch between $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$.

Experimental overview. Experimentally investigating the transition from $P^{(j,\alpha)}(\Delta n)$ to $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ requires photon-number-resolving detectors with a large dynamic range. The standard approach for photon-number-resolved measurements relies on multiplexing binary “click” detectors Hadfield (2009). However, this approach is not practical for our experiment since the required number of click detectors scales prohibitively with dynamic range Fitch et al. (2003). Instead, we use superconducting transition edge sensors (TESs) which can resolve up to $\sim 20$ photons with $>$ 90% efficiency Lita et al. (2008); Humphreys et al. (2015). TESs achieve photon-number resolution by acting as a bolometer, i.e. directly measuring the energy of the absorbed light. More information on TES operation and readout is provided in the Supplemental Material SM .

To prepare our signal $\hat{\rho}_{a}$, we couple femtosecond laser pulses into a periodically-poled potassium titanyl phosphate (ppKTP) waveguide. The waveguide produces an approximately spectrally decorrelated two-mode squeezed vacuum state,

[TABLE]

via type-II parametric down conversion ($\lambda=\tanh{(r)}$ where $r$ is the squeezing parameter) Eckstein et al. (2011). Mode $h$ is sent to a heralding TES detector. By post-selecting on events where the heralding TES detects $j$ photons, we prepare our signal $\hat{\rho}^{(j)}_{a}$ in mode $a$ (see Fig. 2). Due to losses and imperfect detection efficiency ($\eta_{h}=0.395\pm 0.002$), we do not herald a pure photon-number state $\ket{j_{a}}$. Rather, we herald a signal $\hat{\rho}^{(j)}_{a}$ which is a statistical mixture of photon-number states $\geq j$. Despite the losses, $\hat{\rho}^{(j)}_{a}$ is still a nonclassical signal. Namely, we measure that $\hat{\rho}^{(j)}_{a}$ has both sub-Poissonian Leonhardt (1997) and submultinomial Sperling et al. (2017a) photon-number statistics when $j\geq 1$ (see the Supplemental Material SM ).

To obtain a coherent state that is mode-matched to our signal, we prepare $\ket{\alpha_{b}}$ using a second ppKTP waveguide. Unlike in the first waveguide, we stimulate the parametric down-conversion process by seeding it with light from a continuous-wave laser. We use a strong seed to ensure that the spontaneous process is negligibly small relative to the stimulated process Aichele et al. (2002); Liscidini and Sipe (2013). As such, we generate a coherent state of light in the polarization mode orthogonal to the seed, $\ket{\alpha_{b}}$. We then attenuate $\ket{\alpha_{b}}$ to the single-photon level with neutral density filters. Since $\hat{\rho}^{(j)}_{a}$ has no defined phase relative to $\ket{\alpha_{b}}$, $P^{(j,\alpha)}(\Delta n)$ and $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ do not depend on $\theta$ and hence we do not require control of $\theta$.

Both $\hat{\rho}^{(j)}_{a}$ and $\ket{\alpha_{b}}$ are coupled into fibers and temporally overlapped in a balanced fiber beam splitter using a delay stage. The output modes of this beam splitter are then sent to two TESs. We measured the total system efficiency (i.e. coupling, transmission and detection efficiencies combined) in modes $c$ and $d$ to be $\eta_{c}=0.274\pm 0.001$ and $\eta_{d}=0.354\pm 0.002$, respectively. Further details on the experimental setup can be found in the Supplemental Material SM .

Results. In Fig. 3, we show the measured photon-number difference statistics for $j=6$ and three different values of $|\alpha|^{2}$. The measured statistics (black points) can be compared with the statistics expected from the theory model with and without the classical field approximation, i.e. $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ (blue curves) and $P^{(j,\alpha)}(\Delta n)$ (red curves), respectively. Both models include the effects of experimental imperfections such as detection inefficiencies and mode mismatch. The parameters modeling these imperfections are determined from independent measurements (see the Supplemental Material SM ). For $|\alpha|^{2}=0$, $P^{(j,\alpha)}(\Delta n)$ is obtained by projecting $\hat{\rho}^{(6)}_{a}$ onto photon-number states. In this case, $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ is not defined. For $|\alpha|^{2}=6.52$, the data points agree with $P^{(j,\alpha)}(\Delta n)$ but not with $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$, indicating that the classical field approximation is not yet valid. However, for $|\alpha|^{2}=15.41$, there is good agreement between the data and both $P^{(j,\alpha)}(\Delta n)$ and $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$. Surprisingly, it seems that the measurement is already projecting $\hat{\rho}^{(6)}_{a}$ onto quadrature states despite still being far from the regime where $|\alpha|\gg N$ (note that $N=10.4$ for $j=6$ due to loss in the herald mode).

To understand why, we quantify the transition from $\Delta\hat{n}$ to $\Delta\hat{n}_{\mathrm{classical}}$ by computing the sum of the squared residuals,

[TABLE]

where $\nu$ is the number of data points. $\mathcal{S}_{\mathrm{classical}}$ quantifies the discrepancy between the measured data $P^{(j,\alpha)}_{\mathrm{exp}}$ and the classical model $P^{(j,\alpha)}_{\mathrm{classical}}(\Delta n)$. Ideally, $\mathcal{S}_{\mathrm{classical}}$ smoothly converges to zero for increasing $|\alpha|^{2}$ as the validity of the classical field approximation improves. In Fig. 4(a), we plot $\mathcal{S}_{\mathrm{classical}}$ for four different herald outcomes $j$. We found heuristically that an exponential curve adequately models $\mathcal{S}_{\mathrm{classical}}$ for sufficiently large $|\alpha|^{2}$. Thus, for each $j$, we fit $\mathcal{S}_{\mathrm{classical}}$ to $A\exp{\left(-B|\alpha|^{2}\right)}$ from which we determine $|\alpha|^{2}_{\mathrm{min}}$, the coherent state strength required to reach below $\mathcal{S}_{\mathrm{classical}}=6.7\times 10^{-6}$. This threshold corresponds to the sum of the squared residuals obtained with the quantum model [i.e. replacing $P^{(j,\alpha)}_{\mathrm{classical}}(\Delta n)$ with $P^{(j,\alpha)}(\Delta n)$ in Eq. (6)], averaged over all $j$ and $|\alpha|$. We plot $|\alpha|^{2}_{\mathrm{min}}$ for $j=0$ to $6$ in Fig. 4(b). Interestingly, we observe a linear scaling between $|\alpha|^{2}_{\mathrm{min}}$ and the average photon-number in the signal, $N$, instead of a quadratic scaling Tyc and Sanders (2004); Skotiniotis et al. (2017). We believe this relaxed requirement on $|\alpha|^{2}_{\mathrm{min}}$ is due to our detector inefficiency which smooths the fine features in the quadrature distribution of $\hat{\rho}^{(j)}_{a}$. This reasoning agrees with the findings of Ref. Braunstein (1990), which showed that a smaller $|\alpha|^{2}_{\mathrm{min}}$ is required for states with smooth quadrature distributions such as coherent states. Thus, the transition from photon-number to quadrature measurement in Fig. 3 occurred for a weaker $|\alpha|$ than might be expected (i.e. before the regime $|\alpha|\gg N$) since the quadrature distribution of $\hat{\rho}^{(6)}_{a}$ is smooth.

State engineering. So far, we have demonstrated that WFHD can tune between performing photon-number and quadrature measurements. Here, we show it can also be used as a state engineering tool by projecting one part of a photon-number entangled state onto a particular measurement basis, thus steering the possible measurement outcomes on the other part of the entangled state. The concept is shown schematically in Fig. 5(a). As before, mode $a$ of a two-mode squeezed vacuum state $\ket{\Psi_{ha}}$ is sent to the weak-field homodyne detector. The detector projects mode $a$ onto a state $\ket{\Gamma_{a}}$ which depends on the detection outcome $(m,n)$ and $\alpha$ Puentes et al. (2009). This measurement transforms mode $h$ to the state $\ket{\phi_{h}}=\mathcal{N}\braket{\Gamma_{a}}{\Psi_{ha}}$, where $\mathcal{N}$ is a normalization factor. Many different classes of states $\ket{\phi_{h}}$ can be heralded since $\ket{\Gamma_{a}}$ can be continuously tuned between photon-number and quadrature-like states. This versatility makes WFHD a powerful state engineering tool. For example, it was recently shown that $\ket{\phi_{h}}$ can have nearly perfect fidelity with a Schrödinger cat state $\ket{\beta}+\ket{-\beta}$ of arbitrary amplitude $\beta$ Thekkadath et al. (2020).

As a proof of concept, we consider the specific detection outcome $(m,n)=(6,0)$ when $|\alpha|^{2}=15.4$. The measured photon-number distribution $P(k)=\left|\braket{k_{h}}{\phi_{h}}\right|^{2}$ is shown in Fig. 5(b), which we use to calculate the second-order correlation function $g^{(2)}=\sum_{k}(k^{2}-k)P(k)/(\sum_{k}kP(k))^{2}$. When $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$ are not temporally overlapped (using a temporal delay much smaller than the detection window), the two do not interfere and $P(k)$ is a thermal-like distribution (red squares, $g^{(2)}=1.59\pm 0.15$). This is because the average number of photons in mode $b$ (15.4) is much larger than in mode $a$ (1.74), and so the detection outcome provides little information about the number of photons in mode $a$. In contrast, when $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$ are temporally overlapped, $P(k)$ changes drastically and resembles a Poisson distribution (blue circles, $g^{(2)}=1.19\pm 0.11$). This demonstrates that the heralded state is strongly modified by the interference between $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$.

Conclusions. We experimentally demonstrated that weak-field homodyne is a versatile detection scheme that can tune between photon-number and quadrature measurements. We were able to demonstrate this tunability with quantum signals having up to nearly 11 photons. Our work shows that the dynamic range, efficiency, and noise-level of current photon-number-resolving detectors is enough to probe both particle and wave properties of quantum signals using weak-field homodyne detection.

Another important achievement of our work is performing state engineering using weak-field homodyne for the first time. In particular, we showed that the photon-number distribution of a state heralded using weak-field homodyne can be strongly modified by the interference between the local oscillator and the signal. This paves the way towards promising applications in hybrid discrete- and continuous-variable quantum information processing protocols Andersen et al. (2015). For example, weak field homodyne can be used to herald phase-sensitive non-Gaussian states such as two and four-component Schrödinger cat states Thekkadath et al. (2020). An exciting prospect would be to verify the quantum non-Gaussian features of this heralded state using a second instance of the weak-field homodyne detector, which has the added advantage over regular homodyne detection of being robust against losses Royer (1977); Laiho et al. (2010); Lipfert et al. (2015); Harder et al. (2016); Bohmann et al. (2018); Sperling et al. (2020). The same experimental setup could also be used to perform a hybrid Bell test Grangier et al. (1988); Banaszek and Wódkiewicz (1998); Kuzmich et al. (2000); Ji et al. (2010); Donati et al. (2014).

Acknowledgements. The authors are grateful for discussions with J. Sperling, and thank J. Renema for his assistance with the installation of the cryogenic infrastructure. This work was supported by the following: the Natural Sciences and Engineering Research Council of Canada (NSERC); the Networked Quantum Information Technologies Hub (NQIT) as part of the UK National Quantum Technologies Programme Grant EP/N509711/1); the EPSRC Programme Grant ”Building Large Optical Quantum States” (EP/K034480); the Fondation Wiener Anspach; the European Commission project QUCHIP (H2020-FETPROACT-2014) and the ERC grant MOQUACINO (AdG).

Supplemental material

Theory I: Photon-number difference statistics without a classical field approximation

Here, we derive the photon-number difference statistics measured by our detector without invoking a classical field approximation. The general setup is shown in Fig. S6. We begin by deriving the joint photon-number distribution measured at the output of the beam splitter in the idealized case where our signal is in a pure photon-number state (and single spatio-temporal mode), $\hat{\rho}_{a}=\ket{j_{a}}\bra{j_{a}}$, and assume perfect detection efficiency, $\eta_{c}=\eta_{d}=1$. The joint probability of detecting $m$ photons in mode $c$ and $n$ photons in mode $d$ (i.e. the outcome $(m,n)$) at the output of the beam splitter, $P^{(j,\alpha)}(m,n)$, is given by:

[TABLE]

where $\hat{c}$ ($\hat{d}$) is the photon annihilation operator in mode $c$ ($d$). The input and output modes of the beam splitter are related through the following transformation:

[TABLE]

Inserting these expressions into Eq. (S7), we find:

[TABLE]

The non-zero terms in the first sum satisfy $l=m+k-j$. Furthermore, $\hat{b}\ket{\alpha_{b}}=\alpha\ket{\alpha_{b}}$. Using these two facts, Eq. (S9) can be simplified to:

[TABLE]

and $P^{(j,\alpha)}(m,n)=0$ when $m+n<j$.

We now consider the effect of detection efficiency. As shown in in Fig. S6, we model detection efficiency by placing a fictitious beam splitter of transmissivity $\eta_{c}$ ($\eta_{d}$) before just before the detector in mode $c$ ($d$). In this case, the joint probability to measure the outcome $(m,n)$ can be obtained by performing a Bernoulli transformation on Eq. (S10) Leonhardt (1997):

[TABLE]

The intuition for this transformation is as follows. Suppose $m$ photons are detected in mode $c$ (we ignore mode $d$ to make the explanation easier to follow). This detection event could have occurred when there were $x$ photons before the fictitious beam splitter, but $x-m$ photons were reflected, i.e. lost (where $x\geq m$). Let $p(x)$ be the probability for there to have been $x$ photons before the beam splitter. The probability to have transmitted $m$ photons and lost $x-m$ photons is $\eta_{c}^{m}$ and $(1-\eta_{c})^{x-m}$, respectively. Since the detector does not discriminate between the photons, we must include the binomial factor ${x\choose m}$. Thus, the total probability of detecting $m$ photons is the product of each of these probabilities, summed over all possible values of $x$, i.e. $\sum_{x=m}^{\infty}{x\choose m}\eta_{c}^{m}(1-\eta_{c})^{x-m}p(x)$. In Eq. (S11), we apply a similar procedure but consider both modes $c$ and $d$.

We also consider the effect of mode mismatch (e.g. spatial, temporal, spectral, and polarization mismatch) between $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$. In principle, both $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$ can occupy several spatio-temporal modes, thus making a full treatment of mode mismatch quite involved Raymer et al. (1995). However, we found that the simpler approach of decomposing the problem into two effective orthogonal modes is sufficient to model our data. This sort of heuristic model has been used and studied elsewhere, see e.g. Refs. Ourjoumtsev et al. (2006); Tualle-Brouri et al. (2009). Through a Gram-Schmidt process, $\ket{\alpha_{b}}$ can be decomposed in the following way Tichy (2014):

[TABLE]

where $||$ denotes the same mode as the signal, $\perp$ denotes a mode orthogonal to the signal, and $\mathcal{M}\in[0,1]$ is a mode overlap parameter. Since $\ket{\alpha_{b}}$ remains a coherent state when it is split across the modes $||$ and $\perp$, its amplitude is simply scaled by $\sqrt{\mathcal{M}}$ and $\sqrt{1-\mathcal{M}}$, respectively. We wish to determine the joint probability of measuring the outcome $(m,n)$ when $\mathcal{M}\neq 1$, i.e. $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\mathcal{M}}(m,n)$. We treat the detector as mode-insensitive, and so we convolve the joint probabilities of measuring the outcome $(m,n)$ in mode $||$ and mode $\perp$:

[TABLE]

where $*$ denotes a convolution operation. The term $P^{(j,\sqrt{\mathcal{M}}\alpha)}_{\eta_{c},\eta_{d}}(m,n)$ is the joint probability of measuring the outcome $(m,n)$ in the mode of the signal $||$. Similarly, the term $P^{(0,\sqrt{1-\mathcal{M}}\alpha)}_{\eta_{c},\eta_{d}}(m,n)$ is the joint probability of measuring the outcome $(m,n)$ in the mode orthogonal to the signal, i.e. $\perp$ (hence the reason why the signal is vacuum in this mode, i.e. $j=0$).

Finally, we also consider the effect of an imperfect signal preparation. Our signal $\hat{\rho}_{a}$ is prepared by performing a heralding measurement $\ket{j_{h}}\bra{j_{h}}$ on a two-mode squeezed vacuum state $\ket{\Psi_{ha}}=\sqrt{1-|\lambda|^{2}}\sum_{n=0}^{\infty}\lambda^{n}\ket{n_{h},n_{a}}$, where $\lambda$ is a squeezing parameter (see Fig. S7). However, an imperfect detection efficiency $\eta_{h}$ transforms the heralding measurement to $\ket{j_{h}}\bra{j_{h}}\rightarrow\hat{\Pi}_{h}=\sum_{f=j}^{\infty}{f\choose j}\eta_{h}^{j}(1-\eta_{h})^{f-j}\ket{f_{h}}\bra{f_{h}}$ (see text below Eq. (S11) for intuition) Leonhardt (1997). This imperfect heralding measurement $\hat{\Pi}_{h}$ prepares the signal $\hat{\rho}^{(j)}_{a}$ which is given by:

[TABLE]

where $\mathcal{X}$ is a normalization factor such that $\mathrm{Tr}_{a}\left(\hat{\rho}_{a}^{(j)}\right)=1$. We wish to obtain the joint probability of measuring the outcome $(m,n)$ given the imperfect signal $\hat{\rho}_{a}^{(j)}$, i.e. $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\eta_{h},\mathcal{M}}(m,n)$. Eq. (S14) is a statistical mixture of photon number states. Thus, $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\eta_{h},\mathcal{M}}(m,n)$ can be found by computing $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\mathcal{M}}(m,n)$ for each term in the mixture $\hat{\rho}_{a}^{(j)}$ and adding the probabilities with the appropriate weight, that is:

[TABLE]

The photon-number difference statistics, i.e. the probability to measure $\Delta n=m-n$ photons at the output of the beam splitter, can be obtained by performing the following sum:

[TABLE]

Eq. (S16) is used to calculate the red curves in Fig. 3.

Theory II: Photon-number difference statistics with a classical field approximation

Here, we would like to invoke the classical field approximation $\hat{b}\rightarrow|\alpha|e^{i\theta}$ and obtain an expression analogous to $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\eta_{h},\mathcal{M}}(\Delta n)$. We refer to reader again to Fig. S6. As before, we begin in the idealized case where $\hat{\rho}_{a}=\ket{j_{a}}\bra{j_{a}}$ and $\ket{\alpha_{b}}$ are combined on a balanced beam splitter. For perfect detection efficiency, it was shown in Ref. Vogel and Grabow (1993) that the probability to measure $\Delta n$ photons at the output of the beam splitter after invoking the classical field approximation is giving by:

[TABLE]

where $H_{j}$ is a Hermite polynomial of order $j$. Notably, Eq. (S17) is simply the quadrature distribution of a $j$ photon number state scaled by $\alpha$, i.e. $\Delta n\rightarrow\Delta n/\alpha$ Leonhardt (1997). Eq. (S17) is referred to as $P_{\mathrm{classical}}^{(j,\alpha)}(\Delta n)$ in the main text.

The effect of (imbalanced) detection efficiency was also considered in Ref. Vogel and Grabow (1993):

[TABLE]

where * denotes a convolution, $g=\eta\alpha$, $\Delta\tilde{n}=\Delta n-\alpha^{2}(\eta_{d}-\eta_{c})/2$, $\sigma=\sqrt{g^{2}(1-\eta)/\eta}$, and $\eta=(\eta_{c}+\eta_{d})/2$.

We treat the effect of mode mismatch in the same way as in Eq. (S13):

[TABLE]

Similarly, our imperfect signal preparation is treated in the same way as was done previously. Using Eq. (S14), we obtain:

[TABLE]

Eq. (S20) is used to calculate the blue regions in Fig. 3 of the main text.

Theory III: State engineering using weak-field homodyne

We can use Eq. (S15) to calculate the photon-number distribution of $\ket{\phi_{h}}$ since $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\eta_{h},\mathcal{M}}(m,n)$ is the joint probability to measure $j$ photons in mode $h$, $m$ photons in mode $c$, and $n$ photons in mode $d$. That is, we compute $P^{(j,\alpha)}_{\eta_{c},\eta_{d},\eta_{h},\mathcal{M}}(m,n)$ for every $j$ while fixing $(m,n)$ and $\alpha$. The resulting distribution is then normalized, thus yielding $P(j|m,n,\alpha)=\left|\braket{j}{\phi_{h}}\right|^{2}$. This strategy is used to calculate the red and blue bars in Fig. 5 of the main text.

Experimental details

The experimental setup is shown in Fig. S8. The pump is a titanium sapphire oscillator followed by a regenerative amplifier that outputs femtosecond pulses (780 $\pm$ 20 nm [full width at half maximum]) at a rate of 100 kHz. This rate is chosen to accommodate the thermal relaxation time of the transition edge sensor (TES) detectors. The pump pulses are split into two paths, labelled \small{1}⃝ and \small{2}⃝ in Fig. S8.

In path \small{1}⃝, we prepare the signal $\hat{\rho}_{a}$. The pump is filtered to 775 $\pm$ 2 nm using two angle-tuned bandpass (BP) filters. The filtered pump is then coupled into a ppKTP waveguide and generates two-mode squeezed vacuum via type-II spontaneous parametric down-conversion. The pump spectrum is chosen such that the down-converted signal (1554 nm) and idler (1547 nm) modes are approximately spectrally uncorrelated ($g^{(2)}=1.849\pm 0.007$ in idler mode). The pump is then discarded with a longpass (LP) filter and the orthogonally-polarized signal and idler modes are spatially separated with a polarizing beam splitter (PBS). Both modes are sent through a BP filter and coupled in fibers. The idler mode is sent directly to a TES detector, whereas the signal mode is sent to a polarization-maintaining fiber beam splitter (FBS).

In path \small{2}⃝, we prepare the coherent state $\ket{\alpha_{b}}$. The pump is filtered to 783 $\pm$ 2 nm and coupled into a second ppKTP waveguide. In contrast to the previous path, we also couple 1580 nm light from a continuous-wave laser in order to seed the down-conversion process. Through difference frequency generation, light with Poissonian photon-number statistics ($g^{(2)}=1.005\pm 0.002$, Fano factor $=1.015\pm 0.002$) is generated in the polarization mode orthogonal to the seed, i.e. $\ket{\alpha_{b}}$. The seed and pump spectra are chosen to optimize the spectral overlap between the signal $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$ given that there are small differences between the phase matching properties of both waveguides. $\ket{\alpha_{b}}$ is separated from the seed light using a PBS and a BP filter. To further minimize the amount of seed light leakage, we carve 2 ns pulses from the continuous-wave laser using an electro-optic modulator (EOM). $\ket{\alpha_{b}}$ is sent through a half-wave plate (HWP) to match its polarization to that of the signal. The intensity $|\alpha|^{2}$ is adjusted using a neutral-density (ND) filter wheel. Finally, $\ket{\alpha_{b}}$ is coupled into fiber and sent to the FBS.

Model parameters

We now describe how we determined the model parameters from various measurements. Our model (see Eqs. (S13) and  (S20)) requires the following five parameters: the detector efficiencies $\eta_{c}$, $\eta_{d}$, and $\eta_{h}$, the squeezing parameter $\lambda$, and the mode overlap parameter $\mathcal{M}$. The results are summarised in Table I.

We determine $\mathcal{M}$ by measuring an interference signal. In theory, when the signal is a single photon, the probability to detect only one photon at both outputs of the FBS, $P^{(1,\alpha)}(1,1)$, vanishes. In practice, this probability does not vanish due to a number of experimental imperfections, such as background counts, the modal purity of $\hat{\rho}_{a}$, and imperfect mode overlap between $\hat{\rho}_{a}$ and $\ket{\alpha_{b}}$. As such, the visibility $\mathcal{V}$ of this interference signal provides a lower bound on $\mathcal{M}$. We measured $\mathcal{V}=0.800\pm 0.060$, and found that the models agree best with the data by using $\mathcal{M}=0.82$.

The total system efficiencies are obtained via a Klyshko measurement Klyshko (1980). For this measurement, we set $|\alpha|=0$ and decrease the pump power to ensure that $|\lambda|\ll 1$. In this small squeezing limit, $\ket{\Psi_{ha}}$ consists mostly of photon pairs, i.e. higher order photon-number states can be neglected. Then, $\eta_{h}$ is the probability that TES 1 detects a photon given that either TES 2 or TES 3 detected a photon. Next, we determine the probability $\eta_{s}$ that either TES 2 or TES 3 detects a photon given that TES 1 detected a photon. Because of the balanced beam splitter before TES 2 and TES 3, the measured $\eta_{s}$ is the the average of $\eta_{c}$ and $\eta_{d}$, i.e. $\eta_{s}=(\eta_{c}+\eta_{d})/2$. Moreoever, $R=\eta_{c}/\eta_{d}$ can be obtained by looking at the ratio of the number of detected photons at TES 2 and TES 3. As such, we can obtain $\eta_{d}$ and $\eta_{c}$ from $\eta_{d}=2\eta_{s}/(1+R)$ and $\eta_{c}=R\eta_{d}$, respectively. Using this procedure, we find $\eta_{h}=0.395\pm 0.002$, $\eta_{c}=0.274\pm 0.001$, and $\eta_{d}=0.352\pm 0.002$. Note that these are not intrinsic detector efficiencies as they include the efficiencies of fiber coupling, fiber transmission, and waveguide transmission.

Next, we determine the squeezing parameter $|\lambda|$. To do so, we count the number of photons produced by $\ket{\Psi_{ha}}$ per pump pulse when $|\alpha|=0$. We measure $\braket{\hat{n}_{h}}=0.689\pm 0.001$ photons at TES 1. The number of photons before losses is $\braket{\hat{n}_{h}}/\eta_{h}$. Thus, the squeezing parameter is given by $|\lambda|=\tanh[{\mathrm{arcsinh}(\sqrt{{\braket{\hat{n}_{h}}/\eta_{h}}}})]=0.797\pm 0.001$.

Finally, $|\alpha|$ is measured by blocking the pump in path \small{1}⃝ and counting the number of photons arriving at TES 2 and 3, i.e. $|\alpha|=\sqrt{\braket{\hat{n}_{c}}/\eta_{c}+\braket{\hat{n}_{d}}/\eta_{d}}$.

Transition edge sensors

The transition edge sensors (TESs) operate in a dilution refrigerator at a temperature of around 80 mK. Each TES consists of a superconducting tungsten film that is biased with a current in parallel with a 20 m$\Omega$ resistor such that the TES is near its transition to the non-superconducting state. The resistor provides electrothermal feedback which stabilises the TES operating temperature Irwin (1995). When light is absorbed by the film, there is an increase in resistance that depends on the energy (and thus the number of photons) absorbed. This generates an electrical response which is amplified using an array of superconducting quantum interference devices (SQUIDs) coupled to an inductor to provide cryogenic transimpedence gain Welty and Martinis (1991). This signal is then further amplified and filtered at room temperature. The final electrical signal is read by an analogue-to-digital converter (ADC) triggered at the 100 kHz clock of the pump laser. The trigger initiates a 10 $\mu$s acquisition on each TES.

Preliminary to data acquisition, we measure the average electrical response of the TES, which provides us with a matched filter Figueroa-Feliciano et al. (2000). During the data acquisition, we calculate the overlap integral between this average trace and the detector response in real time. The result of this overlap integral is a scalar value corresponding to an estimate of the energy absorbed by the TES. As we only save the matched filter output, we minimize the amount of time spent writing data to the hard drive and thereby minimize the amount of trigger events lost due to data transfer bottlenecks. It should be noted that this method of processing was chosen to maximise the data acquisition rate and that more sophisticated processing techniques exist Humphreys et al. (2015); Levine et al. (2012).

When we plot a histogram of the matched filter output, such as in Fig. S9, we see peaks due to the quantised energy of the optical field. We use a peak finding routine to create bins for our data. We then assign a label to each bin corresponding to its estimated photon number.

Testing the nonclassicality of our signal

Here, we test whether the photon-number statistics of our signal $\hat{\rho}^{(j)}_{a}$ have signatures of nonclassicality. We test for two different signatures, namely sub-Poissonian statistics Leonhardt (1997) and submultinomial statistics Sperling et al. (2017a, b, 2020). For both tests, the coherent state is blocked, i.e. $\ket{\alpha_{b}}=\ket{0_{b}}$. The results are shown in Fig. S10. We find that the state $\hat{\rho}^{(0)}_{a}$ has classical photon-number statistics and hence fails both non-classical tests, as expected. However, for $j>0$, we measure both sub-Poissonian and submultinomial statistics with more than 1 standard deviation of confidence. We now discuss how we performed each test.

To test for submultinomial statistics, we follow the procedure outlined in Refs. Sperling et al. (2017a, b, 2020). Namely, we determine the correlation matrix $M$ from our measured statistics. Let $E^{(j)}(k,l)$ denote the number of events where $j$ photons were detected at TES 1 (i.e. the herald TES), $k$ photons were detected at TES 2, and $l$ photons were detected at TES 3 for $j,k,l\in[0,6]$. Also, let $E_{tot}^{(j)}=\sum_{k,l}E^{(j)}(k,l)$ denote the total number of events for a given herald outcome. Following Eqs. (D1), (D2), and (D3) of Ref. Sperling et al. (2020), the correlation matrix element $M^{(j)}_{x,y}$ for a particular herald outcome $j$ is given by:

[TABLE]

with $\delta$ denoting the Kronecker delta. We compute $M^{(j)}_{x,y}$ for $x,y\in[0,6]$ then calculate the minimum eigenvalue $\mu_{\mathrm{min}}$ of the resulting matrix $M^{(j)}$. $\mu_{\mathrm{min}}<0$ is a signature of submultinomial statistics and hence nonclassicality.

To test for sub-Poissonian statistics, we first determine the photon-number distribution of our signal $P^{(j)}(n)=\mathrm{Tr}(\ket{n_{a}}\bra{n_{a}}\hat{\rho}^{(j)}_{a})$ from our measurements. The quantity $P^{(j)}(n)$ can be obtained by simply combining the statistics of TES 2 and TES 3, that is:

[TABLE]

Using Eq. (S22), we compute the second-order correlation function $g^{(2)}$ for each herald outcome $j$:

[TABLE]

$g^{(2)}<1$ is a signature of sub-Poissonian statistics and hence nonclassicality. We note that, since the value of $g^{(2)(j)}$ is unchanged if loss is applied to $\hat{\rho}^{(j)}_{a}$, the procedure above is not affected by the efficiencies of TES 2 and TES 3.