Parallel tomography of quantum non-demolition measurements in multi-qubit devices

An efficient characterization of QND measurements is an important ingredient towards certifying and improving the performance and scalability of quantum processors. In this work, we introduce a parallel tomography of QND measurements that addresses single- and two-qubit readout on a multi-qubit quantum processor. We provide an experimental demonstration of the tomographic protocol on a 7-qubit IBM-Q device, characterizing the quality of conventional qubit readout as well as generalized measurements such as parity or measurement-and-reset schemes. Our protocol reconstructs the Choi matrices of the measurement processes, extracts relevant quantifiers -- fidelity, QND-ness, destructiveness -- and identifies sources of errors that limit the performance of the device for repeated QND measurements. We also show how to quantify measurement cross-talk and use it to certify the quality of simultaneous readout on multiple qubits.

In this work, we experimentally implement an efficient parallel QND-MT to characterize the most important measurement properties of a 7-qubit IBM-Q quantum computer [19]. The protocol exploits the low correlations between the qubit readout to implement a cheap parallel single-qubit characterization of each measurement, obtaining relevant quantifiers from the Choi operators such as readout fidelity, QND-ness, and destructiveness [47]. We observe that the device is optimized to maximize the fidelity--calibrated at around ∼ 98% for every qubitbut not the QND-ness, which varies more accross the device and it is lower on average ∼ 96.7%. QND-MT also reveals that bit flip errors are the main source of imperfections. Using a two-qubit QND-MT we quantify measurement cross-talk across device. We find a similar correlation strength between local and non-local pairs of qubits, which introduces an error of less than 1% in the simultaneous execution of qubit readout. This validates the application in parallel of the single-and two-qubit tomographic protocol on the IBM-Q device, which can be executed with a constant number of circuits, avoiding the exponential scaling of a full QT. This parallelization is also extended to the post-processing of data on classical computers.
Finally, we demonstrate the generality of QND-MT by reconstructing composite measurement processes relevant to quantum error correction protocols such as parity measurements and measurement-and-reset schemes with classical feedback. Our experiment shows that the parity measurement involves more errors-mainly nondispersive-than a direct QND measurement due to the presence of an entangling gate. In addition, we observe that the measurement-and-reset scheme can enhance the QND nature of the readout.

