Controlled multi-photon subtraction with cascaded Rydberg superatoms as single-photon absorbers

The preparation of light pulses with well-defined quantum properties requires precise control at the individual photon level. Here, we demonstrate exact and controlled multi-photon subtraction from incoming light pulses. We employ a cascaded system of tightly confined cold atom ensembles with strong, collectively enhanced coupling of photons to Rydberg states. The excitation blockade resulting from interactions between Rydberg atoms limits photon absorption to one per ensemble and rapid dephasing of the collective excitation suppresses stimulated re-emission of the photon. We experimentally demonstrate subtraction with up to three absorbers. Furthermore, we present a thorough theoretical analysis of our scheme where we identify weak Raman decay of the long-lived Rydberg state as the main source of infidelity in the subtracted photon number and investigate the performance of the multi-photon subtractor for increasing absorber numbers in the presence of Raman decay.

F uture optical quantum technology relies on precise control over the quantum state of light. Deterministic removal of exactly one, or more generally exactly n sub , photons enables applications such as state-preparation for optical quantum simulation and computing [1][2][3][4][5] or quantum-enhanced metrology 6 . Photon subtraction can also give insight into more fundamental aspects of quantum optics 7 . Heralded single-photon subtraction 8 can be implemented using highly imbalanced beamsplitters 7,9 , but the probabilistic nature limits the scalability of this approach 3,9,10 . Individual absorbers like a single two-level atom in free space seem well-suited for photon subtraction, as they are saturated by just one photon, but this approach is limited by weak atom-photon coupling, stimulated emission and short lifetimes of the saturated state. These problems can be mitigated by enhancing the atom-light coupling using a resonator, and transfer of the absorber to a third, dark state [11][12][13] not coupled to the incoming light as demonstrated with single atoms coupled to a microsphere resonator 14 .
Strong photon-emitter coupling can also be achieved without optical resonators in atomic ensembles with collectively excited and long-lived Rydberg states, also referred to as Rydberg superatoms 15 . Rydberg atoms interact strongly with each other 16 and the resulting excitation blockade 17 can be mapped onto light fields to create strong optical nonlinearities at the single-photon level [18][19][20][21][22] . This has enabled many technical applications such as single-photon sources 23,24 , optical transistors 25,26 , removal of photons from stored light pulses 27,28 , and photon-photon quantum gates 29 with recent efforts to combine these into multinode networks [30][31][32] . Photon subtraction can also be realised using Rydberg superatoms as saturable single-photon absorbers 33,34 combining the blockade, which prevents multi-photon absorption, with rapid dephasing of the superatom into dark collective states to avoid stimulated re-emission.
In this work, we demonstrate a cascaded quantum system of up to three Rydberg superatom absorbers for controlled subtraction of specific photon numbers from a light pulse. In addition to demonstrating controlled multi-photon subtraction, we find that Raman decay is the main source of deviations from the ideal absorber behaviour. This is supported by a detailed theoretical analysis, which also shows that scaling beyond n sub = 3 absorbers with high probability to subtract exactly n sub photons is realistic as long as Raman decay can be suppressed and a sufficiently high single-photon coupling is maintained.

