Observation of Highly Correlated Ultrabright Biphotons Through Increased Atomic Ensemble Density in Spontaneous Four-Wave Mixing

,

Introduction.-Temporallycorrelated biphotons have recently garnered considerable attention in the fields of optical quantum computing and quantum communication, thanks to their exceptional nonclassical properties.Of particular significance is their role as heralded singlephoton sources, which have found applications in diverse domains, including quantum cryptography [1][2][3][4], quantum metrology [5][6][7][8], and quantum imaging [9][10][11][12].Among the various biphoton sources, the spontaneous four-wave mixing (SFWM) mechanism, distinguished by its operation near atomic resonance, stands out for its ability to conveniently manipulate bandwidth and serves as a bridge between different quantum devices, thus attracting significant interest.
In this Letter, we present the thorough investigation of the pairing ratio using the double-Λ SFWM in a cold 87 Rb ensemble.This configuration, chosen for its inherent EIT effect, allows for easy control of the biphoton bandwidth.We achieved an exceptionally high biphoton generation rate of approximately 1.3 × 10 7 per second at an optical depth (OD) of 20, with the pairing ratio measured at 0.61.We also illustrated that this ratio can be enhanced by increasing OD and achieved the highest pairing ratio of 0.89 at an OD of 120.Additionally, at a relatively low biphoton generation rate of 5.0 × 10 4 per second, the signal-to-background ratio for the biphoton wavepacket reached 241, exceeding the Cauchy-Schwarz criterion by approximately 1.5 × 10 4 times.
Experimental setup.-Wetrapped cold 87 Rb atoms using a standard magneto-optical trap.After optically pumping them to the ground state |5S 1/2 , F = 1⟩, as illustrated in Fig. 1, we irradiated the atomic ensemble with a far-detuned driving field characterized by a Rabi frequency Ω d and a nearly resonant coupling field denoted by Ω c .Synchronization between these fields was achieved through injection locking with an external cavity diode laser (not shown).The driving field, with detuning ∆ d , operated on the σ + -transition |5S 1/2 , F = 1⟩ → |5P 3/2 , F = 2⟩, effectively suppressing incoherent fluorescence from one-photon absorption.This allowed for primary emission of Stokes photons via spontaneous Raman scattering [50].Subsequently, the nearly resonant coupling field, with detuning ∆ c , acted on the σ +transition |5S 1/2 , F = 2⟩ → |5P 3/2 , F = 2⟩, inducing the emission of anti-Stokes photons.We used an elongated atomic ensemble to enhance the FWM effect and boosting specific direction scattering probability.
In the SFWM experiment, biphotons were generated using a 10-µs driving pulse in each 2.5-ms cycle.To prevent laser leakage, we adopted a backward configuration where the driving and coupling beams counterpropagated, each with 1/e 2 full widths of 250 and 310 µm.The corresponding optical powers of the 1Γ driving and coupling fields were approximately 7.5 and 11.1 µW.Both fields drove σ + -transitions, resulting in generated photon pairs exhibiting σ + polarization, propagating in opposite directions and passing through respective etalon filter sets (EFS).The intersection angle between the Stokes (anti-Stokes) and driving (coupling) beams was set at 1.7 • for our experiment.Each EFS consisted of two etalons, each with an extinction ratio of roughly 30 dB and an approximate bandwidth of 100 MHz, separated by an optical isolator.The total extinction ratios of the Stokes and anti-Stokes channels were 114 dB and 124 dB, respectively.Biphotons were detected using fiber-coupled single-photon counting modules (SPCM-AQRH-13-FC).Upon detection, an 8-ns pulse was emitted from the SPCM toward the time-of-flight multiscaler (TOF, MCS6A-4T8, not shown).In the coincidence count experiment, we measured the time difference between Stokes and anti-Stokes photons.When generated as a correlated pair, they arrived at the SPCMs within the correlation time, contributing to a nonflat biphoton wavepacket.The TOF generated a histogram of coincident counts based on these data points, providing insight into the source of biphotons.
Theoretical model.-Weused the Heisenberg-Langevin operator approach to analyze the biphoton generation in double-Λ SFWM [49].More details can be found in the Supplemental Material.The photon generation rate, given by R = c L ⟨â † â⟩ and derived from the annihilation operator â, leads to the following expressions: where D jk † ,j ′ k ′ and D jk,j ′ k ′ † are diffusion coefficients, while R s(as) represents the spectrum of Stokes (anti-Stokes) photons.The total photon generation rates, R s and R as , comprise two components: correlated photons governed by the FWM process with coefficients B and C, and uncorrelated photons due to vacuum field fluctuations represented by an integral term with diffusion coefficients.The pairing ratio, r p , denote the ratio of the correlated photons to the total generated photons.Under ideal conditions, R s and R as are nearly equal.However, in experiments, R as is slightly smaller due to phase mismatch and ground state decoherence.Therefore, the biphoton generation rate R B is contingent on R as In the SFWM process within the atomic ensemble, a spontaneously emitted Stokes photon from one atom may interact with nearby atoms, triggering stimulated Raman scattering.Unlike spontaneous Raman scattering, the Stokes photon generated by stimulated Raman scattering shares the same direction as the incident Stokes photon, enhancing directionality.This collective enhancement effect establishes paired correlation with the anti-Stokes photon through the FWM process, reflected in coefficients B and C.However, as SFWM relies on vacuum field fluctuations, the generated photons exhibit isotropic (uncorrelated) nature, represented by the integral term with diffusion coefficients in Eqs. ( 1) and (2).While the paired correlations of these biphotons can be established through the FWM process, it requires a sufficiently high OD within the atomic ensemble.
In biphoton systems, the normalized Glauber secondorder cross-correlation function g (2) s-as (τ ) is often used alongside the photon generation rate.This function is a crucial parameter for evaluating the temporal correlation between biphotons.The derived theoretical expression is as follows: g (2)  s-as The integral term on the right-hand side of Eq. ( 3) reveals the correlation of biphotons.This term is equivalent to the wavepacket of the anti-Stokes single photon, conditioned on the post-selection of a Stokes single photon.This correlation provides valuable information for evaluating the biphoton source.For instance, the peak signal-to-background ratio, denoted as r SB , is defined as the maximum value of [g s-as (τ ) − 1].It serves as a standard metric for assessing the nonclassicality of a biphoton source.In the case of a classical field, the Cauchy-Schwarz inequality universally applies: as-as (0) −1 ≤ 1.The normalized autocorrelation functions of the Stokes and anti-Stokes fields can be derived as g as-as (0) = 2. Nonclassical behavior is observed when r SB > 1.Additional details and initial proofs of the thermal field distributions for both the Stokes and anti-Stokes fields can be found in the Supplemental Material.
Biphoton bandwidth.-Figure 2 presents the experimental biphoton coincidence count rate R C (or biphoton temporal wavepacket) for various coupling field conditions.R C is calculated as R s R as g (2) s-as (τ )∆T + R env (refer to Supplemental Material), where R env accounts for environmental background count rates, arising from laser leakage or SPCM dark counts.For data processing, we used a time bin of ∆T = 1/R s to tally Stokes photons, enabling the post-selection of a single Stokes photon.This ensures that the background and correlated regions of the coincidence count rate correspond to R B +R env and r p , respectively.The time bin for detected anti-Stokes photons, ∆τ = 6.4 ns, aligns with the time interval between experimental data points in Fig. 2.
In Figs.2(a) and 2(b), we set the coupling Rabi frequencies Ω c to 4Γ and 1Γ, respectively.With the OD fixed at 15, both cases yielded a measured biphoton generation rate R B of approximately 3.4 × 10 5 s −1 .The delay time in Fig. 2(b) is significantly longer than that in 2(a).This delay arises from two intrinsic properties of the double-Λ SFWM system: the damped Rabi oscillation with a period denoted as  Figures 2(c) and 2(d) demonstrate how changes in coupling detuning ∆ c affect the biphoton bandwidth.All experimental parameters were consistent with those in Fig. 2(b), except for the ∆ c .A shorter delay time in Fig. 2(c) is observed due to the introduction of ∆ c = 1Γ, which reduces the effective OD and shortens the EIT-induced delay.Conversely, with ∆ c = 3Γ in Fig. 2(d), the tail lengthens again.This behavior is attributed to damped Rabi oscillations, where a larger ∆ c weakens the interaction between the coupling field and the atomic medium, requiring more time to convert spinwave excitations into anti-Stokes photons.The values of R B in Figs.2(c) and 2(d) are 3.4 × 10 5 s −1 and 3.1 × 10 5 s −1 , respectively.This demonstrates that by detuning the coupling field, we can control the biphoton bandwidth without significantly reducing R B .Further discussions can be found in the Supplemental Material.
High-purity biphotons.-Figure3(a) illustrates highpurity biphoton generation achieved with parameters Ω d = 0.5Γ, Ω c = 4Γ, and OD = 10.The theoretical R B is calculated as 5.0 × 10 4 s −1 .In experiments, R B was determined by subtracting the total measured background count rate, R tot = 1.6 × 10 5 s −1 , from the environmental background count rate, R env = 1.1×10 5 s −1 .This yielded an experimental R B of approximately 5.0 × 10 4 s −1 , in close agreement with the theoretical prediction.In Fig. 3(b), we theoretically calculated R B to be 1.9 × 10 6 s −1 at Ω d = 3Γ.The experimentally observed R B , obtained from measurements of R tot = 2.0 × 10 6 s −1 and R env = 1.2 × 10 5 s −1 , also closely matches theoretical prediction.As Ω d increases, both R B and R tot rise significantly.However, this also leads to a notable decrease in r SB , as shown in Fig. 3(c).At Ω d = 0.5Γ, we observed an experimental r SB = 241, surpassing the Cauchy-Schwarz criterion by a factor of approximately 1.5 × 10 4 .If the R env in our experiment could be completely eliminated, it would lead to a more pronounced violation of the Cauchy-Schwarz criterion, exceeding the normal level by a factor of 5.9 × 10 4 .
Figure 3(d) shows the variation of R B and r p with different Ω d values.At Ω d = 0.5Γ and Ω d = 3Γ, the corresponding r p values are 0.63 and 0.59, respectively.These experimental r p values were determined based on the area under the correlated biphoton wavepacket.The r p obtained from the area and those obtained from Eqs.
(1) and ( 2) are equivalent, as detailed in the Supplemental Material.In the SFWM process, atomic ensembles play a crucial role in collectively enhancing the correlation between the Stokes and anti-Stokes fields along the applied light direction.Therefore, with a fixed OD, while increasing Ω d can boost R B , the limited density of atomic ensembles constrains their ability to produce correlated photon pairs, leading to a slight decrease in r p .Highly correlated biphotons.-Figure4(a) showcases the generation of ultrabright biphotons using specific parameters: Ω d = 3Γ, Ω c = 4Γ, ∆ d = 5Γ, and OD = 20, resulting in a remarkable theoretical R B of 1.3 × 10 7 s −1 .This exceeds rates reported in the literature for the double-Λ SFWM scheme.Under these conditions, the experimental total background count rate was R tot = 1.3 × 10 7 s −1 , with environmental background at R env = 2.3 × 10 5 s −1 , accounting for only 1.8% of the total.Thus, in this high R B scenario, the primary source of background count arises from the high photon generation rate rather than environmental factors.Furthermore, in this scenario, the measured r SB was 2.4, surpassing the Cauchy-Schwarz criterion by a factor of 2.9, while r p was only 0.61.Although increasing the coupling power can enhance r SB , as demonstrated in Fig. 2, it does not lead to corresponding improvements in r p .To enhance both r p and r SB , we further increased the OD.In addition to OD = 20, we measured the biphoton wavepacket at OD = 40, 60, 80, 100, and 120.We fine-tuned Ω c to maintain a consistent biphoton bandwidth, while keeping Ω d = 3Γ constant and adjusting ∆ d to maintain a theoretucal R B of 1.3 × 10 7 s −1 .Specific parameters can be found in the Supplemental Material.
In Fig. 4(b), we present the scenario with OD = 120, Ω c = 8.8Γ, and ∆ d = 14.9Γ.Here, the measured r SB at 4.2 exceeds the Cauchy-Schwarz criterion by 6.8 times.An evident positive correlation emerges between increased OD and enhanced r p , resulting in a higher r SB due to augmented coincidence receptions [Fig.4(c)].This enhancement stems from the increased accumulation of biphoton correlations along a specific direction at higher OD values.Photons generated at higher OD levels are more likely to encounter subsequent atoms, amplifying the collective enhancement through the FWM process.Furthermore, while r SB can also be improved by increasing Ω c , this approach does not enhance r p , and therefore cannot improve the generation rate of temporally correlated biphotons.Figure 4(d) illustrates the relationships between R B and r p with OD.At OD = 120, we observed the highest r p of 0.89, indicating a significant improvement in correlated photon pair generation.The experimental R B = 1.3 × 10 7 s −1 signifies the successful generation of approximately 1.2 × 10 7 pairs of correlated photons per second.Additionally, the Fourier transform of (R C − R tot ) reveals a biphoton bandwidth of approximately 24 MHz, resulting in a spectral brightness of the biphoton source at 5.4 × 10 5 s −1 MHz −1 , surpassing the highest achieved by sub-megahertz biphoton sources [52].These results highlight the crucial role of high OD in SFWM-based biphoton sources, allowing for higher values of R B and r p .This enables the generation of a large quantity of high-quality correlated photon pairs for use in various quantum systems.
Conclusion.-Our investigation into the biphoton pairing ratio, utilizing the double-Λ SFWM in cold 87 Rb atoms, revealed a marginal decrease with higher biphoton generation rates.Nonetheless, this trend can be effectively addressed by elevating the atomic ensemble density.The highest pairing ratio observed was 0.89 at an OD of 120, accompanied by an ultrabright biphoton generation rate of up to 1.3 × 10 7 s −1 , surpassing previously reported rates achieved via the double-Λ SFWM scheme.Furthermore, our experiment demonstrated the highest signal-to-background ratio of the biphoton wavepacket at 241, achieved at a low biphoton generation rate of 5.0 × 10 4 s −1 .This outstanding performance exceeded the Cauchy-Schwarz criterion by approximately 1.5×10 4 times.These results underscore the capability of the double-Λ SFWM scheme in advancing biphoton sources for future quantum technologies.
In this research, we used the Heisenberg-Langevin operator approach to characterize the biphotons generated with the spontaneous four-wave mixing (SFWM) process [1].The slowly varying collective atomic operator σjk (z, t) = Nz l=1 σl jk (z, t)/N z was used to model the intricate atomic dynamics during the interaction, where N z is the number of atoms within the infinitesimal one-dimensional (1D) spatial interval dz.Here, σl jk denotes the adiabatic atomic operator of the lth atom that has been transformed by the slowly varying amplitude.The dynamics of slowly varying collective atomic operators are governed by the Heisenberg-Langevin equations (HLEs): ∂ ∂t σjk = i ℏ [ Ĥ, σjk ]+ Rjk + Fjk , where Rjk and Fjk are the adiabatic relaxation term and the Langevin noise, respectively.The interaction Hamiltonian Ĥ of the double-Λ SFWM atomic ensemble is where N and L are the total number of atoms and the atomic medium length, respectively.The applied highintensity coupling and driving fields, characterized by the semiclassical Rabi frequency Ω for the double-Λ SFWM system can be derived as follows: as σ14 e i∆kz − g as âas σ41 e −i∆kz + Γ 31 σ33 + Γ 41 σ44 + F11 , (S2) Here Γ jk is the spontaneous decay rate from the excited state |j⟩ to the ground state |k⟩.γ jk is the decoherence rate corresponding to the states |j⟩ and |k⟩, which is due to the dephasing from the spontaneous emission.Of these, although both the states |1⟩ and |2⟩ are ground states, the term γ 21 is introduced because of some effects such as atomic collisions and magnetic inhomogeneity.In our system, Γ 31 = Γ 32 = Γ 41 = Γ 42 = Γ 2 , where Γ is the spontaneous decay rate of the rubidium-87 D 2 line, and therefore we have γ 31 = γ 32 = γ 41 = γ 42 = Γ and γ 43 = 2Γ.
To simplify the solution of the HLEs, we employed perturbation theory.By neglecting the influence of the quantum field operator âs(as) in the HLEs, we can derive the zeroth-order steady-state collective atomic operators in the form of σ(0) jk = ⟨σ (0) jk ⟩ + mn ϵ mn jk Fmn , with expectation values: ⟨σ 13 , and ⟨σ Despite our assumption that both input fields operate in a single mode described by steady-state zeroth-order HLEs, the solutions of the zeroth-order collective atomic operators still exhibit time-varying terms due to the presence of the vacuum field.These time-varying terms lead to a dynamic equilibrium, which is reflected in the linear combination of Fmn .By substituting the derived zeroth-order adiabatic collective atomic operators back into the HLEs, we can obtain the first-order HLEs as follows: ∂ σ( 1) When substituting the zeroth-order results into the first-order HLEs, we omitted the product of the field ladder operators and the Langevin noises Fmn due to their negligible effects.To investigate biphoton spectra, we applied the Fourier transform on Eqs.(S12)-(S15).We chose the angular frequency shift, represented as ω, relative to the central angular frequency of the Stokes field, ωs , as the frequency-domain variable in the Fourier transform.In other words, a s (z, ω) = 1 2π dtâ s (z, t)e iωt ≡ F{â s (z, t)}.To maintain the factor of the integrand as e iωt , the Fourier transform of the creation operator should be adjusted to a † s (z, −ω) = 1 2π dtâ † s (z, t)e iωt .Similarly, we have a as (z, ω) = F{â as (z, t)}, a † as (z, −ω) = F{â † as (z, t)}, σ jk (z, ω) = F{σ F jk , where {jk} = {21, 23, 41, 43}.We combine these results with the frequency-domain Maxwell-Schrödinger equations used to describe our experimental backward double-Λ SFWM system as follows: We can simplify the expressions by introducing a change of variables, defining a † as,∆k ≡ a † as e i∆kz to remove the exponential phase variation.This allows us to reorganize Eqs.(S16) and (S17) into the following matrix form: The Raman gain coefficient, denoted as g R , is expressed as The solution to Eq. (S18) can be obtained by directly integrating from z = 0 to z = L, as illustrated below: Since the Stokes and anti-Stokes fields are co-propagating outputs in the backward double-Λ SFWM system, we now reorganize Eq. (S19) to describe this biphoton output as follows: A and D represent coefficients for mode preservation, elucidating the influence of the Stokes and anti-Stokes inputs on their respective outputs.On the other hand, B and C represent coefficients for mode conversion, describing the energy conversion between the Stokes and anti-Stokes fields, which is the process of photon mutual conversion in the double-Λ four-wave mixing (FWM) [2,3].

B. Biphoton Generation Rate and Correlation Functions
The biphoton generation rate is related to the photon flux passing through the interaction cross-section area with the atoms.The photon flux can be obtained using the Poynting vector, denoted as ⃗ S = 1 µ0 ⃗ E × ⃗ B * , which represents the energy flux per unit area of an electromagnetic wave.By dividing the surface integral of the Poynting vector by the energy of a single photon, and taking into account the quantized electric field in the form of ⃗ E = ê ℏω 2ϵ0V âe i( ⃗ k•⃗ r−ωt) + H.c., we can determine the quantum photon generation rate as follows: The parameter ê is the unit vector of the polarized direction.The term c 2L on the right-hand side of Eq. (S21) arises from the zero-point energy of the vacuum field and does not contribute to the actual count of generated photons.Thus, the photon generation rate can be expressed as R = c L ⟨â † â⟩.By employing the frequency-domain commutation relation ), we can derive the photon generation rates for the Stokes (R s ) and anti-Stokes (R as ) fields, as given by Eqs. ( 1) and (2) presented in the main text.In our SFWM system, it is worth noting that g p = g s ≡ g.To conform to commonly used experimental parameters, we perform the substitution g 2 N c = ODΓ 4L , where OD denotes the optical depth of the atomic medium.In addition to the photon generation rate, we can also derive the normalized Glauber second-order cross-correlation function, denoted as g (2) s-as (τ ), as depicted in Eq. ( 3) in the main text.Similarly, we can obtain the normalized autocorrelation functions g s-s (τ ) and g (2) as-as (τ ) for the generated Stokes and anti-Stokes photons.

