Entanglement-enhanced dual-comb spectroscopy

Abstract

Dual-comb interferometry harnesses the interference of two laser frequency combs to provide unprecedented capability in spectroscopy applications. In the past decade, the state-of-the-art systems have reached a point where the signal-to-noise ratio per unit acquisition time is fundamentally limited by shot noise from vacuum fluctuations. To address the issue, we propose an entanglement-enhanced dual-comb spectroscopy protocol that leverages quantum resources to significantly improve the signal-to-noise ratio performance. To analyze the performance of real systems, we develop a quantum model of dual-comb spectroscopy that takes practical noises into consideration. Based on this model, we propose quantum combs with side-band entanglement around each comb lines to suppress the shot noise in heterodyne detection. Our results show significant quantum advantages in the uW to mW power range, making this technique particularly attractive for biological and chemical sensing applications. Furthermore, the quantum comb can be engineered using nonlinear optics and promises near-term experimentation.

Introduction

Dual-comb interferometry, a frequency-comb-based precision measurement technique harnessing the interference of two laser frequency combs of slightly different repetition rates in a static device, has emerged to provide unprecedented capability in various applications including spectroscopy^1,2,3, hyperspectral imaging^4,5, and light detection and ranging (LiDAR)^6,7,8,9,10. In combination with on-chip frequency comb generators, dual-comb technique has been demonstrated with various platforms, including quantum cascade lasers^11,12, micro-resonator-based soliton combs^13,14, electro-optic micro-rings¹⁵, and on-chip semiconductor lasers¹⁶.

In terms of spectroscopy, dual-comb interferometry has advantages of (1) rapid effective data acquisitions without mechanical moving parts^17,18,19; (2) broad spectral coverage over the large span of a comb generator^20,21,22,23; (3) spectral resolution reaching the comb line spacing of sub-picometers^21,24; (4) frequency scale calibrated with the accuracy of an atomic clock²⁵; (5) feasibility to on-chip integration^13,14, whose high repetition rates supports rapid measurements. For linear spectroscopy, it acquires thousands of molecular transitions simultaneously^20,21,22, providing rich spectral information for quantitative concentration analysis on a sample^26,27,28. For nonlinear spectroscopy, the high peak power of the comb sources has been utilized for coherent anti-stokes Raman spectroscopy¹⁷ and four-wave-mixing multi-dimensional spectroscopy²⁹ on bio-chemicals, with the potential to improve the measurement speed by several orders of magnitude.

Similar to other broadband spectroscopic techniques, the sensitivity in a dual-comb measurement is inversely proportional to the optical bandwidth as the laser power in photo-detection is constrained due to detector nonlinearity or sample damage^2,30. In this regard, the product of signal-to-noise ratio (SNR) per unit acquisition time and the number of resolved spectral elements is taken as the figure of merit for a dual-comb system². As long-term coherent averaging^21,24,31,32 improves the sensitivity at the cost of acquisition time^20,21,22,23, these techniques hinder real-time sensing and fail to improve the figure of merit. Indeed, one has to suppress various noise sources to achieve a higher figure of merit, including detector noises, laser relative intensity noise (RIN), shot noise and others. While more device engineering efforts can potentially suppress detector noise and RIN, the fundamental shot noise remains. In this regard, state-of-the-art dual-comb sensing systems have shown examples where shot noise is dominant^33,34,35.

To go beyond the shot noise limit, quantum resources such as squeezing and entanglement are necessary. For example, the Laser Interferometer Gravitational-wave Observatory (LIGO)^36,37,38 and the Haloscope At Yale Sensitive To Axion CDM (HAYSTAC) dark matter search³⁹ inject squeezed light to suppress the shot noise. Photonic radar⁴⁰ and optomechanical force sensing⁴¹ adopted entanglement with the distributed sensing paradigm⁴². In terms of spectroscopy, amplitude-squeezing has been demonstrated in nonlinear spectroscopy⁴³ and entangled two-mode squeezed vacuum has been shown to benefit linear absorption spectroscopy⁴⁴. However, none of these advantages directly applies to dual-comb spectroscopy, as its essential component of heterodyne detection presumably precludes the use of squeezing and entanglement.

In this work, we develop a quantum description of dual-comb spectroscopy and then propose an entanglement-enhanced scheme that utilizes quantum combs of light to gain sensitivity enhancement in dual-comb spectroscopy. We first provide a complete quantum model for dual-comb spectroscopy, which recovers the SNR results of ref. ³⁰ in the case of classical source. Furthermore, the quantum model allows us to design a quantum comb composed of pair-wise entanglement around each strong comb line to improve the SNR drastically. The quantum advantage is robust against loss and phase misalignment. We also provide an experimental design to engineer the quantum comb with off-the-shelf components.

Results

Overview of the protocol

Dual-comb spectroscopy employs the interference between the signal comb (shown in red) and the local comb (shown in blue), as shown in Fig. 1a. The protocol’s strength stems from the selection of signal and local combs with slightly different frequency spacings, f_r and f_r + Δf_r, respectively, as illustrated in Fig. 1b. The signal comb interrogates the sample and undergoes loss and phase shift that can be modeled as a bosonic quantum channel, with frequency-dependent transmissivity κ(f) and phase-shift α(f). Meanwhile, the local comb serves as a local oscillator (LO) for the final heterodyne measurement. After mixing the LO and return at a balanced beamsplitter, information about the sample can be extracted from the photocurrent difference ${\hat{N}}(t)={\hat{c}}_{+}^{{\dagger} }(t){\hat{c}}_{+}(t)-{\hat{c}}_{-}^{{\dagger} }(t){\hat{c}}_{-}(t)$, obtained from the photocurrent measurement on both output ports ${\hat{c}}_{\pm }(t)$. While we utilize quantum operator language to describe the measurement to prepare our analyses for quantum combs, our analyses also recover the results obtained from semi-classical analyses for classical combs³⁰.

Fig. 1: Schematic and performance advantage of entanglement-enhanced dual comb spectroscopy.

a Conceptual schematic. The teeth share intermodal entanglement within the signal comb (red beam). LO local osillator, BS beamsplitter, PD photon detector. b Schematic of the quantum comb. Each pair of signal modes beating with the same LO comb tooth (purple line) for the same intermediate frequency is entangled, indicated by a black dashed line connecting a pair of purple circles. c Practical SNR involving NEP-type and RIN-type noises, plotted versus signal power (analog to Fig. 2 of ref. ³⁰), normalized to unit acquisition time T = 1 s. We assume an ideal detector with unity efficiency, and zero loss and noise κ ≈ 1, η = 1. In c, both signal and LO are entangled with equal gain G, which increases from 0 dB (coherent-state) to 30dB in steps of 10 dB, plotted in color from blue to magenta. The NEP/RIN-dictated SNR is presented by green-dot-dashed/black-dot-dashed line, along with the shot noise (SN) limit in blue-dashed. N = 10⁵, λ = 1 μm, RIN = −170 dBc Hz⁻¹⁴⁷ (a more accessible RIN = −150 dBc Hz⁻¹ enforces an earlier saturation shown by the gray dot-dashed line), NEP = 5 × 10⁻¹³ W Hz^−1/2 (NEP = 4.5 × 10⁻¹⁵ W Hz^−1/2 is actually achievable, e.g. by Thorlabs FGA01FC-InGaAs Photodiode), P_LO/P_S = 5. Inset: Quantum advantage in SNR (in decibel unit) versus various values of RIN for total signal power P_S = 10 mW (blue) and P_S = 10 μW (red) at G=20dB. The LO-signal power ratio γ ≡ P_LO/P_S = 5 is fixed—although this figure shows the case of both signal and LO entangled, as we show later, only the signal needs to be entangled under large γ.

In the classical protocol, both combs are classical—the quantum state of each comb line is in a coherent state, obeying the shot-noise-limited standard quantum limit (SQL). We propose to engineer quantum combs to further improve the performance of dual-comb spectroscopy. To begin with, we consider the case of signal comb being quantum engineered. To suppress noise below the SQL, squeezing is commonly adopted in quantum sensing protocols. However, in the case of dual-comb acquisitions, squeezing a single mode alone is inadequate, as heterodyne measurement is necessary to read out quadratures across the entire spectrum. To surpass the SQL, entanglement between different frequency modes is required so that joint quadratures are squeezed. Thus, we propose entangling the side-bands of the signal comb around each local comb line, as indicated in Fig. 1b by the dashed lines. Such an entangled comb with squeezing gain G allows both measured quadratures in heterodyne to be squeezed, resulting in fundamental noise a factor of 1/G below the SQL. Furthermore, as we will detail in later part of the paper, in general, the LO can be engineered to be similarly entangled, which can further improve the performance, especially when the power of the local comb is similar as or lower than that of the signal comb. Note that quantum engineering almost preserves the power of the comb, as the joint-squeezing power is negligible compared to the mean field in dual-comb systems.

To analyze the performance of the proposed quantum dual-comb spectroscopy protocol, we model the signal and all noise involved in the protocol systematically. In this overview, the case of bright LO is considered, where the power P_LO is much larger than the signal comb power P_S, while the full analysis is presented later. The information about the sample is derived from the amplitude decay and phase-shift of the signal comb, ${{\rm {e}}}^{i\alpha }\sqrt{\kappa {P}_{{{{\rm{S}}}}}}$, which is subject to contamination from various sources of noise. For low-loss (κ ~ 1) and low-noise scenarios, typical for room-temperature dilute chemical gas sensing and bio-sensing with thin sample slices at signal wavelength ≲ 5 μm, the fundamental noise stemming from fluctuation properties of the light field ~ 1/2G for the quantum comb. In terms of estimating the transmissivity $\sqrt{\kappa }$, it contributes an inverse-law noise term O(1/P_S). In practice, the detector noises characterized by noise equivalent power (NEP) and the relative intensity noise (RIN) in laser sources also mix in, which dominate the estimation error at the low and high comb power region respectively, as analyzed for the classical protocol in ref. ³⁰. In addition to the inverse-law term due to the fundamental noise, the RIN adds a constant noise term O(1) independent of P_S, denoted as RIN-type noise, and the NEP adds an inverse-square-law noise term $O(1/{P}_{{{{\rm{S}}}}}^{2})$ to the estimation, denoted as NEP-type noise. We ignore detector dynamical range noise, as it can be resolved by engineering the detector array³⁰ and it has similar effects as RIN-type noise.

