Electric-field metrology of a terahertz frequency comb using Rydberg atoms

TL;DR

This paper demonstrates a novel Rydberg atom-based detector for precise, sensitive, and room-temperature calibration of terahertz frequency combs, advancing terahertz spectroscopy into the quantum regime.

Contribution

It introduces a new single-photon detection method using Rydberg atoms for calibrating and detecting terahertz frequency combs with high sensitivity and selectivity.

Findings

Achieved electric field sensitivity down to 45.2 fW/cm²/Hz^0.5

Calibrated comb intensity against atomic properties

Enabled high-precision terahertz spectroscopy in a room-temperature setup

Abstract

Terahertz radiation finds an increasing number of applications, yet efficient generation and detection remain a challenge and an active area of research. In particular, the precise detection of weak and narrowband terahertz signals is notoriously difficult. Here, we employ a novel type of single-photon detector based on Rydberg atoms to both detect and calibrate a terahertz frequency comb over an octave-spanning range, yet with a MHz-level selectivity. We calibrate the intensity of the electric field of the comb against the fundamental atomic properties, while achieving the intensity (power) sensitivity down to 45.2 $\mathrm{fW\ cm^{-2}\ Hz^{-0.5}}$ (1.84 $\mathrm{fW\ Hz^{-0.5}}$) within a single mode of the frequency comb, all in a room-temperature operated setup. Our results elucidate the transition of terahertz frequency combs into the quantum regime, enabling high-precision and high-sensitivity spectroscopy. This breakthrough allows terahertz science to better leverage revolutionary techniques developed for optical frequency combs.

Tables

Table 1. Table 1: Numerical parameters and their uncertainties, used to obtain the frequency comb mode number, tabularized for clarity.

Parameter name	Value and uncertainty
ARC transition frequency	$125.836\pm 0.004$
Repetition rate	$79.4902\pm 0.0020$
Thermal peak center	$4.92\pm 0.05$
Main comb peak center	$-3.21\pm 0.04$
Reference frequency	$125.831\pm 0.004$
Main comb peak frequency	$125.828\pm 0.004$
Mode number	$1582.93\pm 0.06$
Measured transition frequency	$125.8411\pm 0.0032$

Table 2. Table 2: Numerical results with their uncertainties obtained when deriving the sensitivity and noise equivalent power (NEP) of the detector, tabularized for clarity.

Parameter name	Value and uncertainty
Transition frequency	$125.836\pm 0.004$
Transition wavelength	$\mathrm{2.382\hspace{0.17em}41}\pm \mathrm{0.000\hspace{0.17em}08}$
Small antenna aperture	$\mathrm{4.065\hspace{0.17em}04}\pm \mathrm{0.000\hspace{0.17em}26}$
Thermal photon counts	$\left(9.87\pm 0.25\right)\text{×}{10}^3\mathrm{phot}/s$
Thermal field	$36.7\pm 0.7$
Thermal intensity	$1.78\pm 0.07$
Thermal power	$72.5\pm 2.7$
Sensitivity	$45.2\pm 1.7$
NEP	$1.84\pm 0.07$

Equations

$$E=\frac{h}{|d|}\Omega ,$$

$${\mathcal{N}}_{\text{comb}}({\Delta }_{1200})\propto \sum _n{I}_n|\mathcal{G}({\Delta }_{\text{THz}}+n{f}_{\text{rep}},{\Delta }_{1200}){|}^2,$$

$$A=g\frac{{\lambda }^2}{4\pi },$$

$${\mathcal{N}}_{\text{comb}}(t)\propto {\left|\sum _n{A}_n{e}^{-i{\omega }_{n}t}\mathcal{G}({\Delta }_{\text{THz}}+n{f}_{\text{rep}},{\Delta }_{1200})\right|}^2,$$

$${\mathcal{N}}_{\text{beat}}(t)\propto \sum _{n\ne m}{A}_n{A}_m^{\ast }{e}^{-i{\delta }_{nm}t}{\mathcal{G}}_n{\mathcal{G}}_m^{\ast },$$

$${\mathcal{N}}_{\text{beat}}(t)\propto \sum e^{-i2\pi n{f}_{\text{rep}}t}\sum _m{A}_{m+n}{A}_m^{\ast }{\mathcal{G}}_{m+n}{\mathcal{G}}_m^{\ast }.$$

$$\begin{array}{rl}\dot{\rho }(t)=& -\frac{i}{\hslash }\left[H,\rho (t)\right]+{\mathcal{L}}_{\Gamma }\left[\rho (t)\right],\end{array}$$

$$H=-\frac{\hslash }{2}\left(\begin{array}{ccccc}0& {\Omega }_{780}^{\ast }& 0& 0& 0\\ {\Omega }_{780}& i{\Gamma }_1+2{\Delta }_{\text{I}}& {\Omega }_{480}^{\ast }& 0& 0\\ 0& {\Omega }_{480}& i{\Gamma }_2+2{\Delta }_{\text{II}}& {\Omega }_{\text{THz}}^{\ast }& 0\\ 0& 0& {\Omega }_{\text{THz}}& i{\Gamma }_3+2{\Delta }_{\text{III}}& {\Omega }_{1200}^{\ast }\\ 0& 0& 0& {\Omega }_{1200}& i{\Gamma }_4+2{\Delta }_{\text{IV}}\end{array}\right)-\frac{i\hslash {\Gamma }_{\text{tr}}}{2}\mathbb{I},$$