C. Diffusion Coefficients and Einstein Relation
The HLE for the lth atom is given by ∂ ∂t According to the Markovian approximation, the correlation between Langevin noises follows a Dirac delta function as , where D jk,j ′ k ′ represents the diffusion coefficient satisfying the Einstein relation [4].In this case, the Einstein relation can be derived from the following equation: where the term of ⟨ F l jk σl j ′ k ′ ⟩ can be expressed as The term ⟨ F l jk (t)σ l j ′ k ′ (t − ∆t)⟩ yields zero as the Langevin noises at time t are incapable of influencing the atom at time t − ∆t, owing to the nonanticipatory property.Additionally, given that the terms [ Ĥl (t ′ ), σl j ′ k ′ (t ′ )] and Rl j ′ k ′ (t ′ ) evolve at a slower rate compared to the Langevin noises, the integrals of the first two terms in Eq. (S23) both approach to zero.Consequently, Eq. (S23) simplifies to ⟨ F l jk Employing these relations, Eq. (S22) can be further expressed as follows: and this is commonly referred to as the Einstein relation.After establishing these relationships, we can extend our findings from the single-atom scenario to collective situations.In the context of an atomic ensemble system, we define the collective Langevin noise as Fjk =

Nz l=1
F l jk /N z , where N z = N dz/L.With this definition in place, we proceed to formulate By using the relation δ zz ′ = δ(z − z ′ )dz, we can obtain the Einstein relations in both the time and frequency domains for the atomic ensemble system as follows: We employ the reduced density operator approach to derive the density matrices of the Stokes and anti-Stokes fields generated from the double-Λ SFWM system.This method transforms the obtained field operator into the form of a density matrix [5][6][7].Consider an open system with a Stokes field ρ s (t), anti-Stokes field ρ as (t), and vacuum reservoir ρ R (t).Initially, these fields are assumed to be uncorrelated, meaning ρ(0) = ρ s (0) ⊗ ρ as (0) ⊗ ρ R (0).The dynamics of this open system can be described by ρ(t) = U (t)ρ(0)U † (t), where U (t) is the time evolution operator of the total density matrix.The Stokes field ρ s (t) can be obtained by tracing over the degree of freedom of ρ as (t) and ρ R (t) from the total density matrix, i.e. ρ s (t) = Tr as,R {ρ(t)}.The density matrix element of the generated Stokes field in the Fock-state basis can be given by Consider the vacuum input ρ(0) = |0 s ⟩⟨0 s | ⊗ |0 as ⟩⟨0 as | ⊗ ρ R (0) in the double-Λ SFWM system.By rewriting the vacuum density matrix with the annihilation and creation operators |0⟩⟨0| = ∞ l=0 , we can derive the density matrix element of the Stokes field as follows: Using the inverse Fourier transform and combining Eq. (S20), the Stokes field operator in the time domain can be obtained as follows: where F s (ω) is a noise-correlated term that is a linear combination of F 21 , F 23 , F 41 , and F 43 .After substituting Eq. (S30) into Eq.(S29) and considering the initial vacuum input ρ(0), we find that only the condition where m = n contributes to a non-zero value of the density matrix element of the Stokes field.With these results, we can derive the density matrix element of the Stokes field as follows: To further simplify Eq. (S31), we define ζ ≡ dωe −iωt B * a as (L, ω) and η ≡ dωe −iωt F s (ω), resulting in the following expression: where we apply the binomial coefficient , and only when β equals α, there is a non-zero contribution, leading to the following expression: The operator product term ζα ( ζ † ) α can be simplified using the generalized Wick's theorem.Along with the commutation relation [ a as (L, ω), a † as (L, −ω ′ )] = L 2πc δ(ω + ω ′ ), this term can be expressed as: Similarly, the operator product term (η † ) l+n−α ηl+n−α can be simplified as follows: Thus, Eq. (S32) yields where the term L c R s = ⟨â † s (L, t)â s (L, t)⟩ represents the average photon number of the Stokes photons ns .Hence, the density matrix of the Stokes field becomes which shows the characteristic of thermal field distribution.Similarly, we can consider the anti-Stokes field operator in the time domain to have the following form: This identical structure enables us to deduce its density matrix, denoted as ρ as , utilizing the previously mentioned method.Consequently, we acquire the density matrix for the anti-Stokes field as outlined below: also exhibiting the characteristics of a thermal distribution.