With all noises into consideration, we consider the realistic performance of the quantum dual-comb spectroscopy system. For generality, our calculation covers a wide comb power range from 10⁻⁷ to 10⁻¹ W, the same as that in ref. ³⁰. Indeed, in a typical dual-comb application (power ratio γ = 1), the total comb power should be constrained below 50 μW in the near-infrared range²¹ and below 30 μW in the mid-infrared region^25,34 to avoid detector nonlinearity. This is, most of the time, due to the requirement of a high dynamic range to sample the strong center burst in the dual-comb interferogram. A promising route of engineering a high dynamic range detector has been able to extend to operation power up to 2 mW without clear spectral artifacts observed, as demonstrated in ref. ⁴⁵.

We plot the signal-to-noise ratio (SNR) per second versus the total signal power under practical experimental settings in Fig. 1c. In this μW-mW region, the performance of state-of-the-art classical systems (blue solid) is limited by the shot noise (SN) limit (blue dashed). While at low/high signal power limit, the SNR converges to the NEP/RIN-dictated limit (green dot-dashed/black dot-dashed). With the quantum comb, we see that a practical entangled source of 10 dB gain (G = 10, purple solid) yields a quantum advantage up to 4.9 dB over the coherent-state source (blue solid). As the gain increases (from blue to magenta), the quantum advantage improves further until it saturates subject to the limits dictated by NEP-type noise alone (greed dot-dashed) and the RIN-type noise alone (black dot-dashed). In the scenario of Fig. 1c, we observe that the ultimate limit of quantum advantage can go up to 13.4 dB at P_S ≈ 0.1 mW, which is of great interest to bio-sensing applications. Additionally, we provide predictions at the power levels when such saturation happens (dots on the blue-magenta curves), as we detail in the “Methods” section. Our analyses show that for a state-of-art dual-comb system to enjoy quantum advantage, RIN-type noise is often the major constraint: in the inset of Fig. 1, we show that for 10 mW signal power, to enjoy a significant quantum advantage it requires RIN-type noise to be smaller than ~−170 dBc Hz⁻¹, challenging but still possible⁴⁶; For lower power of 10 μW, the requirement is less stringent, RIN = −150 dBc Hz⁻¹ is readily achievable^30,47.

General amplitude and phase detection

In a general spectroscopy sensing process, one is interested in the input–output relation of light for a range of frequencies f. The pattern of the output light reveals information about the composition of the sample under study. The mathematical model for the input–output relation involves a thermal-loss phase-shift channel, which has frequency-dependent transmissivity κ(f) and phase-shift α(f). Given an input light mode described by the annihilation operator ${\hat{a}}_{{\rm {S}}}$ satisfying the canonical commutation relation $[{\hat{a}}_{{\rm {S}}}^{{\dagger} },{\hat{a}}_{{\rm {S}}}]=1$, the output field annihilation operator is given by the linear relation

$${\hat{a}}_{{\rm {R}}}=\sqrt{\kappa }{\rm {{e}}}^{i\alpha }{\hat{a}}_{{\rm {S}}}+\sqrt{1-\kappa }{\hat{a}}_{\rm {{E}}}.$$

(1)

The channel attenuates the mean of input signal mode ${\hat{a}}_{\rm{S}}$ by ${\sqrt{\kappa}}$, shift the phase by α, and mixes in the environment mode ${\hat{a}}_{\rm {E}}$ with mean thermal photon number given by the Bose–Einstein distribution ${{{\mathcal{E}}}}(f)=1/[\exp \left(hf/{k}_{\rm{B}}{T}_{\rm{B}}\right)-1]$. Here h is the reduced Planck constant, k_B is the Boltzmann constant and T_B is the sample environment temperature. Although the thermal noise ${\mathcal{E}}(f){\ll} 1$ is negligible at the frequency of interest, we will keep it in our analyses to tackle the general case.

Although our results work for the simultaneous estimation of phase-shift and transmissivity, as enabled by dual-comb technique, we consider two special scenarios to simplify the SNR analyses. In the first scenario, we are concerned with only the transmissivity κ(f), while the phase-shift is negligible due to phase cancellation via sending both combs to the sample. In the second scenario, the absorption is almost zero (κ(f) ≃ 1), while the phase-shift α(f) (despite also being small) provides the major information. For example, imaging the subtle changes in phase contrast on the order of 0.1 milli-radians allows the study of neural activities at the single neuron level^48,49. Overall, the absorption or phase-shift can be very weak due to the low concentration of the sample, as is the case in atmospheric sensing and human breath analysis²⁸ at parts per billion and in radiocarbon detection at few parts per quadrillion⁵⁰.

In our analyses, we will ignore phase noise, since various dual-comb noise suppression techniques have been developed, based on comb-source engineering^51,52,53, digital phase correction^22,31,45, and active stabilization²¹, with some of them achieving coherent times up to hours^21,22,45.

Quantum model of dual-comb spectroscopy

Now we formulate the quantum theory for dual-comb spectroscopy of N frequency components. In the rotating frame of the carrier frequency ν₀ (ν₀ ≫ f_r), the signal comb is represented by the field operator

$$\hat{A}(t)=\frac{1}{\sqrt{T}}\left[\hat{a}(t)+\mathop{\sum }\limits_{n=1}^{N}{A}_{n}{\rm {{e}}}^{i2\pi n({f}_{{{{\rm{r}}}}}+{{\Delta }}{f}_{{{{\rm{r}}}}})t}\right],$$

(2)

while the local oscillator (LO) comb is represented by

$$\hat{B}(t)=\frac{1}{\sqrt{T}}\left[\hat{b}(t)+\mathop{\sum }\limits_{n=1}^{N}{B}_{n}{\rm {{e}}}^{i2\pi n{f}_{{{{\rm{r}}}}}t}\right].$$

(3)

Here T is the acquisition time, and should be long enough to resolve the frequency of interest. The sum in each comb consists of the strong mean fields of a frequency comb source at discrete frequencies. The light power is mainly contributed by these mean fields. Specifically, the power is ${P}_{{{{\rm{S}}}}}=h{\nu }_{0}\mathop{\sum }\nolimits_{n = 1}^{N}| {A}_{n}{| }^{2}/T$ for the signal, and ${P}_{{{{\rm{LO}}}}}=h{\nu }_{0}\mathop{\sum }\nolimits_{n = 1}^{N}| {B}_{n}{| }^{2}/T$ for the LO, while the additional power due to squeezing is negligible in this paper. The quantum-operator term in each comb describes the noise

$$\begin{array}{r}\hat{z}(t)=\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{\delta =-N}^{N}{\hat{z}}_{n,\delta }{\rm {{e}}}^{i2\pi [n{f}_{{{{\rm{r}}}}}+\delta {{\Delta }}{f}_{{{{\rm{r}}}}}]t}\end{array}$$

(4)

where $\hat{z}\in \{\hat{a},\hat{b}\}$, and we have quantized the frequency modes of the field into the field annihilation operators, satisfying the commutation relation $[{\hat{a}}_{n,\delta },{\hat{a}}_{n,\delta }^{{\dagger} }]=[{\hat{b}}_{n,\delta },{\hat{b}}_{n,\delta }^{{\dagger} }]=1$ and all the other commutators are zero. As indicated in Fig. 2, the double subscripts n ∈ [1, N] and δ ∈ [ −N, N] determines the frequency of the mode nf_r + δΔf_r, with n denoting which comb line the mode is around and δ further specifying which mode around the comb line specified by n. As Δf_r ≪ f_r, the sideband modes around different comb lines will not overlap. Here we have included all frequency modes relevant to the heterodyne measurement. In a classical strategy, the noise property of all modes is vacuum-limited, but in this work, we propose to engineer the noise property via squeezing and entanglement.

Fig. 2: The frequency arrangement of the comb modes.

A red peak B_n in LO beats with the side band modes in the signal (connected by red dashed lines), whose quadrature fluctuation contributes to the overall noise. Therefore, we entangle the side bands (connected by red dashed lines) of the signal to improve the SNR. Similarly, a blue peak A_n in signal beats with the sideband mode pairs in LO connected by blue dashed lines, which are to be entangled to improve the SNR. Around each peak, there are N such pairs. Here we explicitly label the pairs at the edges.

Field propagation through the sample can be formulated by a bosonic quantum channel, as shown in Eq. (1). We are interested in the transmissivity κ(f) and phase α(f) induced by the sample. Note that the non-ideal LO storage also induces a channel of transmissivity η(f) and phase-shift β(f). We assume that the transmissivity and phase-shift spectra are smooth enough such that their values at sidebands of each comb line are identical. For example, κ(nf_r + δΔf_r) = κ_n for all sideband frequencies −N ≤ δ ≤ N. We define α_n and η_n, β_n similarly. Our formalism can be easily generalized to rapidly varying spectra, while the formula will turn much lengthier. After traveling through the sample, channel output fields ${\hat{A}}^{{\prime} }(t)$ for the sample return and ${\hat{B}}^{{\prime} }(t)$ for the LO can be decomposed in the same form as Eqs. (2) and (3). According to Eq. (1), the input–output relation yields ${A}_{n}\to \sqrt{{\kappa }_{n}}{A}_{n}{\rm {{e}}}^{i{\alpha }_{n}}$ and ${B}_{n}\to \sqrt{{\eta }_{n}}{B}_{n}{\rm {{e}}}^{i{\beta }_{n}}$ for mean fields, and