$$\begin{array}{rl}{\Delta }_{\text{I}}& ={\Delta }_{780},\\ {\Delta }_{\text{II}}& ={\Delta }_{780}+{\Delta }_{480},\\ {\Delta }_{\text{III}}& ={\Delta }_{780}+{\Delta }_{480}+{\Delta }_{\text{THz}},\\ {\Delta }_{\text{IV}}& ={\Delta }_{780}+{\Delta }_{480}+{\Delta }_{\text{THz}}-{\Delta }_{1200}.\end{array}$$

$$K=\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ 0& -k_{780}& 0& 0& 0\\ 0& 0& -k_{780}+k_{480}& 0& 0\\ 0& 0& 0& -k_{780}+k_{480}+k_{\text{THz}}& 0\\ 0& 0& 0& 0& -k_{780}+k_{480}+k_{\text{THz}}+k_{1200}\end{array}\right).$$

$$H_v=H+Kv,$$

$$\rho =\int w(v)\rho (v)dv.$$

$$\mathcal{N}\propto {\left|⟨1|\rho |4⟩\right|}^2.$$

$$|\mathcal{G}{|}^2=\frac{\mathcal{N}}{I_{\text{THz}}}.$$

$$\mathcal{N}\propto |\mathcal{G}{|}^2\propto |⟨1|\rho |4⟩{|}^2.$$

Full text

I Introduction

THz (Terahertz) and millimeter-wave radiation occupy the spectral region between the microwave and infrared domains. Historically, this range was one of the least explored parts of the electromagnetic spectrum due to the lack of efficient sources and detectors. However, rapid progress in photonics, semiconductor technologies, and materials science has brought THz technologies to the forefront of modern research. Today, THz radiation is finding growing applications in diverse areas, such as high-speed wireless communications [1], security imaging, non-destructive testing, spectroscopy [2, 3], and advanced material characterization, e.g., study of exciton dynamics and phase transitions in semiconductors [4, 5].

Despite this progress, key challenges remain, particularly in generating and detecting coherent THz signals with high sensitivity and resolution. Pulsed THz sources that generate electromagnetic frequency combs have emerged as powerful tools for addressing these challenges. Frequency combs—well known in the optical domain—enable high-precision spectroscopy and frequency metrology and are now being extended into the THz range [6, 7, 8, 9, 10, 11, 12]. Developments in THz comb generation include techniques based on photomixers driven by continuous-wave lasers, passive mode-locking [13], and referencing to optical frequency combs for improved accuracy. These combs have enabled hyperspectral imaging [14], phase-resolved measurements [15], and studies of phase noise in THz systems [16].

On the detection side, significant progress has been made in fabricating fast, multi-pixel THz detectors, yet fundamental limits in sensitivity—especially at low photon flux—continue to hinder applications in quantum sensing and metrology. This has motivated the exploration of unconventional detection platforms, such as atomic vapor-based sensors.

In particular, Rydberg-atom ensembles have emerged as promising candidates for THz detection due to their extreme sensitivity to electric fields, stemming from their large polarizabilities and strong dipole transitions. These sensors rely on well-understood principles of atomic physics and offer a range of functionalities: fluorescence-based imaging of THz and mmWave fields [17, 18, 19, 20, 21, 22, 23], direct field reception [24, 25], and even THz generation schemes via nonlinear optical processes in Rydberg media [26]. The versatility of these systems has also been demonstrated in demodulation of broadband signals [27], Doppler-free spectroscopy using multi-photon transitions [28], and long-distance sensing in open-air environments [29]. Polarization-insensitive detection schemes [30], access to high-angular-momentum states [31], and applications such as mmWave thermometry in cold-atom and ambient settings [32] have further expanded the scope of Rydberg sensors. Moreover, recent efforts have explored their utility in the calibration of commercial radar systems [33]. A comprehensive review of these developments can be found in [34].

In this paper, we report the observation and analysis of a THz frequency comb using a Rydberg-atom-based detector. Our experimental approach combines two key methods: Autler-Townes splitting measurements, which allow for precise determination of the THz electric field strength in absolute units [35], and a microwave-to-optical conversion technique [36], which extends the sensitivity down to the level of thermal radiation. Using these tools, we characterize the frequency comb by identifying the amplitudes and absolute mode numbers of its individual spectral components.

Furthermore, we explore the interaction of the frequency comb with different Rydberg transitions, providing both experimental measurements and comparisons with theoretical models. The conversion scheme enables us to probe a broad spectral range and observe the structure of the comb around various atomic resonances. We also present direct measurements of beat-note signals between adjacent comb lines, revealing visible high-frequency components and enabling insight into the underlying coherence and stability of the comb structure.

Our results demonstrate the potential of Rydberg-atom detectors in THz frequency metrology and highlight their applicability in emerging quantum-enabled sensing technologies.

II Experimental setup

II.1 Detection

In this work, a detection setup in two configurations is employed: the first one is the standard A-T (Autler-Townes) splitting measurement used for calibration to SI units, while the second is an improvement over microwave-to-optical conversion, enabling ultrasensitive detection of higher RF (radio frequencies).

The first auxiliary method relies upon a 87Rb vapor cell and two optical fields. The laser beams are introduced into the cell from opposing directions and overlaid colinearly. This setup is the typical two-photon excitation scheme utilized in Rydberg atom electrometry, where usually two lasers are used to excite the atoms. The $780\text{\,}$ probe beam tuned to the D2 line is scanned in frequency and, after interaction with atomic vapors is shone upon an avalanche photodiode (Thorlabs APD430A) to obtain the probe transmission spectrum. By the observation of the EIT (electromagnetically induced transparency) effect, introduced by the counter-propagating $483\text{\,}$ coupling laser tuned to the $5^{2}\text{P}_{3/2}(\text{F}=3)\rightarrow 27^{2}\text{D}_{5/2}$ transition, the setup is capable of performing RF electric field measurements by analyzing the Autler-Townes splitting in the domain of the probe laser detuning [35].