II. EXPERIMENTAL DETAILS A. Coincidence Count Rate
The coincidence count rate characterizes the detection rate in the anti-Stokes channel for a specific number of received Stokes photons.This rate is closely related to the normalized Glauber second-order cross-correlation function, defined as s-as (τ )∆T , where ∆T represents the time bin for detecting the Stokes photons.To ensure a fixed value of unity for the post-selected Stokes photon, we chose a time bin of ∆T = 1/R s .This choice leads us to the following expressions: Our experimental approach involves accumulating 2 16 counts in the Stokes channel for each data point.By selecting a total of 2 18 receptions in the Stokes channel, each with a purity P s , we can effectively convert the coincidence count rate R C into coincidence counts N C using the following relations: These relations enable us to quantitatively determine the number of coincidence counts based on the coincidence count rate from the reception characteristics P s in the Stokes channel, the detection efficiency η as in the anti-Stokes channel, and the relevant noise rate R as noise .The first term in Eq. (S41) represents the detection of biphotons, corresponding to the registration of both Stokes and anti-Stokes photons.This term includes the contribution of the background from the biphoton generation rate.The second and last terms are environmental background, which are respectively caused by the pure and impure detections in the Stokes channel.If Stokes photons are detected, the impure registers in anti-Stokes channel results in an undesired background; this is represented by the second term.By contrast, an uncorrelated background is necessarily the case if the Stokes channel receives noised photons; this is indicated by the last term.Finally, the experimental biphoton coincidence count rate can be expressed as N C /(2 18 P s η as ∆τ ) ≡ R C + R env , where R env represents the count rate of the environmental background.