$$\begin{array}{rcl}{\hat{a}}_{n,\delta }^{{\prime} }&=&\sqrt{{\kappa }_{n}}{\rm {{e}}}^{i{\alpha }_{n}}{\hat{a}}_{n,\delta }+\sqrt{1-{\kappa }_{n}}{\hat{e}}_{n,\delta },\\ {\hat{b}}_{n,\delta }^{{\prime} }&=&\sqrt{{\eta }_{n}}{\rm {{e}}}^{i{\beta }_{n}}{\hat{b}}_{n,\delta }+\sqrt{1-{\eta }_{n}}{\hat{f}}_{n,\delta },\end{array}$$

(5)

for noise modes, where ${\hat{e}}_{n,\delta }$’s and ${\hat{f}}_{n,\delta }$’s are environmental noise modes of thermal photon number ${{{{\mathcal{E}}}}}_{n}\equiv {{{\mathcal{E}}}}(n{f}_{{{{\rm{r}}}}})$.

At the receiver, the two combs are combined by a balanced beamsplitter, yielding two output combs, ${\hat{c}}_{\pm }(t)=[{\hat{A}}^{{\prime} }(t)\pm {\hat{B}}^{{\prime} }(t)]/\sqrt{2}$. Then the photon counts of the two output combs are measured, and subtracted from each other $\hat{N}(t)={\hat{c}}_{+}^{{\dagger} }(t){\hat{c}}_{+}(t)-{\hat{c}}_{-}^{{\dagger} }(t){\hat{c}}_{-}(t)={\hat{A}}^{{\prime} {\dagger} }(t){\hat{B}}^{{\prime} }(t)+{\hat{B}}^{{\prime} {\dagger} }(t){\hat{A}}^{{\prime} }(t)$. Taking into account that Δf_r ≪ f_r, we can filter out the direct current (DC) term and the fast-oscillating terms at ∣f∣ ≫ Δf_r. The resulting alternative current (AC) N_AC(t) is a random variable with mean

$$\begin{array}{r}{\overline{N}}_{{{{\rm{AC}}}}}(t)=\frac{1}{T}\left[\mathop{\sum }\limits_{n=1}^{N}\sqrt{{\kappa }_{n}{\eta }_{n}}{{\rm {e}}}^{i({\alpha }_{n}-{\beta }_{n})}{A}_{n}{B}_{n}{\rm {{e}}}^{i2\pi n{{\Delta }}{f}_{{{{\rm{r}}}}}t}+{\rm {c.c.}}\right],\end{array}$$

(6)

where c. c. represents the complex conjugate. One can perform a finite-time-T Fourier transform to obtain the discrete spectrum mean