This measurement scheme has a fundamental limitation on the minimal resolvable field, resulting from the finite EIT linewidth, caused by the movement of atoms at room temperature and the arising Doppler broadening. With the objective of measuring much weaker fields, the A-T setup is then extended by introducing a third laser beam, namely the 1278 nm decoupling field tuned to the $26^{2}\text{F}_{7/2}\rightarrow 5^{2}\text{D}_{5/2}(\text{F}=4)$ transition, realizing the full setup, visible in Fig. 1 a).

This constitutes the second and main measurement method, utilizing parametric RF-to-optical upconversion via six-wave mixing [36]. By transferring the atomic population down from the Rydberg state, it allows for a radiative, 776 nm transition between $5^{2}\text{D}_{5/2}(\text{F}=4)$ and $5^{2}\text{P}_{3/2}(\text{F}=3)$ states to take place, resulting in the emission of optical photons, driven by the incident THz fields. The conversion output beam is then separated using a dichroic mirror and a spectral filter, before being sent to a single photon counter in the form of SPAPD (single photon avalanche photodiode), Excelitas SPCM-AQRH. Finally, the SPAPD photon counts are registered with an externally triggerable IDQ 900 time controller.

Such a detection scheme achieves superior detection sensitivity to the intensity of THz fields, down to the thermal noise limit, with the potential for time-resolved measurements and synchronization with a pulsed source. Thanks to having a strong reference field source, it is possible to calibrate the photon count readout to the electric field amplitude measured from A-T splitting. What is more, as is the case with A-T splitting, we are not limited to resonant RF detection and addressing different frequencies is possible by detuning the coupling or decoupling fields. For simplicity, this work explores tuning just the decoupling laser for mapping the bandwidth of a particular RF transition.

Both detection methods can be carried out in a consecutive manner, with minimal effort required for switching between them. The coupling and decoupling laser frequencies are stabilized to a narrowband frequency-doubled fiber laser (NKT Photonics, $1560\text{\,}$), which is itself stabilized to a Rb cell via a modulation-transfer lock. The conversion scheme requires that the probe laser is also stabilized to the frequency-doubled laser via an optical phase-loop, which also drives the frequency scans in A-T measurements. The beams within the cell are collimated to have the matching Gaussian waists of $250\text{\,}$ with probe, coupling and decoupling beam powers of $140\text{\,}$, $600\text{\,}$, and $126\text{\,}$ respectively. To optimize the conversion process, the vapor cell is heated to $50\text{\,}$ using a custom made heater utilizing hot air to avoid interference with metal elements.

II.2 Signal generation

The THz field is introduced from a high-power source, namely an automotive radar chip TRA_120_045 from Indie Semiconductor. While A-T measurements may be performed with the signal entering the vapor cell from the side, the parametric conversion necessitates phase-matching, so the THz field is combined colinearly with optical fields and focused using a parabolic mirror with a through hole for optical beams. This enables colinear propagation of all fields in the vapor cell, and the setup is further improved using custom-made metamaterials, which allow for beam collimation, attenuation, and polarization control [33]. Both lenses and waveplates are made using an FDM (fused deposition modeling) 3D printer and HIPS (high impact polystyrene) filament, while the polarizers are custom-made flex PCBs (printed circuit boards). All of the designs are fine-tuned for the operating transition frequency using FDTD (finite difference time domain) simulations. As indicated in the Fig. 1, the TRA source is enclosed in an RF absorbing box with a small aperture facing the first lens. The enclosure pre-attenuates the source, minimizes reflections from it and ensures that only the radiation passing through the aperture couples to the atoms.

As such, a sufficiently powerful source, like the aforementioned radar chip, is employed for calibration and connecting the results obtained in both detection methods. Subsequently, we can switch to conversion detection and study lower intensity field sources, in particular, a THz frequency comb. We therefore introduce a PCA (photoconductive antenna – BATOP iPCA-21-05-300-800-h), which, given a bias voltage, emits THz pulses following the intensity characteristic of a femtosecond Spectra-Physics MaiTai laser, which has a repetition rate $f_{\text{rep}}$ of $79.49\text{\,}$. The PCA has a hyperhemispherical silicon lens mounted on top, which partially focuses the output. It is then collimated and sent into the setup with the same path, and a metal sheet is used as a mirror to toggle between the high-power chip and the locally weaker PCA.

All of the fields, including the THz beams, reach the vapor cell with circular polarizations, which drive transitions between the atomic levels of the highest possible total angular momenta, exhibiting larger dipole moments.

III Results

III.1 Calibration

III.1.1 A-T splitting calibration

Due to its low intensity, the THz frequency comb can be observed only via the conversion detection, and as such, its strength can be expressed only in relative units of photon counts per second. However, by performing the absolute calibration of the receiver and tying together A-T splitting and photon counts regimes, the strength of the low intensity field, down to thermal level, can be depicted in terms of electric field amplitude.