B. Pairing Ratio
Pairing ratio is an important quantity that should be considered.Without introducing this concept, it could potentially lead to misinterpretations of the experimental results.In the main text, we discussed two different forms of pairing ratio definitions.It's important to note that these two definitions are equivalent in reality.This core concept arises from the fact that temporal correlation is not inherently present in all generated biphotons; rather, it needs to be acquired by increasing the ensemble density.In the method using R C , we calculate the pairing ratio by integrating R C over delay time τ .When the received Stokes photons are normalized to one, only the proportion of r p exhibits temporal correlation.Consequently, the proportion of anti-Stokes photons correlated in time with this Stokes photons will also be r p , equivalent to the area of R C with background subtracted.However, it's crucial to establish the upper limit of this integral before deriving the pairing ratio from our experimental results.To address this concern, we pinpoint a specific time, denoted as t 0 , based on the theoretical R C .The selection of t 0 ensures that the cross-correlation function g s-as (t 0 ) reaches the background value of 1, and R C (t 0 ) at this time equals the value of R B .For instance, we utilized the same parameters as depicted in Fig. 3(b) of the main text, setting t 0 at 180 ns, as demonstrated in the subsequent figure presenting the theoretically predicted R C .While the choice of t 0 was informed by theoretical curves, our preceding experimental results, including those illustrated in Figs. 2, 3, and 4 of the main text, consistently exhibit a robust agreement between our theoretical model and experimental data.Hence, we deem the selection of t 0 to be well-founded.The same concept is evident in Eqs. ( 1) and ( 2) presented in the main text.Although the biphotons based on SFWM inherently possess frequency anti-correlation, not all generated biphotons exhibit temporal correlation.This characteristic is elucidated in the Stokes photon generation rate [Eq.( 1)] and the anti-Stokes photon generation rate [Eq.( 2)].The generated biphotons are categorized into temporally correlated photons [the first terms of Eqs. ( 1) and ( 2)] and temporally uncorrelated photons [the second terms of Eqs. ( 1) and ( 2)].The pairing ratio is defined as the ratio of temporally correlated photons to the total generated photons.In our system, characterized by a small ground-state dephasing rate, the results obtained from the Stokes and anti-Stokes channels are nearly identical.
In contrast to the pairing ratio, the heralding efficiency addresses the practical application of collection efficiency [9].It focuses on measuring anti-Stokes photons and determining the proportion heralded by Stokes photons, essentially representing the collection efficiency η s of the Stokes channel [10].For instance, assuming a biphoton generation rate of R B and a pairing ratio of one (r p = 1), the detection rates on both sides, R 1 = R B η s and R 2 = R B η as , can be obtained by separately measuring each channel.Through coincidence detection experiments, the temporally correlated biphoton detection rate R 3 = R B η s η as can be obtained.These temporally correlated photons occupy a fraction of the detected photons in the anti-Stokes channel with a proportion of η s .Specifically, when detecting the anti-Stokes channel, only the rate of R 3 is temporally correlated, representing a single-photon Fock state, while the remaining rate of R 2 − R 3 corresponds to a thermal state.

