Towards a microwave single-photon counter for searching axions

The major task of detecting axions or axion-like particles has two challenges. On the one hand, the ultimate sensitivity is required, down to the energy of a single microwave photon of the yoctojoule range. On the other hand, since the detected events are supposed to be rare, the dark count rate of the detector must be extremely low. We show that this trade-off can be approached due to the peculiar switching dynamics of an underdamped Josephson junction in the phase diffusion regime. The detection of a few photons’ energy at 10 GHz with dark count time above 10 s and the efficiency close to unity was demonstrated. Further enhancements require a detailed investigation of the junction switching dynamics.


INTRODUCTION
A long-standing puzzle in quantum chromodynamicswhy charge conjugation parity symmetry seems to be preserved 1,2belongs to the list of unsolved problems in physics. A proposed solution to the problem is novel pseudoscalar particles called axions. Being the imprint of dark matter, axions can also resolve the problem of missing mass in the Universe as well as the anisotropy of the cosmic microwave background [3][4][5] . The theory predicts that in a magnetic field, galactic axions are converted into microwave photons 2,6 , attracting the renewed interest to the axions search [6][7][8][9][10][11] .
The development of quantum technology requires efficient detection of electromagnetic waves at a single-photon level. A practical solution for the microwave frequency range is still challenging. Indeed, the sensitivity of conventional detectors such as transition-edge sensors 12 or superconducting nanowire 13,14 is not sufficient since the microwave photon energy is very low, just 6.6 yoctojoule (6.6 × 10 −24 J) at 10 GHz. Moreover, the sensitivity of the most widespread devices, such as Superconducting Quantum Interference Devices (SQUIDs) or Josephson parametric amplifiers, is inevitably restricted by quantum fluctuations 15 , making them unfavourable for the registration of weak photon fluxes and almost useless for detecting individual photons.
The circuit quantum electrodynamics platform, based on superconducting qubits with microwave-ranged transition frequency, can potentially provide all necessary elements [16][17][18][19][20][21][22][23][24][25][26] to realise a procedure for single microwave photon measurements. Indeed, several protocols, differing in specific purpose, have been demonstrated. However, their detection efficiency and the dark count time are restricted by low coherence times of qubits, which is unacceptable for searching axions.
Specific features of the axions dynamics dictate an optimised choice of the mentioned above detector. The events of the axionphoton conversion are expected to be rare. This means that the low dark count rate of such a detector is required. On the other hand, the destructive nature of the measurements, as well as relatively long initialisation time, is not of practical importance. Consequently, the implementation of the next generation of threshold single-photon counters is needed for the progress in axions search 27,28 .
To realise a microwave Single-Photon Counter (SPC), we use a conventional Al-based tunnel Josephson junction in a currentbiased regime. This junction, coupled to a photon field, caused its switching to the resistive state by the absorbed photons, giving an output voltage of the order of 0.4 mV. While some attempts to exploit Josephson tunnel junctions for SPC applications have been performed earlier [29][30][31][32][33][34] , the following trade-off limits its potential ability. On the one hand, to obtain maximal sensitivity, a junction should be biased close to the photon-assistant switch to the resistive state. On the other hand, to minimise the dark count rate, the bias current should not be close to the transition point.
In this work, we have demonstrated 5-photon detection at 10 GHz frequency with efficiency close to unity with the dark count time of about 10 s. At the same time, the single-photon sensitivity is achieved with an efficiency of 0.02 and with the dark count time above 0.01 s. The dark count times, much larger than predicted from the existing theories, have been reached due to the phase diffusion regime [34][35][36][37][38][39][40][41] . This shows intriguing perspectives for the implementation of a Single-Photon Counter for axion search.