We begin the calibration by measuring the A-T splitting in the probe transmission spectrum. For this measurement, the decoupling laser is covered, and only the effects from the coupling laser and THz field can be seen. The probe transmission spectrum is then measured for attenuations of the THz field between 0 and $1.5\text{\,}$ introduced by two attenuators, each consisting of a half-waveplate and a polarizer. The A-T splitting, created by the THz field, can then be directly observed in the probe transmission spectra as the splitting of the EIT caused by the $483\text{\,}$ laser. Next, we determine the Rabi frequencies $\Omega$ for the corresponding splitting by fitting the numerical solutions of the steady state to the measured probe transmission spectra [37]. Subsequently, knowing the dipole moment $d$ of the considered THz transition, the strength of the incoming field can be expressed in the amplitude of the electric field $E$ by the relation

[TABLE]

where $h$ is Planck’s constant. In the experiment, the dipole moment of the main considered atomic transition, $27^{2}\text{D}_{5/2}\rightarrow 26^{2}\text{F}_{7/2}$, is equal to $d=594a_{0}e$ [38], where $a_{0}$ is the Bohr radius and $e$ is the value of the elemental charge. The calibration of the A-T splitting regime is then determined by fitting a linear function with a fixed slope, equal to $-1$, to the calculated electric field amplitudes as a function of THz field attenuation.

III.1.2 Photon counts calibration

To perform the measurements of converted photon counts, the probe laser is locked at the transition resonant frequency, the decoupling laser is turned on, and both the coupling and decoupling lasers are tuned to maximize the conversion signal coming from the comb tooth closest to the resonance. Then, the photon count measurements of the converted signal are performed for the attenuation range spanning from 0 to 32 dB via the setup of HWPs and polarizers. To further extend the calibration range, the measurements are done for two TRA source power levels, which we will label as: $I_{1}$ used for the A-T splitting measurements, and $I_{2}$, for which the A-T splitting cannot be measured in the probe transmission spectrum. The $I_{1}$ power according to the manufacturer’s documentation is equal to $I_{1}=$-7\text{,}\text{dBm}$$. However, the real power entering the setup is significantly lower due to the pre-attenuation caused by the small aperture in the enclosure. As the conversion signal due to the $I_{1}$ power saturates the SPAPD, for those measurement series, an additional neutral density (ND) filter is added before the photon counter. Finally, we turn off the TRA chip and measure the photon counts corresponding to the frequency comb tooth closest to the atomic resonance. Additionally, we cover both THz sources and measure the conversion signal corresponding to the thermal radiation, and then cover the decoupling laser, which allows for the measurement of the dark counts.

Similarly to the calibration of the A-T splitting regime, the line function with a fixed slope is fitted to the data for both source power levels. Points lying near the thermal noise and saturation level, defined as a point where the response of the receiver becomes nonlinear, are omitted to increase the fit accuracy. Hence, the linear function is fitted only to points lying between $11.5\text{\,}$ (saturation level) and $45\text{\,}$ (near the thermal level) THz signal attenuation. Both of those threshold points have been determined by finding the point for which additional omission of data points does not impact the fit accuracy. As both measurements are performed for the same attenuation range, the ratio $I_{2}/I_{1}$ can be determined precisely using fitted intercept line parameters, with their values being $I_{2}/I_{1}=$25.48\pm 0.07\text{,} which gives $I_{2}=$-32.48\pm 0.07\text{\,}\text{dBm}. This way, a full response curve spanning from atomic saturation to thermal radiation level is created.

Finally, the scaling factor between the electric field amplitude and measured photon counts is calculated by using the fitted intercept line parameters yielding the absolute calibration of the receiver. The measured points, scaled using the aforementioned steps together with the absolute calibration line, the thermal noise, and the comb tooth level, are shown in Fig. 2, with the square and circle plot markers representing the A-T splitting and photon counts measurement accordingly. The colored regions shown in the figure correspond to the measurement regions for both TRA source powers, with vertical axes representing the number of measured photon counts per second in both cases, and horizontal axes representing the THz attenuation values in both measurements. Additionally, the noticeably small amount of measured data points for the A-T splitting measurement can be attributed to the comparably low sensitivity of this measurement technique. Based on the absolute calibration, we then find that the electric field amplitude of the measured THz frequency comb tooth is equal to $E_{\text{c}}=$405\pm 8\text{,} and the amplitude of thermal noise integrated in the whole conversion bandwidth is equal to $E_{\text{noise}}=$36.7\pm 0.7\text{\,}. The calculated values, as well as the entire calibration, are determined without subtracting the dark counts of the detector. It can also be noted based on Fig. 2 that the response to the frequency comb tooth lies in the detector linear response regime, which we define as a linear dependence of the conversion signal on the power of the oncoming THz field. Finally, to extend the context of the measured values, we may derive rough estimates of the total power of the measured fields. Assuming that the THz field ($\lambda=$2.38\text{,}) is focused down to a diffraction-limited area $A=\pi(1.22\lambda N)^{2}\approx 27\ \mathrm{mm}^{2}$ by the parabolic mirror with f-number $N=1$ we get a lower-bounded power for the pre-attenuated TRA (corresponding to $E_{\mathrm{TRA}}=$5\text{\times}{10}^{4}\text{\,}) $P_{\mathrm{TRA}}\approx-30\ \mathrm{dBm}$ and for a single tooth of the comb $P_{c}\approx-72.3\ \mathrm{dBm}$.

III.2 Tunable broadband readout of frequency comb signal