C. Collection Efficiency
In the Stokes channel, a pinhole with a 60% transmission rate is positioned in front of our initial etalon.The first etalon itself has a transmission rate of 60%, while the second etalon has a transmission rate of 75%.An optical isolator with an 85% transmission rate is situated between the two etalons.Subsequent to the second etalon, light is gathered using a multi-mode fiber (MMF) boasting an efficiency of 80%.The MMF is connected to a single-photon counting module (SPCM) with a quantum efficiency of 65%.The coupling efficiency between the MMF and SPCM is 92%.Factoring in the transmission rates of other optical components, the overall collection efficiency is approximately 9%.
The configuration of the anti-Stokes channel mirrors that of the Stokes channel, featuring a 60% transmission rate pinhole, two etalons, and an optical isolator with an 85% transmission rate.The final step involves capturing light using an MMF and coupling it into a SPCM with a quantum efficiency of 65% (coupling efficiency between MMF and SPCM is 92%).Noteworthy differences lie in the efficiencies of the etalons and MMF.These variations stem from the optical path design, resulting in differences in the guiding beam size compared to the Stokes case.Both the first and second etalons maintain a transmission rate of 60%, and the MMF efficiency is 70%.Considering the transmission rates of other optical components, the overall collection efficiency is approximately 6%.
For a collection time T 1 , the total number of received Stokes photons during T 1 is given by R s η s T 1 , where η s denotes the collection efficiency of the Stokes photons.Due to the anti-Stokes collection efficiency η as and the pairing ratio r p , only R s η s T 1 ×η as r p Stokes photons are successfully paired with corresponding anti-Stokes photons.On the other hand, when Stokes and anti-Stokes photons from different pairs are detected, their mutual delay times become stochastic, leading to the appearance of a uniform background signal.Within the collection time T 1 , for each individual detected Stokes photon (total of R s η s T 1 photons), the corresponding anti-Stokes photon detection rate is R as η as .Consequently, the average background count within each time spacing ∆τ can be expressed as R s R as η s η as T 1 ∆τ .We refer this ∆τ to the time bin of the detected anti-Stokes photons.
To ensure sufficient counts at low generation rates, we adopted an alternative approach in our measurements, avoiding the fixed collection time T 1 .Instead, we maintained a constant count of 2 18 in the Stokes channel to ensure an ample number of detections.As a result, the correlated coincidence counts and uncorrelated background counts were adjusted to 2 18 P s η as r p and 2 18 P s R as η as ∆τ , respectively.Here, P s represents the purity of Stokes detection and is calculated as P s = Rsηs Rsηs+R s noise , where R s noise denotes the count rate of the Stokes channel originating from laser leakage, environmental photons, and dark counts from the SPCM.This parameter was determined by conducting single-channel detection on the Stokes channel over a time interval T 2 , resulting in a count of R s η s T 2 + R s noise T 2 .Furthermore, by theoretically calculating the Stokes generation rate of this process, we obtained the Stokes collection efficiency η s .Similarly, the collection efficiency η as for the anti-Stokes channel was determined using the same methodology, yielding a count of R as η as T 2 + R as noise T 2 within the time interval T 2 , where R as noise represents the noise count rate of the anti-Stokes channel.By incorporating these experimentally derived collection efficiencies η s and η as into the correlated coincidence counts and uncorrelated background counts, respectively, we were able to obtain the experimental pairing ratio and biphoton generation rate.
The collection efficiencies obtained from our experiments using guiding beams differ from those acquired by measuring the Stokes and anti-Stokes channels in conjunction with the theoretical biphoton generation rate.Specifically, the collection efficiencies obtained from the guiding beams for the Stokes and anti-Stokes fields were 9% and 6%, respectively, while the corresponding values measured from the biphoton experiments were approximately 2% and 1%, respectively.For convenience, we refer to the method of determining collection efficiency using the guiding beam as the "guiding method."On the other hand, the approach that relies on the measured photon counts from biphoton experiments and the theoretical generation rate is termed the "experimental method."The disparity between these two approaches may arise from differences in the propagation of the guiding beam and the biphoton field, potentially leading to variations in the coupling efficiency with etalons and optical fibers.In this research, we consider the experimental method to be more suitable.This choice is reinforced by the alignment between our experimental results and the theoretical predictions for the signal-to-background ratio and biphoton wavepacket characteristics.Conversely, using the collection efficiencies of 9% and 6% obtained from the guiding method would lead to deviations from the theoretically predicted generation rates.This, in turn, would result in differences between the theoretical signal-to-background ratios and biphoton wavepacket characteristics and their experimental counterparts.

