Direct spectrum analysis using a threshold detector with application to a superconducting circuit

Due to their potential scalability, superconducting circuits provide a promising framework for quantum information processing. However, as solid state systems, they suffer from strong environmental noise sources that limit their quantum coherence. Despite the great improvement in coherence times over the recent years, which is mainly due to clever optimized designs [1–4], and the proof of principle of a correction algorithm on a quantum memory [5], the quantum coherence times are not yet sufficient for implementing non-trivial fault tolerant quantum computations with a realistic number of physical quantum bits per logical qubit [6, 7]. As environmental noise sources cause decoherence their extensive characterization is a key issue in order to identify the origin of the noisy subsystems and improve materials properties, minimize coupling to these sources with better designs, or implement special dynamical decoupling sequences [8]. Until recently, all characterization techniques (except notably [9, 10]) of noise sources in superconducting quantum bit circuits used the decay of coherence functions borrowed from NMR [8, 11–13]. These techniques are the most sensitive and they operate at low temperature ($T<50$ mK) and low probing power, but, through the dependence of decay functions on the control parameters (such as charge or magnetic flux) they give access to the spectrum of the environmental noise only in a form convoluted with a frequency filter which depends on the NMR sequence. Such filters have a frequency width of the order of the inverse of the pulse sequence time (typically the decoherence rate measured by the sequence: tens to hundreds of kHz) which limits the frequency resolution. In contrast, the technique we present here provides a direct access to the noise spectrum with a frequency resolution given, as for any discrete Fourier transform, by the inverse of the total acquisition time (tens of mHz in our case).

In this article, we first discuss the standard operation of a Josephson bifurcation amplifier (JBA), to motivate and introduce the theoretical framework required for using a threshold detector as a spectrum analyzer. Then we apply this technique to the measurement of the frequency fluctuations of a JBA. We demonstrate that this method combines the advantages of state of the art noise measurement techniques in superconducting circuits [14–16] with the advantages of non-dissipative quantum bit readout setups, achieving the four following aims together: a high bandwidth given by half the repetition rate of the measurement, a high sensitivity which may approach the standard quantum limit of a weak continuous measurement, a low temperature of operation ($<15$ mK) due to the absence of on-chip dissipation and a low probing energy (in our case $\approx {{10}^{3}}$ photons in the JBA). This technique can be applied to any detector involving a threshold effect and in particular qubit state measurement setups.

The JBA detector has been extensively studied for the purpose of superconducting quantum bit readout [17–20]. Our device consists of a section of superconducting niobium coplanar waveguide enclosed between input and output capacitors which provide coupling to this resonant Fabry-Perot like structure (see figure 1, more details are given in [21, 22]). An array of superconducting quantum interference devices (SQUIDs) is located in the middle of this structure at anti-nodes of the electric current distribution of odd numbered harmonic modes. This setup provides a strong and tunable non-linearity for these modes. Due to this non-linearity, this system exhibits parametric amplification below a critical number ${{n}_{c}}$ of photons populating the resonator and a bifurcation phenomenon above it (here ${{n}_{c}}\approx 250$ for the third harmonic mode we consider in the following). This bifurcation is a transition between two dynamical states of oscillation: one of small and one of large amplitude, which can be easily detected using commercially available cryogenic amplifiers. The transition rate between these two states depends on experimental parameters (written generically as variable X) which are for instance the resonant frequency of the mode, its quality factor, the frequency and amplitude of the microwave driving (see [23–25] for calculations of this rate). Some of these experimental parameters may be controlled and some others might be the subject of random fluctuations, which can be detected by the JBA.

Zoom In Zoom Out Reset image size

As a first characterization of the sensitivity of this system as a detector, we repeatedly probe the JBA with microwave driving pulses (duration $\tau =35$ μs $\approx \;10/\gamma$ where γ is the linewidth of the mode, in order to damp any transient) at a given rate ${{\nu }_{{\rm{rep}}}}$ (5 kHz and 500 Hz), and we record the state of the resonator at the end of the driving pulse (labeled ${{Y}_{k}}$ at time step ${{t}_{k}}=k/{{\nu }_{{\rm{rep}}}}$) in binary format (${{Y}_{k}}=1$ for the high amplitude state, ${{Y}_{k}}=0$ for the low amplitude state). Counting over $N={{10}^{3}}$ events, one obtains the average switching probability ${{p}_{{\rm exp} }}=1/N\mathop{\sum }_{k=1}^{N}{{Y}_{k}}$ for a given set of parameters X. Actually, as discussed in the following, X might undergo random fluctuations, meaning that the experimentalist can control only the average value $\langle X\rangle$ of X. Recording the switching probability while ramping the control parameter $\langle X\rangle$ across the bifurcation frontier provides the switching probability curve or 'S-curve' whose 10%–90% width $\Delta X$ defines its sensitivity to fluctuations of X, that is, a shift of X by an amount $\Delta X$ can be detected within a single probing pulse with a high level of confidence. Table 1 summarizes all notations used in this paper.