Next, we visualize the THz frequency comb by scanning the decoupling laser, thus tuning the detector to different teeth (frequency modes) of the comb. During the scanning process, both probe and coupling lasers are locked at the optimal frequencies. For the considered transition, we perform a scan over the range of $3\text{\,}$ with a step of $4\text{\,}$ around the atomic resonance. The measurement is performed by detuning the decoupling laser over the frequency range mentioned above and measuring the number of converted photons per second for each of the detunings. As can be seen in Fig. 3 a), the distance between two subsequent peaks is the same as the MaiTai repetition rate $f_{\text{rep}}$ and remains constant across the scan range, with a regime of higher than unity SNR (signal-to-noise ratio) spanning over more than 2 GHz. Additionally, scans without the comb are performed to observe the intensity of converted thermal radiation.

This bandwidth corresponds to a single energy transition. We are, however, able to map different regions of the frequency comb spectrum. We tune both the coupling and the decoupling lasers to different atomic transitions, allowing us to detect the THz comb around different frequencies. Using the ARC library, several candidates have been selected and can be seen in Fig. 3 c)-e). Namely, we detect the THz field at the frequencies of $160\text{\,}$, $182\text{\,}$, and $242\text{\,}$, corresponding to the $25^{2}\text{D}_{5/2}\rightarrow 24^{2}\text{F}_{7/2}$, $24^{2}\text{D}_{5/2}\rightarrow 23^{2}\text{F}_{7/2}$, and $26^{2}\text{D}_{5/2}\rightarrow 28^{2}\text{F}_{7/2}$ atomic transitions respectively. All of them are depicted in Fig. 1 b). Similarly to the previous scan, we perform the measurements by detuning the decoupling laser every $4\text{\,}$ step over the frequency range of $1\text{\,}$, and measure the photon counts per second for each step. Noticeably, the maximum converted signal drops for higher THz frequencies. We attribute this behavior to the strength of the PCA antenna being dependent on the emitted frequency. Using a standard THz detection scheme of time-domain spectroscopy, the spectral shape of an individual pulse is measured to verify the spectral power distribution, visible in Fig. 3 b). While the resolution is limited to $\sim$10\text{,}\mathrm{G}\mathrm{H}\mathrm{z}$$, it confirms an order of magnitude difference in the signal strength. The lower dipole moments also contribute to lower measured photon counts, however the variations are not significant, ranging from $594\text{\,}\mathrm{a}_{\mathrm{0}}\mathrm{e}$ for $125\text{\,}$ to $353\text{\,}\mathrm{a}_{\mathrm{0}}\mathrm{e}$ for $242\text{\,}$ as calculated by ARC [38].

The behavior of the comb frequency scans shown in Fig. 3 can be further understood by comparing them with the numerical predictions. Under the assumption of the linear response, justified by the absolute calibration of the receiver, the expected photon counts as a function of the decoupling laser detuning $\Delta_{1200}$, for the measurement time much longer than $f_{\text{rep}}$ is given as:

[TABLE]

where $I_{n}$ is the intensity of each tooth and $\Delta_{\text{THz}}$ corresponds to the frequency distance of the closest tooth from the transition resonance. The conversion gain $\mathcal{G}$ is the ratio between the amplitude of the THz and converted fields [39, 40], in this case proportional to the atomic coherence between $5^{2}\text{D}_{5/2}$ and $5^{2}\text{P}_{3/2}(\text{F}=3)$ states. As the rest of the lasers remain at the constant frequency, we only consider the dependence of $\mathcal{G}$ on the detuning of the THz field and the decoupling laser. Based on the linear assumption, we then fit the theoretical model to the measured comb frequency scan shown in Fig. 3 a) with the thermal photons subtracted. The theoretical model consists of the steady-state solutions for each of the teeth considered, which are then combined according to Eq. (2) to recover the full scan. The measured comb scan with the fitted theoretical predictions can be seen in Fig. 4. It can be noted that the theoretical model and its envelope match the measurements well, proving the applicability of the detector’s linear response regime. The divergence of the model, visible for the higher decoupling laser detunings, is a result of the influence of other excited states, namely the hyperfine structure of the $5^{2}\text{D}_{5/2}$ state, which has not been taken into account in the simulations.

III.3 Deriving the sensitivity and thermal power

In order to assess the robustness of the THz sensor and perform a quantitative comparison with other platforms, we evaluate its sensitivity and noise equivalent power (NEP). The measured thermal radiation includes other noise sources, such as photon counter dark counts, and as such is considered as the noise floor.

The derivations are performed for the main transition under consideration, i.e., the $27^{2}\text{D}_{5/2}\rightarrow 26^{2}\text{F}_{7/2}$ transition at $125.836\text{\,}$. We measure the system noise with photon counts, using acquisition times of $t=$1\text{,}, and measuring for $30\text{\,}$ in total. We get the mean counts of $\mathcal{N}_{\text{noise}}=$(9.87\pm 0.25)\text{\times}{10}^{3}\text{\,} photons per second. We denote the standard deviation of this statistic as $\sigma_{\text{noise}}=$250\text{,} photons per second. Noticeably, the statistics are super-Poissonian (for Poissonian $\sigma_{\text{noise}}=\sqrt{\mathcal{N}_{\text{noise}}}=$99\text{\,} photons per second should follow), which is expected for thermal light, but for single shot measurements of multimode light at a measurement time far exceeding the bandwidth inverse, the quantum phenomena of thermal bunching should cancel out. What is instead observed are the macroscopic and classical gain fluctuations of the setup, which can be calculated to constitute a relative error of about $1.5\text{\,}$. Shorter measurements are more likely to approach the shot noise (Poissonian) limit, as well as mitigate these macro effects; however, this is not a priority of this work.

