Wide-band Parametric Amplifier Readout and Resolution of Optical Microwave Kinetic Inductance Detectors

TL;DR

This paper demonstrates that a kinetic inductance traveling wave parametric amplifier can be used for MKID readout, achieving near-quantum-limited noise performance and significantly improving photon energy resolution.

Contribution

It introduces a novel application of a kinetic inductance traveling wave parametric amplifier for MKID readout, reaching the quantum limit and enhancing resolution.

Findings

Amplifier achieves ~2.1 quanta of system noise.

Resolving power increased from ~6.7 to 9.3 with the amplifier.

Potential for resolving power up to 25 with detector redesign.

Abstract

The energy resolution of a single photon counting Microwave Kinetic Inductance Detector (MKID) can be degraded by noise coming from the primary low temperature amplifier in the detector's readout system. Until recently, quantum limited amplifiers have been incompatible with these detectors due to dynamic range, power, and bandwidth constraints. However, we show that a kinetic inductance based traveling wave parametric amplifier can be used for this application and reaches the quantum limit. The total system noise for this readout scheme was equal to ~2.1 in units of quanta. For incident photons in the 800 to 1300 nm range, the amplifier increased the average resolving power of the detector from ~6.7 to 9.3 at which point the resolution becomes limited by noise on the pulse height of the signal. Noise measurements suggest that a resolving power of up to 25 is possible if redesigned detectors can remove this additional noise source.

Equations

$$\frac{1}{R_{\text{measured}}^2}=\frac{1}{R_0^2}+\frac{1}{R_{\text{expected}}^2},$$

$$S(\nu )=4hf{Z}_0\left(\frac{1}{e^{\text{nicefrac}hf{k}_{B}T}-1}+\frac{1}{2}\right),$$

$$S_0(\nu )=$$

$$S_1(\nu )=$$

$$S_2(\nu )=$$

$$S_3(\nu )=$$

$$G_{\text{H}}(\nu )=$$

$$S_{\text{H}}(\nu )=$$

$$S_{\text{I}}(\nu )=$$

$$S_{\text{P}}(\nu )=$$

$$A_{\text{I}}(\nu )=$$

$$A_{\text{P}}(\nu )=$$

$$A_{\text{H}}(\nu )=$$

$$A_{\text{sys}}(\nu )\equiv$$

$$\mathcal{L}\left(\bar{X}|X(\nu )\right)=\text{StudentT}\left(\frac{\bar{X}-E_{\nu }[X(\nu )]}{\sqrt{\frac{{\mathrm{Var}}_{\nu }[X(\nu )]}{N_{\nu }}}},N_{\nu }-1\right),$$

