# Quantum Conformance Test

Giuseppe Ortolano1,2 Pauline Boucher1 Ivo Pietro Degiovanni1 Elena Losero1,3
Marco Genovese1 Ivano Ruo-Berchera1 1Quantum metrology and nano technologies
division, INRiM, Strada delle Cacce 91, 10153 Torino, Italy  2DISAT,
Politecnico di Torino, Corso Duca degli Abruzzi 24, 10129 Torino, Italy  3
STI-SB, Ecole Polytechnique Fédérale de Lausanne, Lausanne, 1015, Switzerland

###### Abstract

We introduce a protocol addressing the conformance test problem, which
consists in determining whether a process under test conforms to a reference
one. We consider a process to be characterized by the set of end-product it
produces, which is generated according to a given probability distribution. We
formulate the problem in the context of hypothesis testing and consider the
specific case in which the objects can be modeled as pure loss channels. We
demonstrate theoretically that a simple quantum strategy, using readily
available resources and measurement schemes in the form of two-mode squeezed
vacuum and photon-counting, can outperform any classical strategy. We
experimentally implement this protocol, exploiting optical twin beams,
validating our theoretical results, and demonstrating that, in this task,
there is a quantum advantage in a realistic setting.

##  I Introduction

Significant progress has been made in recent years in the field of quantum
sensing [1, 2] both for continuous parameter estimation [3, 4, 5], and for
discrimination tasks in the case of discrete variables [6, 7, 8]. The use of
quantum states and resources brings an advantage that has been demonstrated in
various specific tasks: both phase [9, 10, 11, 12, 13] and loss [14, 15, 16,
17, 18, 19, 20] estimation, quantum imaging [21, 22, 23], and discrimination
protocols such as target detection [24, 25, 26, 27] and quantum reading [28,
29].

In this class of problems, a parameter of interest is encoded in some quantum
states, or channels, and the values it can take can belong to either discrete
or continuous sets. The estimation of this parameter requires the choice of a
probe, as well as a measurement scheme, i.e. a measurement strategy. Here, we
consider a significant discrimination problem that can find various important
applications, the quantum conformance test. In this problem, one wants to
assess whether a process is conform to a reference, or if it is defective. In
the general case, the process is characterized by a physical parameter,
distributed according to some continuous probability density distribution. We
experimentally have access to the end-products of this process, which we can
perform measurements on. Accordingly, the measurement outcomes belong to a
continuous set of possible values. However, the expected output of the
procedure should be binary: conform or not conform. Such conformity tests
appear frequently in many applications [30], one example being product safety
testing.

Under energy constraints for the probe, quantum mechanical fluctuations set
lower bounds to the probability of error that depend on the decision strategy.
It is important to investigate whether, and to what extent, the use of certain
quantum resources can reduce the probability of error below what is possible
in the classical domain.

In this work, we introduce a formal description of the conformance test
problem in the context of quantum information [31] and consider the
paradigmatic example of bosonic pure loss channels probed with light. We
demonstrate that, under the same energy constraints, i.e. fixing the photon
number, a quantum strategy making use of entangled photons as probe states,
and photon-counting (PC) as of measurement strategy, can deliver better
results than any classical strategy. Finally, we present an experimental
optical implementation of our proposed quantum protocol, which validates our
theoretical model and shows that a genuine quantum advantage persists even in
presence of experimental inefficiencies.

###  I.1 Quantum Conformance test

Figure 1: _Examples of error probabilities for the conformance test._ In
panels A-B, the distributions of the reference process (red),
$\mathcal{P}_{0}$, and defective one (blue), $\mathcal{P}_{1}$, are shown: the
two rows present opposite and archetypal situations. The graphics (C-D-E-F)
are obtained considering a pure loss channel of parameter $\theta=\tau$ as the
SUT, and a PC measurement. In panels C-D, the resulting photon number
distribution $p(n)$, in the case where a classical state is used as a probe,
are shown. The overlaps between the two distributions, highlighted in green
and blue, the two colors distinguishing the cases where the reference or the
defective process are most likely, are visualizations of the conditional error
probabilities, at a given value of $\tau$, while their weighted sum gives an
appreciation of the total error probability. Panels E-F display the case in
which quantum probes are used: quantum correlations enhance the performance of
the discrimination task.

The quantum conformance test (QCT) can be modeled as follows. We define the
binary random variable $x\in\\{0,1\\}$, which corresponds to the process from
which the physical system under test (SUT) was generated, $0$ being the
reference and $1$ the defective one. We consider the monitoring of a physical
process $\mathcal{P}_{x}$, producing a quantum object (the SUT),
$\mathcal{E}_{\theta}$, which depends on a parameter $\theta$. The process
$\mathcal{P}_{x}$ can be described by the ensemble
$\\{g_{x}(\theta),\mathcal{E}_{\theta}\\}$: the parameter $\theta$ is
extracted from a set $\mathcal{A}$, according to the probability distribution
$g_{x}(\theta)$, and it defines the physical object $\mathcal{E}_{\theta}$.
The set $\mathcal{A}$ can either be discrete or continuous. The conformance
test consists in ruling whether an unknown process is conform to a “reference”
process, $\mathcal{P}_{0}$, or if it should be labeled “defective”,
$\mathcal{P}_{1}$, using measurements on the set of objects it produces
$\\{\mathcal{E}_{\theta}\\}$.

A QCT is performed using a probe in a generic quantum state $\rho$. After the
probe has interacted with the SUT, the final state is measured by a positive-
operator-valued measure (POVM) $\Pi$ and, after classical post-processing of
the measurement results, the outcome of the procedure is a final guess on the
nature of the production process, expressed by the binary variable
$y\in\\{0,1\\}$.

The test is successful when $y=x$, i.e. when the guess is correct. On the
other hand, if $y\neq x$, the test fails. Two cases can be distinguished:

  * •

False negative. In this scenario, a SUT produced by a conform process ($x=0$)
is labeled as defective ($y=1$). In an industrial context, this kind of
outcome can be seen as an economic loss for a manufacturer, as a conform
process is considered defective. We will denote the probability of false
negatives as $p_{10}$.

  * •

False positive. A SUT produced by a defective process ($x=1$) is labeled as
conform ($y=0$). This outcome represents a risk since possibly unsafe products
are released. The false-positive probability will be referred to as $p_{01}$.

The analysis of false positives and negatives plays a central part in
conformity testing, and the specific choice, if required, of the tolerance on
either one of those errors vastly depends on the situation considered. In a
general scenario a significant figure of merit – when the energy of the probe
is considered a limited resource and, therefore one compares the results at
fixed energy — that can be used to assess the effectiveness of the test is the
total probability of error, defined as:

$p_{err}=\frac{1}{2}\big{(}p_{01}+p_{10}\big{)}$ (1)

In Fig.(1), we present some examples of error probabilities for the
conformance test. In panels A-B, two possible scenarios for the conform and
defective processes are shown. In panel A, we consider two distributions whose
overlap is negligible, whereas, in panel B, we consider significantly
overlapping distributions. In panels C-D, we show the probability
distributions for the parameter, convoluted with the noise emerging from both
the measurement process and the probe state. From these distributions, the
outcome $y$ must be decided. One can observe that the noise present in the
state which probes the SUT translates into an error in the discrimination.

In the subsequent parts of this work, we focus on a specific class of SUTs:
the bosonic pure loss channels. We show how, in this case, quantum resources
can be used to significantly mitigate the spurious effect of the noise, and
can improve the discrimination performance – evaluated in terms of total
probability of error – over any classical strategy. Subsequently, we analyze
the case in which a constraint is imposed either on the false positives or the
false negatives probabilities.

We note that the quantum reading problem, originally proposed in [28] and
experimentally tested in [29], can be considered to be a special case of the
conformance test (i.e. in the case where $g_{0}$ and $g_{1}$ are both Dirac
delta functions).

##  II Results

###  II.1 Pure loss channels

A bosonic loss channel $\mathcal{E}_{\tau}$ is characterized by its
transmittance $\tau\in\left[0,1\right]$. The input-output relation for such a
channel is $\hat{a}\rightarrow\sqrt{\tau}\,\hat{a}+i\sqrt{1-\tau}\,\hat{v}$,
with $\hat{a}$ the annihilation operator of the input mode [32]. The
annihilation operator of the mode on the second port, $\hat{v}$, acts on the
vacuum, as we are considering pure losses. The configuration for a conformance
test over a pure loss channel is described in Fig.(2). A transmitter
irradiates an optical probe state $\rho$ on the SUT. In a general case, the
state is bipartite, having $M$ signal modes and $L$ idler ones. The state is
measured at the receiver by a joint measurement, and its outputs are processed
to obtain the final outcome.

TransmitterReceiverM (signal)L
(idler)$\mathcal{E}_{\tau}$$\mathcal{I}$$\rho$$\Pi$input statePOVMDP

Figure 2: _Quantum conformance test:_ A probe $\rho$ sends $M$ signal modes
through the loss channel $\mathcal{E}_{\tau}$, while $L$ idler modes directly
reach the receiver – $\mathcal{I}$ represents the identity operator. A
positive-operator-valued measure $\Pi$ is applied to the output state. Using
the result of this measurement and data processing (DP), a decision $y$ is
taken: the process generating $\mathcal{E}_{\tau}$ is identified as conform
($y=0$) or defective ($y=1$).

After the interaction with the SUT, an input state $\rho$ will be mapped into
either $\rho_{0}$ or $\rho_{1}$, with

$\displaystyle\rho_{0}$
$\displaystyle=\mathbb{E}_{\mathcal{P}_{0}}\left[\left(\mathcal{E}_{\tau}\otimes\mathcal{I}\right)\rho\right]\rightarrow\text{reference
process}$ $\displaystyle\rho_{1}$
$\displaystyle=\mathbb{E}_{\mathcal{P}_{1}}\left[\left(\mathcal{E}_{\tau}\otimes\mathcal{I}\right)\rho\right]\rightarrow\text{defective
process}$ (2)

and where $\mathbb{E}_{\mathcal{P}}[\cdot]$ represents the expectation value
over the ensemble $\mathcal{P}$.

In general, the output states of the reference and defective processes,
$\rho_{0}$ or $\rho_{1}$, will overlap. Accordingly, the discrimination
between the processes will be affected by an error probability
$p_{err}^{\rho,\Pi}(\mathcal{P}_{0},\mathcal{P}_{1})$ which is a function of
the processes considered, the input state and the POVM – as well as the
decision procedure applied to the measurement result. The optimization of the
QCT protocol is achieved by a minimization of
$p_{err}^{\rho,\Pi}(\mathcal{P}_{0},\mathcal{P}_{1})$ over all possible input
states $\rho$, and POVMs $\Pi$:

$p_{err}(\mathcal{P}_{0},\mathcal{P}_{1})=\min_{\Pi}\bigg{[}\min_{\rho}\big{[}p_{err}^{\rho,\Pi}(\mathcal{P}_{0},\mathcal{P}_{1})\big{]}\bigg{]}$
(3)

Without constraints on the energy, a trivial strategy is to let the energy of
the probe system go to infinity, which, given a suitable measurement, would
nullify the quantum noise and lead to the minimum possible $p_{err}$ permitted
by the problem, i.e. the overlap between the two initial distributions $g_{0}$
and $g_{1}$ (see Fig.(1)). Nevertheless, in several significant cases, one
cannot dispose of, or use, states of arbitrarily large energy (e.g. in order
not to damage the SUT). In the case we are considering, we fix the total
energy of the probe, and the optimization problem is not easily solved
directly. Following the approaches used, for example in [28, 33], the minimum
probability of error that can be achieved with any classical transmitter can
be derived. Then, if a class of quantum transmitters whose probability of
error is lower than that limit, can be identified, a quantum advantage is
demonstrated. This is the approach we adopt.

###  II.2 Classical limit

First, we turn our attention to the derivation of the minimum error
probability which can be achieved using classical input states. In the context
of quantum optics, classical states can be defined as the states having a
positive P-representation [34]. A generic classical bipartite state can be
represented as:

$\rho^{cla}=\int
d^{\text{\tiny{\hbox{2M}}}}\bm{\alpha}\,d^{\text{\tiny{\hbox{2L}}}}\bm{\beta}\,P\left(\bm{\alpha},\bm{\beta}\right)|\bm{\alpha}\rangle\langle\bm{\alpha}|\otimes|\bm{\beta}\rangle\langle\bm{\beta}|$
(4)

with $|\bm{\alpha}\rangle$ and $|\bm{\beta}\rangle$ respectively $M$\- and
$L$-mode coherent states (corresponding to the signal and idler channels), and
$P\left(\bm{\alpha},\bm{\beta}\right)\geq 0$ a probability density. The energy
constraint for the signal system is expressed, in terms of mean signal photon
number $\bar{n}_{S}$, as:

$\int
d^{\text{\tiny{\hbox{2M}}}}\bm{\alpha}\,d^{\text{\tiny{\hbox{2L}}}}\bm{\beta}\,P\left(\bm{\alpha},\bm{\beta}\right)|\bm{\alpha}|^{2}=\bar{n}_{S}$
(5)

While the choice of imposing the energy constraint only on the signal system
is arbitrary – the alternative being for example a constraint on the total
energy – it is a natural choice for real applications, where the energy
irradiated over the SUT should be limited.

Given an input state $\rho^{cla}$, the output states, after interaction with
the SUT, are $\rho_{0}^{cla}$ and $\rho_{1}^{cla}$, calculated according to
Eq.(2). In this context, the minimum probability of error for the QCT protocol
with classical states is equal to the minimum probability of error in the
discrimination of $\rho_{0}^{cla}$ and $\rho_{1}^{cla}$. The best performance
in this task is achieved by using a POVM measurement, which assumes the
Helstrom projectors as elements [6]. This optimal discrimination procedure
yields a probability of error given by:

$p_{err}=\frac{1}{2}\left(1-D\left(\rho_{0},\rho_{1}\right)\right)$ (6)

where $D\left(\rho_{0},\rho_{1}\right)=||\rho_{0}-\rho_{1}||/2$ is the trace
distance with $||\rho||=\sqrt{\rho^{\dagger}\rho}$. A lower limit for
$p^{cla}_{err}$ can be found by upper bounding
$D\left(\rho_{0},\rho_{1}\right)$. By exploiting the convexity of the trace
distance (see the Supplementary Materials for details), the minimum error
probability for classical states in the QCT protocol is bounded by:

$\mathcal{C}:=p^{cla}_{err}\geq\frac{1-\mathbb{E}_{\mathcal{P}_{0}}\left[\mathbb{E}_{\mathcal{P}_{1}}\left[\sqrt{1-e^{-n_{S}\left(\sqrt{\tau_{0}}-\sqrt{\tau_{1}}\right)^{2}}}\right]\right]}{2}$
(7)

The quantity $\mathcal{C}$ establishes a lower bound for the discrimination
error probability when considering classical resources and an optimal
measurement strategy. We note that this bound is not tight, which means it may
not be reached by any classical receiver.

###  II.3 Quantum strategy

In the following we analyze the particular strategy, involving quantum states,
that is able to surpass the best classical performance $\mathcal{C}$. It uses
a transmitter $\rho$ constituted of $K$ replicas of a two-mode squeezed vacuum
(TMSV) state and a PC receiver [29], whose output is processed by a maximum
likelihood decision.

The TMSV state [35] admits the following expression in the photon number
basis:
$|\psi\rangle=\sum_{n=0}^{\infty}\sqrt{P_{\bar{n}}(n)}|n\rangle_{S}|n\rangle_{I}$
with $P_{\bar{n}}(n)$ a thermal distribution, $\bar{n}$ the mean photon number
and $|n\rangle_{i}$ the $n$-photon state in the $i=S$ signal, or $i=I$ idler,
mode. TMSV states can easily be produced experimentally by parametric down-
conversion [36, 37, 38] or four wave mixing [39, 40, 41]. A multimode TMSV
state is a tensor product of $K$ TMSV states $\otimes^{K}|\psi\rangle$. It
admits a multi-thermal distribution for the total photon number both in the
signal and idler modes, $n_{i}=\sum_{k=1}^{K}n_{i}^{(k)}$, denoted as
$P_{\bar{n},K}(n_{i})$ (the same notation $\bar{n}$ is used for the mean
photon number) and preserves perfect photon number correlation between the
channels i.e.
$P_{\bar{n},K}\left(n_{S},n_{I}\right)=P_{\bar{n},K}\left(n_{I}\right)\delta_{n_{S},n_{I}}$.
The result of the PC measurement is the classical random variable
$\bm{n}=\left(n_{S},n_{I}\right)$, whose distribution
$p\left(\bm{n}|\mathcal{P}_{x}\right)=\langle
n_{S},n_{I}|\rho_{x}|n_{S},n_{I}\rangle$, is conditioned on the nature $x$ of
the process, with $\rho_{x}$ defined in equation 2. Using Bayes theorem and
assuming the defective and reference processes to be equiprobable,
$p\left(\mathcal{P}_{x}=\mathcal{P}_{0}\right)=p\left(\mathcal{P}_{x}=\mathcal{P}_{1}\right)=1/2$,
we can write the a posteriori probability
$p\left(\mathcal{P}_{x}|\bm{n}\right)$ as:

$p\left(\mathcal{P}_{x}|\bm{n}\right)=\frac{p\left(\bm{n}|\mathcal{P}_{x}\right)p\left(\mathcal{P}_{x}\right)}{p\left(\bm{n}\right)}=\frac{p\left(\bm{n}|\mathcal{P}_{x}\right)}{p\left(\bm{n}|\mathcal{P}_{0}\right)+p\left(\bm{n}|\mathcal{P}_{1}\right)}$
(8)

Let’s first consider the action of a pure loss channel on the photon number
distribution of a state. The distribution $p\left(\bm{n}|\tau\right)$ after
the interaction is the composition of the joint distribution
$P_{\bar{n},K}\left(n_{S},n_{I}\right)$ with a binomial distribution:

$p\left(\bm{n}|\tau\right)=\sum_{m=n_{S}}^{\infty}P_{\bar{n},K}\left(m,n_{I}\right)\mathcal{B}\left(n_{S}|m,\tau\right).$
(9)

since the $m$ signal photons can be seen as undergoing a Bernoulli trial each,
with probability of success $\tau$, resulting in the binomial distribution
$\mathcal{B}\left(n_{S}|m,\tau\right)$. Using the linearity of quantum
operations, we can evaluate the effect of a loss channel with transmittance
$\tau$, $p\left(\bm{n}|\mathcal{P}_{x}\right)$, as:

$p\left(\bm{n}|\mathcal{P}_{x}\right)=\int_{\mathcal{A}_{x}}p\left(\bm{n}|\tau\right)g_{x}\left(\tau\right)d\tau.$
(10)

After the measurement, the decision is made by choosing the outcome which
maximizes the conditional probability: $y$ is chosen such that
$y=\operatorname*{arg\,max}_{x\in\\{0,1\\}}p\left(\mathcal{P}_{x}|\bm{n}\right)$.
We note that this condition for the choice of $y$ is equivalent, due the
constant prior, to a maximum likelihood decision, i.e. choosing $y$ such that
$p\left(\bm{n}|\mathcal{P}_{y}\right)\geq
p\left(\bm{n}|\mathcal{P}_{1-y}\right)$. The probability of error for the
quantum strategy, $\mathcal{Q}$, is given by:

$\displaystyle\mathcal{Q}:=p_{err}^{QCT}$
$\displaystyle=\sum_{\bm{n}}\min_{x}p\left(\mathcal{P}_{x}|\bm{n}\right)p\left(\bm{n}\right)$
$\displaystyle=\frac{1}{2}\sum_{\bm{n}}\min_{x}p\left(\bm{n}|\mathcal{P}_{x}\right)$
(11)

###  II.4 Classical states and photon counting

We have derived a limit $\mathcal{C}$ on the performance that can be achieved
with an optimal classical strategy, i.e. optimal classical input states and an
unspecified optimal receiver. We also defined the performance,$\mathcal{Q}$,
of a quantum strategy using TMSV and PC. In this section, we consider the case
in which classical states are paired with a PC receiver. Indeed, the bound
found in Eq.(7) is not tight, meaning that it may not be possible to reach it.
Moreover, in the case where a POVM could be found that saturates the bound,
its implementation may be of difficult practical realization. The analysis of
the best classical performance achievable with the PC receiver will give a
second classical benchmark, whose performance can be experimentally validated.

The analysis is analogous to the one performed for the quantum strategy.
However, since the input states considered are limited to classical ones, the
use of idler modes cannot improve the performance. Classical states are
statistical mixtures of coherent states, that are Poisson distributed in the
photon number: their variance is lower bounded by the Poisson one [34]. In
this scenario, the best performance is achieved using signal states having a
Poisson photon number distribution. Indeed, the error probability $p_{err}$ is
proportional to the overlap of the measurement outcomes, as shown in Fig.(1),
which, in the case of a PC measurement, are the photon number distributions.
Using a state with a narrower photon number distribution will lead to better
discrimination performances. We denote the best performance which can be
achieved using classical states and PC as $\mathcal{C}^{pc}$, and write it in
the form:

$\mathcal{C}^{pc}:=p^{cla,pc}_{err}=\frac{1}{2}\Big{(}1-q_{p}\Big{)}$ (12)

where $0\leq q_{p}\leq 1$. The form of the function $q_{p}$ depends on the
distributions of the considered processes, $g_{0}$ and $g_{1}$. In the
following, we report on the case where both the reference and defective
process have Gaussian distributions $G_{\bar{\tau},\sigma}$ (with $\bar{\tau}$
and $\sigma^{2}$ the mean and the variance). The distributions we consider are
thus $g_{0}=G_{\bar{\tau}_{0},\sigma_{0}}$ and
$g_{1}=G_{\bar{\tau}_{1},\sigma_{1}}$ for the reference and defective
processes respectively. For more general solutions, and a more in depth
analysis, we refer the reader to the Supplementary Materials.

As pointed out in previous sections, the probability of error depends on the
overlap of the two measurement distributions
$p\left(\mathcal{P}_{0}|n_{S}\right)$ and
$p\left(\mathcal{P}_{1}|n_{S}\right)$. Under the assumption of Gaussian
distributions for both processes, and of a large photon number in the initial
probe state, $\bar{n}_{s}\gg 1$, the photon number distribution at the outcome
can be well approximated by a Gaussian:
$p\left(\mathcal{P}_{x}|n\right)\approx
G_{\bar{\tau}_{x}\bar{n}_{s},\sigma_{G^{(x)}}}$, with
$(\sigma_{G}^{(x)})^{2}=\bar{n}_{s}\bar{\tau}_{x}+\bar{n}_{s}^{2}\sigma_{x}^{2}$
(details can be found in Supplementary Materials). An expression for the
overlap can be found, after the determination of the solutions of the equality
$p\left(\mathcal{P}_{0}|n\right)=p\left(\mathcal{P}_{1}|n\right)$. In general,
this equation admits two solutions in the case where the distributions are
Gaussian. If we can assume the variances $\sigma_{G}^{(x)}$ to be of the same
order, and $\bar{\tau}_{0}<\bar{\tau}_{1}$, often only one of the solutions,
labeled $n_{th}$, will be in a range where $p\left(\mathcal{P}_{0}|n\right)$
and $p\left(\mathcal{P}_{1}|n\right)$ are not negligible. Under these
conditions, we can derive a closed form for the function $q_{p}$ in Eq.(12) in
the form of:

$q_{b}=\frac{1}{2}\Bigg{(}\text{erf}\Big{[}\frac{n_{th}-\bar{n}_{s}\bar{\tau}_{0}}{\sqrt{2}\sigma^{(0)}_{G}}\Big{]}-\text{erf}\Big{[}\frac{n_{th}-\bar{n}_{s}\bar{\tau}_{1}}{\sqrt{2}\sigma_{G}^{(1)}}\Big{]}\Bigg{)}.$
(13)

As indicated before, the expression in Eq.(13) holds for certain regimes, and
a more general solution, as well as an explicit expression for $n_{th}$, is
reported in the Supplementary Materials.

###  II.5 Numerical study

Figure 3: _Error probabilities_. The error probabilities as functions of the
mean value of the reference process $\bar{\tau_{0}}$ are displayed, for
different variances of the defective process $\sigma^{2}$ (the numerical
values are reported in the legend), and for the different discrimination
strategies described in the main text: the quantum strategy (solid line), the
classical strategy with PC (dotted line) and the optimal classical strategy
(dash-dotted line). The reference distribution is considered strongly peaked:
$g_{0}\approx\delta\left(\tau-\bar{\tau_{0}}\right)$. In column A, the defect
distribution, $g_{1}(\tau)$, is uniform, in B it is Gaussian. The
distributions $g_{1}(\tau)$ for different values of variance are plotted on
the first row panels. In both bottom panels, the mean number of signal photons
is fixed to $n_{s}=10^{5}$ and the mean value of the defective process is
$\bar{\tau_{1}}=0.997$.

In the strategies described in the previous sections, the performance of the
QCT is a function of the mean photon number $\bar{n}_{S}$, as well as the form
of the distributions of the processes considered, $g_{0}$ and $g_{1}$.