Following the calibration stage, presented in Fig. 2, the converted photon counts can be translated into the equivalent incident electric field of $E_{\text{noise}}=$36.7\pm 0.7\text{,} and further into intensity through $I_{\text{noise}}=E_{\text{noise}}^{2}/2Z_{0}=$1.78\pm 0.07\text{\,}, where $Z_{0}\approx$377\text{,}\text{Ohm}$$ is the impedance of free space. Next, to obtain NEP, the receiver’s effective aperture needs to be considered. For this, we can consider the gain-defined small antenna aperture $A$ as:

[TABLE]

where $g$ is the antenna gain (directivity) and $\lambda=$2.38\text{,} is the wavelength of received field. We obtain the antenna gain $g=9$ from field phase-matching considerations (i.e., the receiver reception pattern) [[36](https://arxiv.org/html/2508.20698v2#bib.bib36)], which leads to the result of $A=$4.07\text{\,}\mathrm{m}\mathrm{m}^{2}. Consequently, the effective registered power corresponding to the thermal radiation equal to $P_{\text{noise}}=I_{\text{noise}}A=$72.5\pm 2.7\text{,}$$ can be retrieved.

Finally, the parameters of interest can be calculated. Power sensitivity $\xi$ is derived as the thermal intensity $I_{\text{noise}}$ scaled by the ratio of thermal photon count fluctuations $\sigma_{\text{noise}}$ to the mean counts $n_{\text{noise}}$, and acquisition time $t$, and calculated as $\xi=I_{\text{noise}}\frac{\sigma_{\text{noise}}}{\mathcal{N}_{\text{noise}}\sqrt{t}}=$45.2\pm 1.7\text{,}. Consequently, NEP equals to this sensitivity multiplied by $A$, yielding a value of $\text{NEP}=\xi A=$1.84\pm 0.07\text{\,}. Note that the aperture used in the above calculation is significantly smaller than the aperture used for estimating the power of the TRA and comb sources, and thus the derived NEP value corresponds to the best case scenario - i.e. perfect matching of the source’s radiation mode to the sensor. The full values and uncertainties considered in this derivation are presented in Tab. 2 in Appendix E.

III.4 Deriving the frequency comb mode number

Energy state levels are usually calculated using the ARC library, which, in the case of the evaluated RF transitions in rubidium-87, relies on experimental data of finite precision. On top of that, there are numerous additional effects, such as Stark or Zeeman shifts, which become even more cumbersome with temperature-dependent optical density variations. Therefore we propose an alternative method to estimate the detected frequency by determining the mode (tooth) number of frequency comb.

Observing particular teeth and performing frequency scans happens through detuning from the resonant case, which happens in relative units. By estimating the mode number, we can move to absolute frequencies, due to the information embedded in the comb peak separation. Knowing that the mode number is an integer, it is possible to infer the transition frequency, which, for some transitions, can be further verified using the radar chip. As it is stabilized using a phase-locked loop, it has a known frequency and can confirm the validity of the derivations.

We have determined the MaiTai repetition rate $f_{\text{rep}}$ to be $79.490\,20\pm 0.000\,01\text{\,}$ using Agilent N9010A spectrum analyzer. In the case of the main transition under consideration, an ARC value of $125.836\pm 0.004\text{\,}$ is calculated, with the uncertainty obtained by redoing the calculations and propagating the uncertainty of quantum defects [41, 42]. The next step involves locating the working point, at which the scan visible in Fig. 3 a) is performed. The thermal profile is equivalent to an integral over the THz frequencies, so by fitting a Lorentzian profile, the resonant $\Delta_{1200}$ value can be found. This particular model shows excellent agreement with the measured thermal photon data, indicating that the resonance is at $\Delta_{1200}=$4.92\pm 0.05\text{,}$$.

The theoretical considerations, presented with the broadband scan in Fig. 4 are then used to yield the position of the tooth closest to the atomic resonance. It is located at the detuning of $-3.213\pm 0.039\text{\,}$, which in absolute units gives $125.828\pm 0.004\text{\,}$. Dividing by $f_{\text{rep}}$ results in $1582.934\pm 0.064\text{\,}$ passing the $3\sigma$ test with respect to the nearest integer. Given the aforementioned lack of consideration for the shifts in energy levels, this result gives a strong indication that the mode number of this particular tooth is equal to 1583. Reverting the above derivations, taking the optical ruler as the ground truth, we arrive at a transition frequency of $125.8411\pm 0.0032\text{\,}$.

Given the small but still relevant deviation from the value provided by ARC, an additional verification is performed using the radar chip source. As it is controlled by a phase-locked loop with a known reference frequency, we are able to assess the frequency of the comb tooth to be located at a target frequency of $125.829\pm 0.001\text{\,}$ corresponding to a mode number estimate of $1582.95\pm 0.04\text{\,}$, confirming the previous result. The parameters used in the derivation of tooth mode number are presented in Tab. 1.

III.5 Beating of adjacent frequency comb peaks

Further evidence of frequency comb detection can be presented by measuring the beat notes between individual peaks. As they are separated by multiples of $f_{\text{rep}}$, the beat note frequencies will also assume values of $nf_{\text{rep}}$. The beating manifests itself in the time domain, so we observe the time-resolved photon conversion response on the timetagger, connected to the single photon counter. The 40 MHz MaiTai internal photodiode readout is used as the external trigger of the timetagger, allowing for measurement of all the potential beat-note frequencies.

As with the frequency scans, detuning the decoupling laser allows adjusting the sensitivity to different frequencies and a value that sits between two comb peaks are selected so as to make the beatings more pronounced in comparison to the resonant case.

The photons are counted for $60\text{\,}$. The observed signals, as seen in Fig. 5 are modulated mainly with the first-order beating contribution, but performing a Fourier transform reveals subsequent multiples of $f_{\text{rep}}$ as well. As detuning the decoupling laser decreases the conversion efficiency, the collected signals are mean-normalized and presented in Fig. 5 a). Positioning our detector frequency sensitivity at a detuning corresponding to $f_{\text{rep}}/2$ corresponds to significantly more pronounced modulation, clearly visible in the temporal signal, while the Fourier transform in Fig. 5 b)-c) reveals the existence of higher-order beat frequencies. It can also be noted that the spectra differ visibly in the contribution of the constant, non-oscillating term, with the spectrum for the detuned case having significantly higher contrast of the oscillating terms. This observation is consistent with the comparison of the signals in the time domain shown in Fig. 5 a).