The natural question which arises now is: Is it possible to infer more information from the array $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$ than just the average switching probability ${{p}_{{\rm exp} }}$? A first improvement is to measure the average switching probability $p_{{\rm exp} }^{n}$ over subsets of n events. This allows the detection of fluctuations of X of order $\Delta X/\sqrt{n}$ still with a high level of confidence but with a much lower bandwidth of the order of ${{\nu }_{{\rm{rep}}}}/n$. Experimentally, we do observe fluctuations of the switching probability ${{p}_{{\rm exp} }}=p_{{\rm exp} }^{N}$ which are well above the expected statistical noise ($\approx 1/\sqrt{N}$ where $N={{10}^{3}}$ events), indicating a low frequency noise present in the experimental parameters. In addition, the measured experimental value for the width of the switching curve as a function of the frequency of the mode: $\Delta \nu \approx 4.5$ kHz is greater by a factor 2 from the theoretical prediction obtained from the Dykman model [26]. Both facts indicate that a non-negligible part of the switching curve width is due to fluctuations of the experimental parameters. We are thus led to consider a 'doubly stochastic process': the outcome of the detector is a random process ${{Y}_{k}}$ depending on a switching probability which is itself a random process (since it depends on a noisy parameter X). After characterizing the noise level of our microwave setup, we can exclude frequency fluctuations of the probing pulses and microwave amplitude fluctuations at the level of the sample. The most likely origin of the noise source is microscopic and on-chip. Such fluctuations may be characterized as inducing fluctuations in the resonant frequency of the cavity.

We will now show that it is possible to extract the spectrum of the frequency fluctuations of the cavity from the binary array $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$. For this purpose, we need to introduce a general model for our detector (which, we note, can be applied to any system involving a threshold effect).