A visualization of the discrimination problem is depicted in Fig.(1), where
the measurement strategy considered at the receiver is photon counting. In
Fig.(1-A) the distributions of the processes to be discriminated are shown,
$g_{0}$ in red and $g_{1}$ in blue. In the first row, the two distributions
barely overlap. Nonetheless, when a probe with finite energy is used to
perform the discrimination, the intrinsic noise of the state on the photon
number distribution results in a significant amplification of the overlap of
the measurement outcome distributions. This is shown in Fig.(1-C), where the
two possible photon number distributions, after the interaction with the SUT,
and for the classical input state discussed in the previous sections, are
plotted. The overlap area, highlighted in green, is a visual representation of
the probability of error: $p_{err}$ is equal to half that area. This overlap
cannot be reduced using classical states without an increment to the total
energy used. However, using the same signal energy, an improvement can be
achieved using quantum correlations. If we consider as input a TMSV state, as
described earlier in the text, we can analyze the photon number distribution
for the signal system, conditioned to having measured a given photon number in
the idler branch. This analysis is shown in Fig.(1-E), where the overlap
between the two distributions is dramatically reduced, leading to much more
efficient discrimination. The nature of the advantage resides in the photon
number correlations, which are used to greatly reduce the photon number noise
of the initial state. In the second row, a situation in which the overlap of
the parameter’s distributions’ is significant is depicted. Indeed, in
Fig.(1-B, the considered distributions have the same mean value, but different
variances. In this case, using a classical discrimination strategy performs
very poorly, as shown in Fig.(1-D), where the two processes are almost
indistinguishable. Fig.(1-F) shows how, once again, quantum correlation can be
used to greatly improve the performance. We note that the best performances
(shown in Fig.(1-E and 1-F) are achieved in the limit of large $K$ so that
$\bar{n}_{S}/K\ll 1$.

From this point on, in order to reduce the number of considered parameters, we
assume the reference process $\mathcal{P}_{0}$ to be strongly peaked around
$\tau_{0}$ so that we can approximate it as
$g_{0}\left(\tau\right)\approx\delta\left(\tau-\bar{\tau_{0}}\right)$, where
$\delta(\tau)$ is the Dirac delta distribution. We also consider two different
noteworthy forms for $g_{1}$: the Gaussian distribution,
$g=G_{\bar{\tau},\sigma}$, and the uniform one,
$g=\mathcal{U}_{\bar{\tau},\delta}$. Using a Gaussian probability distribution
for the transmittance $\tau$ is justified by the wide range of physical
phenomenons it can describe. In this case, the processes are fully
characterized by their mean value $\bar{\tau}$ and variance $\sigma^{2}$. On
the other hand, the uniform distribution, $\mathcal{U}_{\bar{\tau},\delta}$ is
well suited to describe situations in which there is a complete lack of
knowledge of the nature of the process, whose range can be limited by the
physical constraints of the apparatus. Once more, the processes are
characterized by two parameters only: their mean $\bar{\tau}$ and half-width
$\delta$. To make a fair comparison of the error probabilities for Gaussian
and uniform distributions, we choose their parameters such that the resulting
variances are equal. In particular, we will use $\delta=\sqrt{3}\sigma$, since
Var$[\mathcal{U}_{\bar{\tau},\delta}]=\delta^{2}/3$.

In Fig.(3), we show how, for each of the three strategies studied – the
classical optimal, $\mathcal{C}$, the classical strategy with PC,
$\mathcal{C}^{pc}$, and the quantum strategy $\mathcal{Q}$ – the error
probability depends on the reference process’ transmittance, $\tau_{0}$. The
number of signal photons is fixed to $n_{s}=10^{5}$, as well as the mean value
of the transmittance of the defective process $\bar{\tau}_{1}=0.997$. In
column (A), the performances of the different strategies are plotted in the
case of a uniform distributed defect, while in column (B) the case of a
Gaussian distributed one is reported. In both cases, the results for different
values of the variance of the process (denoted by different colors) are shown.
In the range studied, it can be seen that, the quantum strategy (solid lines)
outperforms both the classical lower bound (dashed lines) and the strategy
employing a classical probe and PC (dotted lines). Both classical and quantum
strategies have, in fact, a maximum error probability when
$\tau_{0}=\bar{\tau}$. This shared feature does not depend on the noise but on
the nature of the problem. Indeed, in this configuration, the overlap between
the initial distributions $g$ is maximum and, as a result, the error
probability is maximum as well. In both panels, we see how, in this high
overlap region, the effect of the probe state’s noise is such that, when using
a classical strategy, the processes result are completely ($p_{err}=0.5$), or
almost completely, indistinguishable, while distinguishability is recovered
when quantum states with reduced noise are used.

###  II.6 Experimental results

To validate the presence of a quantum advantage for the conformance test, we
experimentally implemented the protocol described in the section “Quantum
Strategy” – with the substantial difference that we operate with a non-ideal
detection efficiency $\eta$. The detection efficiency encompasses different
processes: the non-unit quantum efficiency of the detectors, the losses due to
the different optical elements, and, in the case of correlated photon sources,
the imperfect efficiency affecting the measurement of the correlated photons.
These losses cannot be distinguished from those caused by a pure loss channel.
Hence, the fit with theoretical curves must be made using the substitution
$\tau\rightarrow\eta\tau$. In general, the effect is equivalent to a reduction
of the probe energy of the same factor $\tau$. For the quantum strategy
however another effect is the reduction of the degree of correlation of the
state. Due to this in general it is expected that the quantum advantage will
be reduced as $\eta$ becomes lower. The setup used is thoroughly described in
the Materials and Methods section.

Figure 4: _Error probabilities for the QCT_. Error probabilities are plotted
as a function of the reference parameter $\tau_{0}$, for different values of
the defective process’ mean value $\bar{\tau}$ and variance $\sigma$:
$A:\\{\sigma=0.001,\bar{\tau}=0.999\\}$,
$B:\\{\sigma=0.002,\bar{\tau}=0.994\\}$ and
$C:\\{\sigma=0.003,\bar{\tau}=0.997\\}$. On the left-hand side, we show the
theoretical curves (quantum and classical probability of error with PC
$\mathcal{Q}$ and $\mathcal{C}^{pc}$, as well as classical bound
$\mathcal{C}$), in the ideal case of unitary detection efficiency
($\eta_{S}=\eta_{I}=1$) for different sets of parameters $\sigma$ and
$\bar{\tau}$. In the right column, we show the experimental error
probabilities obtained using the quantum PC protocol (black dots – Exp Q), and
the classical PC one (black circles – Exp C). The theoretical error
probabilities using the estimated experimental parameters (black and light
green solid lines - $\mathcal{Q}$ and $\mathcal{C}^{pc}$) are plotted for both
protocols, as well as the confidence interval for the quantum case (green
shaded areas – $\pm$ one standard deviation). The classical bound
$\mathcal{C}$ is also displayed as a comparison.

Our results are presented in Fig.(4). Using a sample presenting a varying
spatial transmittance, we acquired data for different values of $\tau$. We
used this experimental dataset to realize different defect distributions: in
Fig.(4), we plot $p_{err}$ (see Eq.(1)) when uniform distributions are
considered with $A:\\{\sigma=0.001,\bar{\tau}=0.999\\}$,
$B:\\{\sigma=0.002,\bar{\tau}=0.994\\}$ and
$C:\\{\sigma=0.003,\bar{\tau}=0.997\\}$. For each acquisition, the mean
numbers of photons $n_{S}$ and $n_{I}$, and the efficiencies $\eta_{S}$ and
$\eta_{I}$ of the signal and idler channels, are estimated, as well as the
electronic noise $\nu$ of the camera. These values are used to draw the
theoretical curves: the quantum and classical error probabilities with PC,
$\mathcal{Q}$ and $\mathcal{C}^{pc}$, as well as the classical optimal bound
$\mathcal{}C$, defined in Eq.7. We display the results in Fig.4 (D-E-F). As a
comparison, we also plot the ideal case $\eta_{S}=\eta_{I}=1$ in Fig.4
(A-B-C). For the quantum error probability, we also report the confidence
region at two standard deviations as a colored band around the curve. The two
sets of experimental data reported in Fig.(4) correspond to the quantum and
classical PC strategies respectively. More details on the elaboration of these
results and data analysis can be found in the section ’Materials and Methods’.

For the three regimes considered, most points fall within the confidence
interval. These results show that, even in the case of degraded detection
efficiency, the quantum strategy always brings an advantage with respect to
the classical one based on PC. However, within the range of the defect
distribution, the optimal classical bound $\mathcal{C}$ on the error
probability becomes smaller than the quantum strategy error probability in
some region. To bring this point into perspective, we point to Fig.4 (A-B-C),
which are constructed with the same experimental parameters but unit detection
efficiency. In this case, the quantum strategy overcomes any classical one. We
note that, as expected, while the classical error probabilities are little
modified by the change in efficiency $\eta$, the quantum one improves
significantly in case of $\eta=1$. This effect stems from the fact that, as
mentioned, spurious losses $\eta$ reduce the photon number correlations
between signal and idler channels.

### Constrained probability of error

Figure 5:  _Cost analysis in QCT_. In panel A, the dependence on the signal
photon number $n_{S}$ of the total error probability ($p_{err}$), false
positive ($p_{01}$) and negative ($p_{10}$) is shown. The maximum likelihood
post-processing used in this case minimizes $p_{err}$. The reference process
is considered strongly peaked around the value $\tau_{0}=0.8$, while the
defective one is chosen uniformly distributed with mean $\bar{\tau}=0.9$ and
half-width $\delta=0.09$. In panel B, a situation similar to that of panel A
is considered, but we use a biased maximum likelihood post-processing with
bias coefficient $b=0.6$ (see main text for details). In panel C, we fix
photon number to $n_{S}=500$ and analyze the dependence of the cost $C$ on
$b$. All the other parameters are equal to those of the previous panels. In
panel D, the optimum value of $b$ as a function of $S$ is shown.

Up until now, the figure of merit considered was the total probability of
error, $p_{err}$ (Eq.1). However, this quantity may not be the most relevant
one to optimize, depending on the nature of the process under test. In some
instances, it can be required to impose a constraint in either one of the
conditional probabilities $p_{01}$ or $p_{10}$, instead of their sum. These
two types of error represent, in fact, different outcomes and the minimization
of one can be deemed more important than that of the other. To proceed
further, we use the notion of cost, which quantifies the fact that in the
mislabelling of a process, the false-positives and false-negatives may not be
equivalent for the operator making the decision: one type of error may be more
“costly” than the other. We introduce the coefficient $0<S<1$, and define the
total cost $C$ of the conformance test as

$C=Sp_{10}+(1-S)p_{01}.$ (14)

From Eq.(14), it follows that if we use the cost $C$ as the figure of merit
for the evaluation of the QCT, the minimum probability of error is achieved
for $S=1/2$. If $S\neq 1/2$, the minimization of the cost will, in general,
not minimize the total error probability $p_{err}$.

In the following, we analyze how the strategies presented in the previous
sections, consisting of PC and a maximum likelihood decision, both for
classical and quantum probes, can be modified to minimize the cost $C$. In the
formalism of the previous sections, we write the general term for the
probabilities $p_{ij}$ as:

$\displaystyle p_{01}$ $\displaystyle=$
$\displaystyle\sum_{\textbf{n}}p(\textbf{n}|\mathcal{P}_{1})\Theta[p(\textbf{n}|\mathcal{P}_{0})-p(\textbf{n}|\mathcal{P}_{1})]$
$\displaystyle p_{10}$ $\displaystyle=$
$\displaystyle\sum_{\textbf{n}}p(\textbf{n}|\mathcal{P}_{0})\Theta[p(\textbf{n}|\mathcal{P}_{1})-p(\textbf{n}|\mathcal{P}_{0})]$
(15)

where $\Theta$ is the step function.

The false-positive and false-negative probabilities are plotted as functions
of the number of signal photons in Fig.(5-A) for the previous section’s
strategies, in the case where $p_{err}$ is the optimized quantity. Both in the
classical (red lines) and quantum (blue lines) cases, the false-negative
probability tends to be smaller than the false-positive one, meaning that the
procedure is more likely to select the reference process. This unbalance comes
from the fact that we consider as a reference process a distribution strongly
peaked around $\tau_{0}$, while, for the defective process, we consider a
uniform distribution. In terms of photon counts probability, this results in
the fact that, when an overlap is present, the most peaked process,
$\mathcal{P}_{0}$, is chosen, while the defective one is selected in a range
that is larger, but also less likely to be measured. This results in a bias
towards $\mathcal{P}_{0}$.

The maximum likelihood post-processing used up until now minimizes the
probability of error by construction. To minimize the cost $C$, a different
post processing is needed. We can modify the maximum likelihood condition used
for the decision, $p\left(\bm{n}|\mathcal{P}_{y}\right)\geq
p\left(\bm{n}|\mathcal{P}_{1-y}\right)$, to make it more likely that one
specific process is selected according to some parameters. We select $y$ when:

$B^{(y)}p\left(\bm{n}|\mathcal{P}_{y}\right)\geq
B^{(1-y)}p\left(\bm{n}|\mathcal{P}_{1-y}\right)$ (16)

where:

$\displaystyle B^{(0)}$ $\displaystyle=\frac{1-b}{2}$ $\displaystyle B^{(1)}$
$\displaystyle=\frac{1+b}{2}$ (17)

with $b\in[-1,1]$ a real constant. We call this new post-processing biased
maximum likelihood. Selecting a positive $b$ results in the defective process
being chosen more often, while a negative $b$ results in a more likely
selection of the reference one. In other words, the value of $b$ can be varied
to shift $p_{01}$ and $p_{10}$, reducing the cost function $C$ at the cost of
the increase of the total error probability. This is shown in Fig.(5.B):
$p_{01}$ and $p_{10}$ have been brought closer to each other and the total
error probability slightly increased, both in the quantum and classical cases.
It is worth highlighting how the quantum advantage that we had with maximum
likelihood processing is well preserved when the biased maximum likelihood one
is used. The optimization of $C$ is performed by varying the coefficient $b$.
As an example, we consider the particular situation in which $S=1/4$, meaning
that each false positive is considered three times as costly as a false
negative. In this situation, it is convenient to use a positive value for $b$,
whose effect is to reduce the number of times the reference process
$\mathcal{P}_{0}$ is selected overall, thus reducing the false positive
probability $p_{01}$. The dependence of $C$ on $b^{(0)}$ is shown in
fig.(5-C). Both in the quantum and classical cases, as expected, there is a
single optimum value for $b$, for which one has a minimum cost. The dependence
of the optimum value of $b$ on $S$, with respect to the cost $C$, is shown in
Fig.(5-D).

The false positive and negative probabilities of the quantum PC strategy can
be evaluated experimentally with the same procedure as the one described in
the Materials and Methods for the total probability of error. The experimental
points, along with the theoretical curves are reported in Fig.(6).

Figure 6: _Theoretical and experimental conditional probabilities of error for
the quantum PC strategy:_ $\bar{\tau}=0.999$ and $\bar{\sigma}=0.001$.

##  III Discussion

In this work, we proposed a protocol addressing the conformance test task,
exploiting specific probe states and measurement strategies. We investigated
the specific scenario in which the systems under test, representing possible
outputs of the unknown process, can be modeled by pure loss channels.

In this context, we found a lower bound on the error probability that can be
achieved in the discrimination task, using any classical source paired with an
optimal measurement. We then showed that a particular class of quantum states,
namely TMSV states, can be used in conjunction with a simple receiver
consisting of a photon-counting measurement and a maximum likelihood decision,
to improve the performances over those of any classical strategy. This
enhancement is due to the high degree of correlation (entanglement) of TMSV
states, whose nature is fundamentally quantum. We also analyzed the particular
case of a classical source paired with photon-counting measurements at the
receiver.

We demonstrated that the quantum advantage persists for a wide range of
parameters, even in presence of losses due to possible experimental
imperfections, and we validated our results by performing an experimental
realization of the protocol. We showed an experimental advantage in a
realistic scenario where experimental losses amounted to more than $20\%$,
highlighting the robustness of the proposed protocol.

We emphasize the fact that these results are particularly significant because
they are achieved using states that are easily produced, as well as a receiver
design of simple implementation, allowing practical applications of the
protocol with present technology. Remarkably, we found an advantage although
the bound $\mathcal{C}$ on the performance of classical states is not tight,
meaning that the actual quantum advantage could effectively be higher. For a
thorough study, one can considered, a scenario in which a set of objects is
tested instead of one. In this scenario, we find an informational limit using
the Holevo bound and we showed how our proposed quantum strategy, relying on
an independent measurement over individual systems, surpasses a more general
classical strategy, where joint measurements over a collection of systems are
allowed. A detailed description of this result will be presented in a
forthcoming paper. The proposed QCT protocol could be used in the foreseeable
future in significant problems concerning the monitoring of production process
of any object probed with quantum states. For example, our results on loss
channel QCT can be used to boost the accuracy in the identifications of issues
in concentration and composition of chemicals production by transmittance
measurement.

##  IV Materials & Methods

###  IV.1 Experimental setup

The experiment is based on spontaneous parametric down-conversion and is
depicted in Fig.(7). A multi-mode TMSV state $\otimes^{K}|\psi\rangle$ is
generated using a 1 cm3 type-II $\beta$-barium borate (BBO) crystal and a
continuous-wave laser at $\lambda_{p}=405$ nm, delivering a power of 100 mW.
The down-converted photons are correlated in momentum. This correlation is
mapped into spatial correlations using a lens in $f$-$f$ configuration: the
“far-field” lens ($f_{\text{FF}}=1$ cm) is positioned at one focal length of
the output plane of the crystal. The absorption sample is positioned in the
conjugated plane, which corresponds to the far-field of the source. It
consists of a coated glass plate, which presents different transmittance
regions, realized with depositions of varying density. A blank coated glass is
inserted in the idler beam’s path, to match the optical paths of both beams.
The sample is imaged using a second “imaging” lens onto a charge-coupled (CCD)
camera (Princeton Instruments PIXIS:400BR eXcellon), working in linear mode,
with high quantum efficiency ($>95\%$ at 810 nm) and low (few electrons per
pixel per frame) electronic noise. Individual pixels of the camera are binned
together in 12x12 macro-pixels, to increase the readout signal-to-noise ratio
and the acquisition speed.

$f_{\textrm{FF}}$$\tau$$f_{\textrm{FF}}$$x$-$x$-$q$$q$n${}_{\textrm{S}}$n${}_{\textrm{I}}$BBOIFCCDS${}_{\textrm{S}}$S${}_{\textrm{I}}$
Figure 7: _Schematic of the experimental setup_. Pumping a BBO crystal with a
laser at $405$ nm, a multi-mode TMSV is generated. Using a lens of focal
$f_{\text{FF}}$, the correlation in momentum is converted into correlation in
position in the sample plane, which is then imaged on the camera using a
second lens. The signal beam passes through the sample of transmittance $\tau$
and is then detected in the $\mathrm{S}_{S}$ region of the CCD. The idler beam
goes directly to $\mathrm{S}_{I}$, without interacting with the sample, though
its optical path is matched with the sample’s one using a non-absorbing glass.
Integrating the signals over the two detection regions, $n_{S}$ and $n_{I}$
are collected.

The photon counts $n_{S}$ and $n_{I}$ of the signal and idler beams are
obtained by integrating over two spatially correlated detection areas
$\mathrm{S}_{S}$ and $\mathrm{S}_{I}$, which are sub-parts of the “full”
correlated areas. The two “full” correlated areas are schematically defined by
the two illuminated regions on the camera, i.e. the spatial extent of the
idler and signal beams. Their precise definition is made using a procedure
described in [22]. The total number of spatial modes collected is $K_{s}\sim
10^{3}$, and that of temporal mode is $K_{t}\sim 10^{10}$ (see [42] for more
details on the estimation). The mean photon number measured in one region
$S_{i}$ is $\bar{n}_{i}\sim 10^{5}$. Hence, the mean occupation number per
mode is very small $\bar{n}_{i}/(K_{s}\times K_{t})\sim 10^{-8}$. Under these
conditions, the multi-thermal marginal photon-number distributions can be well
approximated by Poisson distributions [34]. When the absorptive sample is
replaced with an equivalent coated glass plate without any deposition, we
estimate, for each pair of regions $\mathrm{S}_{S}$ and $\mathrm{S}_{I}$, the
corresponding detection efficiencies $\eta_{S}$ and $\eta_{I}$, by exploiting
the correlations of the SPDC process [43, 44]. Indeed, the phase matching
conditions, the interference filter, the quality of the optical alignment, and
the pixels’ properties on the camera all introduce small differences in the
photon number correlations that need to be accounted for. With the absorptive
sample, for each pair of regions, i.e. for all the different values of $\tau$,
$N_{\mathcal{D}_{\tau}}=2\cdot 10^{4}$ frames are recorded
($\mathcal{D}_{\tau}$ designates the experimental dataset for transmittance
$\tau$).

###  IV.2 Data analysis

Using the experimental data, two error probabilities, $p_{01}$ and $p_{10}$,
must be evaluated to access to the total error probability:

$p_{err}=\frac{1}{2}\left(p_{01}+p_{10}\right).$ (18)

In our analysis, the reference process $\mathcal{P}_{0}$ is strongly peaked
around $\tau_{0}$. Hence, the evaluation of $p_{10}$ with experimental data is
rather straightforward. We use the dataset $\mathcal{D}_{\tau_{0}}$, composed
of $M$ measurement outcome $(n_{s};n_{i})$, when the transmittance $\tau_{0}$
is inserted in the signal branch, to compute the experimental frequency of
error $f_{10}$. In particular, for each data pair
$(n_{s};n_{i})\in\mathcal{D}_{\tau_{0}}$, a label $y\in\\{0,1\\}$ is assigned
according to the maximum likelihood strategy described in the main text. In
the evaluation, experimental parameters, such as the mean photon number of the
source and the detection efficiencies are used, estimated in a calibration
phase. $f_{10}$ is then defined as the number of wrong decisions, which in
this case is simply the number of occurrences of $y=1$ over $M$. $f_{10}$ will
converge to the real probability of error $p_{10}$ as the number of
experimental points in the dataset $M$ becomes large.

On the other hand, the evaluation of $p_{01}$ requires a more careful
approach. In fact, one has to construct an ensemble of experimental data
$\mathcal{D}_{\mathcal{P}_{1}}$, whose data are taken with parameters $\tau$
distributed in a way that is representative of the true probability density
$g_{1}(\tau)$. Experimentally, we acquire $\bar{M}$ photon number pairs for
each of the transmittance $\tau_{i},\;i\in\\{1,...,L\\}$, i.e. $\bar{M}\cdot
L$ data points in total. The transmittance values $\tau_{i}$ are selected in
an interval $[\tau_{min},\tau_{max}]$. They may not be distributed uniformly
in this interval, because of a non perfect experimental control of the
position of the absorption layer, and the uncertainty on the estimation of
$\tau_{i}$. Starting from the experimental dataset, and binning the interval
$[\tau_{min},\tau_{max}]$ in $K$ equal sub-intervals of size $2l$
($2l=(\tau_{max}-\tau_{min})/K$), we define the experimental distribution for
the transmittance, $g$, as the normalized histogram

$g(\tau)=\sum_{k=1}^{K}\frac{d_{k}}{2lL}\Theta_{k}\left(\tau\right),$ (19)

where $d_{k}$ is the number of experimental transmittance values $\tau_{i}$
that fall in the $k$-th bin, so that the relative number of pairs
$(n_{s};n_{i})_{k}$ is $d_{k}\bar{M}$. $\Theta_{k}\left(\tau\right)$ is a step
function, equal to unity in the $k$-th bin and 0 elsewhere. In an ideal case,
the experimental distribution $g$ in Eq.(19), should be a uniform distribution
in the interval $[\tau_{min},\tau_{max}]$, but, in practice, because of the
experimental issues mentioned above, it can present significant discrepancies
with respect to it.

In order to approximate a general probability distribution $f(t)$, we can
multiply the coefficients in Eq. (19), by proper weights $w_{k}$, with
$k\in\\{1,...,K\\}$. Once the weights are determined by a proper optimization
procedure, described in the Supplementary Materials, they are used to modulate
the number of experimental data in each bin according to the substitution
$d_{k}\bar{M}\longmapsto d_{k}w_{k}\bar{M}$. $w_{k}$ represents the fraction
of data kept in bin $k$. The ensemble of data randomly picked from the initial
set, according to the weights, will define the final set
$\mathcal{D}_{\mathcal{P}_{1}}$. We refer to the distribution of
$\mathcal{D}_{\mathcal{P}_{1}}$ after the reweighting procedure, as
$g_{\mathbf{w}}(\tau)$, where $\mathbf{w}=[w_{1},...,w_{K}]$. Since the
experimental set cannot be increased, the constraint $0\leq w_{k}\leq 1$ must
be imposed. The procedure will in general reduce the number of points in the
dataset: an increase of the statistical uncertainty can arise. Thus, it may be
useful to introduce a second constraint on the size of the final dataset,
$M_{T}=\sum_{k=1}^{K}d_{k}w_{k}\bar{M}$, in the optimization procedure to put
a lower bound on how much data is discarded.

For the optimization process, we define the objective function to be maximized
as the Bhattacharyya coefficient between the distributions $f$ and
$g_{\mathbf{w}}$:

$T\left(\mathbf{w}\right)=\int_{\tau_{min}}^{\tau_{max}}\sqrt{f(\tau)g_{\mathbf{w}}(\tau)}d\tau.$
(20)

This quantity measures their similarity and ranges between 0 and 1.
$T\left(\mathbf{w}\right)$ gives a quantitative measure of how close the
experimental dataset can be arranged to resemble the objective distribution.
We define a threshold value $0\leq T_{th}\leq 1$, above which the
approximation of $f\approx g_{\mathbf{w}}$ is deemed “good enough”, i.e.
$T\left(\mathbf{w}\right)\geq T_{th}$. In general, the optimal $\mathbf{w}$
and the corresponding $T\left(\mathbf{w}\right)$, depend on the interval
$[\tau_{min};\tau_{max}]$, which should be chosen so that
$\int_{\tau_{min}}^{\tau_{max}}f(t)$ is close to unity. We note that the
minimum number of data points $M_{T}$ and the initial distribution of the
experimental data can affect the results of the optimization algorithm. For
this reason, in our realization of this algorithm, we sampled the interval
$[\tau_{min},\tau_{max}]$ with a dense array of values $\tau_{i}$, taken
equispaced within the experimental uncertainty, and we took a very large total
number of data points with respect to the target $M_{T}$. The approximated
distribution for the data, reported in the Results section, had a coefficient
$T\sim 1$. A more formal description, as well as further details on the
procedure, are reported in the Supplementary Materials.

## Acknowledgments

### Funding

This work was founded by the EU via “Quantum readout techniques and
technologies” (QUARTET, Grant agreement No 862644).

### Author contributions

GO and IRB proposed the QCT protocol. GO developed the theoretical model and
performed the theoretical and numerical analysis both for the classical limits
and the quantum protocol. PB and EL acquired the experimental data, the
experiment being designed by IRB. PB and GO performed the data analysis, with
contributions from IPD. MG leads the quantum optics group at INRiM and
coordinated the project. All authors contributed to the discussion of the
results and the writing of the manuscript.

### Competing Interests

The authors declare no competing interest.

### Data availability

All data needed to evaluate the conclusions are reported in the paper. Further
data are available under reasonable request to the corresponding author.

## References

  * Pirandola  _et al._ [2018] S. Pirandola, B. R. Bardhan, T. Gehring, C. Weedbrook, and S. Lloyd, Advances in photonic quantum sensing, Nature Photonics 12, 724 (2018).
  * Degen  _et al._ [2017] C. L. Degen, F. Reinhard, and P. Cappellaro, Quantum sensing, Reviews of Modern Physics 89, 035002 (2017).
  * Genovese [2016] M. Genovese, Real applications of quantum imaging, Journal of Optics 18, 073002 (2016).
  * Berchera and Degiovanni [2019] I. R. Berchera and I. P. Degiovanni, Quantum imaging with sub-poissonian light: challenges and perspectives in optical metrology, Metrologia 56, 024001 (2019).
  * Polino  _et al._ [2020] E. Polino, M. Valeri, N. Spagnolo, and F. Sciarrino, Photonic quantum metrology, AVS Quantum Science 2, 024703 (2020).
  * Helstrom [1976] C. Helstrom,  _Quantum detection and estimation theory_ (Academic Press, New York, 1976).
  * Chefles and Barnett [1998] A. Chefles and S. M. Barnett, Quantum state separation, unambiguous discrimination and exact cloning, Journal of Physics A: Mathematical and General 31, 10097 (1998).
  * Chefles [2000] A. Chefles, Quantum state discrimination, Contemporary Physics 41, 401 (2000).
  * Aasi  _et al._ [2013] J. Aasi _et al._ , Enhanced sensitivity of the ligo gravitational wave detector by using squeezed states of light, Nature Photonics 7, 613 EP (2013).
  * Ruo Berchera  _et al._ [2013] I. Ruo Berchera, I. P. Degiovanni, S. Olivares, and M. Genovese, Quantum light in coupled interferometers for quantum gravity tests, Physical Review Letter 110, 213601 (2013).
  * Schäfermeier  _et al._ [2018] C. Schäfermeier, M. Ježek, L. S. Madsen, T. Gehring, and U. L. Andersen, Deterministic phase measurements exhibiting super-sensitivity and super-resolution, Optica 5, 60 (2018).
  * Ortolano  _et al._ [2019] G. Ortolano, I. Ruo-Berchera, and E. Predazzi, Quantum enhanced imaging of nonuniform refractive profiles, International Journal of Quantum Information 17, 1941010 (2019).
  * Pradyumna  _et al._ [2020] S. T. Pradyumna, E. Losero, I. Ruo-Berchera, P. Traina, M. Zucco, C. S. Jacobsen, U. L. Andersen, I. P. Degiovanni, M. Genovese, and T. Gehring, Twin beam quantum-enhanced correlated interferometry for testing fundamental physics, Communications Physics 3, 104 (2020).
  * Tapster  _et al._ [1991] P. R. Tapster, S. F. Seward, and J. G. Rarity, Sub-shot-noise measurement of modulated absorption using parametric down-conversion, Physical Review A 44, 3266 (1991).
  * Monras and Paris [2007] A. Monras and M. G. A. Paris, Optimal quantum estimation of loss in bosonic channels, Physical Review Letter 98, 160401 (2007).
  * Adesso  _et al._ [2009] G. Adesso, F. Dell’Anno, S. De Siena, F. Illuminati, and L. A. M. Souza, Optimal estimation of losses at the ultimate quantum limit with non-gaussian states, Physical Review A 79, 040305 (2009).
  * Matthews  _et al._ [2016] J. C. Matthews, X.-Q. Zhou, H. Cable, P. J. Shadbolt, D. J. Saunders, G. A. Durkin, G. J. Pryde, and J. L. O’Brien, Towards practical quantum metrology with photon counting, npj Quantum Information 2, 16023 (2016).
  * Moreau  _et al._ [2017] P.-A. Moreau, J. Sabines-Chesterking, R. Whittaker, S. K. Joshi, P. M. Birchall, A. McMillan, J. G. Rarity, and J. C. F. Matthews, Demonstrating an absolute quantum advantage in direct absorption measurement, Scientific Reports 7, 6256 (2017).
  * Losero  _et al._ [2018] E. Losero, I. Ruo-Berchera, A. Meda, A. Avella, and M. Genovese, Unbiased estimation of an optical loss at the ultimate quantum limit with twin-beams, Scientific Reports 8, 7431 (2018).
  * Avella  _et al._ [2016] A. Avella, I. Ruo-Berchera, I. P. Degiovanni, G. Brida, and M. Genovese, Absolute calibration of an emccd camera by quantum correlation, linking photon counting to the analog regime, Optics Letter 41, 1841 (2016).
  * Brida  _et al._ [2010] G. Brida, M. Genovese, and I. Ruo Berchera, Experimental realization of sub-shot-noise quantum imaging, Nature Photonics 4, 227 (2010).
  * Samantaray  _et al._ [2017] N. Samantaray, I. Ruo-Berchera, A. Meda, and M. Genovese, Realization of the first sub-shot-noise wide field microscope, Light: Science & Applications 6, e17005 EP (2017).
  * Sabines-Chesterking  _et al._ [2019] J. Sabines-Chesterking, A. R. McMillan, P. A. Moreau, S. K. Joshi, S. Knauer, E. Johnston, J. G. Rarity, and J. C. F. Matthews, Twin-beam sub-shot-noise raster-scanning microscope, Optics Express 27, 30810 (2019).
  * Lloyd [2008] S. Lloyd, Enhanced sensitivity of photodetection via quantum illumination, Science 321, 1463 (2008).
  * Tan  _et al._ [2008] S.-H. Tan, B. I. Erkmen, V. Giovannetti, S. Guha, S. Lloyd, L. Maccone, S. Pirandola, and J. H. Shapiro, Quantum illumination with gaussian states, Physical Review Letter 101, 253601 (2008).
  * Lopaeva  _et al._ [2013] E. D. Lopaeva, I. Ruo Berchera, I. P. Degiovanni, S. Olivares, G. Brida, and M. Genovese, Experimental realization of quantum illumination, Physical Review Letter 110, 153603 (2013).
  * Zhang  _et al._ [2015] Z. Zhang, S. Mouradian, F. N. C. Wong, and J. H. Shapiro, Entanglement-enhanced sensing in a lossy and noisy environment, Physical Review Letter 114, 110506 (2015).
  * Pirandola [2011] S. Pirandola, Quantum reading of a classical digital memory, Physical Review Letter 106, 090504 (2011).
  * Ortolano  _et al._ [2021] G. Ortolano, E. Losero, S. Pirandola, M. Genovese, and I. Ruo-Berchera, Experimental quantum reading with photon counting, Science Advances 7, eabc7796 (2021).
  * Pendrill [2014] L. R. Pendrill, Using measurement uncertainty in decision-making and conformity assessment, Metrologia 51, S206 (2014).
  * Nielsen and Chuang [2011] M. A. Nielsen and I. L. Chuang,  _Quantum Computation and Quantum Information_ , 10th ed. (Cambridge University Press, USA, 2011).
  * Meda  _et al._ [2017] A. Meda, E. Losero, N. Samantaray, F. Scafirimuto, S. Pradyumna, A. Avella, I. Ruo-Berchera, and M. Genovese, Photon-number correlation for quantum enhanced imaging and sensing, Journal of Optics 19, 094002 (2017).
  * Pirandola  _et al._ [2011] S. Pirandola, C. Lupo, V. Giovannetti, S. Mancini, and S. L. Braunstein, Quantum reading capacity, New Journal of Physics 13, 113012 (2011).
  * Mandel and Wolf [1995] L. Mandel and E. Wolf,  _Optical Coherence and Quantum Optics_ (Cambridge University Press, 1995).
  * Weedbrook  _et al._ [2012] C. Weedbrook, S. Pirandola, R. García-Patrón, N. J. Cerf, T. C. Ralph, J. H. Shapiro, and S. Lloyd, Gaussian quantum information, Reviews of Modern Physics 84, 621 (2012).
  * Bondani  _et al._ [2007] M. Bondani, A. Allevi, G. Zambra, M. G. A. Paris, and A. Andreoni, Sub-shot-noise photon-number correlation in a mesoscopic twin beam of light, Physical Review A 76, 013833 (2007).
  * Peřina  _et al._ [2007] J. Peřina, J. Křepelka, J. Peřina, M. Bondani, A. Allevi, and A. Andreoni, Experimental joint signal-idler quasidistributions and photon-number statistics for mesoscopic twin beams, Physical Review A 76, 043806 (2007).
  * Chekhova  _et al._ [2018] M. V. Chekhova, S. Germanskiy, D. B. Horoshko, G. K. Kitaeva, M. I. Kolobov, G. Leuchs, C. R. Phillips, and P. A. Prudkovskii, Broadband bright twin beams and their upconversion, Optics Letter 43, 375 (2018).
  * Glorieux  _et al._ [2011] Q. Glorieux, L. Guidoni, S. Guibal, J.-P. Likforman, and T. Coudreau, Quantum correlations by four-wave mixing in an atomic vapor in a nonamplifying regime: Quantum beam splitter for photons, Physical Review A 84, 053826 (2011).
  * Pooser and Lawrie [2016] R. C. Pooser and B. Lawrie, Plasmonic trace sensing below the photon shot noise limit, ACS Photonics 3, 8 (2016).
  * Wu  _et al._ [2019] M.-C. Wu, B. L. Schmittberger, N. R. Brewer, R. W. Speirs, K. M. Jones, and P. D. Lett, Twin-beam intensity-difference squeezing below 10 hz, Optics Express 27, 4769 (2019).
  * Brida  _et al._ [2010] G. Brida, I. Degiovanni, M. Genovese, M. L. Rastello, and I. Ruo-Berchera, Detection of multimode spatial correlation in pdc and application to the absolute calibration of a ccd camera, Optics Express 18, 20572 (2010).
  * Meda  _et al._ [2014] A. Meda, I. Ruo-Berchera, I. P. Degiovanni, G. Brida, M. L. Rastello, and M. Genovese, Absolute calibration of a charge-coupled device camera with twin beams, Applied Physics Letters 105, 101113 (2014).
  * Brida  _et al._ [2006] G. Brida, M. Genovese, I. Ruo-Berchera, M. Chekhova, and A. Penin, Possibility of absolute calibration of analog detectors by using parametric downconversion: a systematic study, Journal of the Optical Society of America B 23, 2185 (2006).

Supplementary Materials

##  V Derivation of the classical limit

In this section, we show how the bound $\mathcal{C}$ in Eq.(7) of the main
text can be derived. $\mathcal{C}$ is a bound on the minimum probability of
error, in the QCT protocol, achievable using as probe classical states – as
defined in the main text –, $\rho_{cla}$, with fixed signal energy
$\bar{n}_{S}$.

Given an input state $\rho^{cla}$, the output state, after interaction with
the SUT, will be either $\rho_{0}^{cla}$ or $\rho_{1}^{cla}$. The best
performance in discriminating between those states is achieved by using a POVM
measurement, which assumes the Helstrom projectors as elements [6], yielding a
probability of error:

$p_{err}=\frac{1}{2}\left(1-D\left(\rho_{0},\rho_{1}\right)\right)$ (S.1)

where $D\left(\rho_{0},\rho_{1}\right)=||\rho_{0}-\rho_{1}||/2$ is the trace
distance and $||\rho||=\sqrt{\rho^{\dagger}\rho}$. A lower bound on the
probability of error can be found by upper bounding
$D\left(\rho_{0},\rho_{1}\right)$.

###  V.1 Strongly peaked reference process

We consider first a scenario in which only the defective process
$\mathcal{P}_{1}$ is distributed according to some arbitrary distribution
$g_{1}(\tau)$, while the reference process $\mathcal{P}_{0}$ is strongly
peaked, i.e. it can be assumed that all the object produced have the same
trasmittance parameter $\tau_{0}$. In this situation, for any classical input
state $\rho_{cla}$, as defined in the main text, using the strong convexity of
the trace distance, we can write:

$\displaystyle D(\rho_{0},\rho_{1})\leq\int
d^{\text{\tiny$2M$}}\bm{\alpha}\;d^{\text{\tiny$2L$}}\bm{\beta}\;P(\bm{\alpha},\bm{\beta})\times$
$\displaystyle\times
D(\mathcal{E}_{\tau_{0}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)\otimes|\bm{\beta}\rangle\langle\bm{\beta}|,\mathbb{E}_{\mathcal{P}_{1}}[\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)]\otimes|\bm{\beta}\rangle\langle\bm{\beta}|)$
$\displaystyle=\int
d^{\text{\tiny$2M$}}\bm{\alpha}\,P(\bm{\alpha})D(\mathcal{E}_{\tau_{0}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|),\mathbb{E}_{\mathcal{P}_{1}}[\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)]),$
(S.2)