$${\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\sqrt{{\kappa }_{m}{\eta }_{m}}{A}_{m}{B}_{m}{\rm {{e}}}^{i({\alpha }_{m}-{\beta }_{m})},1\le m\le N,$$

(7)

from which we can extract information about the transmissivities and phase shifts across the entire N-line spectrum.

To evaluate the fluctuation of the readout, now we consider the contribution to N_AC(mΔf_r) from noise modes $\hat{a},\hat{b}$. As the amplitudes A_n, B_n ≫ 1, the noise in N_AC(mΔf_r) is

$$\begin{array}{l}{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\simeq \mathop{\sum }\limits_{n=1}^{N}\left[\sqrt{{\eta }_{n}{\kappa }_{n}}{B}_{n}{\hat{X}}_{n,m}+\sqrt{{\eta }_{n}{\kappa }_{n}}{A}_{n}{\hat{Q}}_{n,m}\right.\\ \qquad\qquad\qquad \left.+\,\sqrt{{\eta }_{n}(1-{\kappa }_{n})}{B}_{n}{\hat{X}}_{n,m}^{(e)}+\sqrt{(1-{\eta }_{n}){\kappa }_{n}}{A}_{n}{\hat{Q}}_{n,m}^{(f)}\right].\end{array}$$

(8)

For the full derivation of Eqs. (7) and (8), see the “Methods” section. Here we have adopted the nomenclature widely used in quantum optics^54,55 that defines the joint quadrature operators

$$\begin{array}{l}{\hat{X}}_{n,m}\,\equiv \,{\hat{a}}_{n,m}{{\rm {e}}}^{i({\alpha }_{n}-{\beta }_{n})}+{\hat{a}}_{n,-m}^{{\dagger} }{{\rm {e}}}^{-i({\alpha }_{n}-{\beta }_{n})},\\ {\hat{Q}}_{n,m}\,\equiv \,{\hat{b}}_{n,n+m}{{\rm {e}}}^{-i({\alpha }_{n}-{\beta }_{n})}+{\hat{b}}_{n,n-m}^{{\dagger} }{{\rm {e}}}^{i({\alpha }_{n}-{\beta }_{n})},\end{array}$$

(9)

for the signal (${\hat{a}}_{n,m},{\hat{a}}_{n,-m}$ beat with the strong mean field B_n at frequency nf_r) and for the LO (${\hat{b}}_{n,n+m},{\hat{b}}_{n,n-m}$ beat with A_n at frequency n(f_r + Δf_r)) respectively. Simiarly, we define the quadratures ${\hat{X}}^{(e)},{\hat{Q}}^{(f)}$ for the environment modes $\hat{e}$ and $\hat{f}$ in Eq. (5). Note that these quadratures, along with ${\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}$, are usually not Hermitian (real-valued) observables, thus their variances are $\,{{\mbox{var}}}\,\,\hat{X}\equiv \langle \hat{{X}}^{{\dagger} }\hat{X}\rangle$ for any non-Hermitian complex operator $\hat{X}$.

Dual-comb spectroscopy aims to estimate the transmissivity κ_n, phase-shift α_n or both simultaneously, for all 1 ≤ n ≤ N frequencies of the sample from the photo-current of Eq. (7). We define the amplitude SNR at each comb line as

$$\,{{\mbox{SNR}}}\,=| {\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})| /\sqrt{{{{\rm{var}}}}\left[{N}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\right]}.$$

(10)

The noise, which is defined in Eq. (8), collects the beating modes near all N comb lines. As shown in Methods, it is a good indicator for the minimum mean square error of either the transmissivity or the phase-shift estimation task. Furthermore, a neat figure of merit is the overall quality factor N SNR, which eliminates the dependence on total line number N.

To evaluate the SNR, we make use of the independence between modes around different comb lines and evaluate the variance from Eq. (8),

$${{{\rm{var}}}}\left[{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\right]=\mathop{\sum }\limits_{n=1}^{N}\left[{{{{\mathcal{N}}}}}_{n}+{\eta }_{n}{\kappa }_{n}\left({B}_{n}^{2}{{{\rm{var}}}}{\hat{X}}_{n,m}+{A}_{n}^{2}{{{\rm{var}}}}{\hat{Q}}_{n,m}\right)\right],$$

(11)

where the thermal noise ${{{{\mathcal{N}}}}}_{n}={\eta }_{n}{B}_{n}^{2}(1-{\kappa }_{n})(2{{{{\mathcal{E}}}}}_{n}+1)+{\kappa }_{n}{A}_{n}^{2}(1-{\eta }_{n})(2{{{{\mathcal{E}}}}}_{n}+1)$ is determined by the sample, the LO storage and the environment temperature; the complex quadrature noises ${{{\rm{var}}}}{\hat{X}}_{n,m}$ and ${{{\rm{var}}}}{\hat{Q}}_{n,m}$ are determined by the quantum state of the signal comb source and local comb source.

SNR with entangled quantum comb

From the definition of quadratures in Eq. (9), we see that the noise ${{{\rm{var}}}}{\hat{X}}_{n,m}$ can be suppressed by entangling the modes ${\hat{a}}_{n,\pm m}$ in a two-mode squeezed vacuum state (see the “Methods” section), as indicated in Fig. 2. By such means, the joint quadrature ${\hat{X}}_{n,m}$ is squeezed with suppressed variance (see the “Methods” section)

$${{{\rm{var}}}}{\hat{X}}_{n,m}=\frac{1}{2G}\left[-\left({G}^{2}-1\right)\cos \left(2{\alpha }_{n}-2{\beta }_{n}\right)+({G}^{2}+1)\right],$$

(12)

where squeezing gain G ≥ 1. When phases are perfectly matched as ${\alpha }_{n}-{\beta }_{n}=0,{{{\rm{var}}}}{\hat{X}}_{n,m}$ is minimized to 1/G. Similarly, we can squeeze the joint quadrature ${\hat{Q}}_{n,m}$ of the local comb by gain G_LO. When G = G_LO = 1, the variance reduces to the classical dual-comb spectroscopy. In this case, the variance of complex operator ${{{\rm{var}}}}{\hat{X}}_{n,m}$ is twice of the SQL 1/2, because it is defined as a sum of variances of its real and imaginary parts.

To model the full SNR of the dual-comb spectroscopy system, we involve device and source imperfections. For simplicity, we assume the N comb lines are generated symmetric (${A}_{n}={A}_{{n}^{{\prime} }}$ and ${B}_{n}={B}_{{n}^{{\prime} }}$). In Methods, we derive the full formula of the SNR at intermediate frequency mΔf_r

$$\begin{array}{r}{{{\rm{SNR}}}}^{-2}=\frac{{N}^{2}}{T}\left({a}_{{{\rm {{NEP}}}}}\frac{1}{{P}_{{{{\rm{S}}}}}^{2}}+\frac{{a}_{{{\rm {{quad}}}}}}{{P}_{{{{\rm{S}}}}}}+{a}_{{{\rm {{RIN}}}}}\right),\end{array}$$

(13)

where T is the acquisition time, P_S is the total signal power, the NEP-type noise coefficient a_NEP ≡ NEP²/η_mκ_mγ, the RIN-type noise coefficient a_RIN ≡ RIN/2, and the quadrature noise coefficient

$${a}_{{{\mbox{quad}}}}\equiv \frac{h{\nu }_{0}}{N}\mathop{\sum }\limits_{n=1}^{N}\left[\left({{{\rm{var}}}}{\hat{X}}_{n,m}+\frac{1}{\gamma }{{{\rm{var}}}}{\hat{Q}}_{n,m}\right)+\frac{{{{{\mathcal{N}}}}}_{n}}{{\kappa }_{n}{\eta }_{n}}\right].$$

(14)

Here γ ≡ P_LO/P_S is the LO-to-signal power ratio, hν₀ is the energy per photon. Note the quadrature noises can be suppressed by the entanglement (joint quadrature squeezing) in Eq. (12). The proposed SNR quantifies the performance of both the transmissivity estimation and phase estimation scenarios.

The SNR of Eq. (13) versus total signal power P_S has been evaluated in Fig. 1c for the case of G = G_LO > 1, which highlights the quantum advantages of entanglement. We have taken the case where the phase mismatch α_n − β_n ≪ 1/G are all small. This is the case when one estimates the transmissivity with good phase locking, or estimates small phase-shift caused by weak samples—phase-shift as small as 0.1 milli-radians allows the study of neural activities at the single neuron level^48,49. When there is phase noise, then α_n−β_n cannot be set to zero. Suppose there is a mismatch ϵ ≪ 1, then the leading order of the variance in Eq. (12) ${{{\rm{var}}}}{\hat{X}}_{n,m}=1/G+2{\epsilon }^{2}G$ will mix in the anti-squeezing part. Therefore, phase-locking is crucial for the proposed EA dual-comb spectroscopy, similar to other squeezing-enabled protocols.

Our formula can be regarded as a quantum version of Eq. (2) of ref. ³⁰. We note that our quantum model yields an SNR-γ relation different from the semiclassical model in ref. ³⁰. Specifically, for a fixed signal comb power P_S, we find that the optimum is at γ → ∞ in our quantum model, while the optimum is finite in the semiclassical model. This is because when γ is large, the RIN-type noise increases proportional to γ in the semiclassical model, while the RIN-type noise remains constant in our quantum model. Our result on the RIN-type noise agrees with ref. ⁵⁶.

Strategies for applying quantum combs

In the above analyses, we have allowed both the signal comb and local comb to be quantum-engineered. In general, having quantum entangled combs in both arms might be not necessary and experimentally challenging. Here we address different scenarios of applying the quantum combs.

From Eq. (11), we see that the noises from signal and LO are indeed amplified by the mean field of the other, i.e. the signal noise is amplified by the LO mean and vice versa. Hence, squeezing merely the signal is sufficient to yield a significant advantage when P_LO dominates, while squeezing merely the LO is sufficient when P_S dominates.

Now we evaluate the quantum advantages of various entanglement (or joint squeezing) strategies over the classical coherent-state source in terms of amplitude SNR. In Fig. 3, row 1 shows the ideal advantage without practical noises (note that the ideal advantage depends on the relative power ratio γ = P_LO/P_S between the signal and LO only, not the absolute magnitudes of power). Rows 2 and 3, operating under signal power P_S = 10 mW and P_S = 10 μW, respectively, show the practical advantage with NEP-type and RIN-type noise involved. Along the horizontal axis, all contours are centered at the signal power P_LO = P_S, thus at the left half the signal dominates while at the right half, the LO dominates.

Fig. 3: Quantum advantage in amplitude SNR over the coherent-state source versus various LO-signal power ratios γ ≡ P_LO/P_S and sample transmissivity κ.

We consider the squeezing gain of 10 dB and uniform transmissivity. The quantum advantage is in the decibel unit. Rows: first row, a–c: without practical noise; second row, d–f: with both NEP-type and RIN-type noises, total signal power P_S = 10 mW; third row, g–i: with both NEP-type and RIN-type noises, P_S = 10 μW. Columns: first column, a, d, g: signal squeezed only; second column, b, e, h: local comb (LO) squeezed only; third column, c, f, i: both signal and LO squeezed. We assume an ideal detector with unity efficiency and ideal LO link η_m = 1. N = 10⁵, T = 1 s, λ = 1 μm, NEP = 5 × 10⁻¹³ W Hz^−1/2, RIN = − 170 dBc Hz⁻¹.

In row 1, we verify our results of fundamental limits for the ideal advantages. In Fig. 3a, only the signal is squeezed, we see that the advantage peaks at the LO-dictated region (right half); in Fig. 3b, only the LO is squeezed, the advantage now peaks at the signal-dictated region (left half); finally, in Fig. 3c, both the signal and LO are squeezed, here we enjoy both the advantageous areas in the two squeezing strategies above. It is noteworthy that when LO is squeezed, the quantum advantage survives even when κ → 0 as shown at the bottom-left corners of subplots Fig. 3b and c, which is useful when the sample is lossy and LO power is limited.

In rows 2 and 3, we consider the effect of practical noises. In Fig. 3d–f, the signal power P_S = 10 mW is relatively large. In this scenario, the advantage is mainly constrained by the RIN-type noise which is more significant for large P_LO. We find that the advantages at the LO-dictated region (right half) of all subplots Fig. 3d–f are significantly undermined, which is especially noticeable in subplot Fig. 3d. On the other hand, in Fig. 3g–i, the signal power P_S = 10 μW is extremely small. In this scenario, the advantageous region is mainly affected by the NEP-type noise which is more significant for small P_LO. As expected, the advantages at the left half region of subplots Fig. 3g–i are undermined, which is especially noticeable in subplot Fig. 3h when compared with Fig. 3b. In subplot Fig. 3f or i where both the signal and LO are squeezed, the patterns in the two squeezing strategies shown in subplot Fig. 3d, e or g, h occur simultaneously.

Performance under total power constraints

To begin with, we consider the power dependence of the SNR on P_LO, P_S. Here we explore the scenario where the LO is sent along with the signal, thus η_m = κ_m = κ and the total power exposure P_S + P_LO is to be constrained. Such constraint also naturally appear when detector nonlinearity and saturation is taken into account. Figure 4a, b shows that the total power P_S + P_LO contour and the SNR contour are tangent at P_LO = P_S. This explains that in some applications one tends to use comparable LO and signal (i.e. the power ratio γ = 1) rather than very strong LO to save the total power consumption.

Fig. 4: SNR and quantum advantage versus signal and LO power.

The dependence of a, b the amplitude SNR using the squeezed source and c, d its quantum advantage over the coherent-state source, both in decibel unit, on various power constraints P_LO, P_S. Uniform transmissivity is assumed. Both signal and LO are squeezed at the squeezing gain of 10 dB. a, c κ = 1; b, d κ = 0.5. The black dashed lines are contours of total power P_LO + P_S. We assume the ideal detector with unity efficiency, while the LO is sent along with the signal so that η_m = κ. N = 10⁵, T = 1 s, λ = 1 μm, NEP = 5 × 10⁻¹³ W Hz^−1/2, RIN = − 170 dBc Hz⁻¹.

Now consider the quantum advantage. Note that only the fundamental noise is suppressed by the quantum engineering, we expect to maximize these ratios ${\sigma }_{{{{\rm{quad}}}}}^{2}/{\sigma }_{{{{\rm{NEP}}}}}^{2}$ and ${\sigma }_{{{{\rm{quad}}}}}^{2}/{\sigma }_{{{{\rm{RIN}}}}}^{2}$ to see a significant quantum advantage. In Methods, we show that NEP-type noise ${\sigma }_{{{{\rm{NEP}}}}}^{2} \sim 1/({P}_{{{{\rm{LO}}}}}{P}_{{{{\rm{S}}}}})$, fundamental noise ${\sigma }_{{{{\rm{quad}}}}}^{2} \sim ({P}_{{{{\rm{S}}}}}+{P}_{{{{\rm{LO}}}}})/{P}_{{{{\rm{S}}}}}{P}_{{{{\rm{LO}}}}}$, RIN-type noise ${\sigma }_{{{{\rm{RIN}}}}}^{2} \sim 1$. The ratio ${\sigma }_{{{{\rm{quad}}}}}^{2}/{\sigma }_{{{{\rm{NEP}}}}}^{2}$ is proportional to the total power P_S + P_LO, while ${\sigma }_{{{{\rm{quad}}}}}^{2}/{\sigma }_{{{{\rm{RIN}}}}}^{2}$ is maximized at P_S → 0 or P_LO → 0.

For the NEP-dictated scenario, i.e. the total power is small, the quantum advantage grows with the absolute SNR as total power increases, while it does not depend on γ = P_LO/P_S. In this case one can let γ → 0 or ∞ to make signal or LO dominate so that squeezing on the other is no longer needed, as discussed previously. For the RIN-dictated scenario, i.e. the total power is large, we note that the quantum advantage decreases with the total power and it is minimized at P_S = P_LO given a fixed total power, which is opposite to the absolute SNR case. This is not a preferred scenario for quantum advantage. Figure 4c, d verifies that the total power contour and the quantum advantage contour almost overlap in the region of small total power; they are again tangent at P_S = P_LO in the region of large total power, while the gradient direction of the quantum advantage contour is reversed. Comparing subplot Fig. 4c of lossless sample κ = 1 and subplot Fig. 4d of lossy sample κ = 0.5, we see that the quantum advantage degrades significantly when the sample is lossy now that both signal and LO suffer such loss, while a 1 dB advantage still survives.

Performance limits in biological applications

An important application of the proposed entanglement-enhanced dual-comb spectroscopy system is in sensing fragile bio-tissues. In bio-sensing, the power of signal light shining on the sample is typically between 10 and 100 μW to avoid kill, bleach, or perturb the specimen under analysis^57,58,59. For example, ref. ⁴³ showed an extreme case of power at ~ 10 mW which causes severe sample damage. In another extreme case of retina sensing, the safety standard permits much lower power⁶⁰. The maximum permissible radiant power is a function of the exposure duration, wavelength, and visual angles. In a long exposure time limit, the maximum power can be as low as ≲1 μW (see Fig. 2 of ref. ⁶⁰). In other scenarios, the power can be higher, e.g. ~10−200 μW is adopted in some studies^61,62. To gain better SNR, therefore one cannot simply increase the power, but rather consider suppressing the quantum-limited noise such as in our proposal.

In bio-sensing, a major limitation of the applicable frequency region comes from water absorption. The transmissivity spectrum of water absorption can be derived from the Lambert absorption coefficient spectrum α(f)^63,64,65 via $\kappa (f)=\exp \{-\alpha (f)L\}$, where L is the sample depth. For the optical domain of wavelength λ < 1 μm, the absorption is weak: α ≲ 10⁻⁴ μm⁻¹. We take a typical sample depth of L = 15 μm and evaluate the transmissivity in Fig. 5a assuming the sample absorption is majorly dominated by water, with absorption coefficients taken from refs. ^63,64,65. We see that the absorption is substantial starting around λ ~ 2 μm. In wavelength below 2 μm, the thermal noise described by the Bose–Einstein distribution ${{{\mathcal{E}}}}(f)\ll 1$ is negligible at room temperature. In this ≤5 μm frequency of interest, ${{{\mathcal{E}}}}\le 1{0}^{-4}$ and is negligible at room temperature of 300 K. Even at 10 μm, ${{{\mathcal{E}}}} \sim 0.008$ is still small. With the absorption and noise in hand, we can evaluate the quantum advantage in the absence of any NEP-type or RIN-type noise from Eq. (14). To simplify the evaluation, we assume uniform absorption across all comb lines to get a sense of the quantum advantage. In general, the quantum advantage will be an average across a frequency region analyzed here. In Fig. 5b, we find the advantage indeed appreciable below about 2 μm and increases with the gain G, while above 2 μm water absorption starts to limit the possible advantage. Note that the quantum advantage monotonically increases with G, we expect the contour lines to rise to higher G when transmissivity dips. As expected, we see that the contour lines are almost the reverse of the transmissivity spectrum in Fig. 5a.

Fig. 5: Application in water-dominated tissue slice.

a The absorption spectrum of pure water κ(λ) at room temperature 295 K^63,64,65, evaluated for path length L = 15 μm. b The fundamental limit of quantum advantage in amplitude SNR (in decibel unit) enforced by water absorption. The spectrum data subplot (a) is unknown to the observer. An additive thermal Gaussian noise at room temperature is mixed in. P_LO/P_S = 5.

In real application scenarios, there will be additional loss in experimental implementations, which will be improved as the engineering capability advances. Therefore, in this section, we have focused on the fundamental loss from water to provide a performance limit.

Discussion

Before closing, we discuss the experimental generation of the quantum comb. We believe that the integrated lithium niobate (LN) photonic platform is well suited for the proposed scheme thanks to its large second order and Kerr nonlinearity as well as the capability of achieving quasi-phase matching via electrical poling⁶⁶. An illustration of the LN photonic chip is shown in Fig. 6. Two frequency combs with slightly different comb line spacings can be generated via pumping a continuous wave (CW) laser into a coupled-microresonator device at normal group velocity dispersion (GVD). Such normal GVD-based frequency combs are reported to have a high conversion efficiency and high comb line spacing of ~100 GHz⁶⁷, which is desirable to generate the following quantum frequency comb state. An alternative solution is an electro-optic-modulator-based frequency comb and pulse generator⁶⁸. Among the two combs, one is adopted as the classical signal comb, which will further be combined with entangled sidebands to further generate the quantum signal comb (top side of Fig. 6). The other comb is split into two beams (bottom side of Fig. 6). While one of the beam is stored as the local oscillator comb, the other beam is sent to a periodically poled LN (PPLN) waveguide to generate a comb source at the second harmonic frequencies, which is then used to pump another PPLN waveguide for two-mode squeezing around each of the local comb line via spontaneous parametric down-conversion (SPDC), during which the original local frequency comb serves as the seed. By appropriately controlling the phase between the two frequency combs, one can operate in the parametric de-amplification regime to generate amplitude amplitude-squeezed state at each comb line. The phase matching bandwidth of the SPDC needs to be less than half of the comb line spacing⁶⁹. The sideband entangled outputs are then combined with the classical signal frequency comb (shown in the top side of Fig. 6) at a highly transmissive beamsplitter, forming the quantum frequency comb source in the proposed dual-comb configuration. As compared to the bulk optical system, the integrated photonic platform can achieve the stringent phase-matching bandwidth at each comb line over a broad optical bandwidth since both the GVD and group velocity can be tailored via engineering of the nanophotonic waveguide dimension. In addition, the high repetition rate of the chip-based Kerr OFC sets a higher upper bound of the phase-matching bandwidth. The RIN of Kerr combs can be suppressed with phase-lock-loops and external-cavity acousto-optic or electro-optic actuators⁷⁰, which can be integrated with the thin film LN platform.

Fig. 6: Integrated LN photonic circuits for entanglement-based dual comb source generation.

Two classical comb sources with slightly different comb line spacings are generated by pumping a CW laser into coupled microresonators at normal group velocity dispersion regimes. The comb on the bottom side is split into two beam paths: one beam serves as the local comb oscillator while the other beam is sent to a PPLN waveguide for second harmonic generation, after which the SHG beam is used to generate entangled sidebands around each local comb line using a second dispersion engineered PPLN waveguide. The entangled sideband beam is then combined with the classical signal comb, together serving as the quantum signal comb output. The darker green pulse corresponds to the frequency comb at the fundamental frequency while the red pulse corresponds to the frequency comb as second harmonic frequency and the lighter green pulse illustrates the sideband engineered state at the fundamental frequency. SHG, second harmonic generation; PPLN WG, periodically poled lithium niobate waveguide; CW continuous wave.

Squeezing and two-mode squeezing have been routinely demonstrated in experiments. For example, 3 dB of squeezing has been demonstrated with bulk optics in radio-frequency sensing squeezing⁴⁰. In the on-chip system, 1.1⁷¹, 1.6⁷², and 2.1 dB⁷³ of squeezing have recently been demonstrated, slightly lower than bulk optics due to extra coupling from chip to fiber. From these results, we expect a near-term demonstration of the proposed quantum dual-comb system to provide 3 dB of quantum advantage, with further efforts to improve device imperfections. The major limitation in the measured squeezing in these systems is the loss from the source to the detection, which can be improved with engineering. For example, with ingenious system design and engineering, ref. ⁷⁴ is able to achieve 15 dB of squeezing in a bulk optical platform.

Note that the proposed quantum comb in this work is different from ref. ⁷²: there the comb lines are themselves pair-wise entangled to serve as a resource for quantum computation; while our proposed quantum comb has a side-band of each comb line pair-wise entangled to benefit dual-comb spectroscopy sensing precision.

In this work, we proposed an entanglement-enhanced dual-comb spectroscopy protocol, where both the signal comb and local comb can be quantum-engineered. When the local comb is stronger than the signal comb, the protocol promises signal-to-noise ratio advantages in detecting low-loss samples such as thin slices of bio-tissues and molecular gas. When the local comb is weak compared with a signal comb, quantum engineering of the local comb provides signal-to-noise ratio advantages regardless of the sample loss, making the quantum advantage robust against experimental imperfections.

Dual-comb interferometry is evolving into one of the most powerful tools for broadband laser spectroscopy, ranging, and imaging, and our work extends its advantages beyond the standard quantum limit. Such a potential will enable comb-based spectroscopy and metrology with ultrahigh precision and sensitivity. The boost in SNR will directly lead to orders of magnitude improvement in measurement speeds for real-time sensing and imaging.

Methods

Two-mode squeezed vacuum

In this paper, we are interested in the two-mode squeezed vacuum (TMSV) state, the continuous-variable version of the Einstein–Podolsky–Rosen (EPR) state⁷⁵. Consider two modes ${\hat{a}}_{1}$ and ${\hat{a}}_{2}$, with real and imaginary quadrature operators ${\hat{q}}_{j}\equiv ({\hat{a}}_{j}+{\hat{a}}_{j}^{{\dagger} })/\sqrt{2},{\hat{p}}_{j}\equiv ({\hat{a}}_{j}-{\hat{a}}_{j}^{{\dagger} })/\sqrt{2}i,j=1,2$. The entanglement between the two modes is described by the joint squeezing on modes ${\hat{a}}_{+}=({\hat{a}}_{1}+{\hat{a}}_{2})/\sqrt{2}$ and ${\hat{a}}_{-}=({\hat{a}}_{1}-{\hat{a}}_{2})/\sqrt{2}$, such that the variances of joint quadratures ${\hat{q}}_{+},{\hat{p}}_{-}$ are suppressed to e^−2r/2 ≡ 1/2G. On the other hand, the variances of ${\hat{q}}_{-},{\hat{p}}_{+}$ are amplified to e^2r/2 ≡ G/2. In this case, when r → ∞ it has the ideal EPR correlation of ${\hat{p}}_{1}={\hat{p}}_{2},{\hat{q}}_{1}=-{\hat{q}}_{2}$.

Take the phase-matched case (α_n = β_n) of Eq. (9) and index ${\hat{a}}_{n,m}\to {\hat{a}}_{1}$ and ${\hat{a}}_{n,-m}\to {\hat{a}}_{2}$, the operator ${\hat{X}}_{n,m}\equiv {\hat{a}}_{1}+{\hat{a}}_{2}^{{\dagger} }={\hat{q}}_{+}+i{\hat{p}}_{-}.$ Therefore, the TMSV state enables $\,{{\mbox{var}}}\,{\hat{X}}_{n,m}$ to be suppressed to e^−2r ≡ 1/G. This is the keystone of the SNR improvement proposed in this paper. For the general case, under the same definitions

$$\begin{array}{l}{\hat{X}}_{n,m}\,=\,\cos ({\alpha }_{n}-{\beta }_{n}){\hat{q}}_{+}-\sin ({\alpha }_{n}-{\beta }_{n}){\hat{p}}_{+}\\ \qquad\quad \,+\,i\cos ({\alpha }_{n}-{\beta }_{n}){\hat{p}}_{-}+i\sin ({\alpha }_{n}-{\beta }_{n}){\hat{q}}_{-},\end{array}$$

(15)

from Eq. (9). Then, one can obtain the variance in Eq. (12).

Derivation of the quantum noise in dual-comb spectroscopy

Here we derive the mean and noise formulas Eqs. (7) and (8) of the random readout of the AC field ${\hat{N}}_{{{{\rm{AC}}}}}(t)={\overline{N}}_{{{{\rm{AC}}}}}(t)+{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(t)$.

After the channel described by Eq. (5), the dual-comb channel inputs Eqs. (2) and (3) become

$$\begin{array}{rcl}{\hat{A}}^{{\prime} }(t)&=&\frac{1}{\sqrt{T}}\left[\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\hat{a}}_{n,k}^{{\prime} }{{\rm {e}}}^{i2\pi (n{f}_{{{{\rm{r}}}}}+k{{\Delta }}{f}_{{{{\rm{r}}}}})t}+\mathop{\sum }\limits_{n=1}^{N}\sqrt{{\kappa }_{n}}{A}_{n}{\rm {{e}}}^{i{\alpha }_{n}}{\rm {{e}}}^{i2\pi n({f}_{{{{\rm{r}}}}}+{{\Delta }}{f}_{{{{\rm{r}}}}})t}\right]\equiv \frac{1}{\sqrt{T}}\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\hat{A}}_{n,k}^{{\prime} }{\rm {{e}}}^{i2\pi (n{f}_{{{{\rm{r}}}}}+k{{\Delta }}{f}_{{{{\rm{r}}}}})t},\\ {\hat{B}}^{{\prime} }(t)&=&\frac{1}{\sqrt{T}}\left[\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\hat{b}}_{n,k}^{{\prime} }{\rm {{e}}}^{i2\pi (n{f}_{{{{\rm{r}}}}}+k{{\Delta }}{f}_{{{{\rm{r}}}}})t}+\mathop{\sum }\limits_{n=1}^{N}\sqrt{{\eta }_{n}}{B}_{n}{\rm {{e}}}^{i{\beta }_{n}}{\rm {{e}}}^{i2\pi n{f}_{{{{\rm{r}}}}}t}\right]\equiv \frac{1}{\sqrt{T}}\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\hat{B}}_{n,k}^{{\prime} }{\rm {{e}}}^{i2\pi (n{f}_{{{{\rm{r}}}}}+k{{\Delta }}{f}_{{{{\rm{r}}}}})t},\end{array}$$