Results
Implementation. Figure 1a illustrates the implementation of individual Rydberg superatoms as a saturable single-photon absorber 33,34 . Weak pulses of probe light at λ p ≈ 780 nm, from which photons are to be subtracted, propagate through a small, optically thick ensemble of optically trapped 87 Rb atoms. In combination with a strong, co-propagating control field at λ c ≈ 480 nm, the probe light couples the atomic ground state jgi ¼ j5S (Fig. 1b). As a result of the Raman detuning Δ/2π ≈ 100 MHz, probe photons are only absorbed by the ensemble if the control field is tuned onto Raman resonance with two-photon detuning δ = 0. Strong van der Waals interactions between Rydberg atoms lead to the blockade effect that suppresses multiple Rydberg excitations for atoms separated by r < r B , where r B is the blockade radius that characterises the volume inside which the energy shift V = C 6 /r 6 defined by the van der Waals coefficient C 6 exceeds the excitation linewidth.
If the radial probe beam waist (1/e 2 -waist radius ≈ 6.5 μm) and the extent of the ensemble along the probe axis constrain the excitation volume to a single Rydberg excitation, the superatom is saturated after absorbing just one photon. Consequently, all N atoms in the excitation volume share the excitation in a collective bright state jWi ¼ ∑ j jg 1 g 2 :::r j :::g N i= ffiffiffiffi N p with strongly enhanced collective coupling ffiffiffiffiffiffiffiffi ffi κR in p from the collective ground state jGi ¼ jg 1 g 2 :::g j ::: [35][36][37][38] . Here, g 0 Ω c /(2Δ) is the effective single-probe-photon-single-atom coupling strength between jgi and jri, where g 0 is the single-probephoton-single-atom coupling strength between jgi and jei, R in the incoming probe photon rate and Ω c the control Rabi frequency. Following absorption of a photon, jWi dephases with rate γ D into a manifold of N − 1 collective dark states fjDig that are orthogonal to jWi and no longer couple to the probe such that stimulated emission is suppressed, while maintaining the blockade. Besides dephasing, the excited collective states are also subject to decay of jri with Raman decay Γ ¼ Γ e Ω 2 c =ð2ΔÞ 2 as the dominant contribution, with Γ e being the natural linewidth of jei. Following adiabatic elimination of jei, the superatom dynamics can be described in terms of just jWi, jGi and jDi, a single dark state into which we condense all collective states in fjDig 15,39 . This effective three-level system will form the foundation of our theoretical analysis.
To implement multi-photon subtraction, we place n sub ensembles along the path of the probe and control fields (Fig. 1a)