where $\mathbb{E}_{\mathcal{P}_{1}}$ is the expectation value over the
distribution $g_{1}(\tau)$ and the last equality follows from the fact that
the unaffected reference system doesn’t change the trace distance. Using once
again the convexity of the trace on the second argument of the integrating
term $D$, we have:

$\displaystyle
D(\mathcal{E}_{\tau_{0}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|),\mathbb{E}_{\mathcal{P}_{1}}[\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)])$
$\displaystyle\leq\mathbb{E}_{\mathcal{P}_{1}}[D(\mathcal{E}_{\tau_{0}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|),\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|))]$
$\displaystyle=\mathbb{E}_{\mathcal{P}_{1}}[D(|\sqrt{\tau_{0}}\bm{\alpha}\rangle\langle\sqrt{\tau_{0}}\bm{\alpha}|,|\sqrt{\tau}\bm{\alpha}\rangle\langle\sqrt{\tau}\bm{\alpha}|)]$
(S.3)

To derive this expression, we used the fact that pure loss channels map
coherent states into coherent states, i.e.
$|\alpha\rangle\xrightarrow{\mathcal{E_{\tau}}}|\sqrt{\tau}\alpha\rangle$. The
distance between two pure states is easily calculated:

$\displaystyle
D\left(|\sqrt{\tau_{0}}\bm{\alpha}\rangle,|\sqrt{\tau}\bm{\alpha}\rangle\right))$
$\displaystyle=\sqrt{1-F(|\sqrt{\tau_{0}}\bm{\alpha}\rangle,|\sqrt{\tau}\bm{\alpha}\rangle)^{2}}$
$\displaystyle=\sqrt{1-e^{-|\bm{\alpha}|^{2}(\sqrt{\tau_{0}}-\sqrt{\tau})^{2}}}$
(S.4)