(16)

where we have defined the frequency components at f = nf_r + kΔf_r for the returned signal and LO as

$$\begin{array}{l}{\hat{A}}_{n,k}^{{\prime} }\,=\,{\hat{a}}_{n,k}^{{\prime} }+{\delta }_{n,k}\sqrt{{\kappa }_{n}}{A}_{n}{\rm {{e}}}^{i{\alpha }_{n}},\\ {\hat{B}}_{n,k}^{{\prime} }\,=\,{\hat{b}}_{n,k}^{{\prime} }+{\delta }_{0,k}\sqrt{{\eta }_{n}}{B}_{n}{\rm {{e}}}^{i{\beta }_{n}}.\end{array}$$

(17)

Here δ_n,k = 1 for n = k, zero otherwise and the channel output noises ${\hat{a}}_{n,k}^{{\prime} },{\hat{b}}_{n,k}^{{\prime} }$ are defined by Eq. (5). Note that the noises are zero-mean: $\langle {\hat{a}}_{n,k}^{{\prime} }\rangle =\langle {\hat{b}}_{n,k}^{{\prime} }\rangle =0$.

As we explained in the main text, the receiver collects the two combs together and interferes with the signal comb and the LO comb via a balanced beamsplitter to obtain the fields ${\hat{c}}_{\pm }(t)$. One then performs photodetection to obtain the difference current