D. Biphoton Temporal Profile
We introduce an alternative method to characterize the generated biphoton by directly calculating the biphoton wave function [11].The SFWM is a process of third-order susceptibility χ (3) : z−ω2t .In our calculations, we have supposed k Upon substituting these electric fields into the interaction Hamiltonian, it can be obtained that The first term in Eq. (S43) can be safely neglected, as the vacuum field is inherently undetectable.Moreover, we can disregard the integrand H I (t ′ )H I (t ′′ ) in the last terms, given its minimal impact attributed to χ (3) .Consequently, the expression for the steady-state biphoton wave function |ψ(t → ∞)⟩ can be stated as follows: It is worth noting that in the calculation of this biphoton wave function, a † s and a † as appear together.This implies that this method can only calculate correlated photons participating in FWM and, as a result, cannot directly provide information about the pairing ratio r p .
By utilizing the steady-state biphoton wave function, we can calculate the Glauber second-order correlation function as where ϕ(t) = ks+kas 2 L − ωas τ − (ω s + ωas )t.By employing the semiclassical model and assuming a large driving detuning ∆ d , we derive the forms of susceptibilities as follows: On substituting the χ (3) into G s-as (τ ), we have s-as (τ ) as follows: The coincidence count rate can be determined using the expression R s-as (τ )∆T .Here, we use a time bin of ∆T = 1 second to calculate the correlated biphoton generation rate by integrating R C (τ ), resulting in R B r p (i.e., dτ R C (τ, 1 s) = R B r p ).However, the formula in Eq. (S50) does not account for accidental receptions, which can result from biphotons originating from different pairs or those that remain uncorrelated.
After obtaining Eq. (S50), we proceed to explore the regime of damped Rabi oscillations [11].This regime tends to dominate in cases of low OD or when the coupling field strength Ω c is large.In this scenario, the term of sinc ∆kL 2 + kasLχ * as (−ω) 4 e i kas Lχas (ω) 4 in the integrand of Eq. (S50) can be approximated as sinc(∆kL/2).To further simplify the expression, we assume γ 21 = 0Γ, leading to the following form: .Because the term of sin 2 (ατ /2) represents the oscillating behavior, the decay processes can be described by e − Γ 2 τ (e βτ /2 −e −βτ /2 ) 2 4 . Furthermore, given that e −βτ /2 decays more severely, the overall decay behaviors depend on the term e −( Γ 2 −β)τ , which yields the 1/e decay time of 1/( Γ 2 −β).Keeping other parameters constant, increasing ∆ c extends the tail of the biphoton wavepacket.This is because the heightened ∆ c weakens the coupling between the coupling field and the atoms, making it more challenging for the coupling field to efficiently convert the pre-established spinwave excitation  Here, we aim to analyze the propagation characteristics of the anti-Stokes field.In this context, ω represents the angular frequency denoted as ω as .Consequently, the group delay time can be expressed as: We consider the specific scenario where ∆ c = 0Γ and ω = 0Γ, allowing us to derive the EIT group delay time τ EIT = ΓOD |Ωc| 2 [12].It is important to note that under two-photon resonance conditions (ω = ∆ c ), τ EIT remains constant regardless of variations in ∆ c .Therefore, relying solely on Eq. (S53), we cannot fully elucidate the influence of ∆ c on the delay time.
Figure S2(a) shows τ EIT corresponding to different values of ∆ c and ω.It is evident from the figure that, regardless of the value of ∆ c , the EIT group delay time is equal to Γ OD |Ωc| 2 when the condition ω = ∆ c is met, as mentioned above.Note that for larger ∆ c values near ω = ∆ c , there is a significant variation in τ EIT , which can lead to significant waveform distortion.Conversely, in the vicinity of ω = ∆ c = 0Γ, there is a relatively flat region of τ EIT , allowing waveforms to propagate through the medium without severely distortion.Figure S2(b) illustrates the transmission of a Gaussian pulse through the EIT medium.The Gaussian pulse has a 1/e 2 full width of 400 ns.Although Eq. (S53) states that τ EIT remains constant at the two-photon resonance condition, the pulse behaviors after passing through the EIT medium are noticeably different for various ∆ c values.It is observed that the tail of the EIT pulse decreases as ∆ c increases, especially in the case of ∆ c = 3Γ.This observation further supports the claim made about the significant variation of τ EIT in Fig. S2(a).
By analyzing the damped Rabi oscillation and the EIT group delay, we can gain insights into the behaviors depicted in Figs.The driving Rabi frequency Ω d in Fig. 4(c) and 4(d) of the main text is 3Γ.The other parameters are listed in Table S1.