with $F(\cdot,\cdot)$ is the fidelity between two quantum states. We can then
write:

$D\left(\rho_{0},\rho_{1}\right)\leq\mathbb{E}_{\mathcal{P}_{1}}\Bigg{[}\int
d^{\text{\tiny{\hbox{2M}}}}\bm{\alpha}\,P(\bm{\alpha})\sqrt{1-e^{-|\bm{\alpha}|^{2}(\sqrt{\tau_{0}}-\sqrt{\tau})^{2}}}\Bigg{]}$
(S.5)

Given the condition on the total number of photons on the signal system
$N_{S}$ in Eq.(5) of the main text, we have (see Supplementary Materials of
Ref.[28]):

$\displaystyle\int d^{\text{\tiny$2M$}}\bm{\alpha}\,P(\bm{\alpha})$
$\displaystyle\sqrt{1-e^{-|\bm{\alpha}|^{2}(\sqrt{\tau_{0}}-\sqrt{\tau})^{2}}}$
$\displaystyle\leq\sqrt{1-e^{-N(\sqrt{\tau_{0}}-\sqrt{\tau})^{2}}}$ (S.6)

for any probability distribution $P(\bm{\alpha})$.
Using Eq.(S.1), we can now lower bound the minimum probability of error in the
QCT protocol, when one of the processes is strongly peaked, as:

$p_{err}^{cla}\geq\frac{1-\mathbb{E}_{\mathcal{P}_{1}}\big{[}\sqrt{1-e^{-\bar{n}_{s}(\sqrt{\tau_{0}}-\sqrt{\tau})^{2}}}\big{]}}{2}$
(S.7)

###  V.2 Arbitrary reference process distribution

If we relax the condition of strongly peaked reference process, so that
$\mathcal{P}_{0}$ is distributed as $g_{0}(\tau_{0})$, we have in place of
Eq.(S.3), using the convexity on both arguments:

$\displaystyle
D(\mathbb{E}_{\mathcal{P}_{0}}[\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)],\mathbb{E}_{\mathcal{P}_{1}}[\mathcal{E}_{\tau}(|\bm{\alpha}\rangle\langle\bm{\alpha}|)])\leq$
$\displaystyle\mathbb{E}_{\mathcal{P}_{0}}\mathbb{E}_{\mathcal{P}_{1}}[D(\mathcal{E}_{\tau_{0}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|),\mathcal{E}_{\tau_{1}}(|\bm{\alpha}\rangle\langle\bm{\alpha}|))]=$
$\displaystyle=\mathbb{E}_{\mathcal{P}_{0}}\mathbb{E}_{\mathcal{P}_{1}}[D(|\sqrt{\tau_{0}}\bm{\alpha}\rangle\langle\sqrt{\tau_{0}}\bm{\alpha}|,|\sqrt{\tau_{1}}\bm{\alpha}\rangle\langle\sqrt{\tau_{1}}\bm{\alpha}|)]$
(S.8)

So that repeating the steps of previous section the probability of success
will be given by

$p_{err}^{cla}\geq\frac{1-\mathbb{E}_{\mathcal{P}_{0}}\mathbb{E}_{\mathcal{P}_{1}}\big{[}\sqrt{1-e^{-N(\sqrt{\tau_{0}}-\sqrt{\tau_{1}})^{2}}}\big{]}}{2},$
(S.9)

That is Eq.(7) of the main text

##  VI Classical Limit with photon counting

In this section, we will derive the best performance that can be achieved
using a classical input state, $\rho_{cla}$, paired with a photon counting
measurement. In this scenario, idler modes are not required and the minimum
probability of error is given by an input state whose signal photon number
distribution is Poissonian. The effect of a pure loss channel
$\mathcal{E}_{\tau}$ on the photon number distribution, $P_{\bar{n}_{s}}(n)$,
of a state Poisson distributed with parameter ${\bar{n}_{s}}$ is to map it
into another one having Poisson distribution with parameter
${\bar{n}_{s}}\tau$:

$P_{\bar{n}_{s}}(n)\xrightarrow{\mathcal{E}_{\tau}}P_{{\bar{n}_{s}}\tau}(n)$
(S.10)

The total probability of error, $p_{err}$, of the procedure can be evaluated
as the average over the two possible processes, $\mathcal{P}_{0}$ and
$\mathcal{P}_{1}$, of the conditioned probabilities of success :

$p_{err}=\frac{1}{2}(p_{01}+p_{10})$ (S.11)

where $p_{ij}$ is the probability that the discrimination procedure will give
the outcome $i$ conditioned to the SUT being generated by the process
$\mathcal{P}_{j}$.

A SUT generated by the process $\mathcal{P}_{x}$ will have transmission $\tau$
distributed with known probability density function $g_{x}(\tau)$. The
conditioned photon number distribution, $p(n|\mathcal{P}_{x})$, can be
evaluated by taking the expected value over $g_{x}(\tau)$ of the distribution
conditioned on the specific $\tau$:

$h^{(x)}(n):=p(n|\mathcal{P}_{x})=\mathbb{E}_{g}[p(n|\tau)]=\int_{0}^{1}P_{{\bar{n}_{s}}\tau}(n)g_{x}(\tau)d\tau$
(S.12)

The terms in Eq.(S.11) can be evaluated using the expression in Eq.(S.12) and
the fact that:

$p_{ij}=\sum_{\\{n|y=i\\}}p(n|\mathcal{P}_{j})$ (S.13)

where the sum is over all values of $n$ such that the outcome of the decision
is $y=i$.

To evaluate the integral in Eq.(S.12), we notice that $h(n)$ is a compound
distribution. We have then:

$\displaystyle\mathbb{E}_{h}[n]$
$\displaystyle=\mathbb{E}_{g}[\mathbb{E}_{P_{{\bar{n}_{s}}\tau}}(n|\tau)]={\bar{n}_{s}}\mathbb{E}_{g}[\tau]$
(S.14) $\displaystyle var_{h}[n]$
$\displaystyle=\mathbb{E}_{g}[var_{P_{{\bar{n}_{s}}\tau}}(n|\tau)]+var_{g}[\mathbb{E}_{P_{{\bar{n}_{s}}\tau}}(n|\tau)]=$
$\displaystyle={\bar{n}_{s}}\mathbb{E}_{g}[\tau]+{\bar{n}_{s}}^{2}var_{g}(\tau)$
(S.15)

In the following we will consider two different cases for the distribution
$g(\tau)$, the uniform case and the gaussian one.

### Gaussian distribution

We consider the case in which both the processes distributions have Gaussian
form with mean $\bar{\tau}_{x}$ and variance $\sigma_{x}^{2}$, i.e.
$g_{x}(\tau)=G_{\bar{\tau}_{x},\sigma_{x}}(\tau)$. According to
Eq.(S.14-S.15), we have for the final photon number distribution,
$h^{(x)}_{G}(n)$:

$\displaystyle N_{G}^{(x)}$
$\displaystyle=\mathbb{E}_{h}[n]={\bar{n}_{s}}\bar{\tau}_{x}$ (S.16)
$\displaystyle\sigma_{G}^{(x)}$
$\displaystyle=var_{h}[n]={\bar{n}_{s}}\bar{\tau}_{x}+{\bar{n}_{s}}^{2}\sigma_{x}^{2}$
(S.17)

For ${\bar{n}_{s}}$ large enough, the distribution $h_{G}(n)$ will converge to
a gaussian itself:

$h^{(x)}_{G}(n)\sim G(n,N_{G},\sigma_{G})\qquad\qquad\qquad{\bar{n}_{s}}>>0$
(S.18)

In this approximation the probability of error of the procedure can be
evaluated according to Eq.(S.11) and (S.13) as:

$\displaystyle p_{err}$
$\displaystyle=\frac{1}{2}\Bigg{(}\sum_{\\{n|y=1\\}}p(n|\mathcal{P}_{0})+\sum_{\\{n|y=0\\}}p(n|\mathcal{P}_{1})\Bigg{)}$
$\displaystyle\approx\frac{1}{2}\Bigg{(}\int_{\\{n|y=1\\}}G_{N^{(0)}_{G},\sigma^{(0)}_{G}}(n)+\int_{\\{n|y=0\\}}G_{N^{(1)}_{G},\sigma^{(1)}_{G}}(n)\Bigg{)}$
(S.19)

when used the assumption ${\bar{n}_{s}}>>0$ in Eq.(S.18) to approximate the
the sums with integrals over $n$.

To determine the regions of integration, we impose the condition
$G_{N^{(0)}_{G},\sigma^{(0)}_{G}}(n)=G_{N^{(1)}_{G},\sigma^{(1)}_{G}}(n)$.
This condition defines up to two threshold values, $n_{th}$, for $n$, whose
expression, if $\sigma_{G}^{(0)}\neq\sigma_{G}^{(1)}$, is:

$\displaystyle n_{th}^{(-)}$
$\displaystyle=\frac{N_{G}^{(0)}(\sigma_{G}^{(1)})^{2}-N_{G}^{(1)}(\sigma_{G}^{(0)})^{2}-\sqrt{(\sigma_{G}^{(0)}\sigma_{G}^{(1)})^{2}\left(2\left(\sigma_{G}^{(1)}-\sigma_{G}^{(0)}\right)\left(\sigma_{G}^{(1)}+\sigma_{G}^{(0)}\right)\left(\log\left(\sigma_{G}^{(1)}\right)-\log\left(\sigma_{G}^{(0)}\right)\right)+\left(N_{G}^{(1)}-N_{G}^{(0)}\right)^{2}\right)}}{\left(\sigma_{G}^{(1)}-\sigma_{G}^{(0)}\right)\left(\sigma_{G}^{(1)}+\sigma_{G}^{(0)}\right)}$
$\displaystyle n_{th}^{(+)}$
$\displaystyle=\frac{N_{G}^{(0)}(\sigma_{G}^{(1)})^{2}-N_{G}^{(1)}(\sigma_{G}^{(0)})^{2}+\sqrt{(\sigma_{G}^{(0)}\sigma_{G}^{(1)})^{2}\left(2\left(\sigma_{G}^{(1)}-\sigma_{G}^{(0)}\right)\left(\sigma_{G}^{(1)}+\sigma_{G}^{(0)}\right)\left(\log\left(\sigma_{G}^{(1)}\right)-\log\left(\sigma_{G}^{(0)}\right)\right)+\left(N_{G}^{(1)}-N_{G}^{(0)}\right)^{2}\right)}}{\left(\sigma_{G}^{(1)}-\sigma_{G}^{(0)}\right)\left(\sigma_{G}^{(1)}+\sigma_{G}^{(0)}\right)}$