Ion detection
Photon detection + + + + + + + + + + + + + + + + + + + Fig. 1 Realisation of up to three cascaded single-photon absorbers using Rydberg superatoms. a To create n sub saturable superatom absorbers, we place n sub ensembles of cold 87 Rb atoms in the path of a tightly focussed probe beam. Using an acousto-optical deflector (AOD), we can control the number and position of the optical traps that tightly confine the ensembles below the Rydberg blockade radius r B along the probe direction. b Within r B strong van der Waals interactions restrict each ensemble to a single Rydberg excitation as the probe photons and a control field couple jgi ¼ j5S Raman scheme with detuning Δ/2π ≈ 100 MHz and thus to the absorption of a single photon at a time for a two-photon detuning of δ = 0. The transmitted probe pulses are coupled into a singlemode optical fibre (not shown) and detected on four single-photon counters in a Hanbury-Brown-Twiss configuration. c Representation of the absorber as an effective three-level system in terms of singly excited collective states following adiabatic elimination of jei. Strong dephasing γ D from the bright excited state jWi, with strong coupling ffiffiffiffiffiffiffiffi κR in p from the ground state jGi, into dark excited states jDi prevents stimulated reemission of the absorbed photon and the absorption of further photons until it is subject to Raman decay Γ ≪ γ D , κ.
at distances > 50 μm ≫ r B such that the superatoms do not blockade each other and act independently. We trap the ensembles at the foci of individual optical trapping beams that intersect perpendicularly with the probe and confine each ensemble within r B along the probe axis. We use an acoustooptical deflector (AOD) to control the position and number of ensembles via the number and frequencies of radio-frequency (RF) signals applied [40][41][42] . A reservoir dipole trap (not shown in Fig. 1a) provides additional radial confinement ('Methods'). The ensembles have temperatures of ≈10 μK with an extent < r B along the probe axis and N~10 4 with exact numbers varying with n sub = 3, 2, 1, respectively, due to variations in the trapping and cooling dynamics during ensemble preparation.
Experimental results. In the following, we experimentally demonstrate controlled subtraction of up to three photons by placing the corresponding number of absorbers n sub in the probe path. We measure the transmission for coherent, Tukey-shaped probe pulses with a pulse length of τ = 2.5 μs (FWHM, with 1.0 μs rise/fall time) and 〈n in 〉 ≤ 40, where 〈n in 〉 is the mean incoming photon number per pulse, using four single-photon counters in a Hanbury-Brown-Twiss configuration located behind a singlemode fibre. This configuration allows to analyse second-order photon correlations by calculating the correlation between all distinct pairs and averaging over the results. Without control field, we measure combined optical depths of the ensembles of ≈11, 16 and 20 for n sub = 1, 2 and 3 respectively, and find a probe transmission of >0.99 which is slightly reduced at finite Δ/2π ≈ 100 MHz due to off-resonant scattering of the probe light and which the data below are corrected for.
First, we investigate the difference between 〈n in 〉 and the mean transmitted photon number 〈n out 〉 (Fig. 2a). For 〈n in 〉 > 10, we observe the expected reduction by n sub , while we subtract fewer photons for 〈n in 〉 < 10. This behaviour is expected, as for low R in ∝ 〈n in 〉, the pulse area ffiffiffiffiffiffiffiffi ffi κR in p τ is insufficient to drive the superatom population predominantly into jWi and jDi. This becomes particularly evident in the shape of the transmitted pulses for R in ≈ 1 μs −1 (top row in Fig. 2b) with transmission well below one at their end, whereas we observe the onset of saturation for R in ≈ 5 μs −1 (centre row). Importantly, the duration to reach saturation increases with n sub , because the driving between jGi and jWi reduces alongside the probe intensity following each absorber. For R in ≈ 10 μs −1 (bottom row), saturation sets in even faster, but we observe a slight oscillation in the subsequent transmission, which reflects the superatom dynamics as the probe drives Rabi oscillations between jGi and jWi with strong damping due to γ D 15 . To suppress superradiant re-emission of absorbed photons in the forward direction after the probe pulse 43,44 , γ D has to be sufficiently strong not only compared to 1/τ, but also compared to the coherent dynamics 15,33 . The dephasing is dominated by atomic motion, which is enhanced by the copropagating probe and control beams compared to a counterpropagating configuration, with additional contributions from elastic scattering of the Rydberg electron by ground-state atoms [45][46][47] and the non-uniform AC-Stark shift induced by the trapping light ('Methods'), which can only be compensated for on average. To characterise the system, we determine κ, γ D and Γ by comparing the observed transmission to the predictions of a model of n sub effective three-level atoms strongly coupled to a chiral waveguide ('Methods'), assuming that κ, γ D and Γ are equal for all absorbers. The results of the model are in good agreement with the experiment for both the subtracted photons ( A closer look at the number of subtracted photons (inset in Fig. 2a) reveals that it slightly exceeds n sub at high 〈n in 〉 indicating that a single absorber may subtract multiple photons. This excess cannot be explained by the deexcitation of the absorbers in the Rabi oscillation cycle as the associated spontaneous emission occurs back into the probe mode in forward direction with rate κ. Instead, it can be attributed to the small, but non-zero Raman decay 0 < Γ < κ, γ D , which leads to spontaneous re-emission in random, rather than the forward direction. Its presence, even if weak, however, leads to reduced fidelity to subtract exactly one photon per absorber as we will discuss in detail in our theoretical analysis.
To further demonstrate manipulation of the quantum state of light at the single-photon level, we also investigate the effect of the subtraction on the correlations between transmitted photons. Figure 2c shows the second-order correlation function g (2) (t 1 , t 2 ) = 〈n(t 1 )n(t 2 )〉/(〈n(t 1 )〉〈n(t 2 )〉) for two of the transmitted pulses shown in Fig. 2b and n sub = 3 alongside the theory prediction for three cascaded three-level absorbers. Following initial antibunching, we observe g (2) (t 2 − t 1 = 0) > 1, i.e. bunching of the The inset shows the mean subtracted photons for the same data. All data are corrected for off-resonant absorption of the probe light and background counts. Error bars correspond to one standard error and are smaller than the markers for most data points. b Temporal profiles of probe light following transmission through n sub = 1 (left), 2 (centre) and 3 (right) absorbers for different mean photon numbers. The incoming pulses recorded in the absence of the superatom absorbers are shown in light grey. Besides the experimental data, we also show the results of the three-level model fitted to the experimental data (dark grey). c. Photon correlations following transmission through three superatom absorbers. Second-order correlation functions g (2) (t 1 , t 2 ) for R in ≈ 5.5 (top) and 10 μs −1 (bottom). Besides the experimental data (left), we also show the predictions of the three-level model (right). transmitted light. This bunching is expected as the subtraction operation reduces the mean photon number by n sub = 3, without reducing the width of the photon number distribution compared to the incoming coherent pulses, thus inducing super-Poissonian statistics 15,48 . It also highlights that subtraction of exactly n sub photons is not to be confused with n sub consecutive applications of annihilation operators, which also reduces 〈n in 〉 of the incoming coherent pulses by up to n sub , but maintains Poisson statistics, or probabilistic subtraction 49 . Furthermore, the interaction with the absorber induces correlations within the transmitted pulse for t 2 ≠ t 1 , the timescale of which can be observed on the anti-diagonal axes in Fig. 2c. The correlations disappear towards the end of the pulses due to two effects: as saturation suppresses photon absorption, the incoming coherent light is transmitted without change and, in addition, absorbers are more likely to have undergone random Raman decay, which introduces Poissonian fluctuations in the superatom state that ultimately become reflected in the photon statistics. Except for the early onset and end of the decay of the pulse, where the experimental signal is dominated by noise, we observe good agreement between theory and experiment.
To complement the transmission measurements, we also detect whether the absorbers are in the ground or a collective Rydberg state after the probe pulse by field ionisation of atoms in jri. By resolving the time-of-flight from the atomic clouds to detection on a multi-channel plate (MCP) (Fig. 1a), we can determine from which absorbers the produced ions originate. Figure 3a shows mean detected ions 〈n ion 〉 per pulse and absorber for n sub = 3. The number of detection events from each absorber saturates at the detection efficiency η as expected if no more than one excitation is supported. The slight deviations in 〈n ion 〉 for the individual absorbers result from a slight dependence of η on the superatom position (between 0.18 and 0.25, 'Methods'). We also compare 〈n ion 〉 to the combined populations of jWi and jDi after the probe pulse in the three-level model, again in good agreement with the experimental data following multiplication with the corresponding values of η. To verify that each absorber is saturated by exactly one excitation, and thus represents a single ion source 37,50 , we analyse the counting statistics via the Mandel-Q parameter Q = Var(n ion )/〈n ion 〉 − 1 (Fig. 3b), which gives −η for perfect blockade (imperfect detection leads to a binomial distribution with success probability η), 0 for Poissonian, and >0 for super-Poissonian statistics. Analysing each absorber separately, we find Q ≈ −η for sufficient number of input photons, as expected for saturation at one, while analysis of the combined counts from all absorbers yields Q ≈ −〈n ion 〉/n sub as expected for saturation at n sub excitations.
Finally, Fig. 3c shows the ratio Q/〈n ion 〉 for the three individual absorbers. For large input photon number 〈n in 〉 and perfect blockade, when each absorber contains exactly one excitation, this quantity should give −1. We observe a small deviation as 〈n in 〉 increases which indicates the possibility of additional Rydberg excitations beyond the number of absorbers. To account for these in the model results, we increase the excitation probabilities obtained from the Rydberg populations in the three-level model by a small, photon-number-dependent probability p 2 〈n in 〉, which is independent of the superatoms' states. In addition, we also account for the influence of dark counts in the ion detection by adding an offset independent of 〈n in 〉, based on the experimentally observed dark-count rate of 9 kHz. These fluctuations induce a small Poissonian component in the ion counting statistics and thus shift Q/〈n ion 〉 to values above −1. The model results for 〈n ion 〉 and Var(n ion ) shown in Fig. 3c account for double Rydberg excitations and dark counts and reproduce the experimental results well for p 2 = 3.5, 6.5 and 5.0 × 10 −4 for the first, second and third absorber, respectively. To achieve good agreement for the second absorber, we need to increase the constant noise by a factor 5 compared to the darkcount rate. This is presumably due to the occasional detection of ions originating from the first and third absorber as well as atoms trapped in between the superatoms in the corresponding timeof-flight window attributed to the second absorber, which is expected to be more prone to these events due to its central position.
Our analysis supports the hypothesis that there is a small, 〈n in 〉-dependent probability to create additional Rydberg excitations, which may be caused by several mechanisms which we cannot distinguish in our experiment due to the small magnitude of the effect. First, residual atoms trapped between the superatom ensembles may be excited to jri if 〈n in 〉 is sufficiently high so that power broadening becomes comparable to the AC-Stark shift induced by the tightly confining optical traps. Second, imperfections in the blockade can occur from interaction-induced pairstate resonances on shells within the blockaded volume 51 , similar to the anti-blockade effect 52,53 .
Parameter optimisation and scalability. In the following theoretical analysis, we investigate the parameter space of the superatom photon absorber and discuss the potential for scaling beyond n sub = 3. We base the discussion on the results of the Lindblad master equation for a one-dimensional chain of chirally coupled superatoms described by the three-level model 48,54,55 ('Methods'). On one hand, parameter optimisation is necessary as the Raman decay Γ introduces an uncertainty about the number of absorbed photons and sets an upper limit on the pulse length τ. For short τ on the other hand, γ D must be balanced with the driving strength ffiffiffiffiffiffiffiffi ffi κR in p to yield high absorption probability. We begin our analysis by considering stochastic loss of photons due to incoherent Raman decay from the Rydberg manifold into jGi. It opens a scattering channel into non-observed modes, with a

mean number of lost photons
where P Ryd (t) is the combined population of jWi and jDi. Figure 4a shows the total number of absorbed photons P Ryd (t) + 〈N Raman (t)〉 together with one standard deviation ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi VarðN Raman ðtÞÞ p ('Methods') for Γ = 0, 0.004 and 0.04 μs −1 , for a single superatom at constant driving. The fluctuations in photon absorption increase over time, highlighting that deterministic photon subtraction requires both low Γ and short duration τ.
Consequently, we now analyse the dynamics of a single superatom at τ = 3 μs, similar to the experiment. Figure 4b shows the population of jDi vs. γ D and R in indicating a large parameter regime where photon absorption occurs with high probability. This regime is bounded by three processes with independent timescales, which we indicate by dashed lines. Firstly, the superatom is excited into jWi with rate ffiffiffiffiffiffiffiffi ffi κR in p and the requirement ffiffiffiffiffiffiffiffi ffi κR in p τ ) 1 gives a lower bound for the necessary photon rate. Similarly, jWi decays into jDi with rate γ D and thus the dark state will only be populated for γ D τ ≫ 1. However, at a large γ D , the superatom dynamics enter an overdamped regime in which we can adiabatically eliminate jWi ('Methods') and the effective absorption rate scales asymptotically as γ eff ≃ 4κR in /γ D . Therefore, we also require γ eff τ ≫ 1, limiting the maximal dephasing rate γ D . This analysis is valid until we reach a large ffiffiffiffiffiffiffiffi ffi R in κ p τ where we expect that the blockade mechanism breaks down and the superatom starts to absorb more than one photon.
Lastly, we solve the master equation for a chain of up to n sub driven superatoms. Figure 4c compares the total dark-state population for an ideal system with Γ = 0 μs −1 (solid lines) to Γ = 0.04 μs −1 (dashed lines) and shows the potential to extend our photon subtraction scheme up to n sub = 8. In the simulation, we drive the superatoms with a mean number of 20 photons, which indicates that our set-up can work well even when the number of photons becomes comparable to n sub . The dark-state population never reaches n sub exactly, which is due to the short pulse duration τ considered here and, by decreasing the Raman decay rate, higher absorption probabilities can be achieved.