QND measurement tomography on a multi-qubit device
A generalized quantum measurement of an N -qubit system in state ρ is described by a set of non tracepreserving quantum processes E n , which add up to a trace-preserving one, E = n E n [2]. Each individual process determines a post-measurement state ρ n = E n (ρ)/p(n), conditioned to the measurement outcome occurring with probability p(n) = Tr(E n (ρ)). A representation of quantum processes commonly used in quantum tomography are the Choi matrices [54]. In this representation, a measurement is described by a set of Choi operators {Υ n } whose matrix elements are given by [47] Υ ijkl n = ⟨ij|Υ n |kl⟩ = ⟨i|E n (|k⟩⟨l|)|j⟩ , with {|i⟩} the basis of the measured system with dimension d. In terms of these matrices we can conveniently determine the dynamics of the post-measurement states E n (ρ) = ijkl Υ ijkl n ρ kl |i⟩ ⟨j|, the POVM elements Π n = ijk Υ kjki n |i⟩⟨j|, and the measurement statistics p(n) = Tr{Π n ρ}, where ρ kl = ⟨k| ρ |l⟩ are the components of the density matrix before measurement. Note that Υ n is the transposed of the positive Choi operatorΥ n , whose components are related by ⟨ij| Υ n |kl⟩ = ⟨ik|Υ n |jl⟩ [54].
The Choi matrices Υ ijkl n are a complete description of the quantum processes of a system and from them we can extract all the relevant physical properties of the measurements. We discuss three relevant quantifiers of the measurement: the readout fidelity F , the QND-ness Q, and the destructiveness D [47] (see methods). Comparing them, we can quantify the quality of the mea-surement for particular tasks and discriminate between different types of measurements. The readout fidelity F describes the efficiency of the readout irrespective of the post-measurement state, and it is thus maximal when the POVMs are projectors, Π n = |n⟩⟨n|. Operationally, it is defined as the average probability of successfully detecting a state |n⟩ of the computational basis, after preparing the system in the same state. The QND-ness Q is the fidelity respect to an ideal measurement of an observable O, that is, a measurement that projects the states into the eigenvectors |n⟩ of O and whose Choi matrices are projectors, Υ n = |nn⟩⟨nn|. QND-ness incorporates information of the post-measurement states and can be determined by the average probability that states of the computational basis |n⟩ are preserved in two consecutive measurements. Finally, the destructiveness D quantifies the back-action introduced by the measurement [47]. For D = 0, the measurement is exactly QND which means that it preserves the expected value of the observable O after consecutive measurement ⟨O⟩ = Tr[Oρ] = Tr[OE(ρ)]. For D > 0, the destructiveness signals a deviation from the QND condition, which can occur independent on how ideal the measurement is. Therefore, it is convenient to know the three quantifiers F , Q, and D to provide a more complete analysis of general non-destructive measurements.
Our QND-MT protocol [47] reconstructs the Choi matrices of a QND measurement self-consistently. As shown in Fig. 1(a), it consists of two applications of the measurement interspersed by a unitary gate V i that prepares a complete basis of initial states, and a second gate U j that enables a complete set of measurements. Given sufficient statistics, this protocol provides conditional probability distributions of single QND measurements p(n|i) and of consecutive measurements p(nm|ij). A maximumlikelihood-based classical post-processing [55][56][57] transforms p(n|i) and p(nm|ij) into a set of physically admissible set of POVM elements {Π n } and Choi matrices {Υ n }, requiring to solve a total of N + 1 optimization problems, with N the number of outcomes (see methods). The full characterization of QND detectors with N qubits and N = 2 N possible outcomes demands reconstructing 2 N Choi operators of size 4 N . In the general case, using a strategy based on Pauli observables, QND-MT requires a total of 18 N circuits, corresponding to 6 N initial gates prepared with the tensor product of V i ∈ {I, σ x , e ∓iπσy/4 , e ∓iπσx/4 } and 3 N intermediate unitaries given by the tensor products of U i ∈ {I, e −iπσy/4 , e −iπσx/4 }, with I and σ j the identity and Pauli operators. Figures 1(b)-(d) show the circuits for the particular case of a single-qubit and two-qubits, as explained below.
The exponential scaling in the number of circuits makes quantum tomography-in the QND-MT or in any other form-unfeasible for systems with large numbers of qubits. This scaling may be avoided if the measurements of separate quantum subsystems are shown to be independent. QND measurements in superconducting circuits are implemented via dispersive readout [15][16][17], where each qubit is coupled to an off-resonance cavity, on which one performs homodyne detection to individually extract the outcome of each qubit state. These readouts is thus built to be independent of each other, but imperfections in the device can lead to cross-talk between the qubit measurements [26,58]. Nevertheless, these cor-relations can be characterized by two-qubit QND-MT of each pair of qubits. If the correlations are weak enough, it is possible to execute a highly parallelized and scalable QND-MT for multi-qubit detectors, reducing the number of circuits and the classical post-processing time.
Parallel single-qubit QND measurement tomography Let us first discuss the tomographic reconstruction of every single-qubit measurement of a quantum processor with N qubits. This means using QND-MT to reconstruct 2N single-qubit Choi matrices Υ α n , for qubits α = 1, ..., N and n = 0, 1. Single-qubit QND-MT applies two measurements Υ α n in between single-qubit gates V i and U j , as shown in Fig.1(c), requiring the evaluation of 18 different circuits. For each pair of measurement outcomes the associated Choi matrix Υ α n is estimated as the solution of a maximum likelihood optimization problem.
Initialization and gate errors are accounted for in this optimization by means of single-qubit gate set tomography (GST) [48,[59][60][61]. GST self-consistently characterizes the initial state, the final POVM measurement, and a complete set of linearly independent gates G i of a device. In our case, G i ∈ {I, σ x , e −iπσy/4 , e −iπσx/4 } are the gates used to implement all the V i and U i operations from QND-MT. Note that GST requires 64 circuits, each composed of three gates G i and a measurement, as shown in Fig. 1(b).
The execution of the circuits of the N single-qubit QND-MTs and GST circuits can be efficiently paral-lelized, applying the single-qubit operations simultaneously as sketched in Figs. 1(b-c). This reduces the total number of experiments from O(N ) down to the sum of 18 QMD-MT and 64 GST circuits. With this refinement, we studied the readout of all qubits in the IBM quantum device ibm_perth. This processor has a quantum volume [62] of 32 and CLOPS (Circuit Layer operations per second) of 2.9 × 10 3 [63], and a qubit connectivity graph shown in Fig. 1(e). Each circuit was evaluated with 2 13 shots, resulting in an experiment that takes approximately 2 minutes, with classical post-processing of 30 seconds on a Ryzen-7 5800H processor with 8 cores.
From the reconstructed Choi matrices of each qubit, we derived the three quantifiers-readout fidelity F , the QND-ness Q, and the indestructiveness 1 − D-shown in Figure 2(a). The ibm_perth processor exhibits readout fidelities between 0.969 and 0.992, with an average of F = 0.98. QND-ness varies much more along the device, ranging from 0.951 in the qubit α = 0 to 0.987 in α = 6, with an average ofQ = 0.967. Indestructiveness behaves similarly to fidelity except for qubit α = 1 and ranges between 0.97 and 0.991. The arithmetic mean of F , Q, and 1 − D-see colormap in Fig. 2(a)-characterizes the performance across the device: qubits in the upper sector (α = 0, 1, 2) perform notably worse than those in the lower half of the chip.
The Choi matrices not only provide individual qubit metrics (F, Q, D) but also hint the physical processes behind measurement errors. The Choi matrix element p a→b n = ⟨bb|Υ n |aa⟩ quantifies the probability that the state flips from |a⟩ to |b⟩ when outcome n is detected. Each element informs about deviations from the ideal projective measurement p a→a a = 1, as well as possible origins for those deviations.
Let us first put this into practice using the averaged Choi matricesῩ n = α Υ α n /N , shown in Figure 2(b). Note how the readout of the |0⟩ state (p 0→0 0 = 0.975) is implemented with a better quality than that of |1⟩ (p 1→1 1 = 0.960). Bit flip noise is identified as the main source of errors, dominated by the qubit decay process |1⟩ → |0⟩ (p 1→0 1 = 0.02), and slightly less influenced by the excitation channel |0⟩ → |1⟩ (p 0→1 1 = 0.016). Considering that ibm_perth has a relaxation time T 1 ≈ 100µs and a measurement time T = 700ns [19], we estimate a baseline probability of qubit relaxation p th = 1 − e −T /T1 ≈ 0.007, which accounts for 35% of the observed decay error. The remaining bit-flip error may be due to Purcell-induced decay and other non-dispersive errors that occur during the measurement process itself [47].
This analysis may also be done qubit by qubit. = 0.021 in the outcome |1⟩, and non-dispersive errors given by elements ⟨ab|Υ 0 |00⟩ and ⟨ab|Υ 1 |11⟩, with a ̸ = b. The projection in this outcome is thus not done correctly, which explains the reduction in the QND-ness and indestructiveness.
The parallelized tomography of the qubits has obvious performance advantages, but it could increase the error of the operations [25]. To quantify potential deviations, we have compared the outcome of parallel tomography on ibm_perth with the independent characterization of those qubits, running the O(N ) circuits separately. As shown in Fig 2(e), the differences in the three quantifiers-fidelity |∆F | = |F ind − F par |, QNDness |∆Q| = |Q ind − Q par |, and destructiveness |∆D| = |D ind − D par |-lay below 10 −2 , and are smaller than the non-idealities of those quantifiers (see Fig. 2(a)). Similarly, we have quantified the average distance in diamond norm [64,65] between the Choi operators computed using both strategies |∆Υ| ⋄ = |Υ ind −Υ par | ⋄ ≤ 1, and these lay below 1.4 × 10 −2 (see Fig. 2(e)), validating the use of the parallelized strategy.

Two-qubit QND measurement tomography and cross-talk quantification
The low distinguishability between parallel and independent single-qubit QND-MT suggests that measurement correlations are weak across the device. We can further quantify such correlations comparing the joint measurement process for pairs of qubits (α, β), given by Υ αβ mn for outcomes m, n = 0, 1, with the individual measurement processes Υ α m ⊗ Υ β n . The two-qubit QND-MT requires the evaluation of 324 circuits-two measurement processes interspersed by layers of gates V i and U j (cf. Fig. 1(d)). Since characterizing all N (N − 1)/2 pairs on an N -qubit device is very costly, we first focused on neighboring qubits, which we expect to exhibit the greatest correlations. More precisely, for a device with M physical connections (α, β) ∈ C of M , we aim to reconstruct the 4M two-qubit Choi matrices.
This two-qubit QND-MT can be parallelized by executing similar circuits on non-overlapping pairs of physically connected qubits. This requires dividing the quantum processor into sets of edges that do not share a common qubit. For the 7-qubit ibm_perth and the 65qubit ibm_brooklyn quantum processors, illustrated in Fig. 1(e)-(f), we only need three sets. For a generic planar graph with M vertices, coloring theorems [66] ensure that the number of sets is never larger than 4, setting a bound on the number of circuits 4 × 18 2 that does not grow with the processor's size. Finally, the protocol requires solving 5M optimization problems, a task that can be efficiently parallelized on classical computers. Here, we employ a parallel two-qubit QND-MT to characterize the readout of physically connected qubits on the IBM quantum device ibm_perth. The experiment runs approximately in 38 minutes -using 2 13 shots per circuit-and the post-processing in 3 minutes. Fig. 3(a) shows the experimental results for the quantifiers F , Q, and 1 − D describing the two-qubit measurement. We see an overall decrease in the readout performance of pairs qubits with respect to the single-qubit results in Fig. 2(a), but still we identify the same qualitative behavior: F and Q increase from top to bottom of the device (see inset of Fig. 3(a)) and QND-ness is the worst and most fluctuating quantifier. As discussed below, the two-qubit quantifiers of all connected pairs are very well approximated as products of the single-qubit ones-e.g. F αβ ≈ F α F β and Q αβ ≈ Q α Q β -. This explains the reduction in average two-qubit fidelity and QND-ness between pairs toF = 0.958 andQ = 0.937. Indestructiveness 1 − D is the most stable quantifier, ranging between 0.955 and 0.97, but reduces its average in a similar amount to 1 −D = 0.963. ) 2 ), which suffers more from bit-flip errors. As in the single-qubit case, we estimate the errors introduced by parallelization by comparing the parallelized two-qubit QND-MT with the independent tomography of each pair. Fig. 3(c) shows the error in fidelity, QND-ness, destructiveness, and Choi operators for each pair of physically connected qubits, as defined in the previous section. Parallelization introduces an error below 2 × 10 −2 in all quantifiers and Choi operators.
We quantify the measurement crosstalk by comparing the measurement processes of individual and pairs of qubits. This is done at the level of quantifiers, introducing heuristic measures of separability for the fidelity C[F αβ ] = |F αβ − F α F β | and for the QNDness C[Q αβ ] = |Q αβ − Q α Q β |. It is also done at the level of operators, with estimates of the POVM correlation As hinted above, we observe a good separability of quantifiers. In Fig. 4(a) we see correlations of ibm_perth device below 10 −2 for all pairs, allowing us to estimate the fidelity and QNDness of pairs of qubits as products of the properties of individual qubits. Figure 4(a) shows the POVM and Choi correlations for the ibm_perth device. We certify the presence of measurement cross-talk between all physically connected pairs of qubits: all POVM elements and Choi matrices are non-separable with correlations on the order of 10 −2 , which exceed the statistical error bars from the tomography for most of the qubits. This represents a crosstalk error of about 1%, which is smaller than the physical error found on single-and two-qubit tomography.
QND-MT is not restricted to nearest-neighbor correlations. As example, we have analyzed the correlations between all qubits in ibm_perth and the central qubit α = 3, in 6 sets of separate experiments. This produces pairs at two different distances, the first neighbors (1, 3) and (3,5), and the second neighbors (0, 3), (2, 3), (3,4), and (3, 6). Fig. 4(b) shows the correlations obtained for those pairs. We can see that correlations C[F αβ ] and C[Q αβ ] are of order 10 −3 for all qubits, while C[Π αβ ] and C[Υ αβ ] are approximately 10 −2 . In this small device we do not observe a clear decay of correlations with distance, but we verify that all correlations are smaller than the measurement errors detected for independent single-qubit tomography.

Scaling of QND-MT on larger devices
Parallel single-qubit QND-MT is an efficient technique to characterize large devices, that requires a fixed number of circuits-82 including GST-independently of the device size. Using the execution times obtained in the experiments on the ibm_perth we can extrapolate the performance in larger devices. For the 65-qubit ibm_brooklyn, with a degree 3 connectivity shown in Fig.  1f and a smaller CLOPS number of 1.5×10 3 [63], we estimate 4 minutes for the single-qubit characterization and 5 minutes of post-processing in a Ryzen-7 5800H processor with 8 cores. Notice that all experimental execution times do not depend on the size of the device but they are limited by the number of CLOPS, which are typically lower for larger devices.
We have discussed also three strategies to certify the errors in parallel QND-MT. One strategy is the application of QND-MT of individual qubits in separate, nonparallel experiments. This has a cost that grows linearly O(N ) with the number of sampled qubits, but it is a routine that may be applied with less frequency than the complete calibration. This method enables the development of heat maps of the chip and suggests the order of magnitude of underlying correlations.
The second strategy is the parallelized QND-MT of pairs of neighboring qubits, a method that will provide results consistent with the previous methodology, but also give information about the strength of the cross-talk. In the two-qubit parallelized strategy, our estimate give a total of 1296 independent circuits for any device size, taking 63 minutes for the two-qubit circuit evaluation in a 65-qubit ibm_brooklyn processor, and 30 minutes in a Ryzen-7 5800H processor with 8 cores.
The third and most expensive strategy is to implement a two-qubit QND-MT for all qubit pairs in large devices with O(log(N )) parallel groups [67] and O(N 2 ) optimization problems. In this case, we estimate 2.5 hours for the experiment and a similar amount of post-processing to characterize the ibm_brooklyn device. This is an efficient scaling that enables a very robust calibration of the complete device, to be done only sporadically.
Finally, for larger devices, the execution and postprocessing times could be too long for a complete twoqubit measurement tomography-extending to days for devices with more than 1000 qubits-in which case it makes sense to either randomly sample those pairs, or concentrate the study to specific regions of the chip, that revealed more problematic in the first two methods.

QND measurement tomography of generalized measurements
The QND-MT protocol we introduce can be applied to any kind of generalized measurements [2]. These include synthetic measurements that combine standard detectors with other computing elements, such as local and entangling gates, auxiliary qubits, and resets.
In this work we discuss the application to stabilizer measurements, a relevant example which are widely used in quantum error correction protocols [18,29,30]. Such measurements are usually implemented with controlled operations over an auxiliary qubit, which is finally measured and reset, to discriminate states with different stabilizer value. If we trace over the auxiliary qubit, the generalized measurement is, up to implementation errors, QND, enabling the repetitive monitoring of error syndromes.
As an illustration of how QND-MT works with a generalized measurement, we discuss a single-qubit parity measurement (PM). As shown in Fig. 5(a), this protocol includes an auxiliary qubit, a controlled CNOT operation, and a single-qubit readout and reset. Note that, unlike all higher parity measurements, the single-qubit PM does not entangle multiple system qubits and thus it is not directly applicable to quantum error correction codes. However, it already includes all the underlying operations supporting multi-qubit PM, which can be scaled to characterize multi-qubit measurement errors in practical error correction codes. Here, we study the performance of the single-qubit PM using two fixed qubits of the ibm_perth quantum processor, and we compare it with the performance of the direct measurement (DM) on the same system qubit, as shown in Fig. 5(b). The Choi matrices and the quantifiers obtained for the parity and direct measurements are shown in Figs. 5(d) and (e), respectively.
In this study we observe a decrease of the fidelity and QND-ness of the PM with respect to the DM. The readout fidelity of the parity measurement F PM = 0.958 is close to the product of F DM = 0.973 and the fidelity of the CNOT provided by IBM F CNOT = 0.9897. Therefore, we can conclude that this decrease is mainly due to the CNOT gate as the error from reset is expected to be smaller than 1% [68]. The indestructiveness is the same for parity and direct measurements, 1 − D PM = 1 − D DM = 0.969, which is consistent with the fact that the CNOT and reset operations do not add measurement back-action on the system. In the Choi operators, we can also see the appearance of new bars that describe the noise introduced by the CNOT gate, as well as an increase in the overall error bars and fluctuations.
Another interesting example of generalized measure-ment is the measure-reset-feedback (MRF) operation, shown in Fig. 5(c). It consists of a QND measurement followed by a reset and a classically-conditioned NOT operation that brings the measured qubit exactly to the quantum state selected as measurement outcome-i.e. the qubit is reset to state |0⟩ or |1⟩ when the measurement outcome was deemed n = 0 or 1, respectively. If the reset and NOT operations have high fidelities, measurementand-reset should fix the QND nature of a measurement, bringing the errors 1 − Q and D closer to the measurement infidelity 1 − F . We applied QND-MT to this generalized measurement using a single qubit of the IBM-Q ibm_perth processor. The resulting Choi matrices and quantifiers are shown in Figure 5(f). The MRF scheme has better performance than the DM in the same qubit, having approximately the same fidelity F DM ≈ F MRF ≈ 0.975 and indestructiveness 1 − D DM ≈ 1 − D MRF ≈ 0.969, but with an increase of QNDness from Q DM = 0.954 to Q DM = 0.960. Considering the error bars of the QND-ness and indestructiveness, we find that the worst case MRF provides a QND readout with the same quality as a direct measurement. Moreover, we also witness a reduction in the noisy components of the Choi operator, such as those describing bit-flip errors.

DISCUSSION
In this work we have demonstrated an efficient, highly parallelizable protocol for QND measurement tomography of a state-of-the-art multi-qubit quantum computer, which works with both single-qubit measurement, as well as generalized measurements-e.g. error syndrome measurements, parity measurements, etc. Our method is based on a self-consistent reconstruction of the Choi matrices for single-qubit and two-qubit measurements, which provides information about measurement quality, the QND nature of the measurement and the strength and type of errors.
In the single-qubit scenario, we have developed strategies to massively parallelize the tomography, an approximation that works when multiple measurements can be executed with small crosstalk or correlation. We have applied this protocol in experiments with a 7-qubit IBM quantum computer, obtaining fascinating insight into the performance of the device. First of all, we have found that the chip is well tuned to high-fidelity measurements, with weak and long pulses-much longer than single-or two-qubit gates-that mitigate non-dispersive and discrimination errors, at the expense of increasing incoherent errors, in particular single-qubit bit flip. This limits the QND nature of the measurement which fluctuates along the different qubits of the device.
We have also developed different strategies determine whether single-qubit measurements are independent, and can be parallelized. The most sophisticated strategy involves applying QND-MT to the simultaneous measurement of two qubits, to reconstruct the joint Choi matrices and quantify the degree of correlation. In the setup considered, these correlations lay below 1% and validate the parallelization strategy which, as discussed above, can be efficiently scaled to large multi-qubit processors with an almost fixed cost.
Finally, we have also demonstrated how QND-MT can be generalized to custom measurements, in particular to parity-type measurements relevant to quantum error correcting codes and measurement-and-reset schemes with classical feedback. We used the Choi matrices to identify coherent errors introduced by the CNOT gate in parity measurements and we provided evidence that the reset operation with classical feedback is an appealing way to improve the QND quality of a measurement. This work opens several avenues for further research. The obvious one is to use QND-MT as an input for systematic optimization of the measurement pulses. The goal here is to optimize the driving amplitude and measurement time, minimizing the errors that manifest in the Choi matrices. This would allow us to reduce the decay channels found in the experiment, while keeping other sources of error at bay-e.g. non-dispersive effects [21], discrimination errors [47,69], decoherence [20], leakage to higher levels of the transmon [23], or rotating wave corrections [22]. Another approach is to design alternative schemes for qubit readout that may be more QND [49,50,70], but this would add new error sources that could be similarly identified and characterized with the application of QND-MT.
An additional research avenue is to further understand and mitigate the correlations between simultaneous measurements. In this work we have explored two-qubit correlations, but higher-order correlations, involving 3 or more qubits, could also be analyzed with the help of better tomography methods, such as compressed sensing [71,72]. These methods could also be used to quantify readout and cross-talk errors occurring in multi-qubit stabilizer measurements involving plaquettes of 4 or more qubits as they are required in practical quantum error correction codes such as the surface [18,29,30] or color codes [73].

QND measurement quantifiers
To characterize the most important properties of nondestructive measurements we employ three quantifiers: the readout fidelity, the QNDness, and the destructiveness. Here, we show how to obtain these quantifiers from the reconstructed Choi matrices as introduced in [47].
The fidelity F is the standard quantifier of a detector's readout performance, measured by the probability that an initially prepared eigenstate |n⟩ is successfully identified, The fidelity can be interpreted as the efficiency of the readout as it can be related to the signal-to-noise ratio of the measurement [15,70]. It ignores any information about the post-measurement state and the QND nature of the measurement. The QND-ness Q incorporates information from the post-measurement state and quantifies how close are the Choi matrices with respect to an ideal projective measurement. In quantitative terms, it is the probability that an initially prepared eigenstate |n⟩ is preserved and successfully identified in two consecutive measurements, The destructiveness D [47] asserts precisely the QND nature of generic measurements by verifying the preservation of the expectation value ⟨O⟩ after the measurement.
Operationally, it is defined as the largest change suffered by any observable compatible with O as, where ∥ · ∥ is the Hilbert-Schmidt norm. Unlike F and Q, computing D requires a complete tomographic reconstruction of the measurement process, E † (O c ) = ijkln Υ klij n * O kl c |i⟩⟨j|, but it allows us to quantify the measurement back-action without the bias of Q towards ideal measurements [47]. Note that equation (4) is motivated by the definition of a QND measurement, Tr(Oρ) = Tr(OE(ρ)). Moving into the Heisenberg picture this condition becomes Tr(Oρ) = Tr(E † (O)ρ), where E † is the self-adjoint process of E. Therefore, we can quantify how QND is a measurement by the deviation between O and E † (O), that is ||O − E † (O)||. The last step to obtain Eq. (4) consists in searching for the largest disagreement over the set of all normalized observables compatible with O, so that we ensure that D is an upper bound for the back-action of the measurement.

Quantification of measurement cross-talk and correlations
To quantify the correlations in the measurement of pairs of qubits we introduce heuristic metrics that compare the POVM and Choi matrices derived from twoqubit tomography with the tensor product of the operators obtained from single-qubit tomography. Although this is a comparison between two detector models rather than a intrinsic property of the operators, the outcome provides information about the operation of the device and the effect of including higher order interactions in the QND-MT. We also use these correlations to quantify the distinguishability error of performing the tomography in parallel or independently, as shown below.
First, we define the correlation in two-qubit Choi operators C[Υ αβ ]. Let Υ α n and Υ β m be the process operators of two single qubits α and β, and let Υ αβ nm be the joint measurement of both qubits. We define the Choi matrices correlation as where⊗ is tensor product operation in the superoperator space and || · || ⋄ is the diamond norm [64,65]. This quantity not only evaluates the distance between the processes Υ α n⊗ Υ β m and Υ αβ nm , but is also related to the probability of discriminating the quantum states generated by them. This probability is given by P d = (1 + C[Υ αβ ])/2. Therefore, a small C[Υ αβ ] ≪ 1 means that the post-measurement states are nearly indistinguishable. Notice that the pre-factor in C[Υ αβ ] are chosen normalize the correlation between 0 and 1.
Conveniently, we can use the same definition in Eq. (5) to evaluate the distinguishability of Choi matrices reconstructed in parallel Υ par or independently Υ ind as C[∆Υ] = C[Υ par − Υ ind ], and thereby quantify the error introduced by the parallelization. This is done in Figs. 2(e) and 3(c).
In a similar spirit as done for C[Υ αβ ], we can define the correlation C[Π αβ ] in the POVMs of two-qubit measurements. Let Π α n and Π β m be the POVM elements of two single qubits α and β, and Π αβ nm be the joint POVM element of both qubits. We define the POVM correlation as where || · || 2 is the 2-norm, that is, the largest singular value. This quantity establishes an upper bound for the average error on the probability distribution predicted by the single-qubit reconstruction P S nm = T r(ρ[Π α n ⊗ Π β m ]) compared with the joint measurement P J nm = T r(ρΠ αβ nm ), given by nm |P S nm − P J nm |/4 ≤ C[Π αβ ] for any density matrix ρ. Notice that C[Π αβ ] is normalized between 0 and 1 as C[Υ αβ ].

Maximum Likelihood estimation for QND measurement tomography
Maximum likelihood estimation (MLE) is a statistical inference method widely used in quantum tomography. MLE allows us to recover density matrices, POVMs, or Choi matrices that are meaningful and satisfy all the physical constraints of a measurement. It achieves this goal by optimizing the likelihood function L(θ|f ) of the experimental dataf for a given parametric model M(θ). We employ as a Gaussian distribution as a likelihood function, wheref (i) are the estimated probabilities obtained from the experiment and p(i) are the theoretical probabilities predicted by the model M(θ). We minimize this likelihood function (7) for both the QND-MT and GST. Notice that, for simplicity, the notation of the theoretical probabilities p(i) omits the dependence on the parameters θ, and that the index i may refer to a group of indices as shown below. The QND-MT consists of two steps, first a measurement tomography of the POVM and then a process tomography of each Choi matrix. We reconstruct the POVMs {Π j } by first obtaining the theoretical probabilities, of obtaining the outcome n condition to the application of gate V k . We then minimize the likelihood function of form (7) over the set of feasible matrices {Π n } satisfying Π n ≥ 0 and n Π n = 1 1. Finally, we estimate the Choi matrices Υ n by obtaining the theoretical probabilities, of obtaining the outcome n in the second measurement and the outcome m in the first measurement, condition to the application of gates j and k. We then minimize the corresponding likelihood function of the form (7) over the set of the Choi matrix Υ n that satisfiesΥ n ≥ 0 and the POVM constraint Tr 1 (Υ n ) = Π n . Here, Tr 1 (·) is the partial trace over the first subsystem, andΥ is the transposed Choi matrix, which is a positive matrix with elements ⟨ik|Υ|jl⟩ = ⟨ij|Υ|kl⟩.
To separate experimental errors in gates and state preparation from the measurement errors that we want to characterize, we can apply a GST previously to the QND-MT. The GST gives us an experimental estimate of the set {ρ, Π i , G j }, composed of estimators of the initial state, the POVM elements, and the gates, respectively. Here {G j } are generic trace-preserving processes and not necessarily unitary operations. The theoretical probabilities of obtaining the outcome l are which are condition to the application of gates i, j, k as shown by the circuits in Fig. 1(b). We then minimize (7) by comparing the probabilities (10) with the experimental data to obtain a physically meaningful set {ρ, Π i , G j } that self-consistently accounts for state preparation, gates, and measurement errors. Notice that when we use GST, we can omit the first step of QND-MT as we already have an experimental estimate of the POVMs {Π i }. In addition, the gates {U j } and {V k } needed for the second step of QND-MT must be formed as concatenations of the {G k } processes in order to account for gate errors.
In total, QND measurement tomography of the device requires solving 3N + 5M optimization problems. We solve them using sequential least-squares programming, satisfying the positivity of operators via Cholesky decomposition, and the completeness constraints via Lagrange multipliers.
To quantify the goodness-of-fit of our estimators, we employ the χ 2 -test [74,75]. This is a standard tool for statistical hypothesis testing, that is, for rigorously deciding if there is enough evidence to reject a model. In our work, we apply this test to all single-and two-qubit Choi matrix reconstructions and demonstrate that the fits and models are in agreement with the experimental data within a standard confidence interval of 95%. See Supplementary Methods 1 for a detailed analysis and a description of the method. Let be c(mn|jk) the counts obtained from the QND-MT of measurement process, which are used to obtain an estimator {Υ n }. The goodness-of-fit χ 2 test for this data reads 1. Compute the predicted probabilities, where {ρ,Π n ,F j } is the gate set estimated with GST.
2. Compute the test statistic, where N s is the number of shots used to evaluate the probabilities.
3. Set an error probability q (typically as 0.05) and compute χ 2 q , implicitly defined by where P r is the probability density function of a χ 2 variable with mean r, P r (x) = x (r−2)/2 e −x/2 2 r/2 Γ( r 2 ) .
The mean value r is given by For a N -qubit detector, this mean value is given by Each term from left to right corresponds to: number of circuit (18 N ), independent probabilities per circuit (4 N − 1), number of free parameters of the Choi matrices (2 N × 16 N ), and completeness constraint (4 N ). For the cases of single-and twoqubits detectors, these are r 1 = 26 and r 2 = 3852, respectively.
We apply this procedure to quantify the goodness-offit of our characterization of the ibm_perth device with QND-MT presented. As shown in Suplementary Figure 19, both single-(a) and two-qubit (b) characterizations have a χ 2 value below the threshold χ 2 q (black horizontal line) and therefore are consistent with the experimental data with confidence 95% (q = 0.05).