$$\hat{N}(t)={\hat{A}}^{{\prime} {\dagger} }(t){\hat{B}}^{{\prime} }(t)+{\hat{B}}^{{\prime} {\dagger} }(t){\hat{A}}^{{\prime} }(t),$$

(18)

which consists of both the mean field and the noise. In the difference current $\hat{N}(t)$, we obtain beat notes at frequencies equal to the frequency differences between each frequency component pair in signal and LO including cross-tooth high-frequency terms of frequency f ~ O(f_r). Filtering out the DC term and high-frequency terms of O(f_r), we keep only f ~ Δf_r ≪ f_r frequency components

$${\hat{N}}_{{{{\rm{AC}}}}}(t)=\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum}\limits_{k\ne {k}^{{\prime} }}{\hat{B}}_{n,k}^{{\prime} {\dagger} }{\hat{A}}_{n,{k}^{{\prime} }}^{{\prime} }{\rm {{e}}}^{i2\pi ({k}^{{\prime} }-k){{\Delta }}{f}_{{{{\rm{r}}}}}t}+{{{\rm{h.c.}}}}$$

(19)

By taking the Fourier transform (omitting the delta-function-like envelops that describe the comb linewidth and finite integration time), we obtain the spectrum at mΔf_r for m ≥ 1,

$${\hat{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\hat{B}}_{n,k}^{{\prime} {\dagger} }{\hat{A}}_{n,k+m}^{{\prime} }+{\hat{B}}_{n,k}^{{\prime} }{\hat{A}}_{n,k-m}^{{\prime} {\dagger} }.$$

(20)

Any negative frequency component of m ≤ −1 is fully determined by its positive frequency component as ${\hat{N}}_{{{{\rm{AC}}}}}(t)$ is real and ${\hat{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})={\hat{N}}_{{{{\rm{AC}}}}}^{{\dagger} }(-m{{\Delta }}{f}_{{{{\rm{r}}}}})$, so it is sufficient to measure the m ≥ 1 positive frequency components only.

For the mean, after plugging in Eq. (17) we have

$$\begin{array}{l}\langle {\hat{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\rangle \\ =\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}\left[{\delta }_{0,k}\sqrt{{\eta }_{n}}{B}_{n}{\rm {{e}}}^{-i{\beta }_{n}}{\delta }_{n,k+m}\sqrt{{\kappa }_{n}}{A}_{n}{e}^{i{\alpha }_{n}}\right.\\ \qquad \qquad \quad \left.+\,{\delta }_{0,k}\sqrt{{\eta }_{n}}{B}_{n}{\rm {{e}}}^{i{\beta }_{n}}{\delta }_{n,k-m}\sqrt{{\kappa }_{n}}{A}_{n}{e}^{-i{\alpha }_{n}}\right]\end{array}$$

(21)

$$=\sqrt{{\eta }_{m}}{B}_{m}{\rm {{e}}}^{-i{\beta }_{m}}\sqrt{{\kappa }_{m}}{A}_{m}{\rm {{e}}}^{i{\alpha }_{m}},$$

(22)

where the second term of Eq. (21) vanishes since we collect only m ≥ 1 terms. This recovers Eq. (7).

For the noise, we obtain the leading order contribution by replacing one annihilation operator in the quadratic terms of Eq. (20) with its mean (which contains a strong displacement) and the other with its weak noise part, i.e.,

$$\begin{array}{l}{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\\ \simeq \mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}\langle {\hat{B}}_{n,k}^{{\prime} {\dagger} }\rangle {\hat{a}}_{n,k+m}^{{\prime} }+{\hat{b}}_{n,k}^{{\prime} {\dagger} }\langle {\hat{A}}_{n,k+m}^{{\prime} }\rangle \\ \quad +\,\langle {\hat{B}}_{n,k}^{{\prime} }\rangle {\hat{a}}_{n,k-m}^{{\prime} {\dagger} }+{\hat{b}}_{n,k}^{{\prime} }\langle {\hat{A}}_{n,k-m}^{{\prime} {\dagger} }\rangle \\ =\mathop{\sum }\limits_{n=1}^{N}\mathop{\sum }\limits_{k=-N}^{N}{\delta }_{0,k}\sqrt{{\eta }_{n}}{B}_{n}\left({\rm {{e}}}^{-i{\beta }_{n}}{\hat{a}}_{n,k+m}^{{\prime} }+{\rm {{e}}}^{i{\beta }_{n}}{\hat{a}}_{n,k-m}^{{\prime} {\dagger} }\right)\\ \quad +\,\sqrt{{\kappa }_{n}}{A}_{n}\left({\delta }_{n,k+m}{{\rm {e}}}^{i{\alpha }_{n}}{\hat{b}}_{n,k}^{{\prime} {\dagger} }+{\delta }_{n,k-m}{\rm {{e}}}^{-i{\alpha }_{n}}{\hat{b}}_{n,k}^{{\prime} }\right)\\ =\mathop{\sum }\limits_{n=1}^{N}\sqrt{{\eta }_{n}}{B}_{n}\left({\rm {{e}}}^{-i{\beta }_{n}}{\hat{a}}_{n,m}^{{\prime} }+{\rm {{e}}}^{i{\beta }_{n}}{\hat{a}}_{n,-m}^{{\prime} {\dagger} }\right)\\ \quad +\,\sqrt{{\kappa }_{n}}{A}_{n}\left({\rm {{e}}}^{i{\alpha }_{n}}{\hat{b}}_{n,n-m}^{{\prime} {\dagger} }+{\rm {{e}}}^{-i{\alpha }_{n}}{\hat{b}}_{n,n+m}^{{\prime} }\right)\end{array}$$

(23)

After plugging the noise input–output relation Eq. (5) in, this recovers Eq. (8).

Justification of the SNR definition

Here, we construct the estimator for κ_n, θ_n and show that the estimation error is connected to our definition of SNR in Eq. (10). Consider the real quadratures ${\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\equiv {\rm{Re}}\, {\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})+ {\rm{Re}}\, {\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}}),{\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\equiv {\rm{Im}}\, {\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})+ {\rm{Im}}\, {\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})$ of the complex heterodyne readout N_AC [see Eqs. (7) and (8)]. Their means are given by Eq. (7) as