and, if $\sigma_{G}^{(0)}=\sigma_{G}^{(1)}$:

$n_{th}^{(0)}=\frac{\mu_{1}+\mu_{2}}{2}$ (S.21)

We consider the first case, i.e. $\sigma_{G}^{(0)}\neq\sigma_{G}^{(1)}$,
assuming, without loss of generality, $\sigma_{G}^{(1)}>\sigma_{G}^{(0)}$, we
have $n_{th}^{(+)}>n_{th}^{(-)}$, and the value $y=1$ will be selected if
$n>n_{th}^{(+)}$ or $n<n_{th}^{(-)}$. We have:

$\displaystyle p_{err}\approx\frac{1}{2}\Bigg{(}$
$\displaystyle\int_{-\infty}^{n_{th}^{(-)}}G_{N_{G}^{(0)},\sigma_{G}^{(0)}}(n)+\int^{n_{th}^{(+)}}_{n_{th}^{(-)}}G_{N_{G}^{(1)},\sigma_{G}^{(1)}}(n)+$
$\displaystyle+$
$\displaystyle\int^{\infty}_{n_{th}^{(+)}}G_{N_{G}^{(0)},\sigma_{G}^{(0)}}(n)\Bigg{)}$
(S.22)

when we note how the first integral should run from $0$ rather than $-\infty$
but since we are considering ${\bar{n}_{s}}>>0$ we can approximate the exact
range of integration with the latter ($n_{th}^{(-)}$ can take also negative
values but in that case the integral will not contribute to $p_{err}$).
Defining:

$\displaystyle
E^{x}(n_{th})=\text{erf}\Big{[}\frac{n_{th}-N_{G}^{(x)}}{\sqrt{2}\sigma_{G}^{(x)}}\Big{]}$
(S.23)

Where erf$(x)$ is the error function. We have:

$\displaystyle p_{err}$ $\displaystyle\approx\frac{1}{2}\Big{(}1-q_{G}\Big{)}$
(S.24) $\displaystyle q_{G}$
$\displaystyle:=\frac{1}{2}\Big{(}E^{0}(n_{th}^{(+)})+E^{1}(n_{th}^{(-)})-E^{0}(n_{th}^{(-)})-E^{1}(n_{th}^{(+)})\Big{)}$
(S.25)

To simplify the expression above, we note how, if
$\sigma_{G}^{(0)}\approx\sigma_{G}^{(1)}$ one of the threshold values will be
far from the center of both distributions so that it can be neglected and the
solution can be considered unique. In particular, we can define:

