Detector-Agnostic Phase-Space Distributions

TL;DR

This paper introduces a detector-agnostic method for reconstructing quantum phase-space distributions, allowing characterization of light states without detailed detector models, demonstrated on single- and two-photon states.

Contribution

We develop a measurement scheme that reconstructs phase-space distributions independently of detector specifics, generalizing traditional quasiprobability distributions.

Findings

Successfully implemented with transition-edge sensors

Revealed features of heralded photon states in phase space

Certified nonclassicality with high statistical significance

Abstract

The representation of quantum states via phase-space functions constitutes an intuitive technique to characterize light. However, the reconstruction of such distributions is challenging as it demands specific types of detectors and detailed models thereof to account for their particular properties and imperfections. To overcome these obstacles, we derive and implement a measurement scheme that enables a reconstruction of phase-space distributions for arbitrary states whose functionality does not depend on the knowledge of the detectors, thus defining the notion of detector-agnostic phase-space distributions. Our theory presents a generalization of well-known phase-space quasiprobability distributions, such as the Wigner function. We implement our measurement protocol, using state-of-the-art transition-edge sensors without performing a detector characterization. Based on our approach, we reveal the characteristic features of heralded single- and two-photon states in phase space and certify their nonclassicality with high statistical significance.

Figures

Tables

Table 1. Table 1: For different heralding outcomes, k h subscript 𝑘 ℎ k_{h} , we show the nonclassicality criteria μ min < 0 subscript 𝜇 0 \mu_{\min}<0 and g min < 0 subscript 𝑔 0 g_{\min}<0 . g min subscript 𝑔 g_{\min} is defined in Eq. ( 6 ). μ min subscript 𝜇 min \mu_{\mathrm{min}} is the minimal eigenvalue to the second-order correlation matrix M 𝑀 M defined in Eq. (6) of Ref. SCEMRKNLGVAW17 ; see also the SM supplement . “ − 0 0 -0 ” indicates a slightly negative mean value that rounds to zero.

$k_h$	${\mu }_{\min }$	$g_{\min }$
0	$\left(-0\pm 9\right)\times {10}^{-4}$	$\left(-0\pm 2\right)\times {10}^{-9}$
1	$-0.15\pm 0.03$	$-0.026\pm 0.003$
2	$-0.10\pm 0.03$	$-0.017\pm 0.003$

Table 2. Table 2: For the available heralding bins, k h subscript 𝑘 ℎ k_{h} , the minimal eigenvalues μ min subscript 𝜇 \mu_{\min} of the matrix M 𝑀 M are shown. Significant negativities defines the notion of nonclassical sub-multinomial light SCEMRKNLGVAW17 ; SECMRKNLGWAV17 .

$k_h$		sub-multinomial
0		$-\mathrm{0.000\hspace{0.17em}0}\pm \mathrm{0.000\hspace{0.17em}9}$
1		$-0.15\pm 0.03$
2		$-0.10\pm 0.03$
3		$-0.17\pm 0.07$
4		$-0.3\pm 0.2$

Equations

$$c_{N_0,\dots ,N_K}(\beta )=$$

$$\begin{array}{rl}{g}_{z_0,\dots ,z_K}(\beta )=& \sum _{N_0,\dots ,N_K}{z}_0^{N_0}\cdots z_K^{N_K}{c}_{N_0,\dots ,N_K}(\beta )\\ =& {⟨:{\left(z_0{\hat{\pi }}_0+\cdots +z_K{\hat{\pi }}_K\right)}^N:⟩}_{\hat{\rho }(\beta )},\end{array}$$

$$g_{z_0,\dots ,z_K}(\beta )\stackrel{cl.}{\ge }0.$$

$$G_z(\beta )=$$

$$g_{\min }=\underset{\beta }{\min }\underset{\begin{array}{c}{z}_0,\dots ,z_K:\\ |z_0{|}^2+\cdots +|z_K{|}^2\le 1\end{array}}{\min }{g}_{z_0,\dots ,z_K}(\beta ).$$

$$|{\beta }^{(\mathrm{DI})}|=\sqrt{\sum _{N_0,\dots ,N_K}[0N_0+\cdots +K{N}_K]c_{N_0,\dots ,N_K}^{(\mathrm{vac})}(\beta )},$$

$$G_z(\beta )=\int d^2{\beta }^{\prime }P({\beta }^{\prime };1)G_z^{(\mathrm{vac})}\left(\beta -\frac{t}{r}{\beta }^{\prime }\right).$$

$$\hat{\rho }(\beta )=\int d^2\alpha P(\alpha )|t\alpha -r\beta ⟩⟨t\alpha -r\beta |.$$

$$\begin{array}{rl}{c}_{N_0,\dots ,N_K}(\beta )=& \int d^2\alpha P(\alpha )\frac{N!}{N_0!\cdots N_K!}\\ & \times \prod _{k=0}^K{\underset{\stackrel{def.}{=}{p}_k\left(\frac{t\alpha -r\beta }{\sqrt{N}}\right)}{\underbrace{⟨\frac{t\alpha -r\beta }{\sqrt{N}}|{\hat{\pi }}_k|\frac{t\alpha -r\beta }{\sqrt{N}}⟩}}}^{N_k},\end{array}$$

$$g_{z_0,\dots ,z_K}(\beta )=$$

$$G_z(\beta )=$$

$$\begin{array}{rl}{G}_z(\beta )=& \int d^2\alpha P(\alpha )\exp \left(-[1-z]|t{|}^2\eta {\left|\alpha -\frac{r}{t}\beta \right|}^2\right)\\ =& {⟨:\exp \left[-(1-z)|t{|}^2\eta \hat{n}\left(\frac{r}{t}\beta \right)\right]:⟩}_{\hat{\rho }},\end{array}$$