$$\begin{array}{l}\langle \hat{{q}_{{{{\rm{AC}}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\rangle \,=\,\sqrt{{\eta }_{m}{\kappa }_{m}}{B}_{m}{A}_{m}\cos ({\alpha }_{m}-{\beta }_{m}),\\ \langle {\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})\rangle \,=\,\sqrt{{\eta }_{m}{\kappa }_{m}}{B}_{m}{A}_{m}\sin ({\alpha }_{m}-{\beta }_{m}).\end{array}$$

(24)

The quadrature fluctuations of Eq. (8) are sums of contributions from N comb lines

$$\begin{array}{rcl}\,{{\rm{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})&=&\mathop{\sum }\limits_{n=1}^{N}\left[\frac{{{{{\mathcal{N}}}}}_{n}}{2}+{\eta }_{n}{\kappa }_{n}\left({B}_{n}^{2}{{{\rm{var}}}}\,{\rm{Re}}\, {\hat{X}}_{n,m}+{A}_{n}^{2}{{{\rm{var}}}}\,{\rm{Re}}\, {\hat{Q}}_{n,m}\right)\right],\\ \,{{\rm{var}}}\,\,{\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})&=&\mathop{\sum }\limits_{n=1}^{N}\left[\frac{{{{{\mathcal{N}}}}}_{n}}{2}+{\eta }_{n}{\kappa }_{n}\left({B}_{n}^{2}{{{\rm{var}}}}\,{\rm{Im}}\, {\hat{X}}_{n,m}+{A}_{n}^{2}{{{\rm{var}}}}\,{\rm{Im}}\, {\hat{Q}}_{n,m}\right)\right].\end{array}$$

(25)

The sum of the two noises gives Eq. (11). From Eq. (15) for the TMSV state, we can see that ${\rm{Re}}\, {\hat{X}}_{n,m},{\rm{Im}}\, {\hat{X}}_{n,m},{\rm{Re}}\, {\hat{Q}}_{n,m},{\rm{Im}}\, {\hat{Q}}_{n,m}$ are mutually independent when α_n = β_n.

We begin with the estimation of κ_n’s, assuming perfect phase matching, α_n = β_n. In this case, two-mode squeezing has $\,{{\mbox{var}}}\,\,{\rm{Re}}\, {\hat{X}}_{n,m}=\,{{\mbox{var}}}\,\,{\rm{Im}}\, {\hat{X}}_{n,m}=\frac{1}{2}\,{{\mbox{var}}}\,\,{\hat{X}}_{n,m},\,{{\mbox{var}}}\,\,{\rm{Re}}\, {\hat{Q}}_{n,m}=\,{{\mbox{var}}}\,\,{\rm{Im}}\, {\hat{Q}}_{n,m}=\frac{1}{2}\,{{\mbox{var}}}\,\,{\hat{Q}}_{n,m}$, thus

$$\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\,{{\mbox{var}}}\,\,{\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\frac{1}{2}\,{{\mbox{var}}}\,\,{N}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}}).$$

(26)

Also, note that ${\hat{q}}_{{{{\rm{AC}}}}}$ and ${\hat{p}}_{{{{\rm{AC}}}}}$ commute, and indeed they are mutually independent Gaussian variables [which can be verified from Eq. (15)]. Thus, we can define the distribution of the readouts q, p for ${\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}}),{\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})$ as P_q(q)P_p(p). Then the minimum mean square error (MMSE) for unknown parameters κ_n is given by the Cramér–Rao lower bound of Gaussian distribution:

$$\begin{array}{l}{[{{\mbox{MMSE}}}\sqrt{{\tilde{\kappa }}_{m}}]}^{-1}\,\equiv \,{\left[\int{{\mbox{}}}{\rm {d}}q{\left(\frac{\partial \log {P}_{q}(q)}{\partial \sqrt{{\kappa }_{m}}}\right)}^{2}{P}_{q}(q)\right]}^{-1}\\ \qquad\qquad\qquad\qquad +\,{\left[\int{{\mbox{}}}{\rm {d}}p{\left(\frac{\partial \log {P}_{p}(p)}{\partial \sqrt{{\kappa }_{m}}}\right)}^{2}{P}_{p}(p)\right]}^{-1}\\ \quad\qquad\qquad\qquad\,=\,\frac{{\left\vert \frac{d{\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}{d\sqrt{{\kappa }_{m}}}\right\vert }^{2}}{\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}+2\cdot \frac{{\left\vert \frac{d{{\rm{var}}}{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}{d\sqrt{{\kappa }_{m}}}\right\vert }^{2}}{2{[{{\mbox{var}}}\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})]}^{2}}\\ \quad\qquad\qquad\qquad\,\simeq \,\frac{{\left\vert \frac{d{\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}{d\sqrt{{\kappa }_{m}}}\right\vert }^{2}}{\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}=\frac{{\eta }_{m}{B}_{m}^{2}{A}_{m}^{2}}{\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}.\end{array}$$

(27)

In the last (approximate) equality, we have assumed that the modulation on the readout variance is negligible, which is true due to squeezing power much lower than the comb power. Thus, in comparison with Eq. (10),

$${[{{\mbox{MMSE}}}\sqrt{{\tilde{\kappa }}_{m}}]}^{-1}=2/{\kappa }_{m}\cdot {{{\mbox{SNR}}}}^{2}.$$

(28)

Now we estimate the phase mismatch θ_m ≡ α_m−β_m. Note that ${\left\vert \frac{{\rm {d}}{\overline{N}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}{{\rm {d}}{\theta }_{m}}\right\vert }^{2}=2\cdot \kappa \eta {A}^{2}{B}^{2}$. Similar to the transmissivity estimation above, we assume a TMSV input state such that $\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\,{{\mbox{var}}}\,\,{\hat{p}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})$ and the independence between ${\hat{q}}_{{{{\rm{AC}}}}}$ and ${\hat{p}}_{{{{\rm{AC}}}}}$ still holds. Similarly, the Cramér–Rao lower bound gives

$$\begin{array}{r}{[{{\mbox{MMSE}}}\,{\theta }_{m}]}^{-1}=\frac{{\eta }_{m}{\kappa }_{m}{B}_{m}^{2}{A}_{m}^{2}}{\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})}\end{array}$$

(29)

Note that $\,{{\mbox{var}}}\,\,{\hat{q}}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})=\frac{1}{2}\,{{\mbox{var}}}\,\,{N}_{{{{\rm{AC}}}}}(m{{\Delta }}{f}_{{{{\rm{r}}}}})$. Thus, in comparison with Eq. (10),

$${[{{\mbox{MMSE}}}\,{\theta }_{m}]}^{-1}=2\cdot {{{\mbox{SNR}}}}^{2}.$$

(30)

Full formula of estimation error

Here we derive Eq. (13), including the inverse-square-law, the inverse-law, and constant noise terms with respect to source power P_S. We denote the total power of the signal or the LO as P_S or P_LO, and define their ratio γ ≡ P_LO/P_S. We can identify ${{TP}}_{\rm{S}}=h{\nu}_{0}\mathop{\sum }\nolimits_{m = 1}^{N}| {A}_{m}{| }^{2},{{TP}}_{{{{\rm{LO}}}}}=h{\nu }_{0}\mathop{\sum }\nolimits_{m = 1}^{N}| {B}_{m}{| }^{2}$, where hν₀ is the energy per photon. To simplify the formulas we assume symmetric comb lines A_m = A, B_m = B for any 1 ≤ m ≤ N. To connect to the SNR of Eq. (10), we normalize each noise by the power of mean field Eq. (7), ${\overline{N}}_{{{{\rm{AC}}}}}{(m{{\Delta }}{f}_{{{{\rm{r}}}}})}^{2}={\eta }_{m}{\kappa }_{m}{A}^{2}{B}^{2}={\eta }_{m}{\kappa }_{m}\cdot ({P}_{{{{\rm{S}}}}}T/Nh{\nu }_{0})\cdot ({P}_{{{{\rm{LO}}}}}T/Nh{\nu }_{0})$.

The inverse-square-law noise comes from the detector noise. It is modeled as a fluctuation of constant noise equivalent power (NEP) on the photon current readout, including dark current, Johnson noise, amplifier noise figure, etc. Note that NEP-type noise is circular symmetric in the phase space, it affects both real and imaginary parts of N_AC. Since we defined var N_AC as the sum of the two quadrature noises, the physical NEP-type noise is 2NEP²Δf in power. In photon number, the noise is $2{({{\mbox{NEP}}}\cdot T/Nh{\nu }_{0})}^{2}{{\Delta }}f$. Here the bandwidth is defined as Δf = 1/2T for single-sided NEP spectral density. According to Eq. (10), the detector noise results in the normalized NEP-type noise at intermediate frequency mΔf_r

$${\sigma }_{{{{\rm{NEP}}}}}^{2}=2\frac{{({{\mbox{NEP}}}\cdot T/Nh{\nu }_{0})}^{2}{{\Delta }}f}{{\eta }_{m}{\kappa }_{m}({P}_{{{{\rm{LO}}}}}T/Nh{\nu }_{0})({P}_{{{{\rm{S}}}}}T/Nh{\nu }_{0})}=\frac{{N}^{2}}{T}\frac{{{{\mbox{NEP}}}}^{2}}{{\eta }_{m}{\kappa }_{m}\gamma {P}_{{{{\rm{S}}}}}^{2}}.$$

(31)

Here NEP has the unit of W Hz^−1/2. Now we see that the NEP-type noise-power relation is ${\sigma }_{{{{\rm{NEP}}}}}^{2}\propto \frac{1}{{P}_{{{{\rm{S}}}}}\cdot {P}_{{{{\rm{LO}}}}}}$, which is an inverse-square-law term $\sim O(\frac{1}{{P}_{{{{\rm{S}}}}}^{2}})$ when γ is fixed.

The inverse-law noise comes from the intrinsic quantum noise Eq. (8). For the case where phase noise is negligible, the quadrature fluctuation leads to the normalized quadrature noise