RESULTS
The aluminium superconductor-insulator-superconductor (SIS) tunnel junction with the area of 0.8 μm 2 has been measured at the temperatures of 10-300 mK in a low noise environment 42,43 . The scheme of the experimental setup is shown in Fig. 1. The experiment was performed in three steps. The first one was the lifetime measurements without an itinerant microwave field to get the SPC dark count time. The second step was the procedure of calibration to find out the actual power absorbed by the SIS junction. And the third step was to test the SPC performance as a function of a microwave power P 0 and bias current I.
For SIS junctions with critical currents of the order of 1 μA and above, it has been shown theoretically in 44,45 that it is possible to tune the parameters to get the single-photon sensitivity along with a sufficiently long lifetime. In this work, we choose for the experimental study a SIS junction with a small critical current of only 23 nA, expecting high single-photon sensitivity due to the lower threshold, with a trade-off of decreasing the dark count time.
The conventional theory of escape, both in the thermal [46][47][48][49] , and the quantum 47,48 limits, predicts a lifetime for such junctions of the order of nanoseconds. However, the dark count time can be seriously increased due to the phase diffusion regime [34][35][36][37][38][39][40][41] at the millikelvin temperatures. The supercurrent across the SIS junction is defined by the phase difference φ(t) between its superconducting electrodes. The dynamics of this phase difference can be modelled by an effective particle, moving in the washboard potential 50 . Basically, the escape of the particle from the potential well leads to its continuous motion down the potential profile, i.e., to switch into the resistive state with a finite voltage. However, if the dissipation is high enough, the escape of the particle from the well does not always lead to the appearance of a constant voltage because it stops in the adjacent well, which corresponds to a return into the superconducting state. Such SIS junction dynamics with phase slips is called the phase diffusion regime [34][35][36][37][38][39][40][41] , and the corresponding parameter space plots, where this regime is realised, are given in 38,41 .
To confirm the findings above, we have calculated the junction lifetime through numerical simulations of the junction switching in the frame of the resistively-capacitively shunted junction (RCSJ) model with thermal noise 50 , as well as by the Kramers' formula for the thermal escape of a particle from a well 51 : Here k is the Boltzmann constant and Φ 0 = h/(2e) is the flux quantum, with the Planck constant h and the electron charge e. The numerical simulations are performed with the following parameters: critical current I c = 23 nA, capacitance C = 46 fF, normal state resistance R N = 4200 Ω, temperature T = 300 mK. The calculations with the Kramers' theory (1) are straightforward, while by the numerical simulations we reconstruct two dependencies, φ(t) and dφ(t)/dt. The latter, due to the second Josephson law, is the voltage V(t) across the junction. The phase dynamics is computed by numerical solution of the second-order differential Langevin equation with noise in the frame of the RCSJ (a mechanical pendulum) model using the Heun scheme 52,53 . Here, we can compute the lifetime either as the escape time of the potential barrier by the phase variable, or by the first passage time of the value 0.2 mV by the voltage variable. In the latter case, we automatically consider multiple phase slips before the switching to the finite voltage state.
The obtained results are presented in Fig. 2a. Considering the approximate nature of the formula (1) in the considered parameter range of moderate damping 50 , the obtained agreement between the Kramers' theory (green diamonds) and the first phase-slip event in the RCSJ model (violet dashed line) is rather good. Moreover, the experimental data (blue dots with error bars) are perfectly fitted by the voltage simulations (red dots), considering multiple phase slips before approaching a finite voltage, see the upper curves in Fig. 2a. Since we are restricted by switching times above 0.01 s due to experimental constraints and simulation times grow exponentially, it is impossible to cover the full experimental parameter range; that is why we show only partial overlap between curves. Interestingly, the voltage simulation was performed with the only slight fitting of capacitance due to possible variation of oxidation barrier thickness, but with other fixed parameters defined independently above. Thus, the phase diffusion regime is fully confirmed and indeed, the junction lifetime is drastically increased; see the voltage traces in the inset of Fig. 2a.
At low temperatures, the SIS junction represents the phase qubit with discrete energy levels in a potential well. The main escape mechanism here is the tunnelling under the potential barrier. It is intuitively clear that the change of the escape origin (tunnelling through the barrier instead of thermal escape) would not seriously modify the phase dynamics 38 , expecting the quantum phase diffusion regime. Since in the case of tunnelling under the barrier, the potential energy will be smaller than in the case of thermal escape above the barrier, the quantum phase diffusion will develop in a more broad parameter range than the classical one. In this case, we cannot any more recourse to the numerical modelling of the classical Langevin equation, since it only works for higher temperatures well above the quantum crossover. However, Fig. 1 The scheme of the measurement setup with thermal mounting and two filtering stages. The sample was mounted in an RF-tight box with superconducting shielding on the 10 mK cryostat plate. The DC bias wires were filtered with feed-through capacitors at room temperature and RC filters at the 10 mK plate, minimising the effect of low-frequency noise. A pulse of 10 GHz signal after a two-stage attenuation at room temperature was further attenuated by 30 dB by the twisted pair and sent to the loop antenna, weakly coupled to the sample. The coupling efficiency and the power received by the SPC was measured in situ. In such a setup, a classical pulse (violet rectangle) was attenuated until it was absorbed by the SPC as separate photons during the response time δt, see the text.
the comparison between the experiment and the quantum tunnelling theory 32 , indeed demonstrates our expectations: the calculated lifetime is again much smaller than the experimental data, see Fig. 2b.
The advantage of the presented setup for single-photon counting is the ability to calibrate the absorbed power in situ, using the photon-assisted tunnelling (PAT) steps at the currentvoltage characteristics of the SIS junction 54,55 , see Methods. The 70 fW power level is easily seen with the help of PAT steps and further attenuation from this power level by 10-20 dB can be done by a remotely controlled room temperature attenuator without any other change in the route of the microwave signal, see Fig. 1.
To characterise the performance of the photon counter, a photon source is needed. However, there are no available microwave single photon sources on demand. Therefore, in the microwave range, one can only use classical sources with subsequent strong signal attenuation. Since the monochromatic beam can be presented as a flux of photons obeying Poissonian statistics 56 , the probability of appearing of single photons, photon pairs, triples and so on, in the source output scales with the beam power. For an incident photon flux with the power P 0 at the frequency ν, the average number of photons N within the SPC response time δt can be expressed as N = P 0 × δt/(hν). If N is small, the Poisson distribution (see Eq. (6) in Methods) dictates that the switching probability scales as log p sw % j log N. The latter sets the slope of the switching probability p sw curve in logarithmic axes by the number of photons j. The same methodology has been used before for infrared single-photon detectors 57 .
In the experiment, we used the pulses of the duration t pulse = 50 ms, which greatly exceeds the SPC response time δt. The SPC was first initialised in the zero-voltage state and after that was biased to a working point by quasi-adiabatic ramping up current 58 during 50 ms. If, by incoming pulse, the SPC switched into the resistive state, the run was counted as 1, otherwise, as 0. The obtained results are presented in Fig. 3. Each experimental curve here is the average of 10 3 −10 4 runs.
The obtained experimental data were fitted by making use of Eqs. 6-8 (see Methods), which can be reduced to where q[0] is the probability of the false detection without a photon, and q [1], q [2], q [3], etc., are the detection efficiencies of 1, 2, 3, etc., photons. All the fitting curves (see Fig. 3) are obtained using the same values of δt = 0.2 ns and the number of attempts M = 150, and only differ by the number of the detected photons j and the detection efficiency q. These curves are plotted along the "mean photon number N" axis, which is determined by the value of the controlled attenuator voltage. We also give for reference the power axis in Watts under the N axis, which is obtained by calculating the absorbed power at PAT steps on the resistive branch of the IV curve, see Methods. Even though the PAT calibration method is a bit rough in terms of absolute input power, we don't need this to be exactly correct, since the number of photons, initiating the switching, come from fitting to the slopes of the switching curves. Thus, the relative power along the x-axis of Fig. 3 is more important than small offset in the absolute power.
The obtained results are summarised in Table 1. Here the detection efficiencies q[i] versus bias currents, corresponding to the detection of 1, 2, 3, etc. photons, are shown (see the red curves in Fig. 3).