The dynamics of the beating of adjacent comb teeth can be understood by once again investigating the function describing the expected converted signal, which we write this time as:

[TABLE]

where $A_{n}$ is the amplitude of the $n$-th frequency component, and $\omega_{n}=2\pi nf_{\text{rep}}$. For the convenience of notation, we denote the conversion gain $\mathcal{G}(\Delta_{\text{THz}}+nf_{\text{rep}},\Delta_{1200})$ as $\mathcal{G}_{n}$. The equation above can be split into two parts, namely the non-oscillating part, described by Eq. (2), and the time-dependent part given as:

[TABLE]

with $\delta_{mn}=2\pi(n-m)f_{\text{rep}}$ being the beating frequency between $m$-th and $n$-th comb components. It is worth noting that for the signal averaged over a time longer than the period of the slowest oscillating term, all the oscillating terms vanish, leading to the expected converted signal given by Eq. (2). The terms visible can be further grouped by the value of $\delta_{nm}$:

[TABLE]

As can be seen in the equation above, the beating only occurs at frequencies that are multiples of $f_{\text{rep}}$. Moreover, as all of the terms in the mentioned equation are time-dependent, they do not contribute to the constant term of the measured beat note signal, and thus the difference in the beating contrast can be attributed mostly to the part of the conversion signal described by Eq. (2).

IV Discussion

In this paper, we have demonstrated the detection of a THz frequency comb using a microwave-to-optical converter based around a warm Rydberg atomic ensemble. We further characterized the frequency comb by deriving its electric field amplitude and performing wide frequency scans of the conversion signals. Finally, we measured the beating between adjacent frequency components. We have also verified the linear behavior of the converter in the wide frequency range by comparing the frequency scans and beat note signal to the theoretical model.

The obtained THz intensity detection NEP values show an improvement of 2 to 4 orders of magnitude over state of the art room temperature alternatives [43, 44, 45, 46, 47, 48, 49], showcasing great potential for Rydberg atoms as an alternative platform in this field. Similar conclusions can be drawn from other attempts at using a Rydberg atom-based converter for THz detection [50], where we also achieved an improvement of about 4 orders of magnitude in sensitivity. Direct comparison of NEP is difficult in this case, due to different assumptions in aperture derivation.

The are multiple issues to be solved that are beyond the scope of this work. Namely, due to the THz frequency comb being weak, determining its electric field amplitude requires an additional, stronger source. This issue, in theory, can be overcome by increasing the sensitivity of the A-T splitting detection, which can be achieved by considering the Doppler canceling scheme [28], leading to a narrower EIT peak and the possibility of direct measurement of the A-T splitting caused by the frequency comb.

We believe that our results pave the road for future applications of Rydberg atom-based detectors in the rapidly growing field of THz electrometry. The unique possibility of determining the strength of weak fields in absolute units, as well as overall high sensitivity, makes Rydberg atom-based setups a great candidate for ultra-sensitive and easy-to-calibrate THz detectors. Additionally, the setup can be easily extended to perform multiple tasks, such as THz imaging or spectroscopy.

Acknowledgments

We thank M. Lipka for his initial contribution to the construction of the pulsed THz source. This research was funded in whole or in part by the National Science Centre, Poland, grants No. 2021/43/D/ST2/03114 and No. 2024/53/N/ST7/02730. The "Quantum Optical Technologies" (FENG.02.01-IP.05-0017/23) project is carried out within the Measure 2.1 International Research Agendas programme of the Foundation for Polish Science, co-financed by the European Union under the European Funds for Smart Economy 2021-2027 (FENG). S.B. is a recipient of the Foundation for Polish Science START 2025 scholarship.

Disclosures

The authors declare no conflict of interest.

Data availability

Data underlying the results presented in this paper are available in the Ref. [51].

Code availability

The codes used for the numerical simulation are available from W.K. upon request.

Appendix A: Calibration of the half-waveplates’ transmission

To determine the attenuation imposed by the HWP combined with the polarizer, the transmission of such a pair needs to be calibrated as a function of the HWP rotation angle, for which we use the photon counts measurement as a function of angles of both HWP for the $I_{2}$ TRA source intensity. For the sake of calibration of each of the individual waveplates, we choose the cross-sections of measured photon counts through the rotation angle of the other waveplate corresponding to its measured maximum transmission. As the HWPs are not fully rigid and could move by about $\pm 0.5^{\circ}$ on their own, the uncertainty of the transmission functions is taken as a difference between the interpolated function and the interpolated function shifted by $0.5^{\circ}$. The interpolated transmission functions together with the data points and uncertainties can be seen in Fig. 6. Additionally, based on the found transmission functions, maximum attenuations of both pairs can be estimated as $A_{1}=$19.6\pm 1.7\text{,} and $A_{2}=$23.6\pm 0.4\text{\,} for the first and second pair, respectively.