$${\sigma }_{{{{\rm{quad}}}}}^{2}=\frac{\,{{\mbox{var}}}\,\,{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}}{{\eta }_{m}{\kappa }_{m}{B}^{2}({P}_{{{{\rm{S}}}}}T/Nh{\nu }_{0})}=\frac{{N}^{2}}{T}{c}_{\gamma }\frac{4h{\nu }_{0}}{{P}_{{{{\rm{S}}}}}},$$

(32)

where

$${c}_{\gamma }\equiv \frac{\,{{\mbox{var}}}\,\,{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}}{4{\eta }_{m}{\kappa }_{m}{B}^{2}\cdot N}$$

(33)

with $\,{{\mbox{var}}}\,\,{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}$ defined in Eq. (11). For vacuum input G = 1, G_LO = 1, unit transmissivities η_n = κ_n = 1 and zero noise ${\tilde{{{{\mathcal{E}}}}}}_{n}=0$ for any 1 ≤ n ≤ N, we have ${c}_{\gamma }=\frac{1}{4}(1+\frac{1}{\gamma })$. In this case, $\,{{\mbox{var}}}\,\,{\hat{{{\Sigma }}}}_{{{{\rm{AC}}}}}\propto {P}_{{{{\rm{S}}}}}+{P}_{{{{\rm{LO}}}}}$, we see that the fundamental noise–power relation is ${\sigma }_{{{{\rm{quad}}}}}^{2}\propto \frac{{P}_{{{{\rm{S}}}}}+{P}_{{{{\rm{LO}}}}}}{{P}_{{{{\rm{S}}}}}\cdot {P}_{{{{\rm{LO}}}}}}$, which is an inverse-law term ~O(1/P_S) when γ is fixed.

The constant noise results from the relative intensity noise (RIN). Consider the power of each comb tooth, hν₀A²/T for signal and hν₀B²/T for LO. RIN is modeled as a single-sided white noise $\,{{\mbox{var}}}\,\,(h{\nu }_{0}{A}^{2}/T)=(\,{{\mbox{RIN}}}\,{{\Delta }}f){P}_{{{{\rm{S}}}}}^{2},\,{{\mbox{var}}}\,\,(h{\nu }_{0}{B}^{2}/T)=(\,{{\mbox{RIN}}}\,{{\Delta }}f){P}_{{{{\rm{LO}}}}}^{2}$ on the power of field amplitudes A, B, generated from additive amplified spontaneous emission (ASE) from the laser or from any subsequent optical amplification. Similar to NEP-type noise, the bandwidth is defined as Δf = 1/2T for single-sided RIN spectral density. RIN-type noise is also circular symmetric in the phase space, which affects both real and imaginary parts of N_AC and yields a factor of 2 in the complex-observable variance var N_AC. Consequentially, the physical noise is

$$\begin{array}{l}2{\eta }_{m}{\kappa }_{m}[\,{{\mbox{var}}}\,\,(A)\cdot {B}^{2}+{A}^{2}\cdot \,{{\mbox{var}}}\,\,(B)]\\ \qquad\,=\,2{\eta }_{m}{\kappa }_{m}\left[\frac{\,{{\mbox{var}}}\,\,({A}^{2})}{4{P}_{{{{\rm{S}}}}}}\cdot {P}_{{{{\rm{LO}}}}}+{P}_{{{{\rm{S}}}}}\cdot \frac{\,{{\mbox{var}}}\,\,({B}^{2})}{4{P}_{{{{\rm{LO}}}}}}\right]\\ \qquad\,=\,2{\eta }_{m}{\kappa }_{m}{\left(\frac{T}{h{\nu }_{0}}\right)}^{2}\cdot \,{{\mbox{RIN}}}\,{{\Delta }}f\cdot {P}_{{{{\rm{S}}}}}^{2}\left(\frac{1}{4}\cdot \gamma +\frac{{\gamma }^{2}}{4\gamma }\right)\\ \qquad\,=\,{\eta }_{m}{\kappa }_{m}{\left(\frac{T}{h{\nu }_{0}}\right)}^{2}\cdot \,{{\mbox{RIN}}}\,{{\Delta }}f\cdot \gamma {P}_{{{{\rm{S}}}}}^{2}\end{array}$$

(34)

Here RIN is in unit 1 Hz⁻¹, and in the first equality we have used $\,{{\mbox{var}}}\,\,({A}^{2})={(\frac{\partial }{\partial A}{A}^{2})}^{2}\,{{\mbox{var}}}\,\,(A)=4{A}^{2}\,{{\mbox{var}}}\,\,(A),\,{{\mbox{var}}}\,\,({B}^{2})=4{B}^{2}\,{{\mbox{var}}}\,\,(B)$, assuming A² ≫ var A and B² ≫ var B. It is valid to consider fluctuations on amplitudes A, B instead because the RIN physically comes from the amplified spontaneous emission, which is modeled as a Gaussian noise on the field quadratures in quantum optics. The amplitude fluctuations result in the normalized RIN-type noise

$${\sigma }_{{{{\rm{RIN}}}}}^{2}=\frac{{\eta }_{m}{\kappa }_{m}\,{{\mbox{RIN}}}\,{{\Delta }}f\cdot \gamma {({P}_{{{{\rm{S}}}}}T/h{\nu }_{0})}^{2}}{{\eta }_{m}{\kappa }_{m}({P}_{{{{\rm{LO}}}}}T/Nh{\nu }_{0})({P}_{{{{\rm{S}}}}}T/Nh{\nu }_{0})}=\frac{{N}^{2}}{T}2{c}_{{\gamma }^{2}}\,{{\mbox{RIN}}}\,.$$

(35)

Here the coefficient ${c}_{{\gamma }^{2}}\equiv 1/4$. We immediately see that the RIN-type noise-power relation is ${\sigma }_{{{{\rm{RIN}}}}}^{2} \sim O(1)$, which does not depend on the signal or LO power.

Overall, we can write the full formula of the SNR Eq. (10) at intermediate frequency mΔf_r as ${{{\mbox{SNR}}}}^{-2}={\sigma }_{{{{\rm{NEP}}}}}^{2}+{\sigma }_{{{{\rm{quad}}}}}^{2}+{\sigma }_{{{{\rm{RIN}}}}}^{2}$, which gives Eq. (13) in the main text.

To recover the classical results in ref. ³⁰, we further assume that the frequency spectra of all parameters are almost uniform, such that κ_m ≈ κ, η_m ≈ η for any 1 ≤ m ≤ N. When η, κ → 1 and classical source are used (G = G_LO = 1), our result recovers the formula of σ_H [Eq. (2)] in ref. ³⁰ mostly [H(f) is the transfer function of the electrical field, equivalent to $\sqrt{\kappa (f)}$]. Note that in the quantum model, the shot noise term a_shot is instead formulated as the quadrature fluctuation a_quad, while the quadrature fluctuation incidentally gives a similar result ${c}_{\gamma }=\frac{1}{4}(1+\frac{1}{\gamma })$, up to an extra 1/2 factor. Also, ${c}_{{\gamma }^{2}}=\frac{1}{4}\cdot \frac{2\gamma }{2\gamma }$ is different from (1 + γ²)/2γ in Eq. (2) of ref. ³⁰, which is the result of our derivation of the physical RIN-type noise Eq. (34) in the balanced detection. The independence of ${c}_{{\gamma }^{2}}$ and thereby of ${\sigma }_{{{{\rm{RIN}}}}}^{2}$ in our Eq. (13) agrees with the recent results of Eqs. (58) and (59) in ref. ⁵⁶ assuming white noise spectrum. Meanwhile, our derivation recovers the result ${c}_{{\gamma }^{2}}^{{{{\rm{unbal}}}}}=\frac{1+{\gamma }^{2}}{2\gamma }$ in ref. ³⁰ for the unbalanced detection case where the physical noise is ~var(A²) + var(B²) instead. To summarize, our result when applied to the classical dual-comb SNR can be obtained by letting ${c}_{\gamma }=\frac{1}{4}(1+\frac{1}{\gamma }),b=1,{c}_{{\gamma }^{2}}=1/4$ in Eq. (2) of ref. ³⁰.

We note that our quantum model yields an SNR-gamma relation different from the semiclassical model in ref. ³⁰. For example, at a fixed P_S there is a finite-optimum value of γ to maximize the SNR in the semiclassical model while the optimum is at γ → ∞ in our quantum model. This is because when γ is large, the RIN-type noise increases with γ in the semiclassical model, while RIN-type noise remains a constant in our derivation which agrees with ref. ⁵⁶.

Finally, we address the saturation in the SNR with respect to the squeezing gain G, due to NEP or RIN noise. For simplicity, we consider κ_m ≈ κ, η_m ≈ η for any 1 ≤ m ≤ N. The saturation of SNR begins when the squeezing gain G is large enough such that NEP-type or RIN-type noise overwhelms the fundamental noise. By solving ${\sigma }_{{{{\rm{NEP}}}}}^{2}\, > \,{\sigma }_{{{{\rm{quad}}}}}^{2}$ or ${\sigma }_{{{{\rm{RIN}}}}}^{2}\, >\, {\sigma }_{{{{\rm{quad}}}}}^{2}$, we derive the saturation threshold as ${P}_{{{{\rm{S,sat}}}}}^{{{{\rm{NEP}}}}}=G\cdot {{{\mbox{NEP}}}}^{2}/\{h{\nu }_{0}(\gamma [G(1-\kappa )+\kappa ]+\kappa )\}$ for the NEP-type and ${P}_{{{{\rm{S,sat}}}}}^{{{{\rm{RIN}}}}}=h{\nu }_{0}(\gamma (G(1-\kappa )+\kappa )+\kappa )/\{G\gamma \kappa \cdot \,{{\mbox{RIN}}}\,\}$ for the RIN-type. The above thresholds are indicated in Fig. 1 as the dots on the SNR curves, showing a good agreement with numerical results.