$$\displaystyle c_{N_{0},\ldots,N_{K}}=\left\{\begin{array}[]{ll}c_{0,\ldots,0,N_{i}=2,0,\ldots,0}&\text{ for }0\leq i\leq K,\\ c_{0,\ldots,0,N_{i}=1,0,\ldots,0,N_{j}=1,0,\ldots,0}&\text{ for }0\leq i<j\leq K,\\ 0&\text{ otherwise},\end{array}\right.$$

$$\overline{f}=$$

$$=$$

$$+\sum _{0\le i<j\le K}\frac{E(i,j)+E(j,i)}{E}\underset{\stackrel{\mathrm{def}}{=}f(i,j)=f(j,i)}{\underbrace{f_{0,\dots ,0,N_i=1,0,\dots ,0,N_j=1,0,\dots ,0}}}$$

$$=$$

$$E(k_1,k_2)=\stackrel{\text{(symmetric part)}}{\overbrace{\frac{E(k_1,k_2)+E(k_2,k_1)}{2}}}+\stackrel{\text{(asymmetric part)}}{\overbrace{\frac{E(k_1,k_2)-E(k_2,k_1)}{2}}}.$$

$$\left|\overline{f}-{\overline{f}}^{(\mathrm{sym})}\right|\le \sum _{k_1,k_2}|f(k_1,k_2)|\left|\frac{E(k_1,k_2)-E(k_2,k_1)}{2E}\right|={ϵ}_f,$$

$$\begin{array}{rl}{M}_{i,j}=& N\overline{N_i(N_j+{\delta }_{i,j})}-(N-1)\overline{N_i}\overline{N_j}\\ =& N^2(N-1)\left(⟨:{\hat{\pi }}_i{\hat{\pi }}_j:⟩-⟨:{\hat{\pi }}_i:⟩⟨:{\hat{\pi }}_j:⟩\right),\end{array}$$

$$\overline{N_i}=$$

$$\overline{N_i(N_j+{\delta }_{i,j})}=$$

$$=$$

$$I(\beta )=$$

$$=$$

$${\left|{\beta }^{(\mathrm{DI})}\right|}^2=\sum _{i=0}^{K}i\overline{N_i}\in [0,NK],$$

$$G_z(\beta )=\sum _{N_0+\cdots +N_K=N}{z}^{0N_0+\cdots +K{N}_K}{c}_{N_0,\dots ,N_K}(\beta ).$$

$$\mathrm{Prob}(|X-X^{\prime }|\le \delta )=\underset{-\delta }{\overset{+\delta }{\int }}du\underset{-\infty }{\overset{+\infty }{\int }}dzp(X=z)p(X^{\prime }=z+u),$$

$$\mathrm{Prob}\left(G_z^{(k_h)}\ne G_z^{(k_h^{\prime })}\right)=1-\mathrm{Prob}(\forall n:|X_n-X_n^{\prime }|\le {\delta }_n)$$

$$=$$

$${\left[\mathrm{Prob}\left(G_z^{(k_h)}\ne G_z^{(k_h^{\prime })}\right)\right]}_{k_h,k_h^{\prime }=0,1,2}=\left[\begin{array}{ccc}7.5\%& 100.\%& 100.\%\\ 100.\%& 7.5\%& 98.9\%\\ 100.\%& 98.9\%& 7.5\%\end{array}\right],$$

$$G_z^{(m)}(\beta )=\sum _{j=0}^m\left(\genfrac{}{}{0ex}{}{m}{j}\right)\frac{1}{j!}{\left[\frac{|t{|}^2}{|r{|}^2}\right]}^j{\partial }_{\beta }^j{\partial }_{{\beta }^{\ast }}^j{G}_z^{(0)}(\beta ),$$

$$g_{z_0,\dots ,z_K}(\beta )=\int d^2\alpha P(\alpha )g_{z_0,\dots ,z_K}^{(0)}\left(\beta -\frac{t}{r}\alpha \right),$$

Full text

Introduction.—

The characterization of quantum light is a main challenge one encounters when implementing classically infeasible tasks, such as quantum communication protocols KMSUZ16 ; BFV09 ; RL09 . On a more fundamental level, studying the peculiarities of quantized radiation fields leads to a profound understanding of the role of quantum physics in nature in general, and how it is distinct from classical wave theories in particular. As in classical systems, quantum-optical phase-space distributions offer a versatile instrument to directly visualize unique features of nonclassical light, such as demonstrated for squeezing W83 ; SBRF93 ; BSM97 . Moreover, negativities in certain phase-space functions directly point at quantum properties of light; see, e.g., Refs. LHABMS01 ; KVPZB08 ; LKCGS10 ; DEWKDKCM13 ; HSRHMSS16 ; BFS17 . For the above reasons, the representation of quantum light in phase space is one of the most frequently applied methods to characterize nonclassical light.

However, the estimation of phase-space distributions from experimental data is a cumbersome task. Consequently, this reconstruction problem inspired a wide range of research WVO99 ; S07 ; LR09 , leading to sophisticated analytical tools, such as solving inversion problems T97 ; SSG09 , employing diverging pattern functions R96 ; LMKRR96 , performing maximum-likelihood estimations H97 ; L04 ; KWR04 , and using data pattern recognition RMH10 ; MIMSRH13 . In addition, each family of detection devices has to be equipped with its own precise model to reliably extract information about phase-space functions WVO99 ; S07 ; LR09 . This treatment comprises a comprehensive analysis that assesses (i) how a detector responds to incident light KK64 ; FW91 , including, e.g., nonlinear detection responses JA69 ; AML11 , and (ii) how the light absorption is influenced by a number of possible imperfections, e.g., efficiencies PM07 ; WCLMSPTW09 . Moreover, applying these methods can also require universally applicable, yet rather demanding theoretical and experimental techniques in practice, such as performing detector tomography and calibration K80 ; LS99 ; F01 ; AMP04 ; FLCEPW09 ; CLPFSMSEPW09 ; RFZMGLDFE12 ; PHMH12 ; BKSSV17 .

Despite these challenges, phase-space distributions constitute a highly successful approach to revealing nonclassical properties of light LHABMS01 ; KVPZB08 ; LKCGS10 ; DEWKDKCM13 ; HSRHMSS16 ; BFS17 . For example, $s$-parametrized quasiprobabilities CG69 ; AW70 , as well as their non-Gaussian generalizations K66 ; KV10 , can exhibit negativities that are incompatible with classical light. Even if a phase-space function does not exhibit negativities, observable patterns render it possible to identify quantum features, for instance, via the nonnegative Husimi function H40 ; LGC18 ; MS04 or through marginal distributions A93 ; PLLSZZZKN17 . Because of its success, the concept of phase-space functions has been further extended to other physical scenarios; see Refs. SW18 ; SV19 . To name a few, atomic ensembles A81 ; DAS94 ; MZHCV15 ; LMKMIW96 and entanglement STV12 ; SV09 ; SMBBS18 have been successfully characterized using quasiprobability distributions. Nevertheless, there remains a dependency on well-defined detection schemes and reconstruction algorithms.

In this contribution, we circumvent the reconstruction problem by devising a measurement protocol that results in detector-agnostic phase-space (DAPS) distributions, which can be directly estimated, encompass known quasiprobabilities, and apply to arbitrary quantum states of light. We demonstrate our scheme with transition-edge sensors (TESs), which have sophisticated physics underlying their operation, and analyze our data without relying on any specific detector models. Our DAPS functions reveal nonclassical features expected from our heralded multiphoton states with high statistical significance. Moreover, the measurement of vacuum alone enables us to predict the unique structures of DAPS distributions as demonstrated for our experimentally generated states.

Theory framework.—

Our measurement scheme is a combination of unbalanced homodyning WV96 and a multiplexed detection layout PTKJ96 ; see Fig. 1. A signal light field, $\hat{\rho}$, is mixed with a local oscillator (LO), $|\beta\rangle$, on a beam splitter. One of the output states, represented through $\hat{\rho}(\beta)$, is injected into a multiplexing scheme that consists of $S$ steps. In each step, light is split into output fields with the same intensity, which then can be split again. The finally obtained $N=2^{S}$ output beams are individually measured with unknown detectors, which are not specified but assumed to operate in the same manner. Each detector returns one of the possible outcomes $\mathcal{K}=\{0,\ldots,K\}$. In Refs. SCEMRKNLGVAW17 ; SECMRKNLGWAV17 , we have shown for the multiplexing part that independently of the detector response, the probability to simultaneously measure $N_{k}$ times the outcome $k$ ($\forall k\in\mathcal{K}$) follows a quantum version of the multinomial distribution; its generalization to $\hat{\rho}(\beta)$ reads

[TABLE]

where ${:}\cdots{:}$ denotes the normal ordering and $\{\hat{\pi}\}_{k\in\mathcal{K}}$ is the unknown positive operator-valued measure of the detectors.

The only assumptions made are a balanced splitting in the multiplexing and identical response functions for the $N$ detectors, including all imperfections. We can account for deviations from both assumptions by including a systematic error, directly estimated from asymmetries in the measured data; see the Supplemental Material (SM) for details supplement .

A probability distribution is entirely characterized through its generation function, which can be expressed as

[TABLE]

for $z_{0},\ldots,z_{K}\in\mathbb{R}$. The second line is a result of the multinomial form of the statistics in Eq. (1). One salient feature is that classical light fields have a nonnegative generation function $g_{z_{0},\ldots,z_{K}}$. To see this, first recall that a classical light field is described through a nonnegative Glauber-Sudarshan distribution G63 ; S63 , which is not affected by displacements and describes a state as a statistical mixture of coherent states. Furthermore, for all even $N$, we can define the operator $\hat{f}=\hat{f}^{\dagger}=\left(z_{0}\hat{\pi}_{0}+\cdots+z_{K}\hat{\pi}_{K}\right)^{N/2}$. Since for any nonnegative Glauber-Sudarshan function $\langle{:}\hat{f}^{\dagger}\hat{f}{:}\rangle\geq 0$ holds true TG65 ; M86 ; VW06 ; A12 , we conclude

[TABLE]

A violation of this inequality certifies the nonclassicality of the signal light, $\hat{\rho}$. We can also define a special case of this generating function,

[TABLE]

Similarly to the expression in Eq. (2), $G_{z}$ is straightforwardly estimated from the measured detector outcomes $c_{N_{0},\ldots,N_{K}}(\beta)$ by setting $z_{k}=z^{k}$, and $G_{z}$ is nonnegative for classical light.

As an example, we may consider photocounting KK64 . Although this model is not required for our approach and does not apply to our experiment (TESs have a finite photon-number resolution, a non-unit detection efficiency, and a nonlinear response function SECMRKNLGWAV17 ), it demonstrates how $G_{z}$ generalizes the concept of well-known phase-space distributions. For photocounting, we find supplement ; WM64

[TABLE]

with $\eta$ and $\hat{n}$ being the efficiency and the photon-number operator, respectively, and $s=1-2/[\eta|t|^{2}(1-z)]$ commentsParameter . Thus, $G_{z}(\beta)$ resembles the $s$-parametrized distributions $P(r\beta/t;s)$ CG69 ; AW70 . Beyond photoelectric detectors, we refer to $g_{z_{0},\ldots,z_{K}}$ and $G_{z}$ as DAPS distributions as Eqs. (2) and (4) apply without any knowledge of the measurement operators $\{\hat{\pi}_{k}\}_{k\in\mathcal{K}}$. In this context, it is worth emphasizing that the first line in Eq. (2) enables the estimation of our DAPS distributions as a result of the measured coincidence statistics $c_{N_{0},\ldots,N_{K}}(\beta)$ alone.

Implementation.—

By implementing a single multiplexing step, $N=2^{S}=2$ for $S=1$, we demonstrate how to apply our theoretical framework of DAPS distributions. To realize our protocol in Fig. 1, we produce heralded photon states $\hat{\rho}$ and different LO amplitudes $\beta$. The detectors used for the multiplexing measurement and the heralding are TESs, which count photons up to a maximal number $K$. In the following, we describe the experimental setup supplement .

Femtosecond pulses with a $100\,\textrm{kHz}$ repetition rate from a titanium sapphire laser are coupled into two separate, periodically poled potassium titanyl phosphate (ppKTP) waveguides. Both pulses are filtered to a full-width at half-maximum of $\pm 2\,\textrm{nm}$ using angle-tuned bandpass filters. With the first ppKTP waveguide, prepare the signal $\hat{\rho}$. When filtering the pump at $775\,\textrm{nm}$ the waveguide produces two-mode squeezed vacuum in approximately a single spatio-temporal mode via type-II parametric down conversion (PDC) ECMS11 . The signal mode at $1554\,\textrm{nm}$ and the herald mode at $1547\,\textrm{nm}$ are separated with a polarizing beam splitter, then filtered and coupled into single-mode optical fibers. The herald mode is then sent to a single TES detector. With the second ppKTP waveguide, we prepare the LO. In contrast to the signal state generation, we filter the pump at $783\,\textrm{nm}$ and stimulate the PDC process by seeding it with $2\,\mathrm{ns}$ pulses carved with an electro-optic modulator from a $1580\,\textrm{nm}$ continuous-wave laser. Because of the strong seed signal, this nonlinear mixing generates coherent light to an excellent approximation in the polarization mode orthogonal to the seed LS13 . We separate the LO from the seed with a polarizing beam splitter, then pass through a bandpass filter at $1554\,\textrm{nm}$. By pumping the two waveguides at different wavelengths, we are able to create an LO that is well mode-matched to the signal using a seed laser that is detuned from the heralding mode. This avoids a potential source of noise due to the seed laser passing through the filters for the heralding TES. The generated LO is attenuated to the single-photon level and coupled into single-mode optical fiber. Crucially, this process prepares an LO with Poissonian photon statistics with a measured second-order correlation function $g^{(2)}(0)$ of $1.005\pm 0.002$.

Finally, the LO $|\beta\rangle$ and signal $\hat{\rho}$ are combined on a $90:10$ fiber beam splitter. We consider the port that uses $|r|^{2}=10\%$ of the LO and transmits $|t|^{2}=90\%$ of the signal. The light from this port, $\hat{\rho}(\beta)$, is then impinged on a $50:50$ fiber beam splitter for realizing a multiplexing step; both outputs are then sent to two separate TESs. See Fig. 1.

Our experiment uses three TES detectors that can have efficiencies above $\eta=90\%$ LMN08 . TESs are superconducting photon-number-resolving detectors that we operate in a dilution refrigerator at a temperature of around $80\,\mathrm{mK}$. Their response is amplified using an array of superconducting quantum interference devices WM91 , followed by further amplification and filtering at room temperature. This electrical signal is read by an analogue-to-digital converter and processed using a matched filter technique FCMPSW00 , which outputs a single value when triggered by a clock signal from the laser. We bin these values to assign a photon number. It should be noted that it is possible to extract slightly more accurate estimates of photon number, however, using more sophisticated signal processing techniques, yet without affecting the applicability of the DAPS distribution approach SECMRKNLGWAV17 ; HMGHLNNDKW15 .

We record the binned outcome at all three TESs for various LO amplitudes ($|\beta|^{2}$ from [math] to ${\sim}28$ in steps of ${\sim}1$). The amplitude is controlled by varying the seed laser power. To obtain data for a specific heralded state $\hat{\rho}$, we consider the subset of trials with the appropriate detection outcome (i.e., heralding bin $k_{h}$) at the herald TES.

Verification of nonclassicality.—

In a first step, we apply our DAPS distribution to uncover nonclassical features of our prepared states through the violation of condition (3). The optimal negativity we obtain from the DAPS function [Eq. (2)] is given by the minimum

[TABLE]

To assess the quality of this approach, we compared our verification of nonclassicality with other methods. In Ref. SCEMRKNLGVAW17 , we demonstrated that a correlation matrix, $M$, obtained from the measured statistics in Eq. (1), is positive semidefinite for classical light, described through a nonnegative minimal eigenvalue $\mu_{\min}$ of $M$. The resulting notion of sub-multinomial light, $\mu_{\min}<0$, was shown to be a better figure of merit than other means of verifying nonclassicality SECMRKNLGWAV17 , such as sub-Poisson light M79 ; SM83 and sub-binomial light SVA12 ; BDJDBW13 .

The comparison of $g_{\min}$ and $\mu_{\min}$ for our data is shown in Table 1 for different heralding bins $k_{h}$. For the heralded one-photon (two-photon) states, we confirm $g_{\min}<0$ with $9$ ($6$) standard deviations, while the sub-multinomial behavior is less significant, $5$ ($3$) standard deviations. For the vacuum state, i.e., $k_{h}=0$, both measures are consistent with the classical expectation, $g_{\min}=0=\mu_{\min}$.

Reconstructed distributions.—

From the data, we can directly estimate our DAPS distributions. The results of our extended analysis are shown in Fig. 2. The estimation procedure is this: we run the experiment twice, once with the signal blocked and once with the signal unblocked.

To have full detector-agnostic approach, we first define a detector-independent coherent amplitude,

[TABLE]

which is given by the statistics $c^{\mathrm{(vac)}}_{N_{0},\ldots,N_{K}}(\beta)$ measured by blocking the signal supplement . In case of photocounting, this gives $|\beta^{\mathrm{(DI)}}|=\sqrt{\eta}|r||\beta|$. As we do not record a phase, we consider full phase randomization. This does not affect the DAPS distributions of our heralded photon states. In Fig. 2, our DAPS distributions $G_{z}$ are shown as a function of the amplitude in Eq. (7), determined by means of the vacuum measurement.

The same measurement renders it possible to theoretically predict the DAPS distribution of arbitrary states. Namely, a general DAPS distribution can be described as a convolution of the measured vacuum distribution $G_{z}^{\mathrm{(vac)}}$ and the Glauber-Sudarshan distribution $P(\beta^{\prime};1)$ of the state under study supplement ,

[TABLE]

In our case, the Gaussian shape of $G_{z}^{\mathrm{(vac)}}$ implies that heralded single-photon (two-photon) states should follow a Gaussian distribution multiplied with a first-order (second-order) polynomial in $|\beta^{\mathrm{(DI)}}|^{2}$. In Fig. 2, this prediction (dashed lines) is confirmed as it correctly represents the DAPS distributions of the measured heralded photon states. The heralding to $k_{h}=1$ gives a characteristic dip at the origin $|\beta^{\mathrm{(DI)}}|=0$, and the two-photon case, $k_{h}=2$, leads to additional oscillations together with the appearance of a peak at the origin. We emphasize that the functional behavior $\beta\mapsto G_{z}(\beta)$ depends on the measurement operators, but the estimation of $G_{z}(\beta)$ is done without any specification of the detector operators, according to the first line in Eq. (2). Moreover, we are able to characterize defining features of other states without any other prior knowledge about the detectors from the data obtained using the vacuum state input [Eq. (8)].

Based on our reconstruction, we were able to determine a number of other properties of the experimentally produced states supplement . For instance, we can determine how well the DAPS distributions enable us to perform a quantum state discrimination task. The single- and two-photon states [plots (b) and (c) in Fig. 2] can be distinguished from each other with more than $98\%$ certainty. Furthermore, we found that, for $z<-2.4$, the central dip of $G_{z}$ becomes negative, similar to the behavior of other phase-space quasiprobabilities. The negativity has the highest statistical significance for $z=-4.85$, where $G_{z}(0)=-0.51\pm 0.08$ is more than 6 standard deviations below the classical threshold of zero.

Summary and discussion.—

We have developed a theory and realized an experiment to characterize quantum light by means of phase space that functions for any type of detector and without performing a prior detector characterization. Our framework is based on the generating function derived from the properties of a balanced linear optical network, enabling our DAPS distribution, to be directly estimated from measured correlations and that applies to arbitrary states.

To demonstrate this concept, we showed that a single multiplexing step is already sufficient for applying our method. This renders it possible to verify the nonclassicality of multiphoton states based on DAPS distributions, which results in greater statistical significance than obtained with earlier approaches SCEMRKNLGVAW17 ; SECMRKNLGWAV17 , which themselves already outperformed previous quantifiers of nonclassicality. Moreover, our approach encompasses prominent phase-space quasiprobabilities and straightforwardly generalizes to multimode light.

Our general theory also includes more recent phase-space functions based on on-off detectors LSV15 , constituting the special case $K=1$ and being applicable to off-the-shelf detectors (e.g., avalanche photodiodes in Geiger mode and single-photon nanowire detectors); to prove this, see the corresponding experiment with $S=1,2,3$ multiplexing steps BTBSSV18 . Furthermore, recent advancements in detector technology (see, e.g., Refs. DetNew1 ; DetNew2 ; DetNew3 ; DetNew4 ) offer new photon counters to which our detector-agnostic framework is also readily applicable.

With our approach, we are further able to predict defining phase-space features of any states by measuring vacuum as a reference. Thus, we can compare a target state with the actually reconstructed DAPS distribution, thus enabling us to estimate other quantum properties as well. As a practical example, a state discrimination task based on our DAPS distributions resulted in distinguishing one- and two-photon states with almost unit certainty, despite high losses in our setup. It is also worth noting that our DAPS distribution includes the full quantum information of the state that is accessible with the detectors used and does not require demanding reconstruction algorithms and detection models.

Our experiment comprises state-of-the-art detectors combined with an advantageous method to create coherent states, well mode matched to our signal. As our method is detector agnostic, the detector efficiency need not be specified or even known; the number of data points merely has to be sufficient to produce statistically meaningful results. Also, our approach is not restricted to any specific states; currently, we are mainly limited by the available sources of nonclassical light.

In the future, recording the LO’s phase would be beneficial for applying our scheme to phase-sensitive nonclassical states as well. Furthermore, generalizing other interferometric measurement schemes in a detector-agnostic manner is feasible, e.g., as done for on-off detectors SVA15 . In addition, we encounter the imperfections stemming from imbalances by assigning systematic errors. It may be possible to avoid this by using more sophisticated strategies LFPR16 .

Generalized phase-space distributions are becoming increasingly important in identifying vastly different notions of quantumness; see Refs. SW18 ; SV19 for thorough overviews. To date, however, such universally applicable techniques are also highly dependent on the particular response of the detectors. Our DAPS approach, however, sets a precedence that such limitations can be overcome in theory and experiment.

In conclusion, our detector-agnostic framework provides a universally applicable approach to the robust characterization of quantum light in phase space under conditions where detailed knowledge of the measurement apparatus is not available, and forms a basis for future research.

Acknowledgments.—

The authors are grateful to William R. Clements for helpful discussions and Jelmer J. Renema for his assistance with the installation of the cryogenic infrastructure. The authors also thank Scott Glancy, Arik Avagyan, and Tim Bartley for valuable comments. The Integrated Quantum Optics group acknowledges financial support from the Gottfried Wilhelm Leibniz-Preis (Grant No. SI1115/3-1). This work received funding through the Networked Quantum Information Technologies (NQIT) hub (part of the UK National Quantum Technologies Programme) under Grant No. EP/N509711/1. G. S. T. acknowledges financial support from the Natural Sciences and Engineering Research Council of Canada and the Oxford Basil Reeve Graduate Scholarship. A. E. is supported by EPSRC (project EP/K034480/1 BLOQS). T. A. W. W. is supported by Fondation Wiener - Anspach. J. L. thanks the European Commission (H2020-FETPROACT-2014 grant QUCHIP). I. A. W. acknowledges ERC (Advanced Grant MOQUACINO). This work was supported by the Quantum Information Science Initiative (QISI).

Supplemental Material

Appendix A Details on the experiment

The experimental setup is shown in Fig. 3. Our pump laser is a titanium sapphire (Ti:Saph) regenerative amplifier that generates femtosecond pulses (center wavelength $780\,\mathrm{nm}$, full width at half maximum [FWHM] $15\,\mathrm{nm}$) at a rate of $100\,\mathrm{kHz}$. This rate is chosen to accommodate the thermal relaxation time (${\sim}10\,\mathrm{\mu s}$) of the transition-edge sensors (TESs). We split the pulses into two paths, each pumping a periodically poled potassium titanyl phosphate (ppKTP) waveguide.

We pump the first waveguide (right in Fig. 3) using filtered (center $775\,\mathrm{nm}$, FWHM $2\,\mathrm{nm}$) pulses from the Ti:Saph. Pumping the ppKTP waveguide generates two-mode squeezed vacuum via type-II parametric down-conversion. The pump is then discarded using a longpass filter. The two down-converted modes (signal $1554\,\mathrm{nm}$, idler $1547\,\mathrm{nm}$) are orthogonally polarized and separated by a polarisation beam splitter. Each mode is sent through a bandpass filter (FWHM $10\,\mathrm{nm}$). The idler mode is sent to a TES detector to herald photon-number states in the signal mode (heralding efficiency ${\sim}40\%$) by postselecting to a specific outcome $k_{h}$.

In the second waveguide (left in Fig. 3), we prepare the local oscillator (LO). Since the first and second waveguides have slightly different phase-matching properties, a different pump spectrum (center $783\,\mathrm{nm}$, FWHM $2\,\mathrm{nm}$) is used. This pump spectrum is chosen to maximize the spectral overlap between the LO and signal (SI). We also carve $2\,\mathrm{ns}$ square seed pulses from a continuous-wave laser (center $1580\,\mathrm{nm}$) using an electro-optic modulator. The pump and seed pulses are temporally overlapped and coupled into the second ppKTP waveguide. Through difference frequency generation, the LO is generated in the polarization orthogonal to the seed. The LO is separated from the seed using a polarisation beam splitter. As before, we discard the pump by a longpass filter. The LO’s polarization is adjusted with a half-wave plate to match the SI’s polarization. Then, the LO is sent through a bandpass filter (FWHM $10\,\mathrm{nm}$) to further eliminate seed light as well as increase the spectral overlap with the signal. Finally, neutral-density filters attenuate the LO to the single-photon level.

The SI and LO are combined on a 90:10 beam splitter. The resulting light field of one output then enters the multiplexing step (50:50 beam splitter) and the then resulting beams are measured with two TESs. The recorded coincidences give the detection events $E(k_{1},k_{2})$, which we use for our analysis.

In addition, we characterized the mode overlap of the signal and LO by combining the two on a 50:50 beam splitter. We consider the specific case of a single photon ($k_{h}=1$) and a weak LO ($|\beta|\ll 1$). By scanning the delay between the SI and LO, we expect to measure a Hong-Ou-Mandel-type dip in two-fold coincidences at the output of the beam splitter. We measured a dip of ${\sim}80\%$ visibility, suggesting that the mode overlap is at least $80\%$. By blocking the signal, this setup constitutes a Hanbury Brown-Twiss interferometer that allows us to measure the LO’s $g^{(2)}(0)$. We measured $g^{(2)}(0)=1.005\pm 0.002$, which is consistent with the expected Poisson distribution for the LO’s photon statistics.

Appendix B Details on the theory

B.1 General approach

Let us formulate some additional details on the theory. As we can expand any SI state in the Glauber-Sudarshan decomposition, $\hat{\rho}=\int d^{2}\alpha\,P(\alpha)|\alpha\rangle\langle\alpha|$, it is sufficient to consider the propagation of coherent states $|\alpha\rangle$. Our detection scheme consists of a combination of the SI with the LO state $|\beta\rangle$ on a beam splitter, the multiplexing, and the detection.

Applying a beam splitter transformation, we map an input, consisting of SI and LO, as follows: $|\alpha\rangle\otimes|\beta\rangle\mapsto|t\alpha-r\beta\rangle\otimes|r^{\ast}\alpha+t^{\ast}\beta\rangle$, where $t$ and $r$ define the trasmissivity and reflectivity ($|t|^{2}+|r|^{2}=1$). When tracing over the second mode, we obtain the state that enters the multiplexing stage,

[TABLE]

Further, the multiplexing distributes the coherent state components in Eq. (9) among the $N=2^{S}$ output beams, where $S$ is the depth of the multiplexing scheme, $|\gamma\rangle\mapsto|\gamma/\sqrt{N}\rangle^{\otimes N}$, resulting in $\int d^{2}\alpha\,P(\alpha)[|(t\alpha-r\beta)/\sqrt{N}\rangle\langle(t\alpha-r\beta)/\sqrt{N}|]^{\otimes N}$. Using some combinatorics (see Ref. SECMRKNLGWAV17 for details), we find

[TABLE]

where $\{\hat{\pi}_{k}\}_{k=0,\ldots,K}$ is an unknown positive operator-valued measure (POVM) that describes the detector.

Consequently, the two types of generating functions under consideration read

[TABLE]

using the Glauber-Sudarshan $P$ function. As long as $N$ is even and $P\geq 0$, both expressions are necessarily nonnegative.

B.2 Photoelectric counting with loss

Let us analyze our scheme for the special case of photoelectric counting. A simple photoelectric detection is described through POVM elements $\hat{\pi}_{k}={:}e^{-\eta\hat{n}}(\eta\hat{n})^{k}/k!{:}$ for $k=0,1,\ldots$ ($K=\infty$), where $\eta$ is the quantum efficiency and $\hat{n}$ is the photon-number operator. In this scenario, the generating function in Eq. (12) can be further evaluated WM64 and reads

[TABLE]

where $\hat{n}(\gamma)$ is the displaced photon-number operator. Since we have $P(\gamma;s)=2(\pi[1-s])^{-1}\langle{:}e^{-2\hat{n}(\gamma)/[1-s]}{:}\rangle_{\hat{\rho}}$ WV96 , the above expression can be related to $s$-paramterized distributions.

Also note that according to the characterization performed in Ref. SECMRKNLGWAV17 (Sec. III), our TESs are more precisely described through POVMs of the form $\hat{\pi}_{k}={:}e^{-\Gamma(\hat{n})}\Gamma(\hat{n})^{k}/k!{:}$ for $k=0,\ldots,K-1$ and $\hat{\pi}_{K}=\hat{1}-\sum_{k=0}^{K-1}\hat{\pi}_{k}$. Therein, $K<\infty$ reflects the finite photon-number resolution, and the response function $\Gamma$ has a nonlinear form, $\Gamma(\hat{n})\approx\eta\hat{n}+\eta^{(2)}\hat{n}^{2}$, where the quantum efficiency $\eta$ is not one ($\eta<1$) and the nonlinear contribution does not vanish ($\eta^{(2)}\neq 0$). As $\eta^{(2)}$ is small, the nonlinear behavior mainly affects higher LO and SI intensities.

Appendix C Coincidences and systematic errors

A single multiplexing step, $N=2$, was implemented. Thus, it is convenient to formulate the data processing in terms of measured coincidences. For this purpose, we denote with $E(k_{1},k_{2})$ the number of coincidence events for the measurement outcomes $k_{1}$ and $k_{2}$ ($k_{1},k_{2}\in\{0,\ldots,K\}$), resembling the detection bins of the TESs $1$ and $2$, respectively. $E=\sum_{k_{1},k_{2}}E(k_{1},k_{2})$ defines the total number of events.

The coincidences are directly related to the desired quantum version of a multinomial distribution, $c_{N_{0},\ldots,N_{K}}$, cf. Eq. (1) in the main text. Since $N_{0}+\cdots+N_{K}=N=2$, we have

[TABLE]

where we can identify $c_{0,\ldots,0,N_{i}=2,0,\ldots,0}=E(i,i)/E$ and $c_{0,\ldots,0,N_{i}=1,0,\ldots,0,N_{j}=1,0,\ldots,0}=(E(i,j)+E(j,i))/E$. To estimate the value $\overline{f}$ of a function $f_{N_{0},\ldots,N_{K}}$, we can recast the standard sampling formula as follows:

[TABLE]

See the Supplemental Material to Ref. SCEMRKNLGVAW17 for the generalization to $N>2$.

In order to apply the multinomial framework described in the main text, one has to satisfy the premise that the coincidence statistics is symmetric, $E(k_{1},k_{2})=E(k_{2},k_{1})$. However, in reality, this is only true to a limited extent since the beam splitters in the multiplexing might not be perfectly balanced, and the detectors after the multiplexing might have slightly different responses. In Ref. SCEMRKNLGVAW17 , we provided a rough systematic error estimate to account for such imperfections, which is further refined in the following.

The premises mentioned above state that the coincidences are symmetric. The actual measurements $E(k_{1},k_{2})$ naturally exhibit a certain amount of asymmetry; if not, no systematic error needs to be assigned. We can decompose the coincidences as follows:

[TABLE]

Furthermore, assume we estimate a function $f(k_{1},k_{2})$ to obtain the mean $\overline{f}=\sum_{k_{1},k_{2}}f(k_{1},k_{2})E(k_{1},k_{2})/E$. Inserting the above decomposition and denoting with $\overline{f}^{\mathrm{(sym)}}$ the value obtained from the symmetric part in Eq. (19), we apply the triangle inequality and find

[TABLE]

which is the systematic error resulting from the asymmetry in the measured data.

As we use the typical quadratic error propagation—rather than the linear form used for the above derivation—, we replace the right-hand-side expression in Eq. (20) with $\epsilon_{f}^{2}=\sum_{k_{1},k_{2}}|f(k_{1},k_{2})|^{2}\left|[E(k_{1},k_{2})-E(k_{2},k_{1})]/[2E]\right|^{2}$. Recall that the general relation between linear and quadratic error expansion for a function $F(x_{1},x_{2},\ldots)$ is given by $\Delta^{\mathrm{(lin)}}F=\sum_{j}|\partial F/\partial x_{j}|\Delta x_{j}$ and $\Delta^{\mathrm{(quad)}}F=(\sum_{j}|\partial F/\partial x_{j}|^{2}[\Delta x_{j}]^{2})^{1/2}$. In addition, let us remind ourselves that the random error reads $\sigma_{f}=[(\overline{f^{2}}-\overline{f}^{2})/(E-1)]^{1/2}$, which is combined with the systematic error to give the overall uncertainty, $\Delta f=\sqrt{\epsilon_{f}^{2}+\sigma_{f}^{2}}$.

Appendix D Sub-multinomial light

We assess our results with our previously derived nonclassicality criteria SCEMRKNLGVAW17 ; SECMRKNLGWAV17 . Let us briefly recapitulate this approach and its implementation for a self-consistent reading. The previously devised method is based on the observation that a correlation matrix, $M$, is positive semidefinite for classical light, i.e., $M=(M_{i,j})_{i,j=0,\ldots,K}\geq 0$, where

[TABLE]

with $\delta$ denoting the Kronecker symbol. It was shown that the required first- and second-order moments can be obtained from coincidences as SCEMRKNLGVAW17

[TABLE]

Finally, the minimal eigenvalue $\mu_{\min}$ of the correspondingly reconstructed matrix $M$ is computed to probe for positive semidefiniteness.

The heralding with a TES enables us to generate higher-order photon-number states. In Table 2, we listed the nonclassicality in terms of the criteria $\mu_{\min}<0$ for data with a coherent amplitude zero. The observed nonclassicality in Table 2 for heralding bins larger than two is no longer significant within a three-standard-deviation error margin as the total number of events $E$ is too small in those cases. For this reason, we restrict our considerations to heralding bins $k_{h}\in\{0,1,2\}$.

Appendix E Local oscillator amplitudes

From the first derivative of Eq. (13) for the photoelectric model, we can infer the dimensionless and displaced intensity,

[TABLE]

where the latter expression results from the definition $G_{z}(\beta)=\sum_{N_{0},N_{1},\ldots}z^{0N_{0}+1N_{1}+\cdots}c_{N_{0},\ldots,N_{K}}(\beta)$. Most importantly, in the case that the SI is vacuum, we find $\langle{:}\hat{n}(\gamma){:}\rangle_{|0\rangle\langle 0|}=|\gamma|^{2}$. We abstract this observation and define for the general detection scenario a detector-independent amplitude $|\beta^{\mathrm{(DI)}}|=\sqrt{I(\beta)}=\sqrt{\partial G_{z}(\beta)/\partial z|_{z=1}}$ for general POVMs and the SI $\hat{\rho}=|0\rangle\langle 0|$. This also explains the following Eq. (25) as well as Eq. (7) in the main text.

We block the SI to infer the amplitude of the LO. The intensity, here represented through the dimensionless quantity

[TABLE]

has been varied, realized via $n=0,\ldots,28$ power settings of the seed laser. These settings have been also applied when the SI is not blocked. The settings are chosen such that an equidistant intensity grid is generated. This is confirmed through the linear fit in Fig. 4.

Appendix F Additional analysis and results

F.1 Typical data and error estimates

In order to further assess the impact of uncertainties, let us consider a typical example. In Fig. 5, we depict a typical data set, there for a single photon (i.e., heralding to bin $k_{h}=1$) and a vanishing LO amplitude (i.e., the setting $n=0$).

The resulting (systematic and random) observational errors for the generating function are depicted in Fig. 6. The expected trend of monotonicity and the diverging behavior of the uncertainties as a function of $|z|$ are clearly visible. As $G_{z=1}(\beta)=1$ corresponds to the total probability, which is not subject to fluctuations, the random error vanishes for $z=1$ (Fig. 6, top-left panel). By construction, the systematic errors are symmetric with respect to $z=0$ (top-right plot in Fig. 6). See also the following discussion in Sec. F.2.

For probing the nonclassicality through the generating function, obtained as $g_{z_{0},\ldots,z_{K}}(\beta)=\sum_{k_{1},k_{2}}z_{k_{1}}z_{k_{2}}E(k_{1},k_{2})/E$, we compute the eigenvector $Z=[z_{0},\ldots,z_{K}]$ to the normalized coincidence matrix $[E(k_{1},k_{2})/E]_{k_{1},k_{2}=0,\ldots,K}$ that corresponds to the minimal eigenvalue. For the example under study, we get $Z=[0.228,-0.973,-0.033,-0.001,-0.000]$. With this information, we can now estimate the general generating function and get $g_{z_{0},\ldots,z_{K}}(\beta)=-0.026\pm 0.003$ as the optimal negativity, here for $\beta=0$. The minimum over all measured LO amplitudes then yields $g_{\min}$, cf. Eq. (6) in the main text.

F.2 Optimal quasiprobability distribution and error estimates

We found that the parameter $z=-4.85$ is optimal in the sense that $G_{z}(\beta)$ has the most statistically significant negativity at the origin (see Fig. 7), i.e., $-G_{z}(0)/\Delta G_{z}(0)$ is maximized. In the plotted scenario, the error estimates $\Delta G_{z}(\beta)$ are rapidly increasing for increasing $|\beta|$, which we discuss in the following based on the sampling formula

[TABLE]

For increasing LO amplitudes, the components of $c_{N_{0},\ldots,N_{K}}(\beta)$ that relate to a higher power of $z$ have a higher contribution to the estimate of this function. Similarly, $|z|>1$ also leads to a most relevant term that corresponds to a higher power of $z$. Recall that the exponent $0N_{0}+\cdots+KN_{K}$ relates to the overall intensity [cf. Eq. (25)]. Consequently, both a large LO amplitude and $|z|>1$ result in the fact that the contribution for $z^{p}$ for larger $p$ becomes the most relevant one. Standard error propagation then implies a relative error scaled by the large factor $p$, for increasing LO amplitudes and increasing $|z|$ values, which also explains why the increase of observational uncertainties in those scenarios is expected.

Appendix G Quantum state discrimination

As our distributions can, in principle, take arbitrary forms for arbitrary detectors, let us formulate the statistical model to discriminate states based on the reconstructed phase-space functions alone. The probability that two distributions, described through multivariate random variables $X$ and $X^{\prime}$, are indistinguishable within a $\delta$-uncertainty can be expressed as

[TABLE]

where $p(X)$ and $p(X^{\prime})$ are the probability densities of the uncertainties for the two random variables. We identify $X=[G^{(k_{h})}_{z}(\beta_{0}),\ldots,G^{(k_{h})}_{z}(\beta_{28})]$ and $X^{\prime}=[G^{(k^{\prime}_{h})}_{z}(\beta_{n})]_{n=0,\ldots,28}$ for different heralded states and LO settings and use a Gaussian error model with a mean and variance that corresponds to the reconstructed distributions for each measured setting $n$. Consequently, we get from Eq. (27) the following probability for the discrimination:

[TABLE]

where we set the vector $\delta$ to correspond to three combined standard deviations for each setting, $\delta_{n}=3\Delta_{n}$ and $\Delta_{n}=[\Delta G^{(k_{h})}_{z}(\beta_{n})^{2}+\Delta G^{(k_{h})}_{z}(\beta_{n})^{2}]^{1/2}$, and using the mean values $\mu_{n}=\overline{G^{(k_{h})}_{z}(\beta_{n})}$ and $\mu^{\prime}_{n}=\overline{G^{(k^{\prime}_{h})}_{z}(\beta_{n})}$. Note that $\mathrm{Err}[z]=\int^{z}d\xi e^{-\xi^{2}/2}/\sqrt{2\pi}$ denotes the error function.

For instance, we find for our measured data that the likelihood to discriminate the phase-space distributions for $k_{h}$ from the one for $k^{\prime}_{h}$ for $z=-1.5$ is given by the matrix

[TABLE]

where “$100.\%$” corresponds to a value which is $100\%$ within the used numerical precision. Note that the diagonal elements are nonzero as identical distributions could still represent different states when considering a finite error margin.

Appendix H Fit model from vacuum measurements

From the measurement in which the SI is blocked (i.e., vacuum SI), we can extrapolate the general shape of the phase-space distribution for photon states, without relying on any particular detector model. This approach is also used to fit the reconstructed distributions for the heralding to $k_{h}$.

Using the data where the signal is blocked, we find that a Gaussian distribution describes the reconstructed phase-space distribution for vacuum quite well; see Fig. 8. This information can be used to predict the phase-space distributions of $m$-photon states as well. Because of Eq. (12) and the known representation $P^{(m)}(\alpha)=\sum_{j=0}^{m}\binom{m}{j}j!^{-1}\partial_{\alpha}^{j}\partial_{\alpha^{\ast}}^{j}P^{(0)}(\alpha)$ VW06 , where $P^{(0)}(\alpha)$ describes the delta distribution centered at the origin, we find that the $m$-th photon state is given by

[TABLE]

where $G^{(0)}_{z}(\beta)$ is experimentally obtained by blocking the signal (Fig. 8) and which is determined without relying on any detection models. For deriving Eq. (30), note that the argument of the vacuum function implies that a partial integration of Eq. (12) with derivatives of delta distributions results in $\left.\partial^{j}_{\alpha}\partial^{j}_{\alpha^{\ast}}f(\beta-t\alpha/r)\right|_{\alpha=0}=(-t/r)^{j}(-t^{\ast}/r^{\ast})^{j}\partial^{j}_{\beta}\partial^{j}_{\beta^{\ast}}f(\beta)$.

From Eq. (30) and the fit obtained from $G^{(0)}_{z}(\beta)$, we can therefore predict the phase-space distribution of an $m$-photon state. In our case, this means that $G^{(m)}_{z}(\beta)$ is a Gaussian function multiplied with a fixed $m$th-order polynomial in $|\beta|^{2}$. In this context, also recall the linear relation between the actual intensity (via the setting number $n$) and the detector-independent intensity in Fig. 4. In addition, it is known (see, e.g., Ref. SVA14 ) that imperfect heralding for the kind of photon source used leads to additional noise contributions. For this reason, our fit for an heralding to the $k_{h}$th bin is described through $G^{(k_{h})}_{z}(\beta)=\sum_{j=0}^{k_{h}}f_{j}|\beta^{\mathrm{(DI)}}|^{2j}\exp[-b|\beta^{\mathrm{(DI)}}|^{2}]$, which constitutes the generalized fit function used in the main text [Figs. 2(a)–(c)] and is determined from the vacuum measurements alone and without relying on any detection models.

As a final remark, it is worth mentioning that the above treatment can be straightforwardly generalized to predict $G_{z}(\beta)$ for arbitrary states, resulting in Eq. (8) in the main text. This is based on the fact that the $P$ function of an arbitrary state can be written as a convolution, $P(\alpha)=\int d^{2}\alpha^{\prime}\,P(\alpha-\alpha^{\prime})P^{(0)}(\alpha^{\prime})$, recalling that the vacuum state is described by a delta distribution, $P^{(0)}$. Thus, Eq. (12) implies that $G_{z}$ of an arbitrary state, represented through the Glauber-Sudarshan distribution $P$, is predicted to resemble the convolution of the already measured vacuum state’s $G_{z}^{(0)}$ and $P$. Even more generally, we can write

[TABLE]

where $g^{\mathrm{(0)}}_{z_{0},\ldots,z_{K}}(\beta)=\left[\sum_{k=0}^{K}z_{k}p_{k}\left(-r\beta/\sqrt{N}\right)\right]^{N}$, cf. Eq. (11), to predict the phase-space distribution for a state, theoretically described through $P(\alpha)$, via the measured $g^{\mathrm{(0)}}_{z_{0},\ldots,z_{K}}$.