$\displaystyle q_{G}^{(+)}({\bar{n}_{s}},\tau_{0},\bar{\tau},\sigma):=$
$\displaystyle\frac{1}{2}\Bigg{(}E^{0}(n_{th}^{(+)})-E^{1}(n_{th}^{(+)}\Bigg{)}$
(S.26) $\displaystyle q_{G}^{(-)}({\bar{n}_{s}},\tau_{0},\bar{\tau},\sigma):=$
$\displaystyle\frac{1}{2}\Bigg{(}E^{1}(n_{th}^{(-)})-E^{0}(n_{th}^{(-)})\Bigg{)}$
(S.27)

and in the approximation of similar variances, we have:

$q_{G}\approx\begin{cases}q^{(+)}_{G}&\text{if}\hskip
5.0pt\bar{\tau}_{0}\leq\bar{\tau}_{1}\\\ q^{(-)}_{G}&\text{if}\hskip
5.0pt\bar{\tau}_{0}>\bar{\tau}_{1}\end{cases}$ (S.28)

The function $q_{G}$ it’s monotonically increasing in the number of photons
${\bar{n}_{s}}$. The approximation in Eq.(S.28) is the result reported in the
main text.

### Uniform distribution

If we consider a uniform distribution having mean $\bar{\tau}$ and half width
$\delta$, $\mathcal{U}_{\bar{\tau},\delta}(\tau)$ the integral in Eq.(S.12)
can be evaluated analytically. We have:

$\displaystyle h_{U}(n)=p(n|\mathcal{P}_{1})$
$\displaystyle=\frac{1}{2\delta}\int_{\bar{\tau}-\delta}^{\bar{\tau}+\delta}P_{{\bar{n}_{s}}\tau}(n)d\tau=$
$\displaystyle=\frac{\Gamma(n+1,{\bar{n}_{s}}(\bar{\tau}-\delta))-\Gamma(n+1,{\bar{n}_{s}}(\bar{\tau}+\delta))}{2\delta{\bar{n}_{s}}n!}$
(S.29)

where $\Gamma(a,b)$ is the incomplete gamma function. Using Eq.(S.14)-(S.15)
we have:

$\displaystyle N_{U}$
$\displaystyle=\mathbb{E}_{h}[n]=\bar{\tau}{\bar{n}_{s}}$ (S.30)
$\displaystyle\sigma^{2}_{U}$
$\displaystyle=var_{h}[n]=\bar{\tau}{\bar{n}_{s}}+\frac{1}{3}\delta^{2}{\bar{n}_{s}}^{2}$
(S.31)

If we assume, as in the previous section, ${\bar{n}_{s}}>>0$, the distribution
in Eq.(S.29) has two noteworthy extreme cases:

  1. 1.

if $\sigma_{U}>{\bar{n}_{s}}\delta$, $h_{U}(n)$ can be approximated with a
Gaussian one, $h_{U}(n)\sim G(n,N_{U},\sigma_{U})$.

  2. 2.

if $\sigma_{U}<<{\bar{n}_{s}}\delta$, $h_{U}(n)$ can be approximated with a
uniform distribution,
$h_{U}(n)\sim\mathcal{U}_{N_{U},{\bar{n}_{s}}\delta}(\tau)$.

The two cases are showed in Fig(S.1).

Figure S.1:  Limiting cases for the distribution $h_{U}(n)$. In figure are
shown the two limiting cases for the distribution $h_{U}(n)$, described in the
main text. In the first panel the parameters are fixed to
${\bar{n}_{s}}=1.5*10^{5}$, $\bar{\tau}$ and $\delta=0.002$, so that
$\sigma_{U}>{\bar{n}_{s}}\delta$, the situation described in case $(1)$. In
the second panel the parameters are ${\bar{n}_{s}}=1.5*10^{5}$,
$\bar{\tau}=0.92$ and $\delta=0.7$, and we have the situation of case (2).

If we consider case (1), that is the situation we have in our experimental
range of parameters, since we are dealing with a Gaussian distribution the
probability of success can be found using the analysis of the previous
section.

### Efficiency $\eta$

In the classical case, inefficiencies resulting in a loss of photons can be
accounted for by introducing the efficiency of the channel $\eta<1$. The
effect of the efficiency can evaluated using the fact that the action of a
loss $\eta$ and the one of the SUT, commute. The final effect of $\eta$ is to
reduce the energy available for the discrimination and the classical bounds
found can be evaluated accounting for the inefficiency performing the
substitution ${\bar{n}_{s}}\rightarrow\eta{\bar{n}_{s}}$ in the formulas
referring to the perfect efficiency case.

##  VII Experimental data analysis

The probability of error of the QCT procedure, $p_{err}$, can be divided in
two separate contributions, namely the probability of error conditioned to the
object being produced by the process $\mathcal{P}_{0}$, $p_{10}$, and one
referring to the process $\mathcal{P}_{1}$, $p_{01}$:

$p_{err}=\frac{1}{2}\Big{(}p_{10}+p_{01}\Big{)}$ (S.32)

The process $\mathcal{P}_{0}$ is supposed to be strongly peaked around a value
of transmission $\tau_{0}$, so to evaluate experimentally the contribution
$p_{10}$, one can use an object having transmission $\tau_{0}$ and repeat the
discrimination procedure $N_{\mathcal{D}}$ times. The frequency of error
$f_{\mathcal{D}}(\mathcal{P}_{1}|\mathcal{P}_{0})$ on the dataset
$\mathcal{D}$ will converge to the probability of error $p_{10}$ as
$N_{\mathcal{D}}\to\infty$.

The process $\mathcal{P}_{1}$, on the other hand, will have an arbitrary, but
known, probability density $f(\tau)$. To estimate the probability $p_{01}$, a
dataset $\mathcal{D}$, representative of the density $f(\tau)$, must be
constructed. A convenient way to do this is to create experimentally a
collection of $L$ dataset $\mathcal{D}_{i}$ each generated from a different
value of trasmission $\tau_{i}$, $1\leq i\leq L,$ in a range
$[\tau_{min},\tau_{max}]$ determined by the distribution $f(\tau)$ one wants
to approximate. Suppose for example the target distribution is a uniform one
centered in $\bar{\tau}$ and having half-width $\sigma$. To approximate this
distribution, measurements are taken with $L$ different transmission $\tau$
equispaced in the interval
$[\tau_{min}=\bar{\tau}-\sigma,\tau_{max}=\bar{\tau}+\sigma]$. For a Gaussian
distribution, having mean $\bar{\tau}$ and variance $\sigma^{2}$, measurements
can be taken in the interval $[\bar{\tau}-k\sigma,\bar{\tau}+k\sigma]$, where
$k$ can be chosen depending on the accuracy required for the approximation. In
general, for an arbitrary density $f(\tau)$, one can chose a value $0\leq
r\leq 1$ and determine $\tau_{min}$ and $\tau_{max}$ such that:

$\int_{\tau_{min}}^{\tau_{max}}f(\tau)d\tau=r$ (S.33)

Each dataset $\mathcal{D}_{i}$ will have $N_{i}$ points. The final dataset
used to estimate the gain can be constructed in different ways. The most
straightforward one is to take a number of measurements $N_{i}$, such that the
union dataset $\mathcal{D}_{U}=\bigcup_{i=1}^{L}\mathcal{D}_{i}$ is directly
distributed as $f(\tau)$.

A more flexible approach is to fix the number of measurements taken for each
$\tau$ to a given number $\bar{N}$, i.e. $N_{i}=\bar{N}$ $\forall i$.
Following this procedure, the union dataset $\mathcal{D}_{U}$ will be composed
of $N_{U}=L\bar{N}$ measurements extracted from an ideally uniform
distribution in $[\tau_{min},\tau_{max}]$. The final dataset can then be
rearranged with a procedure of statistical weighting, described in detail in
the following, consisting in discarding a certain amount of measurement
generated from selected $\tau_{i}$. The downside of this approach is that the
number of measurements to be performed is higher then the number of
measurements that will be used in the final dataset approximating $f(\tau)$,
$\mathcal{D}_{f}$. There are however different practical advantages. In
practical scenarios in fact the control in selecting an exact value for the
transmission $\tau$ is limited, so that it may be difficult to have exactly
equispaced values. The situation is showed in Fig.(S.2).

Figure S.2:  Sampling of a trasmissivity interval. In panel A an ideal
equispaced sampling for the trasmissivity $\tau$ is reported in the interval
$[0.5,1]$. In the $x$ axis are reported the labels $i$(ascending ordered) of
the values of trasmissivity $\tau_{i}$ reported in the $y$ axes. In panel B
the same plot is reported in a noisy case, i.e. where the $\tau_{i}$ are taken
equispaced but a value can be assigned up to a certain precision. In panel C
and D the histograms referring to the situation of panel A and B respectively
are reported.

In panel A and B the plots of an equispaced sampling of the interval
$[\tau_{min}=0.5,\tau_{max}=1]$ are reported in the ideal and noisy cases
respectively. As the $\tau_{i}$ are affected by random noise their
distribution is affected as well as showed in panels C and D. Panel C refers
to the ideal case where the empirical distribution of $\tau$ is uniform as
expected. In this case, acting before the measurement only on the coefficients
$N_{i}$ could be an effective way to achieve a good approximation of the
target distribution. In the noisy case however, whose distribution is reported
in panel D, the distribution is not only non uniform but, in most cases,
unknown before the experiment is performed. For this reason, acting only on
the coefficient $N_{i}$, supposing an uniform distribution of the sampling
could result in a bad approximation of $f(\tau)$.

The amount of noise varies depending on experimental conditions so in cases in
which the deviation from an uniform distribution of the empirical one is
small, or in cases in which it can be well characterized beforehand, acting on
the coefficient $N_{i}$ remains a viable solution. Another advantage of the a
posteriori statistical weighting procedure is the fact that while some
measurements are discarded to construct the final dataset, the procedure can
be repeated with the same initial data with a different target distribution,
meaning that the gain over different distributions can be evaluated without
having to perform the experiment again.

In the following will be given a formal description of the statistical
weighting procedure used to create a dataset $\mathcal{D}_{f}$ approximating
one generated by a probability density $f(\tau)$.

### Approximation using statistical weighting

Using the notation of the previous section, the initial experimental dataset
$\mathcal{D}_{U}$ will be composed by $N_{U}=L\bar{N}$ measurement results
where each group of $\bar{N}$ measurements has been performed on a different
transmission $\tau_{i}$, $0\leq i\leq L$, where the values $\tau_{i}$ are a
sampling of the interval $[\tau_{min},\tau_{max}]$ determined using Eq.(S.33).

Our aim is to generate a dataset $\mathcal{D}_{f}$ to be a good approximation
of an ideal dataset generated with measurement with objects distributed as
$f(\tau)$.

To do this let’s define the map $\mathcal{H}_{\rho}^{\Pi}$, mapping a
distribution of $\tau$ onto the distribution $F(\textbf{X})$ of the
experimental outcomes X (e.g. the distribution of photon pairs):

$F(\textbf{X})=\mathcal{H}_{\rho}^{\Pi}[f(\tau)]$ (S.34)

The subscript $\rho$ and superscript $\Pi$ denote the dependence of the map on
the initial state used to probe $\tau$ and the measurement performed and will
be omitted from here on forward.

Our aim is to operate a transformation $W$ on our initial dataset,
$\mathcal{D}_{U}\xrightarrow{W}\mathcal{D}_{f}$, such that the point in
$\mathcal{D}_{f}$ are distributed as $h(\textbf{X})$, i.e.
$\textbf{X}_{f}=\textbf{X}\in\mathcal{D}_{f}\sim F(\textbf{X})$. Let us
introduce the weight vector $\textbf{w}=[w_{1},...,w_{L}]$. The initial
dataset $\mathcal{D}_{U}$ is composed of $L$ subsets $\mathcal{D}_{i}$, each
composed of $\bar{N}$ experimental points,
$\mathcal{D}_{i}=[\textbf{X}_{i}^{1},...,\textbf{X}_{i}^{\bar{N}}]$. We define
the transformation $W_{i}$ as the one acting on each $\mathcal{D}_{i}$,
$\mathcal{D}_{i}\xrightarrow{W_{i}}\mathcal{D}^{\prime}_{i}$, extracting at
random a fraction $N^{\prime}_{i}=\lfloor w_{i}\bar{N}\rfloor$ of experimental
points from $\mathcal{D}_{i}$. The final dataset will be composed as the union
set of all the $\mathcal{D}^{\prime}_{i}$ and the transformation W can be
defined as:

$\mathcal{D}_{U}\xrightarrow{W}\mathcal{D}_{f}=\bigcup_{i=1}^{L}\mathcal{D^{\prime}}_{i}$
(S.35)

Therefore a successful approximation is dependent on the optimization of the
vector w.

To perform this optimization, we have to introduce a figure of merit. To do so
let’s define a probability density function, $g(\tau)$, approximating the
empirical distribution in $\mathcal{D}_{U}$. This can be done using a
normalized histogram. We fist divide the interval $[\tau_{max},\tau_{min}]$ in
$K$ sub intervals with equal width, where $K$ is selected via an adequate
algorithm. We can define the categories vector
$\textbf{$\tau$}^{C}=[\tau_{1}^{C},...,\tau_{K}^{C}]$, taking $\tau_{i}^{C}$
as the middle point between the edges of the $i$-th subinterval. The empirical
distribution $g(\tau)$ can then be defined as:

$g(\tau)=\sum_{i=1}^{K}\frac{c_{i}}{2lL}\Theta_{i}(\tau)$ (S.36)

where $c_{i}$ are the number of occurences on the $i$-th sub interval, $L$ is
the number of different values of $\tau$ and $\Theta_{i}(\tau)$ is a step
function and $l$ is the half width of each sub interval:

$\displaystyle\Theta_{i}(\tau)$
$\displaystyle=\begin{cases}1&\tau_{i}-l\leq\tau\leq\tau_{i}+l\\\
0&\text{otherwise}\end{cases}$ (S.37) $\displaystyle l$
$\displaystyle=\frac{\tau_{max}-\tau_{min}}{2K}$ (S.38)

The probability distribution in Eq.(S.36) defines a distribution for the
experimental data, $G(\textbf{X})$ according to Eq.(S.34):

$G(\textbf{X})=\mathcal{H}[g(\tau)]$ (S.39)

Since the map $\mathcal{H}$ does not depend on the probability distribution we
can act on the applying a transformation on the distribution of the dataset we
can apply a transformation of the distribution $g(\tau)$:

$\displaystyle g(\tau)\xrightarrow{W}g^{\prime}(\tau)$
$\displaystyle\longrightarrow\mathcal{H}[g(\tau)]=G(\textbf{X})\rightarrow\mathcal{H}[g^{\prime}(\tau)]=$
$\displaystyle=G^{\prime}(\textbf{X})\sim F(\textbf{X})=\mathcal{H}[f(t)]$
(S.40)

The transformation $W$ must then act on $g(\tau)$, transforming it in
$g^{\prime}(\tau)\sim f(\tau)$. We introduce than the parameterized function
$g_{\textbf{w}}(\tau)$:

$g_{\textbf{w}}(\tau)=\sum_{i=1}^{K}\frac{c_{i}w_{i}}{2lL}\Theta_{i}(\tau)$
(S.41)

depending on the weight vector $\textbf{w}=[w_{1},...,w_{K}]$ under the
normalization constraint:

$\int
g_{\textbf{w}}(\tau)d\tau=1\longrightarrow\sum_{i=1}^{K}\frac{c_{i}w_{i}}{L}=1$
(S.42)

Note how in this configuration, the transformation $W$ acts separately on each
of the $K$ sub interval rather than on each of the $L$ values of $\tau$
measured, so that the vector w has dimensionalty $K\leq L$.

It is clear that the initial distribution is $g(\tau)=g_{\textbf{w}_{0}}$,
where $\textbf{w}_{0}=[1,...,1]$ and the target function is
$g^{\prime}(\tau)=g_{\textbf{w'}}$, where w’ is the optimal weight such that
we have the best approximation $g^{\prime}(\tau)\sim f(\tau)$. In a formal
way, we can define the objective function to be maximized:

$T(\textbf{w})=\int\sqrt{f(\tau)g_{\textbf{w}}(\tau)}dt$ (S.43)

and the optimal vector w’:

$\textbf{w'}=\text{argmax}_{\textbf{w}}T(\textbf{w})$ (S.44)

$T(\textbf{w})$ is the Bhattacharyya coefficient (BC) between the distribution
$g_{\textbf{w}}(\tau)$ and $f(\tau)$, a measure of their similarity, and
ranges between $0$ and $1$. $T(\textbf{w}^{\prime})$ gives a quantitative
measure of how close our experimental dataset can be arranged to one produced
by the distribution $f(\tau)$. We can define a threshold value $0\leq
T_{th}\leq 1$ such that we accept the approximation if
$T(\textbf{w}^{\prime})\geq T_{th}$ and reject it otherwise. Before performing
the optimization, some other constraint must be imposed.

The statistical weight represents the fraction of data that we select from
each dataset, and considering that the data are considered fixed, i.e. new
data cannot be added after the initial dataset has been created, we must
impose for the vector w the condition:

$0\leq w_{i}\leq 1\qquad\forall i$ (S.45)

where the first inequality follows from the fact that the statistical weight
has to be non negative. The fixed data condition in Eq.(S.45) is not
compatible with the normalization constraint in Eq.(S.42), but this issue is
easily solved taking into account the number of experimental points in each
dataset in the formulation of our problem.

The number of experimental data in the original union dataset is
$N_{\mathcal{U}}=L\bar{N}$. Since the transformation $W$ consists exclusively
in discarding data, the number of points in the final dataset will be
$N_{T}\leq N_{\mathcal{U}}$. It is clear that $N_{T}$ will influence the
variance of the estimated gain, so it makes sense to fix it before the
optimization process rather then leaving it as a free parameter. We can then
redefine $g_{\textbf{w}}(\tau)$ as:

$g_{\textbf{w}}(\tau)=\sum_{i=1}^{K}\frac{N_{\mathcal{U}}}{N_{T}}\frac{c_{i}w_{i}}{2lL}\Theta_{i}(\tau)$
(S.46)

The new normalization condition is then:

$\sum_{i=1}^{K}\frac{N_{\mathcal{U}}}{N_{T}}\frac{c_{i}w_{i}}{L}=1$ (S.47)

The new definition for $g_{\textbf{w}}(\tau)$ shifts the normalization
condition in Eq.(S.42), imposed on the number $L$ of different $\tau$, to the
new one in Eq.(S.47) where it is imposed on the number of experimental data in
the final dataset. With this new definition, the initial distribution will be
$g(\tau)=g_{\textbf{w}_{0}}(\tau)$, where in this case
$\textbf{w}_{0}=[w_{0},...,w_{0}]$ and $w_{0}=N_{T}/N_{\mathcal{U}}$. In this
formulation, w is exactly a statistical weight vector, i.e. defines the
fraction of experimental data taken from each sub interval dataset, and we can
find the optimal vector $\textbf{w}^{\prime}$, maximizing $T(\textbf{w})$
under the fixed data constraint of Eq.(S.45) and the normalization one in
Eq.(S.47), that are now compatible conditions. In general, the optimization of
eq.(S.43) under the $K+1$ constraint described can be solved by numerical
methods. This procedure will yield a value $\textbf{w}^{\prime}$ and its value
$T(\textbf{w}^{\prime})$. It is clear that the choice of sub interval number
$K$ influences the outcome of the procedure and could be included in the
optimization process. We decided however to non include it to avoid
unnecessary numerical complications, as good approximation can be reached in a
wide variety of cases performing the binning before the optimization, using
well known algorithms. In our case the binning was performed using Sturges
rule.

Figure S.3:  Distribution approximation with different initial dataset. In row
A are showed the initial dataset and the best approximation for different
values, $\bar{N}=150$ and $\bar{N}=500$, of experimental data taken (see main
text for more details). In row B the same situation is reported in case of an
initial dataset taken in a smaller interval $[\tau_{min},\tau_{max}]$ w.r.t.
row A. In row C is considered a faulty initial dataset, having few
measurements around the peak of the target distribution.

In Fig.(S.3) are reported some examples of approximation of a distribution
from simulated experimental dataset. The distribution to be approximated is a
Gaussian one with mean $\bar{\tau}=0.65$ and standard deviation $\sigma=0.1$
and is reported in each panel as a red curve. In the first row, Fig.(S.3).A,
the experimental values of $\tau$ are taken from an uniform distribution in
$[\bar{\tau}-3.5\sigma,\bar{\tau}+3.5\sigma]$, meaning that $\tau_{min}$ and
$\tau_{max}$ are selected imposing $r=0.9995$ in Eq.(S.33). Starting from the
left, the first panel represent the initial experimental distribution. The
second and third panels show the best approximation in the situation in which
are considered measurement with $L=100$ different values of $\tau$ and, for
each $\tau$, the experimental points taken are $\bar{N}=150$(second) and
$\bar{N}=500$(third). In both cases, the number of point of the final dataset
is fixed to $N_{T}=10000$. As can be seen, if not enough point can be
discarded, as it is the case in the second panel, a good approximation is
difficult to reach. When a lot of experimental points are taken, as in the
situation of the third panel, as can be expected, a good approximation of most
the distributions can be reached. In quantitative terms, we have for the
second panel a value of $T^{\bar{N}=150}(\textbf{w}^{\prime})\approx 0.9$
while for the third $T^{\bar{N}=500}(\textbf{w}^{\prime})\approx 0.99$. In
general, the BC tends to have high values in the situation we are considering,
for reference consider that the unaltered initial distribution has a value
$T^{0}(\textbf{w}^{\prime})\approx 0.8$. In Fig.(S.3).B are reported the same
scenario considered in row A but considering an initial dataset taken in a
smaller interval, such that $r=0.997$ (half width of $[\tau_{min},\tau_{max}]$
equal $3\sigma$). While a direct comparison may be unfair given that the
dataset are randomly generated, it can seen that reducing the interval may be
useful since more data will be taken at the peak of the target distribution.
In the second panel, referring again to the situation of $\bar{N}=150$ we can
see a better approximation with respect to the same panel in row A. Finally,
in row C the same situation is reported in case of a faulty initial dataset,
i.e. one having very few values around the peak of the target distribution. In
this case even with $\bar{N}=500$ experimental data for each $\tau$ (third
panel) a good approximation cannot be reached, indicating that in those
situation it may be worth to perform the measurements again.