DISCUSSION
We have developed a photon counter prototype, which can detect an energy equivalent of 5 photons at frequency of 10 GHz with a dark count time above 10 s and the efficiency close to unity (see the curve for 12 nA in Fig. 3 and compare with the lifetime in Fig. 2b). While smaller number of photons is also detected, the efficiency and dark count time have to be improved. Such a threshold detector exploits a current biased Al SIS junction coupled to an incident photon flux. Due to the phase diffusion regime, the achieved dark count times are much higher than previously thought, enhancing the prospects for using an SIS junction as a single-photon counter. Looking forward, the counter performance can be further improved by altering junction parameters, aiming to an optimisation of a relationship between the counter sensitivity and its dark count rate.

Sample fabrication
The samples studied in this work were fabricated at Chalmers University of Technology. All layers, except for the SIS junctions themselves, were formed by the method of lift-off lithography, performed using a laserwriter, and the subsequent deposition of thin films by an electron beam.
For the fabrication of aluminium tunnel junction, an electronic lithograph and a shadow evaporation technique at three different angles were used, which made it possible to deposit tunnel junctions without breaking the vacuum.

Methodology of single-photon detection
A photon absorption in an SPC is a complex process involving both active and reactive part of the junction impedance 44 . If the Josephson plasma frequency matches the incoming signal frequency, the signal causes the oscillations in the reactive branches with the highest amplitude. The active part of the impedance equals to R qp 50 , while the parallel resonant RLC circuit increases the current through the Josephson inductance L J by a factor of The photon energy hν is divided between the energy E s stored in the tank circuit, and the energy E d dissipated in the resistor R qp , so that E s + E d = hν.
The amplitude I J of current oscillations induced in a Josephson inductance is related to E s by Since E s /E d = Q/(2π), we find an expression for I J by combining (3) and (4): If the current reaches a value close to the critical current, the probability of the switching into the resistive state increases significantly.
We analyse the switching statistics under the assumption that several photons can participate in an SPC switching, but only if they hit the SPC synchronously and simultaneously. If the photons are absorbed sequentially one after another, the energy of the first photons is dissipated over the times set by the Q-factor of the resonant circuit of the SPC, which is about 10 in our case.
Consider the photon flow with the power P 0 and the frequency ν, accepted by the SPC. The average number of photons within a time interval δt is N = P 0 × δt/(hν). Then the probability p[j] to find j photons during δt is given by Poisson distribution 56 : The power P 0 can be measured in situ by PAT steps, see below, that links the experimental probabilities with Eq. (6).
For the counter, the Poissonian statistics of the incoming signal means that the switching may happen within δt due to the arrival of one photon, and also of 2, 3, 4 and so on, photons at a given bias current. The probability of the SPC switching during δt at a given bias current can be fully described by an array of the detection efficiencies q[j], where q[0] is the probability of false switching if no photons are absorbed (dark count), q [1] is the probability to switch if one photon is absorbed, and so on. For an ideal single-photon detector, the array q[j] will look like [0, 1, 1, . . . , 1], which means absence of dark counts and 100% probability of switching due to one and consequently more than one photon.
The expression for the total switching probability p δt during δt in an experiment, where we do not know the exact number of the absorbed photons, has the form of a sum over the products of the probability p[j] to   Fig. 3 The switching probability of the SPC (blue dots with error bars computed as standard errors) versus the normalised power of the signal for different bias currents I from 16 to 12 nA, corresponding to detection of 1-5 photons, respectively. The fitting red curves are plotted against the middle x-axis "mean photon number". The lower x-axis is the power in Watts, obtained from the PAT steps (the correspondence between the normalised power and the power in Watts is given in Methods). A green dashed line is given to guide eyes, indicating the beginning of the dark count floors.
absorb j photons and the probability q[j] to switch due to j photons: p½j q½j: Here the probability q[j] is the detection efficiency. Let us discuss the meaning of the SPC switching time δt. We assume that δt is determined as the time needed to dissipate the energy of absorbed photons. For a threshold detector, such as an SPC, it is determined from the non-switching event (the photon is absorbed, but the jump to the resistive state does not happen) using the dynamics of the Josephson phase. If the next photon or a group of photons comes later than δt, the SPC already forgets that there has been a previous photon (or several photons). In this case, we can say that the experiment starts from the beginning, and this will be the second attempt to switch. Now, if we open the shutter of the photon source for a longer time, there will be more attempts for the SPC to switch unless it happened, for example, during the previous steps. Thus, the switching probability p sw will be given by the following expression: where M is the number of "elementary" experiments (attempts). This expression is used to fit the experimental switching probabilities (the simplified expression (2) is given in the main text) with three fitting parameters: the time interval δt, the number of attempts M, and the detection efficiency q. It should be noted that all p sw curves can be fitted by the same values of δt and M, while their difference is described by the only detection efficiency q.