We consider the JBA as a generic bistable system: its state labelled Y is a random variable which can take two different values (Y = 1 or 0) with probabilities dependent on whether some parameter X is above or below a threshold value ${{x}_{0}}$. By offsetting X, we set ${{x}_{0}}$ to 0 in the following. Considering first the ideal case where thermal and quantum noises are absent from this detector, we have a 'sharp' threshold: Y = 0 with certainty if $X<0$ and Y = 1 with certainty if $X>0$. In this case, this detector is completely analogous to a 1-bit analog to digital converter (see figure 2(a)): $Y=Q\left[ X \right]$ where Q is a digitizer function ($Q\left[ X \right]=1$ if $X>0$ and $Q\left[ X \right]=0$ if $X<0$). Now, consider a time varying ${{X}_{t}}$ named an 'input' signal, which is sampled by this threshold detector at regular time intervals ${{t}_{k}}=k/{{\nu }_{{\rm{rep}}}}$ ($k\in \left\{ 1,..,N \right\}$ to obtain a binary array of outcomes $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$. What can be inferred about the temporal variations of ${{X}_{t}}$ from the $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$ array? Provided that the signal has some frequency components only for $f<{{\nu }_{{\rm{rep}}}}/2$ (the Shannon criterion) and that the mean value $\left\langle X\right\rangle$ of the signal is within a fraction of its standard deviation ${{\sigma }_{X}}$ of the threshold, then the answer is: only some crude information about the largest Fourier component of ${{X}_{t}}$ can be inferred. For instance, a sinusoidal input oscillation of X around the threshold will be converted into a rectangular output, thus corrupting the spectrum of ${{Y}_{k}}$ with harmonic generation (see figures 3(a) and (b)).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

However, in analog to digital conversion systems, it is common practice to add a small amount of noise to the input signal ${{X}_{t}}$ prior to digitization (see figure 2(c)) in order to control or 'shape' the associated distortion. For instance, this technique is implemented in fast oscilloscopes where ${{X}_{t}}$ is an input voltage and the noise is generated by the pre-amplifier stage of the scope. A careful engineering of this input noise can increase the effective resolution of the converter at the expense of the sampling rate: this is the so-called 'oversampling' technique [27].

We will now consider an analogous situation for our threshold detector and focus on its implication for the spectrum of the digitized signal. A random variable ${{D}_{k}}$ is added to the input signal prior to thresholding (see figure 2(c)), such that the output signal is now ${{Y}_{k}}=Q\left[ {{X}_{k}}+{{D}_{k}} \right]$, where ${{X}_{k}}={{X}_{{{t}_{k}}}}$ is the input sampled at time ${{t}_{k}}$. We assume that the ${{D}_{k}}$ variables are independent and generated from a stationary random process with zero mean and probability density ${{P}_{D}}\left( d \right)$. We focus first on the statistics of order one of this probabilistic model.

3.1. First order statistics of the model

For a single sampling of the input signal ${{X}_{k}}$ at a given time ${{t}_{k}}$ having the value x, the outcome of the digitization process ${{Y}_{k}}$ takes the value 1 with probability p(x) and the value 0 with probability $1-p\left( x \right)$ where p(x) is the conditional probability

$$\begin{eqnarray}&&p(x)=\mathbb{P}\left( {{Y}_{k}}=1/{{X}_{k}}=x \right)\end{eqnarray} \tag{ 1 }$$

to observe ${{Y}_{k}}=1$ knowing that ${{X}_{k}}=x$. This probability can be related to the cumulative distribution of the ${{D}_{k}}$ variables, ${{P}_{D}}\left( d \right)$:

$$\begin{eqnarray*} p\left( x \right) & = & \mathbb{P}\left( {{X}_{k}}+{{D}_{k}}>0/{{X}_{k}}=x \right) \\ {} & = & \mathbb{P}\left( {{D}_{k}}>-x \right)=\mathop{\int }_{-x}^{\infty }{{P}_{D}}\left( u \right){\rm{d}}u \\ \end{eqnarray*}$$

As a result of the D fluctuations, the 'sharp' threshold of the ideal quantizer is broadened by the distribution ${{P}_{D}}$ (see figure 2(d)). However, experimentally ${{X}_{t}}$ can fluctuate over time so we cannot access directly p(x). Instead, one is measuring an average probability ${{p}_{{\rm exp} }}$ calculated over many probing pulses $({\rm{at}}\;{{t}_{1}},{{t}_{2}},\ldots ,{{t}_{N}})$:

$$\begin{eqnarray}&&{{p}_{{\rm exp} }}=\frac{1}{N}\mathop{\sum }_{k=1}^{N}{{Y}_{k}}\end{eqnarray} \tag{ 2 }$$

which, in the limit of the law of large numbers, can be approximated by:

$$\begin{eqnarray}&&\approx \frac{1}{N}\mathop{\sum }\limits_{N}^{k=1}\langle {{Y}_{k}}\rangle \end{eqnarray} \tag{ 3 }$$

where $\langle {{Y}_{k}}\rangle =\mathbb{P}\left( {{Y}_{k}}=1 \right)$. So what is the value of ${{p}_{{\rm exp} }}$ knowing that ${{X}_{t}}$ can fluctuate over time? To answer this question, we need to assume two more hypotheses on the process ${{X}_{t}}$: first, ${{X}_{t}}$ undergoes a random stationary process with a distribution probability ${{P}_{X}}$ centered around an average value $\langle X\rangle$ which can be controlled experimentally (such as an average magnetic flux or an average gate voltage). Then we need to assume a quasi-static approximation: the fluctuations of ${{X}_{t}}$ should be slower than the sampling time, (i.e. the duration of a single microwave probing pulse in the case of the JBA). With these two hypotheses, $\mathbb{P}\left( {{Y}_{k}}=1 \right)$ does not depend on k and is the average of p(x) weighted by the distribution of X:

$$\begin{eqnarray*} {{p}_{{\rm exp} }} & \approx & \mathbb{P}\left( {{Y}_{k}}=1 \right)=\mathop{\int }_{-\infty }^{+\infty }\mathbb{P}\left( {{Y}_{k}}=1/{{X}_{k}}=x \right){{P}_{X}}\left( x \right){\rm{d}}x \\ {} & = & \mathop{\int }_{-\infty }^{+\infty }p\left( x \right){{P}_{X}}\left( x \right){\rm{d}}x. \\ \end{eqnarray*}$$

Setting ${{X}_{t}}=\delta {{X}_{t}}+\langle X\rangle$, we can rewrite ${{p}_{{\rm exp} }}$ as a function of $\langle X\rangle$:

$$\begin{eqnarray}&&{{p}_{{\rm exp} }}(\langle X\rangle )\approx \mathop{\int }_{-\infty }^{+\infty }p(x){{P}_{\delta X}}(x-\langle X\rangle ){\rm{d}}x=p\otimes {{P}_{\delta X}}(\langle X\rangle )\end{eqnarray} \tag{ 4 }$$

The experimental probability ${{p}_{{\rm exp} }}$ of detection considered as a function of the control parameter $\langle X\rangle$ is thus the convolution of the response of the detector p(x) with the probability distribution of X. The p(x) response curve of the detector, already broadened by the fluctuations of the ${{D}_{k}}$, is further broadened by the fluctuations of the ${{X}_{t}}$ process (see figure 2(f)). The threshold of the digitization process is no longer 'sharp', it has an 'S like' shape with a 10%–90% width (defined as $\Delta X$), which can be related to the standard deviations of X and D: $\Delta X\approx 2.56\sqrt{\sigma _{D}^{2}+\sigma _{X}^{2}}$ in the case of gaussian distributions for ${{D}_{k}}$ and ${{X}_{k}}$. We will see that this relation is useful for calibrating our detector. Having studied the first order statistics of our detection model, we now focus on the second order statistics and demonstrate that the spectrum of the X parameter can be extracted from the experimental binary array $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$.

3.2. Second order statistics of the model: autocorrelation and spectral density

We consider the autocorrelation of the $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$ array and show here that it can be related to the autocorrelation of the ${{X}_{t}}$ process. We define first the fluctuation $\delta {{Y}_{k}}={{Y}_{k}}-{{p}_{{\rm exp} }}$. We know that $\langle {{Y}_{k}}\rangle ={{p}_{{\rm exp} }}$ and $\langle Y_{k}^{2}\rangle =\langle {{Y}_{k}}\rangle$ (since ${{Y}_{k}}$ takes only two values 0 or 1). As a consequence the variance $\sigma _{Y}^{2}$ of the ${{Y}_{k}}$ process is:

$$\begin{eqnarray}&&\sigma _{Y}^{2}=\langle \delta Y_{k}^{2}\rangle ={{p}_{{\rm exp} }}(1-{{p}_{{\rm exp} }})\end{eqnarray} \tag{ 5 }$$

Then for $q\ne 0$, we have:

$$\begin{eqnarray} \left\langle {{Y}_{k+q}}{{Y}_{k}}\right\rangle & = & \mathbb{P}\left( \left\{ {{Y}_{k+q}}=1 \right\}\;{\rm{AND}}\;\left\{ {{Y}_{k}}=1 \right\} \right) \\ {} & = & \mathbb{P}\left( \left\{ {{X}_{k+q}}+{{D}_{k+q}}>0 \right\}\;{\rm{AND}}\;\left\{ {{X}_{k}}+{{D}_{k}}>0 \right\} \right) \\ {} & = & \mathop{\int }_{0}^{\infty }\mathop{\int }_{0}^{\infty }{{P}_{U,U}}\left( {{u}_{1}},{{u}_{2}} \right){\rm{d}}{{u}_{1}}{\rm{d}}{{u}_{2}} \\ \end{eqnarray} \tag{ 6 }$$

where $\mathbb{P}\left( \left\{ {{Y}_{k+q}}=1 \right\}\;{\rm{AND}}\;\left\{ {{Y}_{k}}=1 \right\} \right)$ is the joint probability to have the events ${{Y}_{k+q}}=1$ and ${{Y}_{k}}=1$. ${{P}_{U,U}}$ is the joint probability density of ${{U}_{1}}={{X}_{k}}+{{D}_{k}}$ and ${{U}_{2}}={{X}_{k+q}}+{{D}_{k+q}}$. Such a probability density is the double convolution of the joint probability of X: ${{P}_{X,X}}$, with the joint probability of D: ${{P}_{D,D}}\left( {{d}_{1}},{{d}_{2}} \right)$:

$$\begin{eqnarray}&&{{P}_{U,U}}({{u}_{1}},{{u}_{2}})=\left( {{P}_{D,D}}\otimes {{P}_{X,X}} \right)({{u}_{1}},{{u}_{2}})\end{eqnarray} \tag{ 7 }$$

Because the two random variables ${{D}_{1}}$ and ${{D}_{2}}$ are assumed to be independent and identically distributed, we have that ${{P}_{D,D}}\left( {{d}_{1}},{{d}_{2}} \right)={{P}_{D}}\left( {{d}_{1}} \right)\;{{P}_{D}}\left( {{d}_{2}} \right)$. To go further we need to make more assumptions about the statistics of the ${{X}_{t}}$ process.

A Gaussian hypothesis for ${{X}_{t}}$ is physically reasonable since we are dealing with a condensed matter system where sources of noises involve a priori large numbers of uncorrelated fluctuating subsystems. In addition, since we are interested only in the second order statistics, such a gaussian hypothesis for ${{X}_{t}}$ is all that is necessary. Finally, this hypothesis will provide analytical formulas for the relation between the autocorrelations of ${{Y}_{k}}$ and ${{X}_{k}}$. We thus assume that the ${{X}_{t}}$ process is stationary and has the following joint probability density:

$$\begin{eqnarray}&&{{P}_{X,X}}\left( {{x}_{1}},{{x}_{2}};t \right)=\frac{1}{2\pi \sigma _{X}^{2}\sqrt{1-\rho _{X}^{2}}}{{e}^{-\frac{\left( x_{1}^{2}+x_{2}^{2}-2{{\rho }_{X}}{{x}_{1}}{{x}_{2}} \right)}{2\left( 1-\rho _{X}^{2} \right)\sigma _{X}^{2}}}}\end{eqnarray} \tag{ 8 }$$

where ${{\rho }_{X}}\left( t \right)=\left\langle \delta {{X}_{t}}\delta {{X}_{0}}\right\rangle /\sigma _{X}^{2}$ is the normalized autocorrelation (in equation (8), the dependence of ${{\rho }_{X}}$ on time is omitted for brevity).

The double convolution of equation (7) and the integration in equation (6) can then be calculated analytically and gives the main result of this paper: a direct relationship between the normalized autocorrelations of ${{X}_{k}}$ and ${{Y}_{k}}$: ${{\rho }_{X}}$ and ${{\rho }_{Y}}$,

$$\begin{eqnarray}&&{{\rho }_{X}}\left( {{t}_{q}} \right)=\frac{\left\langle \delta {{X}_{k+q}}\delta {{X}_{k}}\right\rangle }{\sigma _{X}^{2}},\qquad {{\rho }_{Y}}\left( {{t}_{q}} \right)=\frac{\langle \delta {{Y}_{k+q}}\delta {{Y}_{k}}\rangle }{\sigma _{Y}^{2}}\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{{\rho }_{Y}}\left( {{t}_{q}} \right)=\frac{2}{\pi }{\rm arctan} \left( \frac{{{\rho }_{X}}\left( {{t}_{q}} \right)}{\sqrt{{{\left( 1+\frac{\sigma _{D}^{2}}{\sigma _{X}^{2}} \right)}^{2}}-{{\rho }_{X}}{{\left( {{t}_{q}} \right)}^{2}}}} \right).\end{eqnarray} \tag{ 10 }$$

Note that equation (10) is valid for $q>0$ only: two different samples have to be considered (when q = 0, there is no relation between $\langle \delta Y_{k}^{2}\rangle ={{p}_{{\rm exp} }}\left( 1-{{p}_{{\rm exp} }} \right)$ and $\langle \delta {{X}^{2}}\rangle$). This point will be important when considering the Fourier transform of the autocorrelation to obtain the spectrum. It is then useful to define a transfer function ${{\rho }_{Y}}/{{\rho }_{X}}$ which is plotted as a function of ${{\rho }_{X}}$ in figure 4 for different values of ${{\sigma }_{D}}/{{\sigma }_{X}}$. We consider the two limiting cases:

• When ${{\sigma }_{D}}\to 0$, the additive noise D disappears, and equation (9) simplifies to:
  $$\begin{eqnarray}&&{{\rho }_{Y}}({{t}_{q}})=\frac{\langle \delta {{Y}_{k+q}}\delta {{Y}_{k}}\rangle }{\sigma _{Y}^{2}}=\frac{2}{\pi }{\rm arcsin} ({{\rho }_{X}}({{t}_{q}})).\end{eqnarray} \tag{ 11 }$$
  This is a strong non-linear relationship between $\langle \delta {{Y}_{k+q}}\delta {{Y}_{k}}\rangle$ and $\langle \delta {{X}_{k+q}}\delta {{X}_{k}}\rangle$, which accounts for harmonic generation.
• The other limit case is the more interesting one: when ${{\sigma }_{X}}\lesssim {{\sigma }_{D}}$, one finds a quasi-linear relation between the autocorrelation of X and the autocorrelation of Y:
  $$\begin{eqnarray}&&\langle \delta {{Y}_{k+q}}\delta {{Y}_{k}}\rangle \approx \frac{2}{\pi }\frac{\sigma _{Y}^{2}}{\sigma _{D}^{2}+\sigma _{X}^{2}}\langle \delta {{X}_{k+q}}\delta {{X}_{k}}\rangle \end{eqnarray} \tag{ 12 }$$
  valid for $q>0$.

Zoom In Zoom Out Reset image size

The linearity of equation (12) gives a direct access to the autocorrelation of ${{X}_{k}}$ from the experimentally measured autocorrelation of the $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$ array. The harmonic distortion due to the thresholding is suppressed by the addition of the noise D to the input signal. This is done at the expense of the 'reduced gain' ${{\rho }_{Y}}/{{\rho }_{X}}$ which decreases as ${{\sigma }_{D}}$ increases. There is thus a tradeoff between linearity and gain. For instance, when ${{\sigma }_{D}}={{\sigma }_{X}}$, ${{\rho }_{Y}}/{{\rho }_{X}}\approx 0.32$ and is constant within $\pm 1\%$.

The gain $\frac{2\sigma _{Y}^{2}}{\pi \left( \sigma _{D}^{2}+\sigma _{X}^{2} \right)}\approx \frac{2}{\pi }{{\left( \frac{2.56}{\Delta X} \right)}^{2}}\sigma _{Y}^{2}$ may be calibrated. The variance ${{\sigma }_{Y}}$ of ${{Y}_{k}}$ may be obtained directly from experiment via $\sigma _{Y}^{2}={{p}_{{\rm exp} }}\left( 1-{{p}_{{\rm exp} }} \right)$. The width of the switching curve $\Delta X$ is also easily accessible experimentally. Further, when the repetition rate increases, the correlation between two successive outcomes $\langle \delta {{Y}_{k+1}}\delta {{Y}_{k}}\rangle$ converges to:

$$\begin{eqnarray}&&\mathop{{\rm lim} }\limits_{{{\nu }_{{\rm{rep}}}}\to \infty }\langle \delta {{Y}_{k+1}}\delta {{Y}_{k}}\rangle =\frac{2}{\pi }\frac{\sigma _{Y}^{2}}{\sigma _{D}^{2}+\sigma _{X}^{2}}\sigma _{X}^{2}\end{eqnarray} \tag{ 13 }$$

These relations may be used, in principle, to provide estimates of the values of ${{\sigma }_{X}}$ and ${{\sigma }_{D}}$.

Equation (12) can be rewritten in the frequency domain: by Fourier transforming equation (12) and using the Wiener–Khinchin theorem, one obtains the relation between the spectral densities of ${{Y}_{k}}$ and ${{X}_{k}}$ $({\rm{respectively}}\;{{S}_{Y}}\;{\rm{and}}\;{{S}_{X}})$:

$$\begin{eqnarray}&&{{S}_{Y}}\left( \nu \right)=\frac{2}{\pi }\frac{\sigma _{Y}^{2}}{\sigma _{D}^{2}+\sigma _{X}^{2}}{{S}_{X}}\left( \nu \right)+2\frac{\sigma _{Y}^{2}}{{{\nu }_{{\rm{rep}}}}}\end{eqnarray} \tag{ 14 }$$

which shows that the digitization noise $\sigma _{Y}^{2}$ (coming from the autocorrelation of ${{Y}_{k}}$ at q=0) is spread as a white background over the acquisition bandwidth ${{\nu }_{{\rm{rep}}}}/2$. This constant background gives the sensitivity at which ${{S}_{X}}$ can be measured. It is important to stress that this noise level can be squeezed down just by increasing the sampling rate.

As an example, we consider again the case of a sinusoidal ${{X}_{k}}$ such that ${{\sigma }_{D}}={{\sigma }_{X}}$ (see figure 3(c)). In this case, the harmonic distortion of the noiseless 1-bit $A/D$ converter is suppressed and replaced by a white background in the spectrum of ${{Y}_{k}}$ (figure 3(d)).

The point we want to make now is that this 1-bit analog to digital conversion with an additive random noise is analogous to our JBA with thermal and quantum noises taken into account. The effects of thermal and quantum noises should be seen as the addition of a gaussian random variable ${{D}_{k}}$ to ${{X}_{k}}$ prior to thresholding. The ${{D}_{k}}$ are assumed to be independent since the sampling interval is much larger than the reset time of the detection process. The probability distribution of the noise ${{D}_{k}}$ is directly related to the switching curve (see figure 2(f)), and can be obtained from the Dykman model in both the thermal and quantum regimes [23–25].

We set the working point of our experiment to ${{p}_{{\rm exp} }}\approx 1/2$ and then record the outcome of the JBA as a binary array over a time $\tau \approx 3$ min at a temperature of 11 mK $\pm 2$ mK (measured with a PdFe magnetic susceptibility thermometer) and for two different repetitions rates: 500 Hz and 5 kHz. We first note that we do not see any dependence of the switching curves on the repetition rate, which allows us to exclude heating effects as a source of correlations (see figure 5(b)). We note from figure 5(c) that ${{\rho }_{Y}}\leqslant 0.16$ which places this data set in the linear regime whatever the value of ${{\sigma }_{D}}$. Hence equation (14) is valid. The spectral density ${{S}_{Y}}$ is then computed from the array $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$ (see figure 5(a)) using a fast Fourier transform routine. From the width of the switching curve $\Delta X\approx 2.56\sqrt{\sigma _{X}^{2}+\sigma _{D}^{2}}\approx 4.5$ kHz and the variance $\sigma _{Y}^{2}={{p}_{{\rm exp} }}\left( 1-{{p}_{{\rm exp} }} \right)\approx 0.25$ we obtain, with the use of equation (14), the 'gain' between the autocorrelations of ${{X}_{k}}$ and ${{Y}_{k}}$: $\frac{2}{\pi }\frac{\sigma _{Y}^{2}}{\sigma _{X}^{2}+\sigma _{D}^{2}}\approx 5\times {{10}^{-8}}$ Hz$^{-2}$ and the conversion factor from ${{S}_{Y}}$ to the fractional frequency noise spectrum ${{S}_{\delta \nu /\nu }}$. This is displayed as the right-hand scale of figure 5(a).

Zoom In Zoom Out Reset image size

As an aside, we note from figure 5(c) that our experimental repetition rate was not fast enough for $\left\langle \delta {{Y}_{k+1}}\delta {{Y}_{k}}\right\rangle$ to saturate so we are unable here to use equation (13) to extract values for ${{\sigma }_{X}}$ and ${{\sigma }_{D}}$, though we are able to put bounds on these values as follows. Assuming that ${{\rho }_{Y}}$ increases as the repetition rate increases we may put a lower bound on the value of ${{\sigma }_{X}}$. From equation (13),

$$\begin{eqnarray*}&&\frac{2}{\pi }\frac{\sigma _{Y}^{2}}{\sigma _{D}^{2}+\sigma _{X}^{2}}\sigma _{X}^{2}\geqslant 0.16*\sigma _{Y}^{2}\end{eqnarray*}$$

where we have used the maximal experimental value of ${{\rho }_{Y}}$ from figure 5(c). Hence $\frac{2}{\pi }{{\left( \frac{2.56}{\Delta X} \right)}^{2}}\sigma _{X}^{2}\geqslant 0.16$ such that ${{\sigma }_{X}}\geqslant 0.9$ kHz. This implies a maximum possible value of ${{\sigma }_{D}}\leqslant 1.5$ kHz given the experimental width of the switching curve $\Delta X\approx 4.5$ kHz. For the same reason the maximum possible value of ${{\sigma }_{X}}$ is 1.8 kHz for ${{\sigma }_{D}}=0$. For this experimental data set, we conclude that 0.9 kHz $\leqslant \;{{\sigma }_{X}}\leqslant 1.8$ kHz and $0\leqslant {{\sigma }_{D}}\leqslant 1.5$ kHz or that $0\leqslant \frac{{{\sigma }_{D}}}{{{\sigma }_{X}}}\leqslant 2$. Dykman predicts ${{\sigma }_{D}}\approx 0.9$ kHz [21]. We note that the upper bound on the variance ${{\sigma }_{X}}\leqslant 1.8$ kHz $\approx \;0.4$ ppm of the frequency of the resonant mode is comparable to the state of the art in superconducting quantum bits achieved with 3D cavities [4]. Figure 5(a) shows as expected that the white background noise corresponding to the digitization noise is present, its level agrees well with the prediction which is plotted on figure 5(a) for the two repetition rates 500 Hz and 5 kHz. This white background gives us the sensitivity of the spectrum measurement and can be squeezed down by increasing the repetition rate .

Table 1.  Notations

Symbol	Definition/result
${{D}_{k}}$	Added noise prior to thresholding at time step ${{t}_{k}}$ considered as a random discrete variable.
${{n}_{{\rm{ph}}}}$	Photon number in the third harmonic mode of the superconducting cavity.
${{\nu }_{{\rm{rep}}}}$	Repetition rate of the acquisition process (typically up to a few kHz).
${{P}_{X}}\left( x \right),{{P}_{D}}\left( d \right)$	Probability density of the random variables X and D.
${{P}_{X,X}}\left( {{x}_{1}},{{x}_{2}};t \right)$	Joint probability density of the random stationary process ${{X}_{t}}$.
$\mathbb{P}\left( {{\omega }_{1}}/{{\omega }_{2}} \right)$	Conditional probability for event ${{\omega }_{1}}$ to happen knowing that event ${{\omega }_{2}}$ has happened.
p(x)	Shorter notation for the conditional probability : $\mathbb{P}\left( {{Y}_{k}}=1/{{X}_{k}}=x \right)$.
${{p}_{{\rm exp} }}$	Experimental switching probability, obtained by counting bifurcation events over $\approx {{10}^{3}}$ sampling pulses.
Q	Quantizer function: $Q\left[ x \right]=1$ if $x>0$ and $Q\left[ x \right]=0$ if $x<0$.
${{\rho }_{X}}\left( t \right)$	Normalized autocorrelation of the ${{X}_{t}}$ process: $\langle \delta {{X}_{t}}\delta {{X}_{0}}\rangle /\sigma _{X}^{2}$
${{\rho }_{Y}}\left( {{t}_{q}} \right)$	Normalized autocorrelation of Y: $\langle \delta {{Y}_{k+q}}\delta {{Y}_{k}}\rangle /\sigma _{Y}^{2}$
${{\sigma }_{X}},{{\sigma }_{Y}},{{\sigma }_{D}}$	Standard deviations of the random variables $X,Y,D$.
${{S}_{Y}}\left( \nu \right)$	Spectral density of the binary array $\left\{ {{Y}_{1}},{{Y}_{2}},\ldots ,{{Y}_{N}} \right\}$.
${{S}_{X}}\left( \nu \right)$	Spectral density of the random process ${{X}_{t}}$.
${{t}_{k}}=k/{{\nu }_{{\rm{rep}}}}$	$k{\rm{th}}$ sampling time.
$\langle X\rangle$	Average of the random variable X.
$\delta X=X-\langle X\rangle$	Fluctuation of X.
$\Delta X$	10%–90% width of the switching curve: ${{p}_{{\rm exp} }}$ as a function of $\langle X\rangle$.
${{X}_{t}}$	Input of the detector considered as a time dependent random process.
X	Shorter notation for ${{X}_{t}}$ when the time dependence can be omitted.
${{X}_{k}}$	Shorter notation for ${{X}_{{{t}_{k}}}}$, the input sampled at time step ${{t}_{k}}$ considered as a discrete random variable.
$\langle \delta {{X}_{k+q}}.\delta {{X}_{k}}\rangle$	Autocorrelation of the input signal of the detector.
${{Y}_{k}}$	Output of the detector at time step ${{t}_{k}}$ considered as a discrete random variable.
$\langle \delta {{Y}_{k+q}}.\delta {{Y}_{k}}\rangle$	Autocorrelation of the output signal of the detector.

It is interesting to compare this sensitivity to a fundamental scale which is the standard quantum limit of a weak continuous measurement of the frequency of a resonator [28] in comparable experimental conditions: average photon number in the cavity (here $\bar{n}\approx 2{{n}_{c}}\approx 500$) leaking at rate $2\pi \gamma$ (here $\approx 2$ MHz). The frequency fluctuations of the resonator equivalent to the shot noise of the driven coherent state are given by $S_{\delta \nu /\nu }^{sn}=2\nu _{0}^{2}/\left( {{\gamma }^{2}}\dot{n} \right)$ where $\dot{n}=2\pi \gamma \bar{n}$. Remarkably, for the maximal theoretical repetition rate of this detector (${{\nu }_{{\rm{rep}}}}\approx \gamma /5\approx 50$ kHz) the theoretical prediction for the sensitivity of the bifurcation as a noise spectrum analyzer would be comparable to the standard quantum limit. Experimentally, we used a maximum repetition rate of 5 kHz, giving a sensitivity within an order of magnitude of the standard quantum limit.

In addition to the digitization noise, a significant $A/f$ frequency noise is present in our sample with $A\approx {{10}^{-15}}$. From flux modulation measurements [21], we can put an upper bound on the contribution of flux noise at the optimal working point ($\phi =0$ where sensitivity to flux noise is only second order), and show that it has negligible contribution. In addition, because of the small value of the participation ratio ${{L}_{{\rm{SQUIDs}}}}/{{L}_{{\rm{tot}}}}\approx 2.5\%$ (where ${{L}_{{\rm{SQUIDs}}}}$ is the total inductance of the SQUID array, and ${{L}_{{\rm{tot}}}}$ the total inductance of the cavity), critical current noise has also negligible contribution. Finally, we note that the noise amplitude and frequency dependence observed is comparable to previous observations made in Kinetic Inductance Detectors [14, 29, 30]. We conclude that dielectric noise is probably the source for the observed $1/f$ noise in this device.

We have presented a model that provides a deeper insight into threshold detectors. This model allows direct access to the spectral density of any noise source coupled to such detectors and is reminiscent of noise shaping with 'dithering' in analog to digital conversion. It was applied to measure the frequency fluctuations of a JBA demonstrating the presence of a $1/f$ noise whose amplitude and spectrum is compatible with previous observations of dielectric noise in superconducting resonators. The main advantage of this technique as an on-chip detector, is its dispersive nature which avoids the dissipation and backaction associated with the voltage state of a SQUID amplifier or switched hysteretic junction. This allows a lower thermalization temperature of the degrees of freedom considered. The sensitivity of this technique as a noise spectrometer is potentially of the order of the standard quantum limit of a weak continuous measurement. The potential of this technique for the extensive characterization of decoherence sources in superconducting quantum bits circuits is thus high. It could provide in situ measurement of noises of any origin, including magnetic, charge, critical current, dielectric, kinetic inductance noises. They can be measured most effectively if the coupling is tunable. In addition, the detection bandwidth of this method is half the repetition rate which is in our case limited by the reset time of the bifurcation detector. As the sensitivity of the bifurcation depends weakly on the quality factor [21], a lower Q than that used in our experiment could allow repetition rates of order of several hundreds of MHz. Obtaining the noise spectrum over this frequency range with a lower digitization noise would be of great interest. Finally, we note that only partial information on a random process is provided by the second order statistics. As a consequence, it would be interesting to generalize this method to higher order correlators. Apart from qubit diagnostics, the technique may have important applications for the measurement of the full counting statistics of a quantum conductor [31].

We wish to thank D Estève and all at the Quantronics Group at CEA Saclay for their support over many years, especially P Bertet and A Palacios-Laloy who fabricated the sample, A Tzalenchuk and T Lindström (NPL) for helpful discussions and the loan of equipment, and John Taylor and Howard Moore for technical help. G Ithier acknowledges financial support from the Leverhulme Trust (Early Career Fellowship SRF-40311) and P J Meeson acknowledges financial support from the EPSRC (grants EP/D001048/1 and EP/F041128/1) and the EXL03 Microphoton project of the European Metrology Research Program. The EMRP is jointly funded by the EMRP participating countries within EURAMET and the European Union.