Discussion
The theoretical analysis of our experimental results reveals that the main contributions to imperfections in the subtracted photon number are two-fold and are not necessarily unique to our scheme. First, a finite lifetime of the saturated state leads to excess absorption and the random nature of a decay process like Raman decay in our system introduces a probabilistic component into an initially deterministic scheme. This applies to any scheme which employs excited or metastable states and the severity of the impact depends on the decay strength compared to the absorber-photon coupling and the pulse duration. While the loss through Raman decay may initially seem as a disadvantage of our implementation, it should be noted that insertion loss into waveguides and cavities can lead to similar probabilistic fluctuations. Second, the slight deviation from n sub for coherent input pulses of finite duration is more general and affects all subtraction schemes relying on irreversible transfer into a dark state or separate optical modes irrespective of the absorber nature, also including hybrid systems of waveguide-and cavity-coupled single quantum emitters 14,56 .
In our scheme, Raman decay could be further suppressed by reducing either Ω c or increasing Δ, with the latter also reducing residual absorption on the probe transition. To compensate the associated reduction in κ, one can increase the number of atoms N per superatom via the ensemble density or increase the probe waist with r B as upper constraint. Meanwhile, for increasing κ combined with fine-tuning of γ D , the dark-state population converges towards n sub shifting the curves in Fig. 4c upwards. In this context, performance limitations will ultimately occur for high R in as power broadening causes a breakdown of the blockade.
While high-fidelity preparation of quantum states of light may require more substantial performance improvements, limitations are less stringent for other applications of our set-up. A more readily implementable application is number-resolved detection of up to n sub photons based on the number of absorbers in a Rydberg state. Currently, performance would be limited by the low-efficiency η to detect the superatom state via field ionisation, but this could be significantly improved by replacing the MCP by another model or using optical detection [57][58][59] . By increasing n sub well beyond the expected photon number, a weak photon-absorber coupling κ could also be compensated. Meanwhile, it is still important to minimise Raman decay as it reduces the detection efficiency for each absorber.
In summary, we have experimentally demonstrated controlled multi-photon subtraction from weak coherent probe pulses using a cascaded chain of saturable Rydberg superatom absorbers in free space. Our theoretical analysis has identified both technical and fundamental sources of imperfections, including the introduction of probabilistic fluctuations through decay into nonobserved modes resulting in a probability slightly below unity to transfer an absorber into its saturated state for coherent input pulses of finite duration.
Obvious next steps include improving the subtraction fidelity via the measures discussed above, changes to the optical trapping of the atomic ensembles to further increase n sub , and the implementation of optical readout for number-resolved photon detection. More generally, our system of cascaded superatoms is also well suited to study the behaviour of emitters coupled to a chiral waveguide 54,60,61 , as the superatoms not only introduce photon correlations, but the photons also coherently mediate interactions between the superatoms, which should become evident in the limit γ D ≪ κ and with increasing n sub .