Experimental setup
The radio frequency (RF) signal from an external microwave synthesiser was attenuated using constant attenuators from 1 dB to 30 dB and voltage controlled room-temperature attenuator, calibrated with a commercial spectrum analyser. The measurements of the power passed through the variable attenuator versus voltage on the attenuator are shown in Fig. 4a, where the decrease in the power by >3 orders of magnitude is seen.
To study the dynamics of an Al/AlO x /Al SIS tunnel junction, we thermally anchored the sample to the mixing chamber of a He 3 /He 4 Oxford dilution refrigerator Triton 200. The reverse branch of the current-voltage characteristic of the investigated SIS junction with parameters I c = 23 nA, R N = 4.2 kΩ, C = 46 fF is presented in Fig. 4b.
For a high-frequency experiment, a microwave signal was fed into the cryostat via a phosphor bronze twisted-pair with attenuation of −15 dB per meter at 10 GHz and with a loop antenna at the end near the Josephson junction. The level of the RF signal passing through two metres of a twisted pair inside the cryostat was also measured as a function of the control voltage on the attenuator by a spectrum analyser. In total, the initial signal with a power of the order of mW was attenuated to a level of several nW passing through several constant attenuators and two metres of the twisted pair. In turn, the absorbed signal (described further) appears to be smaller by several orders of magnitude (~pW) due to poor coupling between the loop antenna and the SIS junction. For switching statistics measurements, this signal was reduced even further to a few femtowatts using a controlled attenuator.
The absorbed power was measured using the well-pronounced photon assisted tunnelling (PAT) steps at the IV curve, see Fig. 4b. Its dependence on the attenuator voltage follows the same curve as the calibration curve, see Fig. 4a. That's how we find the correspondence between relative attenuated power and power in Watts. The lower power, for which the PAT steps are not visible yet, can be probed only with the help of the switching histograms and by the decrease of the superconducting state lifetime.
For the switching probability measurements, to bias the junction at the initial state, the bias current of the junction was ramped up quasiadiabatically, with slow start and stop parts to decrease the disturbance of the quantum system 58 . The voltage was measured using a low noise roomtemperature differential amplifier AD745 and further digitised by a highspeed NI-DAQmx card.