Appendix B: Theoretical model

The experimental setup consists of a five-level atomic open quantum system, which is described by the GKSL (Gorini–Kossakowski–Sudarshan–Lindblad) equation:

[TABLE]

where $H$ is the Hamiltonian of the system, $\rho$ is the density matrix, and $\mathcal{L}_{\Gamma}$ is a Lindblad superoperator, which accounts for decoherence mechanisms. We model a five-level atomic ladder in the RWA (rotating-wave approximation), which produces the following form of the Hamiltonian in dipolar approximation:

[TABLE]

where $\Omega_{ij}$ are Rabi frequencies driving transitions between states $\ket{i}\rightarrow\ket{j}$, $\Gamma_{i}$ is the decay rate of state $\ket{i}$ and $\Gamma_{\text{tr}}$ is a parameter modeling the effective broadening coming from limited interaction time of atoms with thin laser beams. Finally, the detunings of the incident electric fields from the atomic resonances are grouped into common terms:

[TABLE]

Room temperature atomic vapors are significantly influenced by the Doppler effect, stemming from the atoms moving along the electric field propagation axis. Factoring in the broadening of the resonance can be simply implemented as additional diagonal values of the Hamiltonian. For that, we define the $K$ matrix, whose diagonal elements represent cumulative sums of the wave vectors $k$ of the fields required to excite the system to the corresponding energy levels, with each subsequent level involving an additional driving field. The signs are picked to reflect the counter propagation of the probe beam, yielding:

[TABLE]

Then we can write our full Hamiltonian as:

[TABLE]

where $v$ denotes the velocity class of atoms, and obtain the final solution by averaging over the Maxwell-Boltzmann velocity distribution, represented by a statistical weight of $w(v)$:

[TABLE]

The obtained Hamiltonian is then evaluated numerically.

Modeling of Autler-Townes splitting

The Autler–Townes effect is observed as a splitting and broadening of the EIT in the absorption spectrum of the probe laser. Probe laser absorption is proportional to expected values of $\bra{1}\rho(\Delta_{780})\ket{0}$. In the experiment, the measurement of absorption happens after some time from the moment the interaction started, which allows the system to relax into a steady state. In the macroscopic picture, the steady state is time independent, $\dot{\rho}=0$, which simplifies the calculations significantly.

Modelling of terahertz frequency comb

On the other hand, for the analysis of the conversion signal caused by the THz frequency comb, we focused mostly on the conversion gain $\mathcal{G}$ and the conversion signal $\mathcal{N}$. In light of the model described above, the conversion signal is given as:

[TABLE]

To better understand the conversion gain, it is worth first considering the conversion of the monochromatic THz signal. Under the assumption of the detector’s linear response, the conversion gain $\mathcal{G}$ can be interpreted as the ratio between the emitted and THz fields, usually expressed in terms of their intensities as:

[TABLE]

It can be noted that due to the THz field being monochromatic, the value of $I_{\text{THz}}$ is just a scaling factor. Thus, the conversion gain becomes proportional to the conversion signal, and the following relation can be established:

[TABLE]

Appendix C: Determining Rabi frequency from A-T splitting

As mentioned in the main text, for the sake of the calibration of the receiver, the Rabi frequency of the THz field is determined based on the measured A-T splitting visible in the probe transmission spectra. The values of the Rabi frequencies are found by numerically fitting the steady-state solution for the theoretical model described in Appendix B to the probe transmission signal. In such a fit, the Rabi frequency of the THz field is treated as one of the parameters of the fitted function. Exemplary probe transmission spectra, together with the fitted theoretical model, can be seen in Fig. 7. In order to perform the fitting properly, both the experimental and numerical data are min-max normalized. To further increase the accuracy of the numerical fits, the common parameters, which remained constant during the measurement, such as the Rabi frequencies of laser fields, are fitted simultaneously to all the measured spectra.

Appendix D: Selective observation of comb teeth

In order to get a better intuition of the non-standard measurements performed in this paper, we have used the theoretical framework and experimental data to produce a plot of the conversion gain $\mathcal{G}\left(\Delta_{\text{THz}},\Delta_{1200}\right)$ values. It is important to note that the derivation of this value is phenomenological in its nature, as it relies on the min-max normalization of the experimental data and outputs of the theoretical model. The fitting procedure on normalized detuning scans is then performed to retrieve the parameters found in the system’s Hamiltonian, which is subsequently used to calculate the value of $\mathcal{G}$ for arbitrary detunings.

The calculated conversion gain can be seen in Fig. 8, showing two resonances that give a good working intuition for performing comb spectroscopy with this setup. A single photon count measurement, given a frequency comb and a fixed value of $\Delta_{1200}$, can be understood as a discrete sum of conversion gain values along a vertical line of the map, as expressed in Eq. (2). The frequency comb teeth remain static when tuning the decoupling laser, meaning that we adjust the conversion gain to the positions of the teeth. The resonances observed in Fig. 3 are spaced at $f_{\text{rep}}$, because we count two teeth at this interval. On the other hand, the thermal photon measurements correspond to taking an integral across the $y$-axis.

The zoomed-in map in Fig. 9 showcases the asymmetric profile of the conversion gain that could be observed in the performed scans.

Appendix E: Sensor’s thermal power parameters

The sensor’s parameters, used in the considerations of thermal noise level, as well as the sensor’s sensitivity and NEP, are presented in Tab. 2.