FIG. 3 .
FIG. 3. High-purity biphotons.The red lines represent the theoretical curves, while the black circles indicate the experimental data points.The time bin for detecting the anti-Stokes photons is ∆τ = 1.6 ns.The remaining parameters are set to OD = 10, Ωc = 4Γ, ∆ d = 10Γ, γ21 = 0.001Γ, ∆kL = 0.37π, with (a) Ω d = 0.5Γ and (b) Ω d = 3Γ.(c) The peak signal-tobackground ratio rSB versus Ω d .The black squares represent the experimental data, and the black line is the curve fitted to these experimental data points.(d) The biphoton generation rate RB and pairing ratio rp as a function of Ω d .The experimental data points for RB and rp are represented by the unfilled blue and solid magenta circles, respectively.The theoretical curves for RB and rp are depicted by the blue and magenta lines, respectively.

FIG. 4 .
FIG. 4. Highly correlated ultrabright biphotons.The red lines represent the theoretical curves, while the black circles indicate the experimental data points.The time bin for detecting the anti-Stokes photons is ∆τ = 1.6 ns.The remaining parameters are set to Ω d = 3Γ, γ21 = 0.001Γ, ∆kL = 0.37π, with (a) OD = 20, Ωc = 4Γ, ∆ d = 5Γ, and (b) OD = 120, Ωc = 8.8Γ, ∆ d = 14.9Γ.(c) The peak signal-to-background ratio rSB versus OD.The black squares represent the experimental data, and the black line is the curve fitted to these experimental data points.(d) The biphoton generation rate RB and pairing ratio rp as a function of OD.The experimental data points for RB and rp are represented by the unfilled blue and solid magenta circles, respectively.The theoretical curves for RB and rp are depicted by the blue and magenta lines, respectively.
c(d) , induce the generation of Stokes and anti-Stokes photons, which are treated as quantum fields.Here, g s = d 32 ωs /2ℏϵ 0 V is the coupling constant between the Stokes field and atoms, where d 32 indicates the dipole moment of the |2⟩ ↔ |3⟩ transition, ωs is the center angular frequency of Stokes field, ϵ 0 is the vacuum permittivity, and V is the interaction volume.The coupling constant g as = d 41 ωas /2ℏϵ 0 V between the anti-Stokes field and atoms is described in the same manner.The operator âs(as) is the slowly varying annihilation operator of the Stokes (anti-Stokes) field.The phase shift ∆kz ≈ ( ⃗ k d − ⃗ k s + ⃗ k c − ⃗ k as ) • ⃗ z represents the phase-mismatching value in the 1D approximation.All relevant HLEs