Methods
Ensemble preparation. We start from a cigar-shaped ensemble of 87 Rb atoms in a crossed optical dipole trap (wavelength 1064 nm, 1/e 2 -waist ≈ 55 μm, intersection angle 30°) loaded from a magneto-optical trap (MOT). Following a final compression of the MOT, the atoms are evaporatively cooled as we reduce the trap light intensity in two stages. For additional cooling and to reduce atom loss, we employ Raman sideband cooling for 16 ms during each of the linear evaporation ramps. To create multiple ensembles for multiple superatom absorbers, we generate multiple, tightly focused optical traps with an elliptical cross-section (wavelength 805 nm, 1/ e 2 -waists ≈9 μm along and ≈29 μm perpendicular to the probe) that intersect perpendicularly with the cigar-shaped ensemble as well as the probe and control beams by feeding several RF signals into an AOD (as shown in Fig. 1 and discussed in 'Results'). An objective system translates the resulting differences in diffraction angle into different trap positions that can be tuned via the signals' frequencies over a range of order 100 μm, which is limited by the axial extend of the crossed region of the reservoir trap. In our experiments the separation between the ensemble centres is ≈75 μm for n sub = 2 and ≈50 μm for n sub = 3.
Before experiments, we ramp the crossed dipole trap intensity to zero to release atoms confined between the dimples before increasing it again to provide confinement in the radial probe direction for the superatom ensembles. In combination with the 1/e 2 -waist radius of the probe (≈6.5 μm), the dimple confinement restricts the excitation volume below the blockade range. The focus of the control beam is larger (≈14 μm) to limit the variation of Ω c across the excitation volume.
Experimental sequence. Following the preparation outlined above, we turn the crossed dipole trap off for 14 μs every 100 μs, while the dimple traps are left on to maintain confinement along the probe direction. The resulting AC-Stark shift is compensated by adjusting the probe frequency accordingly and we ensure that all superatom absorbers have the same resonance frequency by individual fine-tuning via the power for each RF signal applied to the AOD. Besides the axial confinement, the AC-Stark shift also helps to suppress Rydberg excitation of residual atoms trapped in between the dimple potentials. Following each single experimental shot, we field-ionise any Rydberg atoms to gain information about the absorber state and avoid the presence of residual Rydberg excitations during the next iteration of the superatom excitation. The ions are detected on a MCP. In total, we repeat the cycle described above 500 times before releasing the atoms to obtain reference pulses of the probe light in the same manner in the absence of any atoms and preparing new atomic ensembles.
Site-resolved ion detection. In order to attribute the ions detected following each experimental shot to Rydberg excitations in different ensembles, we use a time-offlight method. We find that the pulses generated by the MCP detector following detection of an ion occur in a time window with a width of ≈30 ns and a typical separation of several 10 ns between the arrival times for two superatom ensembles separated by 50 μm. Combined with the 3 ns time resolution of our data acquisition, this allows us to attribute a detection windows of 75 ns to the location of each individual superatom (with a 40 ns gap between the windows for the superatoms furthest from the detector). In Fig. 3, the detection efficiency varies between η ≈ 0.18 and 0.25 depending on the position of a superatom and is generally highest in the single absorber case. The variation is caused by a grid of steel wires, which is placed in front of the detector to shield the atoms in the experimental region from the strong electric field produced by the MCP front plate and partially obstructs the ion trajectories. For a single superatom, the applied ionisation and steering fields can be adapted to minimise the influence of the grid, but for multiple absorbers we cannot avoid that a fraction of ions is blocked by the wires, which depends on the location of their origin.
Theoretical description. We describe each superatom i as an effective three-level atom whose ground state jG i i is coupled to a collective excited state, the bright state, jW i i by the coherent probe field α(t) (with R in = |α(t)| 2 ). The photon absorber relies on shelving excitations into a non-radiating dark state jD i i, which we model by an incoherent decay of the bright state with rate γ D . Assuming the dipole and rotating wave approximation and no Raman decay, the superatoms obey the master equation 54 The master equation consists of the action of the probe field a hopping term due to the exchange of virtual photons and the dissipative decay terms D ∑ N i¼1 σ À W i h i and ∑ N i¼1 D½σ þ D i σ À W i , describing the collective decay of the superatoms and dephasing of each bright state into the respective dark state. We use the notation σ À A i jG i ihA i j, σ þ A i jA i ihG i j and D½σρ σρσ y À fσ y σ; ρg=2.
To understand non-deterministic effects in our photon absorber set-up, it is essential to include Raman decay in our model. In the simplest description, Raman decay enters as an additional decay term ΓðD½σ À W i ρ þ D½σ À D i ρÞ for each superatom. To gain access to the number statistics of the emitted Raman photons, we further introduce a virtual spin chain to our model and modify the Raman decay terms so that each Raman decay yields an excitation in the spin chain. This allows us to calculate the standard deviation of the number of emitted Raman photons, as depicted in Fig. 4a.
At a large γ D , the superatoms enter an overdamped regime and further increasing γ D negatively impacts the photon absorption rate. In the overdamped dynamics, we may adiabatically eliminate the bright state from the master equation, which we achieve by setting ∂ t ρ WW ¼ 0 ¼ ∂ t ρ WG ¼ ð∂ t ρ GW Þ Ã , with the shorthand notation 〈A|ρ|B〉 = ρ AB . Under the adiabatic elimination, the master equation for a single superatom reduces to a classical rate equation for the groundstate and dark-state population where the effective decay rate reads as

Data availability
The data generated in this study have been deposited in the Zenodo database under accession code https://doi.org/10.5281/zenodo.4984099 Received: 25 March 2021; Accepted: 23 June 2021; Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/ licenses/by/4.0/.