Absorbed power
Even if we know the losses in each element of the microwave path including the antenna pattern, antenna-receiver matching, etc., it is still useful to establish an independent way to measure the power, absorbed by the SIS junction. The PAT steps 55 give us such an opportunity because the height of the steps is unambiguously related to the power absorbed by the quasiparticles of the superconducting electrodes of the SIS junction.
In the presence of a high-frequency signal with a frequency exceeding the scale of the gap smearing, the I-V characteristic (IVC) of the SIS junction takes on a characteristic stepped shape, see Fig. 4b. In this case, the IVC can be represented as the sum of the IVCs without a signal I 0 (V), shifted relative to each other along the voltage axis by an integer hν/e: Here J n (x) is the Bessel function whose argument x = eV rf /(hν) is proportional to the voltage V rf through the junction at the frequency of the external signal ν. Approximation of the experimental IVC by (9) with the adjustable parameter x makes it possible to determine V rf . Further, the absorbed power is expressed as follows: where σ is the conductivity of the SIS junction at the frequency ν. In Ref. 55 , it is proposed to determine σ through the difference of the shifted IVCs since the IVC changes so sharply with the signal that the usual expression dI/dV is not applicable: σ ¼ I 0 ðV þ hν=eÞ À I 0 ðV À hν=eÞ 2hν=e : Fig. 4 Determination of absorbed power from PAT steps. a The output power through the variable attenuator versus the attenuator voltage. b An example of obtaining a PAT step from the IVC using Eq. (9) and fitting an experimental I-V characteristics. c The power absorbed in the SIS junction for the attenuator voltages from 1.66 V to 3 V.
Let us find the absorbed power using the described method for the junction with the critical current of 23 nA. Figure 4b shows the reverse branch of the IVC without the external signal, with a fitted analytical function (red). A good agreement with the reverse branch is observed both below and above the gap. The only difference from the experiment is the absence of the heating effect with an increasing current. To estimate the power of the RF signal, we will analyse the IVCs with only the first PAT step. In this case, Eq. (9) can be restricted to the first three terms.
To find V rf , we use Eq. (9), where we substitute I 0 (V) with the found analytical function. The argument x of the Bessel function is the only fitting parameter. Figure 4b shows the terms from the sum of Eq. (9) without coefficients (dashed red) and the result of the sum in the form of an IVC with a PAT step (dashed blue). Furthermore, Fig. 4b shows real IVCs without (black dots) and with a PAT step (grey dots). As it can be seen, the experimental data and the fitting results are in a good agreement.
The absorbed power is obtained from Eq. (10). Figure 4c shows the power dependence versus the junction voltage. The maximal power is absorbed at the step. When the attenuator voltage is 3 V, the maximum absorbed power is 6.8 × 10 −14 W. This means that at least this power always comes to the junction, but whether it is absorbed or not depends on the voltage at the junction.
To estimate the power absorbed when SIS is biased at supercurrent, we use the power level at voltages above the gap, since it is believed that at zero-voltage across the SIS junction, its quasiparticle resistance is equal to R N 59 .

DATA AVAILABILITY
The data that support the findings of this study are available from A.V. Gordeeva upon reasonable request.