|gs| 2 c
ϵ 23 + iω c .The first component, iN |gs| 2 c ϵ 23 , characterizes the Raman process, while the second component, iω c , represents the phase term arising from free evolution.The electromagnetically induced transparency (EIT) profile coefficient, denoted as Γ as , takes the form: Γ as = iN |gs| 2 c η 41 + i∆k − iω c .In this expression, the first component, iN |gs| 2 c η 41 , describes the EIT characteristics, and second term, i∆k, represents the phase-mismatching effect.Furthermore, we have the Stokes (anti-Stokes) coupling coefficients κ s and κ as , defined as: κ s = iN g * s g * as c η 23 and κ as = iN gsgas c ϵ 41 .The coefficient κ s (κ as ) quantifies how the anti-Stokes (Stokes) field affects the Stokes (anti-Stokes) field.
) D. Proof of Thermal States for Stokes and Anti-Stokes Fields

Figures 2 (
a)-2(d), 3(a), 3(b), 4(a), and 4(b) in the main text present results accumulated over four data points, corresponding to 2 18 receptions in the Stokes channel.The error bars represent the standard deviation across these four data points.
FIG. S1.The coincidence count rate corresponds to the parameters outlined in Fig. 3(b) of the main text.The red solid curve depicts the theoretically predicted coincidence count rate, while the black dashed line represents the uncorrelated background, which is equivalent to the biphoton generation rate.
2(c) and 2(d) of the main text.In Fig. 2(c) of the main text, we introduce ∆ c = 1Γ, which leads to a reduction in τ EIT , as indicated by the red curve in Fig. S2(b).Increasing ∆ c to 3Γ, as demonstrated in Fig. 2(d) of the main text, further decreases τ EIT due to the reduced effective OD.However, introducing ∆ c = 3Γ simultaneously extends the 1/e decay time 1/( Γ 2 − β) of the damped Rabi oscillation.As a result, the elongated tail in Fig. 2(d) of the main text is attributed to the regime of damped Rabi oscillation.E. Experimental Parameters of Figure 4

9 TABLE S1 .
Parameters in Fig. 4(c) and 4(d) of the main text.