$$f(x)=\{\begin{array}{ll}\frac{1}{\beta -1}{e}^{-\frac{x-1}{\beta -1}}& x\ge 1\\ 0& x<1\end{array}$$

$$S_0(\nu )=$$

$$S_1(\nu )=$$

$$S_2(\nu )=$$

$$S_3(\nu )=$$

$$S_{21}(f)=\left(g_0+g_{1}f+g_2{f}^2\right)e^{i{\varphi }_0-2\pi i\tau f}\left(1-\frac{Q{\hat{Q}}_c^{-1}}{1+2iQx}\right)$$

$$Qx=Q\frac{f-f_r}{f_r}+\frac{a}{1+4Q^2{x}^2}$$

$$\mathbf{m}(t;E,t_0)=\left(\begin{array}{c}{A}_{\theta }(E)s_{\theta }(t-t_0,E)\\ A_d(E)s_d(t-t_0,E)\end{array}\right).$$

$${\chi }^2={\int }_{-\infty }^{\infty }df(\tilde{\mathbf{d}}-\tilde{\mathbf{m}}{)}^{†}{\mathbf{S}}^{-1}(\tilde{\mathbf{d}}-\tilde{\mathbf{m}}),$$

$${\sigma }_{\alpha }^2=$$

$$=$$

$${\delta }_{\theta }$$

$${\delta }_d$$

$${\sigma }_{{\delta }_{\theta }}^2={\sigma }_{{\delta }_{\text{TLS}}}^2+{\sigma }_{{\delta }_{\text{sys}}}^2.$$

$$\mathbf{d}(t;E)=\left(\begin{array}{c}{A}_{\theta }(E)s_{\theta }(t-{\delta }_{t_0},E)+{\delta }_{\theta }\\ A_d(E)s_d(t-{\delta }_{t_0},E)+{\delta }_d\end{array}\right)+\mathbf{n}(t).$$

$$\mathbf{d}(t;E)=\left(\begin{array}{c}\begin{array}{rlr}\theta (f(& A_{\theta }(E)s_{\theta }(t-{\delta }_{t_0},E))+f({\delta }_{\text{TLS}}))& +{\delta }_{\text{sys}}\\ & A_d(E)s_d(t-{\delta }_{t_0},E)& +{\delta }_{\text{sys}}\end{array}\end{array}\right)+\mathbf{n}(t).$$

$$R\equiv\frac{E}{\Delta E}=\frac{$1$}{$2.355$}\sqrt{\frac{E}{(F+J)\epsilon}},$$

$$J=\eta_{t}\frac{\Omega_{D}}{\epsilon}\frac{l_{ph}}{d}\frac{$12$($1$+\lambda)}{$11$($1$+\lambda)+$3$}\>g_{1}\!\!\left(\frac{\Omega_{D}}{\Omega_{1}}\right).$$

$$l_{ph}=\frac{\hbar v_{s}}{$8$\pi N_{0}\left<\alpha^{2}\right>_{\text{avg}}\Omega_{D}},$$

$$\Omega_{1}=\Omega_{D}\sqrt{\frac{$2$}{$3$}($1$+\lambda)N_{0}\Omega_{D}},$$

$$\lambda=\frac{$1.04$+\mu^{*}\ln\left(\nicefrac{{\Omega_{D}}}{{$1.45$k_{B}T_{c}}}\right)}{(1-$0.62$\mu^{*})\ln\left(\nicefrac{{\Omega_{D}}}{{$1.45$k_{B}T_{c}}}\right)-$1.04$}.$$

Full text

Optical MKIDs Szypryt et al. (2017) are superconducting, single photon counting, energy resolving sensors which are sensitive to radiation in the ultraviolet to near infrared range. Advantages over semiconductor detectors in this wavelength band include the absence of false counts (read noise, dark current, and cosmic rays), intrinsic spectral resolution, high speed, and radiation hardness. Other superconducting detectors have shown promise at these wavelengths, Martin et al. (2006); Burney et al. (2006) but they are difficult to chain together into large arrays. Optical MKIDs, however, are naturally frequency domain multiplexed, which has enabled full-scale instruments at the Palomar observatory Mazin et al. (2013); Meeker et al. (2018) and Subaru telescope. Walter et al. (2018) In the future, these detectors will be included in a balloon borne mission. Cook et al. (2015)

Photon counting MKIDs operate differently than MKIDs designed for longer wavelength detection in the bolometric regime. Instead of measuring a constant flux of photons, they record individual photon events similarly to an X-ray calorimeter. To achieve a measurable detector response for a single photon event they tend to be smaller and able to handle less signal power than their longer wavelength bolometric counterparts. In these conditions, amplifier noise can be comparable in magnitude to the detector phase noise that originates from microscopic two-level system (TLS) states on the surface or between layers of the device. Gao et al. (2008) The TLS noise can be mitigated through careful sample preparation, de Graaf et al. (2018) fabrication, Megrant et al. (2012); Calusine et al. (2018) and device design Noroozian et al. (2009); Mignot et al. (2018); Vissers et al. (2012) while the effect of amplifier noise can be addressed by designing detectors that can handle higher signal powers. Zobrist et al. (2018); Beldi et al. (2019) These routes are actively pursued, but, for optical MKIDs, improving the main readout amplifier’s noise floor offers an additional path to lowering the total system noise.

Quantum mechanics imposes an uncertainty relationship between the two quadratures of an electromagnetic signal. Caves (1982) This relationship results in a lower limit to the noise that a high gain, phase-preserving, linear amplifier adds to its input signal, equal to that of the electromagnetic zero-point fluctuations $\left(A=\nicefrac{{1}}{{2}}\right)$. To readout an optical MKID array, the primary amplifier must have a moderately high saturation power, significant dynamic range, and a large bandwidth. High electron mobility transistor (HEMT) amplifiers are often used to satisfy these requirements. Wadefalk et al. (2003); Schleeh et al. (2013) State of the art commercial HEMT amplifiers operating over $48\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ at $5\text{\,}\mathrm{K}$ typically reach added noise numbers, in units of quanta, as low as $A=$8.0$$ ($2.3\text{\,}\mathrm{K}$ noise temperature at $6\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$). Schleeh et al. (2012); LNF (2019) A quantum limited amplifier that could operate in these conditions would reduce the added amplifier noise by a factor of $16$.

Parametric amplifiers have been shown to perform at or near to the quantum limit and have been used in several experiments with superconducting resonators. A $\sim 1\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$, bandwidth format was used to investigate a superconducting resonator’s noise properties in the dissipation quadrature. Gao et al. (2011) However, only one resonator could be measured at a time, and a carrier suppression tone was needed, which made the amplifier incompatible with single photon measurements. At low signal powers, larger bandwidth formats have been able to read out up to $20$ superconducting qubits. Macklin et al. (2015) This readout system, too, does not have enough dynamic range for an optical MKID array as the amplifier’s saturation power is on the order of the required signal tone’s power.

Traveling wave parametric amplifiers (TWPA) based on a superconductor’s nonlinear kinetic inductance were designed to handle this wide-band, high power and dynamic range case. Zmuidzinas and Day (2014); Chaudhuri et al. (2017); Eom et al. (2012) In this paper we show how one of these parametric amplifiers can be integrated with a large MKID array and demonstrate its ability to measure single photon events with quantum limited amplifier noise.

The full readout is divided between three temperature stages with most of the large electronics at room temperature. The rest of the components are cooled to either $3.3\text{\,}\mathrm{K}$ or $100\text{\,}\mathrm{m}\mathrm{K}$ in a Leiden Cryogenics CF-200 dilution refrigerator with the MKID array and parametric amplifier inside of a Amumetal 4K magnetic shield developed by Amuneal. Outside of the cryostat we employ a homodyne readout system with the signals digitized at a sample rate of $2\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$ after being low pass filtered at $1\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$ to prevent aliasing. Figure 1 shows a schematic of the setup used for this experiment.

The detector that we tested was designed for the MEC instrument at the Subaru telescope on Mauna Kea in Hawaii. Walter et al. (2018) It has ten niobium coplanar waveguide feedlines and a $20,440$ pixel platinum silicide MKID array multiplexed over $48\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$. Szypryt et al. (2017) Each pixel is a lumped element resonator, capacitively coupled to one of the feedlines. A microscope image of a pixel is shown in the inset of figure 1.

For the photon measurements, we used five single mode laser diodes with wavelengths from $8081310\text{\,}\mathrm{n}\mathrm{m}$. The $1120\text{\,}\mathrm{n}\mathrm{m}$ diode was purchased from Eagle Yard, and the others were purchased from Thor Labs. The diode light was coupled into a fiber using an integrating sphere and sent into the fridge. Near the detector, the fiber output was collimated and directed toward a microlens array on top of the device. The microlens array focused the light onto the inductor of each pixel and was purchased from Advanced Microoptic Systems. It is made from 1 mm thick STIH53 glass with $140\text{\times}146$ lenslets at a $150\text{\,}\mathrm{\SIUnitSymbolMicro}\mathrm{m}$ pitch.

We chose a resonator from the array at $5.67446\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ with a coupling and internal quality factors of $15,100$ and $190,000$ respectively for this test. All data taken off resonance was captured at a frequency of $5.675\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$. Near the resonance frequency the parametric amplifier had a constant gain of $13.7\text{\,}\mathrm{d}\mathrm{B}$, and the resonator was driven almost to saturation at about $-106\text{\,}\mathrm{d}\mathrm{B}\mathrm{m}$. When the parametric amplifier is turned on, we see changes in the coupling quality factor of some resonators on the order of $\sim 10\text{\,}\mathrm{\char 37\relax}$. More isolation between the parametric amplifier and MKID array could remove this effect, but since this difference is smaller than the variation intrinsic to the design, modifications were unnecessary for this test.

The parametric amplifier itself is a wide-band, traveling-wave, kinetic inductance amplifier of the general type first described by Eom et al., Eom et al. (2012) but is an updated version and differs in several respects. The coplanar wave guide transmission line structure uses finer features (e.g. the center line width and gaps are $320\text{\,}\mathrm{n}\mathrm{m}$). In order to lower the characteristic impedance to $50\text{\,}\mathrm{\SIUnitSymbolOhm}$, added capacitance is provided with a interdigitated structure in a similar manner to the amplifier described by Chaudhuri et al. Chaudhuri et al. (2017) Additionally, the amplifier is operated in a three-wave mixing mode by applying a DC bias current using the technique shown by Vissers et al. Vissers et al. (2016) Further details about this amplifier are in preparation and will be published separately.

Operation of the amplifier using three-wave mixing requires the use of both a pump tone and a DC bias current. The pump tone (from synthesizer B in figure 1) at $14.765\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ is attenuated by 20 dB at $3.3\text{\,}\mathrm{K}$ and filtered using a $14.517\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ Marki FB-1575 bandpass filter to ensure that phase noise from the pump generator does not leak into the signal band. The pump tone power at the input of the amplifier is about $-23\text{\,}\mathrm{d}\mathrm{B}\mathrm{m}$. The tone is then combined with a DC current of approximately $0.7\text{\,}\mathrm{m}\mathrm{A}$ to produce about $15\text{\,}\mathrm{d}\mathrm{B}$ of gain from $510\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$. The two Anritsu K250 bias tees and the two Marki DPX-1114 diplexers isolate the DC current and pump tone from the other components, respectively. The diplexers have over $60\text{\,}\mathrm{d}\mathrm{B}$ of isolation at the pump frequency, and the pump tone is terminated on the cryostat ground at $3.3\text{\,}\mathrm{K}$ where there is more cooling power. The parametric amplifier is then protected from reflections off of the warmer components with a $311\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ Pamtech CTH1365K10 isolator.

\ctable

[ caption=Maximum likelihood estimates for the input, parametric amplifier, HEMT, and total system noise in units of quanta at $5.675\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ (from left to right). $1\sigma$ statistical errors then systematic errors are reported to the right of the estimates. The HEMT noise level corresponds to a noise temperature of about $5.3\text{\,}\mathrm{K}$ and is within the manufacturer’s specifications for the amplifier. The effective noise temperature of the combined parametric and HEMT amplifiers is $0.43\text{\,}\mathrm{K}$. The measurement used to determine these numbers and errors is detailed in appendices A and B., label=tab:noise, width=pos=t, botcap ] C C C C

$A_{\text{I}}$ $A_{\text{P}}$ $A_{\text{H}}$ $A_{\text{sys}}$

0.71\stackanchor+0.03-0.03 \stackanchor+0.2-0.01 0.59\stackanchor+0.02-0.02 \stackanchor+0.2-0.07 19.5\stackanchor+0.1-0.1 \stackanchor+2-2 2.13\stackanchor+0.01-0.01 \stackanchor+0.8-0.1

Input signal saturation occurs when the amplified signal power at the output of the amplifier is approximately $15\text{\,}\mathrm{d}\mathrm{B}$ below the pump power. At that point, the pump amplitude becomes depleted and the operating point of the amplifier is altered. For $15\text{\,}\mathrm{d}\mathrm{B}$ gain, this results in an input signal saturation power of around $-53\text{\,}\mathrm{d}\mathrm{B}\mathrm{m}$. The large saturation power obviates the need for any carrier suppression tones. For this array, even if all of the signal tone power from a simultaneous measurement of a full feedline of resonators reached the amplifier input, it would take $-86\text{\,}\mathrm{d}\mathrm{B}\mathrm{m}$ of power per tone before reaching saturation, well above the typical operating point of an optical MKID. In practice, most of the signal tone power is reflected by the MKIDs before reaching the parametric amplifier, so higher signal powers might be usable.

The second low temperature amplification stage is a $412\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$ CIT412 HEMT amplifier from Cosmic Microwave Technology with about $32\text{\,}\mathrm{d}\mathrm{B}$ of gain and thermalized to the $3.3\text{\,}\mathrm{K}$ temperature stage with a copper heat strap. A low pass filter is included before the amplifier input to ensure that the pump tone does not leak past the diplexers and saturate the HEMT. This amplifier is required because the parametric amplifier does not have enough gain to boost the signal above a standard room temperature amplifier’s noise floor.

Since the transmission through the parametric amplifier is near unity when unpowered, we can perform an extended Y-factor measurement of the system. This method allows us to accurately determine the noise components and assign them to different elements in the setup. Table LABEL:tab:noise shows this breakdown for a signal tone slightly detuned from the MKID resonance frequency. The details of the procedure used to collect this data are laid out in appendices A and B.

While the parametric amplifier is performing near optimally, some aspects of the system may be improved. The observed system noise is about two times the achievable limit. Almost $75\text{\,}\mathrm{\char 37\relax}$ of this excess comes from the HEMT amplifier. This effect could be mitigated by either increasing the parametric amplifier’s gain or using a lower noise secondary amplifier. Additionally, the input noise is about $0.2\text{\,}\mathrm{q}\mathrm{u}\mathrm{a}\mathrm{n}\mathrm{t}\mathrm{a}$ larger than one would expect from a $100\text{\,}\mathrm{m}\mathrm{K}$ termination. It could be brought closer to its lower bound by increasing the lowest temperature attenuation to block more of the thermal noise from the $3.3\text{\,}\mathrm{K}$ and $300\text{\,}\mathrm{K}$ stages.

Figure 2 shows the noise in the bandwidth of our signal with the signal tone set to the MKID resonance frequency. This state is used for all of the single photon measurements. The flat component of the noise decreased by a factor of 9.6 when the parametric amplifier was powered, with the noise level in the dissipation response quadrature matching the off resonance noise floor for the system. In the phase response quadrature, there is significant low frequency noise which we attribute to TLS noise in the detector.

When the detector is illuminated, we measure a two dimensional pulse record for each photon event. From a record, we can calculate a maximum likelihood estimate of the photon’s energy and arrival time. More details of this estimation are given in appendices D and E. Since the pulse decay time is about $35\text{\,}\mathrm{\SIUnitSymbolMicro}\mathrm{s}$, care was taken to use low count rates ($<200\text{\,}\mathrm{H}\mathrm{z}$) and to exclude photons which arrived within $2\text{\,}\mathrm{m}\mathrm{s}$ of each other to ensure clean event records.

The resolving power of the detector, $\nicefrac{{E}}{{\Delta E}}$, is determined by measuring many photon events of known energy and evaluating the full-width half-max (FWHM) of the resulting fitted-energy distribution. Figure 3 shows this distribution for $808\text{\,}\mathrm{n}\mathrm{m}$ photons with the parametric amplifier’s pump tone and DC current both on and off, along with the kernel density estimation of the distribution that was used to compute the FWHM. We note that even though the dissipation noise is significantly smaller than the phase noise when the parametric amplifier is on, the smaller dissipation signal for this detector means that roughly $90\text{\,}\mathrm{\char 37\relax}$ of the expected resolving power can be achieved with the phase signal alone.

Using the parametric amplifier clearly improves the resolving power, but the results do not match the expected values for the noise level in the system. Table LABEL:tab:R shows this discrepancy over the five measured photon energies. Our noise model does not account for the small, non-stationary decrease in phase noise during a photon event, but we still see this discrepancy when estimating the energy using only the dissipation signal, which has stationary noise. Additionally, since the actual phase noise for a given event is smaller than a stationary noise model predicts, non-stationary noise is unlikely to explain the measured discrepancy in resolving power.

The resolving power is approximately the inverse of a scaled standard deviation, so we can consider additional noise on the pulse height, corresponding to a resolving power of $R_{0}$, using

[TABLE]

where $R_{\text{measured}}$ and $R_{\text{expected}}$ are the resolving powers from the columns in table LABEL:tab:R. Since $R_{\text{expected}}$ represents a lower bound, $R_{0}$ calculated from equation 1 can be interpreted as a upper bound. For this data, $R_{0}$ ranges from $7.311$ with the parametric amplifier off and from $9.610$ with the amplifier on. These results are in rough agreement with a constant detector-related energy uncertainty of $\sim 10$.

$R_{0}$ may be attributable to more than one source. In appendix H, we detail how the non-uniform current density in the resonators can contribute to an uncertainty in pulse height and account for the skewed distribution in figure 3. This effect is strongly dependent on the diffusion constant for platinum silicide, which is unknown for our films. However, the diffusion constants of similar films Baturina et al. (2005) suggest $R_{0}\sim$2040$$. Pulse shape variations may also contribute to $R_{0}$, but we expect them to be small since the quasiparticle distribution averages out the current non-uniformities at time scales on the order of the quasiparticle recombination rate. In our data, we do not see significant variations in the pulse shape for photons of the same energy with different pulse heights.

Interactions between the quasiparticles and phonons in the superconductor could also introduce fluctuations in the pulse height. The standard Fano-limit Fano (1947) is well above the measured resolving powers, but hot phonon escape from the superconductor to the substrate during the initial energy down-conversion may be an important factor. Kozorezov et al. (2007, 2008) For an $800\text{\,}\mathrm{n}\mathrm{m}$ photon absorbed in platinum silicide on a sapphire substrate, this effect could contribute anywhere from $R_{0}\sim$1330$$—with the uncertainty being dominated by the unknown electron-phonon interaction energy in platinum silicide. See appendix I for more details on this calculation.

Together, the worst case estimates of $R_{0}$ for the current non-uniformity and the hot phonon escape effects may account for the measured value of $R_{0}\sim 10$. More detailed measurements of the superconducting properties of platinum silicide are needed to make a more accurate estimation. Nevertheless, the lowered amplifier noise reveals that the detector’s design is the limiting factor in its resolving power. With the help of a parametric amplifier, future designs addressing this issue could reach resolving powers of up to $25$ without requiring lower TLS noise.

In summary, we have demonstrated a straightforward readout scheme for a full MKID array which incorporates a quantum limited amplifier and presented measurements showing that this amplifier is compatible with single photon detection. Detailed calculations of the noise in this system show that it significantly outperforms the standard HEMT amplifier readout scheme at the cost of a modest increase in complexity and exposes an additional noise source on the pulse height for these detectors.

N.Z. was supported throughout this work by a NASA Space Technology Research Fellowship. The MKID arrays used were developed under NASA grant NNX16AE98G. This research was carried out in part at the Jet Propulsion Laboratory and California Institute of Technology, under a contract with the National Aeronautics and Space Administration.

\ctable

[ caption=The measured increase in the resolving power is shown by comparing data with the parametric amplifier unpowered and powered for different photon energies. A lower bound on the expected resolving power of the estimate, from equation S12, is also specified (see appendix E). The detector response starts to saturate, as designed, at the highest energies, so the resolving power begins to decrease. Results are shown for data reduced using both the phase and dissipation, just the phase, and just the dissipation signals., label=tab:R, pos=b, width=botcap ] @Xd3.3d3.2d3.3d3.2d3.3d3.3

Energy [$\mathrm{e}\mathrm{V}$] Resolving Power $\left[\nicefrac{{E}}{{\Delta E}}\right]$

 Phase and Dissipation    Phase    Dissipation  

  

Measured

Expected

Measured

Expected

Measured

Expected

$1.53$ ($808\text{\,}\mathrm{n}\mathrm{m}$) $5.8$ , $8.9$ $9.5$ , $23$ $5.4$ , $8.8$ $8.7$ , $22$ $1.4$ , $5.9$ $1.8$ , $8.7$

$1.35$ ($920\text{\,}\mathrm{n}\mathrm{m}$) $7.4$ , $9.4$ $10$ , $24$ $6.6$ , $9.1$ $9.7$ , $23$ $1.4$ , $5.8$ $1.7$ , $7.8$

$1.27$ ($980\text{\,}\mathrm{n}\mathrm{m}$) $7.5$ , $9.6$ $11$ , $25$ $6.6$ , $9.3$ $10$ , $23$ $1.5$ , $6.8$ $1.6$ , $8.8$

$1.11$ ($1120\text{\,}\mathrm{n}\mathrm{m}$) $6.6$ , $9.6$ $9.3$ , $24$ $6.1$ , $9.2$ $8.8$ , $22$ $1.9$ , $6.9$ $1.8$ , $9.9$

$0.946$ ($1310\text{\,}\mathrm{n}\mathrm{m}$) $6.0$ , $9.2$ $8.7$ , $23$ $5.5$ , $8.9$ $8.3$ , $20$ $1.9$ , $6.1$ $1.9$ , $9.9$

Appendix A Noise Measurement

The noise figures presented in table LABEL:tab:noise are determined by the following measurement, which consists of taking four noise data sets using a microwave switch to change the input to the HEMT. From each noise data set, the power spectrum is computed for each quadrature and summed. Superconducting coaxes are used between the HEMT and the switch and between the switch and the terminations to ensure an accurate calibration.

1. $\bullet$ $S_{0}(\nu)$:

The input of the HEMT is terminated at {T_{\text{hot}}\sim3.3\text{,}\mathrm{K}}. 2. $\bullet$ $S_{1}(\nu)$:

The input of the HEMT is terminated at {T_{\text{cold}}\sim100\text{,}\mathrm{m}\mathrm{K}}. 3. $\bullet$ $S_{2}(\nu)$:

The HEMT is connected to the parametric amplifier and MKID with the pump tone and DC current off. 4. $\bullet$ $S_{3}(\nu)$:

The HEMT is connected to the parametric amplifier and MKID with the pump tone and DC current on.

The sum of the single sided power spectral densities for each quadrature around a frequency $f$, on a terminated transmission line with impedance $Z_{0}$, and at an equilibrium temperature $T$ is given by

[TABLE]

where the first component comes from the thermal Johnson noise and the second from the zero point fluctuations. Caves (1982) Equation S1 holds for $\nu\ll f$ and shows that the noise is independent of spectral frequency.

The noise added by the two different temperature terminations is calculated using equation S1 and is labeled $S_{\text{hot}}(\nu)$ and $S_{\text{cold}}(\nu)$ for $T_{\text{hot}}$ and $T_{\text{cold}}$ respectively. The four measurements can then be written in terms of their components, where $S_{\text{I}}(\nu)$ is the noise at the input of the parametric amplifier, $S_{\text{P}}(\nu)$ is the noise added by the parametric amplifier, and $S_{\text{H}}(\nu)$ is the noise added by the HEMT amplifier.

[TABLE]

$G_{\text{H}}(\nu)$ is the total gain of the system excluding the parametric amplifier and must be calculated. $G_{\text{P}}$ is the gain of the parametric amplifier and can be accurately determined by measuring the amplitude of a probe tone tuned off resonance with the parametric amplifier turned on and off. This system of equations can be solved for the unknown parameters.

[TABLE]

The added noise in units of photon quanta for the components of the system are given in equation S4, where we have defined $A_{\text{sys}}(\nu)$ to be the total noise of the system off resonance.

[TABLE]

Appendix B Noise Measurement Statistical and Systematic Errors

As written in appendix A, $G_{H}(\nu)$, $S_{\text{H}}(\nu)$, $A_{\text{H}}(\nu)$, $S_{\text{P}}(\nu)$, $A_{\text{P}}(\nu)$, $S_{\text{I}}(\nu)$, $A_{\text{I}}(\nu)$, $A_{\text{sys}}(\nu)$ are all functions of $\nu$. However, they should be frequency independent. We can look for this characteristic to check for any unexpected properties of the data. It also allows us to easily determine the statistical uncertainties from the marginal likelihood distribution for the mean.

[TABLE]

where $X(\nu)$ is one of the above parameters, $\operatorname{E}_{\nu}[\cdot]$ is the mean of the values over frequency, $\operatorname{Var}_{\nu}[\cdot]$ is the variance of the values over frequency, and $N_{\nu}$ is the number of frequency bins. For the amplifier added noise numbers, an additional prior can be included for only allowing physical values, $\bar{X}\geq\nicefrac{{1}}{{2}}$.

We also account for any systematic errors which affect the noise measurement. Error terms are estimated and a Monte Carlo simulation is done to compute a posterior distribution for the noise numbers.

Losses at the $100\text{\,}\mathrm{m}\mathrm{K}$ stage are small but they may exist. We introduce two loss parameters. $l_{\text{P}}$ represents losses between the parametric amplifier and switch, and $l_{\text{H}}$ represents losses between the switch and HEMT. We expect these losses to be not much greater than $1\text{\,}\mathrm{d}\mathrm{B}$, so we draw random samples from the following exponential distribution:

[TABLE]

with $\beta=$1\text{,}\mathrm{d}\mathrm{B}$\approx 1.26$.

We also consider thermometer calibration errors which we do not expect to be more than $10\text{\,}\mathrm{\char 37\relax}$. The temperatures for $S_{\text{cold}}(\nu)$ and $S_{\text{hot}}(\nu)$ are drawn from the normal distribution $\mathcal{N}\left(T,\nicefrac{{T}}{{10}}\right)$ for their respective temperatures, where the second value is the standard deviation.

The resulting equations representing our measurement are modified from equation S2 to

[TABLE]

$S^{\prime}_{\text{hot}}(\nu)$ and $S^{\prime}_{\text{cold}}(\nu)$ are from the calibration terminations if their temperatures were different than our expectation. The solution to equation S7 for many different samples of the loss and temperature distributions give the second set of confidence intervals presented in table LABEL:tab:noise.

Appendix C Resonance Fitting

Each MKID is a notch resonator which has near perfect transmission far from the resonance frequency and a minimum in transmission near the resonance frequency. These resonances can be well modeled by an asymmetric Lorentzian function.

[TABLE]

The $g_{0}$, $g_{1}$, $g_{2}$, $\phi_{0}$, and $\tau$ terms represent a quadratic gain background, constant phase offset and cable delay for the system. The loop asymmetry is set by the complex coupling quality factor, $\hat{Q}_{c}$ and is related to the internal and total quality factors of the resonator using $Q^{-1}=Q_{i}^{-1}+Q_{c}^{-1}=Q_{i}^{-1}+\text{Re}\left[\hat{Q}_{c}^{-1}\right]$. Khalil et al. (2012)

We use a nonlinear kinetic inductance model to describe the frequency detuning parameter, $x$, since the resonator described in this paper is driven close to its bifurcation power. Swenson et al. (2013)

[TABLE]

The parameter, $a$, is linearly related to the generator power and ranges from 0 to about 0.77 before the onset of bifurcation. $f_{r}$ is the resonance frequency when $a=$0. A fit to the phase and magnitude of equation [S8](#E8) is shown in figure [S1](#A3.F1), which gives the parameter values $Q_{i}=$190,000, $Q_{c}=$15,100$$, and $a=0.28$.

Appendix D Detector Response

When an optical MKID is illuminated at photon count rates much lower than the characteristic decay rate of the sensor, we can measure individual photon energies and arrival times. The quasi-particle density resulting from a single-photon absorption can be measured by monitoring either the resonator’s phase or dissipation quadratures. Each signal resembles that of a standard x-ray calorimeter response, and much of the formalism used to analyze this data is adapted from corresponding techniques generalized to this two-dimensional case. Moseley et al. (1988); Szymkowiak et al. (1993); Fowler et al. (2017)

In an ideal detector for a fixed energy, the detector response should not vary. We model this response as an energy dependent pulse shape, $s_{i}(t,E)$, and amplitude, $A_{i}(E)$, where $i\in\{\theta,d\}$ corresponds to either the phase or dissipation response. The amplitude and shape are separated because for our detector the pulse shape is only weakly dependent on energy. $A_{i}(E)$ is also a good indicator of detector linearity. In practice, it is not possible to know these functions over a continuous range of energies, but they can be approximated by averaging many pulse records together for a discrete set of known photon energies. As long as the number of pulse records is large, this averaging procedure gives us noise free estimates of the detector response.

To determine the pulse shapes and amplitudes at all energies, we estimate the functions by interpolating between the known energies (five are used in our case). For the amplitudes, we use second order splines with “not-a-knot” boundary conditions, fixing the zero point. The amplitude calibration is shown in the left of figure S2. The detector is linear up to the highest energies tested where saturation effects begin to flatten the response. At these energies, the detector has the highest resolving power. For the pulse shapes, we linearly interpolate over each frequency in their Fourier transforms. We find that this procedure preserves the pulse shape normalization to a very good approximation as long as enough known energies are used. The pulse decay time is defined as the integrated area under the normalized pulse shape and is shown as a function of energy on the right of figure S2. This definition is adopted for it’s simplicity and because the pulse shapes do not fit well to an exponential function. The pulse shapes for each energy are plotted in figure S3. Both of these interpolations introduce some systematic calibration error for photon measurements of significantly different energy than that of the calibration points. In this paper, however, the measured energies are the same as the calibration energies, so this effect is negligible.

The pulse shapes change very little except for the decay time in the phase quadrature which decreases with energy. This effect might be explained by a higher quasi-particle recombination rate at higher quasi-particle densities. However, the same tendency is not seen in the dissipation quadrature. The decay time in the dissipation quadrature is additionally much lower than in the phase quadrature. This behavior is a characteristic of the superconducting PtSi from which the resonator patterned. Similar results have been seen in TiN superconducting resonators Gao et al. (2012) and may be attributable to a small population of sub-gap states.

Appendix E Pulse Fitting

After the detector response amplitude and shape has been accurately measured, determining an absorbed photon’s energy and arrival time becomes an optimization problem. We model the pulse records with the following equation:

[TABLE]

Equation S10 represents a two-dimensional time dependent vector that we expect to measure for a given photon energy, $E$, and arrival time, $t_{0}$. The actual data contains noise, so we compute a maximum likelihood estimate for the photon energy and arrival time by minimizing a $\chi^{2}$ function with respect to $t_{0}$ and $E$. For a given data record, $\bm{d}(t)$,

[TABLE]

where a tilde denotes a Fourier transform with respect to time and $\bm{S}$ is the spectral noise covariance matrix determined by analyzing data with no photon events.

We compute $\chi^{2}$ in the frequency domain to simplify the computation, and because by not including the zero frequency in the integration, any small DC offset in the data is ignored. Correcting for the DC offset is important since our detectors have significant low-frequency noise below the measurement bandwidth. Because $\bm{m}(t;E,t_{0})$ is a nonlinear function of it’s parameters, no closed form solution can be used to find the best estimates of $E$ and $t_{0}$. Instead, we use a standard least-squares minimization routine to compute the results.

This formalism assumes following:

The photon records are completely isolated in time. 2. 2.

Each photon of the same energy creates the same detector response up to an additive noise. 3. 3.

The added noise is Gaussian and stationary (i.e. does not change as quasi-particles are broken in the detector).

When all of these assumptions are satisfied, the expected variance of each estimator can be calculated from equation S11 giving

[TABLE]

where $\alpha\in\{E,t_{0}\}$, $\hat{\alpha}$ represents the maximum likelihood value, and $\partial_{\alpha}$ represents a partial derivative with respect to the estimated parameter.

Equation S12 is a good estimate of the expected variance of the energy and arrival time model parameters only if the above assumptions are satisfied. We address the first assumption by ensuring that the photon absorption count rate is much slower than the characteristic decay rate of the detector and by excluding photons with nearby arrival times from the analysis. The second assumption is a detector property and may or may not be true depending on the detector design. The third assumption, however, is more difficult to assess. While the noise is very nearly Gaussian, the detector’s phase noise decreases as the resonance frequency moves off of the signal tone during a photon absorption event. Adequately accounting for this effect is still a debated problem and beyond the scope of this paper. Fixsen et al. (2014); Shank et al. (2014); Fixsen et al. (2004) However, since the phase noise is smaller at the peak of the pulse record, Monte Carlo simulations of our data, presented in appendix G, suggest that equation S12 is a lower bound on the true variance of the estimate as long as assumption 2 is satisfied. This intuitively makes sense, since the average noise in the pulse record is less than that used to make the estimation.

Appendix F Pulse Measurement Baseline

Noise below the measurement bandwidth contributes a small DC offset to each recorded photon event. By measuring the pre-trigger baseline of each pulse, we can see, in figure S4, how this baseline varies over the measurement duration. On the $\sim 10\text{\,}\mathrm{min}$ time scale, the baseline stays constant to within the measurement error.

The spread of measured baselines in the phase quadrature does not change much when the parametric amplifier is turned on because the low frequency noise is dominated by the two-level system noise in the detector. The two-level system noise corresponds to real resonance frequency shifts which could change the pulse height or shape through the nonlinear phase to frequency relationship shown in figure S1. We will find it useful to decompose the baseline into these two effects in order to simulate our detector. The baselines can be modeled as normal random variables:

[TABLE]

where the variances are related by the usual equation,

[TABLE]

Appendix G Pulse Fitting Monte Carlo Simulations

To be confident in the results of equation S12 and in our fitting routine, we simulate our data assuming a uniform pulse height for each energy and compute the resulting fitted energies. Our fake data follows the general form

[TABLE]

Our pulse trigger is accurate to the level of the sample spacing, $0.5\text{\,}\mathrm{\SIUnitSymbolMicro s}$, but there is some jitter. We model this effect with the random variable $\delta_{t_{0}}$ drawn from a normal distribution with a standard deviation of $0.25\text{\,}\mathrm{\SIUnitSymbolMicro s}$. A covariance matrix for the noise is estimated using the noise data taken directly before the measurement. $\bm{n}(t)$ is then drawn from the multivariate normal distribution with the corresponding covariance.

As discussed in appendix F, baseline drifts due to two-level systems could cause extra variations in the pulse height and shape. We cannot model this effect perfectly because the phase and dissipation trajectory of our pulses is not well understood. However, if we ignore the dissipation response, we can estimate the magnitude of this effect using $\theta(f)$ and $f(\theta)$ from figure S1. Equation S15 becomes

[TABLE]

As a final check, we can also account for the non-stationary noise in our simulated pulse records. We use the $2000$ pulse records with estimated energies closest to the laser energy and have the same arrival time as the model. The noise is recovered by subtracting off the model from the pulse records, and a covariance matrix is estimated similarly to before. This procedure produces noise that looks realistic, but we expect there to be problems near the rise of the pulse. Small variations in the arrival time between the model and data can cause significant errors in the subtraction in this region. For this reason and because the exact way to account for baseline effects is unknown, the Monte Carlo simulations are not used for the expected values in table LABEL:tab:R. The results are included only as a check that the non-stationary noise does not degrade the energy resolution.

At $808\text{\,}\mathrm{n}\mathrm{m}$ the device calibration is the most nonlinear, so in figure S5 we present the Monte Carlo results for this energy using the simulated data described above. We find good agreement between equation S12 and the simulations for the stationary noise case. Baseline drifts slightly degrade the resolving power, but the decrease is mostly negligible. For non-stationary noise, we recover a higher resolving power as expected, and two-level system baseline effects significantly reduce the simulated resolving power. These results support the claim that the resolving power given by equation S12 represents a lower bound on the expected resolving power for this data.

Appendix H Current Distribution

The data presented in table LABEL:tab:R show a discrepancy between the measured and expected resolving powers. However, the pulse fitting algorithm performs as expected on simulated data computed using the enumerated assumptions discussed in appendix E. The discrepancy between the real and simulated data, therefore, must be either in the noise model or response model for the detector. Here we investigate one possible explanation for the detector response not being fixed for a given photon energy.

The response of a MKID induced by a uniform distribution of quasiparticles differs from the response from a localized quasiparticle distribution. In the localized case, the detector response is proportional to the square of the current density multiplied by the quasiparticle distribution integrated over each position in the resonator. Gao (2008) In an ideal lumped element resonator, all of the current is in the inductor and none is in the capacitor. When designing these detectors, then, it is important to ensure a uniform current distribution in the inductor.

The current distribution was computed in a Sonnet simulation for the resonator used in this paper and is shown in figure S6. It takes the resonator about $1\text{\,}\mathrm{\SIUnitSymbolMicro}\mathrm{s}$ to respond to a photon absorption event, and in that time, the quasiparticles will have diffused an amount determined by the superconductor’s diffusion constant. This time scale roughly sets the amplitude of a photon pulse since at longer time scales the quasiparticles begin to recombine (see figure S3). In similar films the diffusion constant has been measured to be $\sim 28\text{\,}\mathrm{c}\mathrm{m}^{2}\mathrm{/}\mathrm{s}$, Baturina et al. (2005) so we use this range for our analysis.

To model the effect of the current non-uniformity, we randomly selected photon absorption locations on the inductor, solved the diffusion equation for each one, and computed the resulting response. Photon strikes in the capacitor have normalized responses $\ll 0.5$ and are excluded from the calculation to increase the throughput. The distribution of responses for several different diffusion constants was computed to estimate its effect on the resolving power. The results are shown in figure S7.

Appendix I Hot Phonon Escape

After a photon is absorbed by the superconducting PtSi in our detector, the resulting quasiparticle distributes its energy among more broken Cooper pairs and phonons. This energetic cascade ultimately results in an equilibrium distribution of quasiparticles on a roughly $1\text{\,}\mathrm{n}\mathrm{s}$ timescales. Kozorezov et al. (2007) Phonons created in the intermediary stages of the cascade can escape the superconductor to the substrate, effectively removing energy from the measurement. The statistics of this process offer another potential explanation for the resolving power discrepancy presented in table LABEL:tab:R.

The resolving power of our detector absorbing a photon with energy, $E$, is given by

[TABLE]

considering only Fano fluctuations and fluctuations in energy loss due to the escape of athermal phonons. $\epsilon$ is the mean energy required to generate one quasiparticle in the superconductor. Kozorezov et al. (2007) $F$ is the Fano factor, typically assumed to be about 0.4 for most materials. Fano (1947) The $J$ factor accounts for the hot phonon escape and can be calculated using [][.Note:equation~\ref{eq:J}containsanextrafactorof4whichwasmissingfromtheoriginalmanuscript.Thiswasdeterminedthroughprivatecommunicationwiththeauthors.]Kozorezov2008

[TABLE]

$\Omega_{D}$ and $\Omega_{1}$ are the Debye energy and and the lower boundary energy for the phonon controlled down-conversion. $l_{ph}$ is the mean free path for a phonon at the Debye frequency. $d$ is the film thickness. $\lambda$ is the effective electron-phonon coupling constant. $g_{1}(\cdot)$ is a dimensionless function defined by Kozorezov et al.

\ctable

[ caption=Parameters used in the evaluation of equation S17 are tabulated here. The top half corresponds to known or estimated values. Parameters computed from those given above are listed on the bottom half., label=tab:parameters, width=pos=t, botcap ] XXp2.8inl

Parameter Value Description Comment

F $0.4$ Fano factor typical value for superconductors Fano (1947)

$T_{c}$ $930\text{\,}\mathrm{m}\mathrm{K}$ superconducting transition temperature typical value for a MEC array Szypryt et al. (2017)

$N_{0}$ $1\text{\,}\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{t}\mathrm{e}\mathrm{/}\mathrm{e}\mathrm{V}\mathrm{/}\mathrm{i}\mathrm{o}\mathrm{n}$ density of states at the Fermi energy theoretical Koc, Deligöz, and Mamedov (2011)

$v_{s}$ $3600\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$ speed of sound in PtSi theoretical Koc, Deligöz, and Mamedov (2011)

$\Omega_{D}$ $31\text{\,}\mathrm{m}\mathrm{e}\mathrm{V}$ Debye energy theoretical Koc, Deligöz, and Mamedov (2011)

$\left<\alpha^{2}\right>_{\text{avg}}$ $1.5\text{\,}\mathrm{m}\mathrm{e}\mathrm{V}$ averaged squared electron-phonon interaction typical value for superconductors Kaplan et al. (1976)

$\mu^{*}$ $0.13$ Coulomb pseudo-potential typical value for superconductors McMillan (1968)

$\eta_{t}$ 1 phonon transmission coefficient assumed for PtSi-sapphire boundary

$d$ $50\text{\,}\mathrm{n}\mathrm{m}$ PtSi film thickness typical value for a MEC array Szypryt et al. (2017)

$\Delta_{0}$ $0.14\text{\,}\mathrm{m}\mathrm{e}\mathrm{V}$ zero temperature superconducting gap energy $$1.764$\;k_{B}\;T_{c}$

$\epsilon$ $0.25\text{\,}\mathrm{m}\mathrm{e}\mathrm{V}$ mean quasiparticle creation energy $$\sim 1.75$\;\Delta_{0}$ Kozorezov et al. (2000)

$l_{ph}$ $2.0\text{\,}\mathrm{n}\mathrm{m}$ phonon mean free path at $\Omega_{D}$ equation S19

$\lambda$ $0.43$ effective electron-phonon coupling constant equation S21

$\Omega_{1}$ $5.3\text{\,}\mathrm{m}\mathrm{e}\mathrm{V}$ threshold energy for the phonon cascade equation S20

$g_{1}\!\left(\nicefrac{{\Omega_{D}}}{{\Omega_{1}}}\right)$ $0.55$ special function numerically evaluated Kozorezov et al. (2008)

$J$ $2.5$ phonon noise factor equation S18

$\eta_{t}$ is the phonon transmission coefficient from PtSi to sapphire. While sapphire is significantly harder than PtSi, PtSi is more dense. This relationship leads to the acoustic impedance of the two materials being fairly well matched, so we will take $\eta_{t}\sim 1$.

We can write some of the parameters in equation S18 in terms of more tractable material parameters. Kaplan et al. (1976); Kozorezov et al. (2008)

[TABLE]

and

[TABLE]

where $N_{0}$ is the single spin density of states at the Fermi energy, $v_{s}$ is the speed of sound in the superconductor, and $\left<\alpha^{2}\right>_{\text{avg}}$ is the averaged squared electron-phonon interaction energy.

For the case of PtSi, we don’t have a good measurement of $\lambda$, but it can be estimated using McMillan (1968)

[TABLE]

$T_{c}$ is the superconducting transition temperature. $\mu^{*}$ is the Coulomb pseudo-potential which is also unknown, but generally takes a very small range of values for different materials. We make the approximation $\mu^{*}\sim$0.13$$.

The values used to compute the resolving power from equation S17 are listed in table LABEL:tab:parameters, resulting in $R\sim$20 for $808\text{\,}\mathrm{n}\mathrm{m}$ photons and $R\sim$15 for $1310\text{\,}\mathrm{n}\mathrm{m}$ photons. While a few of the parameters needed for this calculation are well known for our array, some had to be taken from theoretical calculations and could be inaccurate. Moreover, the values of $\epsilon$, $F$, $\mu^{*}$, and $\left<\alpha^{2}\right>_{\text{avg}}$ are completely unknown for PtSi and had to be estimated from values for other superconductors. The largest error likely comes from our choice in $\left<\alpha^{2}\right>_{\text{avg}}$. Considering the extremes of the values presented by Kaplan et al., we expand these estimates to $R\sim$1330 for $808\text{\,}\mathrm{n}\mathrm{m}$ photons and $R\sim$1024 for $1310\text{\,}\mathrm{n}\mathrm{m}$ photons