Modulation spectroscopy of Rydberg atoms in an optical lattice

We develop and study quantum and semi-classical models of Rydberg-atom spectroscopy in amplitude-modulated optical lattices. Both initialand target-state Rydberg atoms are trapped in the lattice. Unlike in any other spectroscopic scheme, the modulation-induced ponderomotive coupling between the Rydberg states is spatially periodic and perfectly phase-locked to the lattice trapping potentials. This leads to a novel type of sub-Doppler mechanism, which we explain in detail. In our exact quantum model, we solve the time-dependent Schrödinger equation in the product space of center-of-mass (COM) momentum states and the internal-state space. We also develop a perturbative model based on the band structure in the lattice and Fermi’s golden rule, as well as a semi-classical trajectory model in which the COM is treated classically and the internal-state dynamics quantum-mechanically. In all models we obtain the spectrum of the target Rydbergstate population versus the lattice modulation frequency, averaged over the initial thermal COM momentum distribution of the atoms. We investigate the quantum-classical correspondence of the problem in several parameter regimes and exhibit spectral features that arise from vibrational COM coherences and rotary-echo effects. Applications in Rydberg-atom spectroscopy are discussed.


I. INTRODUCTION
The interaction of an electron with an electromagnetic field consists of a termê A ·p/m and a term e 2Â2 /(2m), with electron mass m, elementary charge e, the field's vector potential A(r), and electron position and momentum operatorsr andp [1,2]. Under certain conditions, periodically modulated and inhomogeneous fields can drive electronic transitions via a ponderomotive interaction, e 2Â2 /(2m). These transitions can occur in Rydberg atoms immersed in light fields, a case in which the field frequency (hundreds of THz) exceeds the Rydberg atom's evolution frequency (tens to hundreds of GHz) by several orders on magnitude. Due to the quasi-free nature of the Rydberg electron on optical time-and energyscales, the Rydberg-atom ponderomotive effect is related to free-electron Kapitza-Dirac scattering [3][4][5]. Further, ponderomotive level shifts in high-intensity laser fields were observed in atoms [6,7] and in molecules [8], as well as in high-intensity multi-photon ionization [9], in highintensity zero-kinetic-energy photoelectron (ZEKE) spectroscopy [10], and in optical [11,12] and microwave [13,14] above-threshold ionization. Ponderomotive forces on atomic electrons are important in atom dynamics in highintensity laser pulses [15][16][17]. Ponderomotive effects are also known from Paul ion traps, where the ponderomotive force drives the secular motion, while the micromotion occurs at the trap's radio-frequency drive [18].
In Rydberg atoms, ponderomotively driven transitions (transitions driven by the e 2Â2 /(2m) operator) are free of multipole selection rules [19] that govern traditional methods of laser and microwave spectroscopy (which are * Corresponding author: graithel@umich.edu based on the properties of theê A ·p/m operator in first or higher orders). For ponderomotive spectroscopy to be effective, the field must be modulated at a (sub-)harmonic of the transition frequency, and the field intensity must vary within the extent of the electron wavefunction. These conditions are quite naturally satisfied by Rydberg-atom transitions in amplitude-or phasemodulated optical lattices [19], because Rydberg atoms have sizes on the order of typical optical-lattice periods [20]. Further, microwave amplitude and phase modulators for optical-lattice lasers are readily available. In addition to driving microwave transitions, the ponderomotive interaction can serve as a tool to trap the Rydberg atoms in a ponderomotive optical lattice (POL) [21][22][23][24][25]. Hence, the modulated POL can satisfy two functionalities at once: it can trap the Rydberg atoms and, at the same time, serve as a spectroscopic probe for a wide variety of Rydberg transitions [26]. POLs also offer great flexibility in designing Rydberg-state-mixing properties [27][28][29] and magic-transition traps, where two or more states have (near-) identical trapping potentials [30]. In these applications perturbing effects from Rydberg-atom photoionization [31][32][33][34][35] are typically irrelevant, in particular for S-type Rydberg states of rubidium and for highangular-momentum Rydberg states of any species.
The combination of the aforementioned features enables high-precision spectroscopy on long-lived circularstate Rydberg atoms trapped in optical lattices [36], utilizing a scheme in which the transition frequency between the circular Rydberg levels is measured via resonant POL-modulation at microwave frequencies. This experimental platform may also be useful for quantum simulators [37] that are based on circular-state Rydbergatom arrays. For such applications, the effects of latticepotential-induced level shifts and interaction-time broadening on the achievable spectroscopic resolution have to arXiv:1911.05671v1 [quant-ph] 13 Nov 2019 be reduced. To that end, it is important to quantitatively model vibrational sidebands, anharmonic corrections, and vibrational quantization and tunneling in POL modulation spectroscopy. This necessitates a fully quantized description of the center-of-mass (COM) and internalstate dynamics of the atoms. In the present paper we develop such models, investigate quantum-classical correspondence in POL modulation spectroscopy, and exhibit the quantum features in the spectra. (right) in identical, modulated, ponderomotive optical-lattice potentials. The spectrum consists of a Doppler-free central peak (change in vibrational quantum number ∆ν = 0) and vibrational sidebands (∆ν = ±2). In the depicted case, the anharmonicity-induced substructure of the sidebands is spectroscopically resolved.
In Fig. 1 we illustrate several quantum aspects of POL modulation spectroscopy. The curves with the light-gray drop areas visualize the sinusoidal COM lattice potentials versus position for a pair of Rydberg states |1 and |2 ; the figure shows the case of a magic POL transition for a lattice depth of 500 kHz. The depicted band structure of the COM dynamics is for 85 Rb atoms in a POL formed from two counter-propagating beams of 1064 nm wavelength. Amplitude modulation (AM) of the lattice drives electronic transitions the Rabi frequencies of which have previously been calculated [19]. While POL modulation spectroscopy has no selection rules for l − l (l and l are the angular-momentum quantum numbers of Rydberg states |1 and |2 , respectively), for even-parity transitions (|l − l| = 0, 2, . . .) the vibrational quantum number ν of the COM motion can only be changed by even numbers, whereas for odd-parity transitions (|l − l| = 1, 3, . . .) it can only be changed by odd numbers. The case depicted in Fig. 1 is for even-parity transitions. The spectrum, sketched on the right, shows the transition probability from |1 into |2 as a function of the detuning of the POL modulation frequency, δ, from the atomic transition frequency. The spectrum has red-shifted (∆ν = ν − ν = −2), unshifted (∆ν = 0) and blue-shifted (∆ν = 2) spectral components, which arise from transitions of the lattice-trapped atoms. The substructure of the ∆ν = ±2 components is a COM quantum effect that results from the anharmonicity of the POL potential and that requires sufficient spectral resolution to be observed. For odd-parity transitions, the spectrum would have two major components, corresponding to changes of the vibrational quantum number by ∆ν = ±1.
The theory leading to spectra as sketched in Fig. 1 is developed and discussed in Secs. III-V. In our models both ground and excited wave functions evolve on one-dimensional, sinusoidal optical-lattice trapping potentials [21] with generally different depths but fixed relative spatial phase (the experimentally most relevant case). The ponderomotive coupling between the electronic states that arises from lattice modulation is described by an effective Rabi frequency, Ω(z), that depends on the atom's COM position z transverse to the lattice planes. The novelties described in this work depend critically on the fact that Ω(z) has a sine-like dependence on z, with the same spatial period as the lattice itself [19]. The spatial phase between the ponderomotive coupling and the lattice depends on whether the modulation-driven transition is between Rydberg states with same or opposite parity. The phase of the ponderomotive coupling Ω(z) exhibits a quite peculiar behavior, as it proceeds in discrete steps of π as a function of the coordinate z. This behavior differs radically from the optical phase of typical plane-wave or Raman couplings, ∆k · r, which is a continuous function of position ( ∆k is the photon momentum transfer). In Sec. VI B we explain why the peculiar phase behavior of Ω(z) in modulated POLs leads into a new paradigm of sub-Doppler spectroscopy.
To obtain the spectrum of the excited-state population as a function of lattice modulation frequency, we numerically solve the time-dependent Schrödinger equation (TDSE) in momentum representation of the COM. The results of these numerical solutions are averaged over the thermal COM momentum distribution of the Rydberg atoms. In addition, we analyze the band structure of the problem. The band structure is employed to model the POL modulation spectra with transition rates between Bloch states, averaged over an initial thermal momentum distribution of atoms loaded into the lattice. The solutions presented in this work account for quantum features such as band structure, vibrational quantization, band curvature and tunneling, and rotary-echo effects. In addition, we discuss the convergence of quantum and semi-classical results in the appropriate limits. In applications, the advanced modeling afforded by our work will enable a reduction of systematic errors caused by latticeinduced shifts of the Rydberg-atom transition frequency.

A. Position representation
We consider a Rydberg atom moving in a onedimensional POL formed by two counter-propagating laser beams of equal polarization, wavenumber k, and wavelength λ = 2π/k. The atom is initially prepared in state |1 . The optical lattice is amplitude-modulated at a frequency that effects ponderomotively driven transitions into Rydberg state |2 . Lattice potentials and Rabi frequencies for this scheme have been derived in [19,21]; here we recite relevant, previously-proven findings. The lattice potentials for the two Rydberg states can be written in the form V i,0 + V i cos(2kz), with i = 1, 2 and constants V i,0 ≥ V i . Here, it is sufficient to consider Rydberg levels |n, l, j, m j without lattice-induced state mixing [19,21]. Assuming azimuthal symmetry and choosing the quantization axis along the direction of the laser beams, amplitude modulation of the lattice can generally drive, in first order, any transition |n, l, j, m j ↔ |n , l , j , m j with m j − m j = ∆m j = 0, and with no other applicable selection rules [19]. To drive the transitions efficiently, the lattice modulation frequency must be close to the atomic transition frequency [26] or a sub-harmonic [30]. The Rabi frequencies for ponderomotive transitions driven by amplitudemodulated lattices [19] take the form Ω cos(2kz) for evenparity (|l − l| = 0, 2, . . .) and Ω sin(2kz) for odd-parity (|l − l| = 1, 3, . . . ) Rydberg transitions.
We define spinor wavefunctions Ψ i (z) via where the base kets |1, z are in the product space of the internal (Rydberg) state space {|1 , |2 } and the position space {|z } of the z-COM degree of freedom. We further define an effective atom-field detuning δ = δ 0 +V 2,0 −V 1,0 , where V i,0 are the lattice-potential offsets defined above, and δ 0 is the lattice-free atomic transition frequency, ω A , minus the optical-lattice modulation frequency, ω L , or its relevant overtone, p ω L (p = 2, 3, . . .). In the most generic case, considered here, only the states |1 and |2 are close to resonance, while for all other transitions Ω δ. Hence, in a dressed-atom picture the near-resonant atomfield states are |1, N 0 and |2, N 0 − p , where N 0 is a constant number of spectator photons/phonons in the modulator crystal, and p is the modulation order that drives the transition. In a holistic picture, p can be interpreted as the number of light-modulator energy quanta (photons/phonons) that are absorbed in the POL-modulation transition. For the cases studied here, we neglect natural and black-body Rydberg-atom decay.
With these definitions, the Schrödinger equation for even-parity (|l − l| = 0, 2, . . . ) Rydberg transitions is where the spinor wavefunctions Ψ i (z) are now in the dressed-atom picture (rotating-frame). We write detunings and couplings in units rad/s. For odd-parity (|l − l| = 1, 3, . . . ) Rydberg transitions it is (3) The potential depths V i are generally different. There exist cases in which V 1 = V 2 ; these are "magic" lattices that lead to particularly narrow spectral lines. Magic lattices are well known from optical clocks (see, for instance, [38,39]). The same concept translates to spectroscopy in modulated POLs for Rydberg atoms. Magic transitions in modulated POLs require a suitable combination of lattice period and Rydberg levels; such transitions have already been demonstrated [30]. Magic lattices are particularly useful for high-precision spectroscopy because they minimize lattice-induced shifts of the transition frequency to be measured. Several examples discussed below are for magic transitions, where V 1 = V 2 . We also consider a generic case in which the lattices for |1 and |2 have different depths, For even-parity transitions, the Rabi frequency and the lattice potentials all share the same spatial modulation ∝ cos(2kz), while for odd-parity transitions the lattice potentials are ∝ cos(2kz) and the Rabi frequency is ∝ sin(2kz). In either case, the drive term Ω(z) is real and alternates between positives and negatives, amounting to discrete phase jumps of π at every π/(2k) step in z. This uncommon behavior greatly differs from the Rabifrequency behavior of first-and higher-order multipole transitions effected by theê A ·p/m term, which (in onedimensional cases) typically is Ω m exp(i∆k z), where ∆k denotes the recoil momentum and Ω m the Rabi frequency of the transition. We see that in the modulated-POL case the spatial phase of the drive follows a (real-valued) staircase function and the magnitude of the coupling varies in z as |Ω cos(2kz)| or |Ω sin(2kz)|, whereas in the latter case the spatial phase is a linear function with fixed slope ∆k and the magnitude of the coupling, |Ω m |, is fixed. In Sec. VI B we show that these facts enable a novel type of sub-Doppler method that is realized automatically in modulated-POL spectroscopy.
The Hamiltonian in the above Schrödinger equation for the even-parity case can be conveniently written in matrix formĤ where V ± = (V 2 ± V 1 )/2,σ i are the Pauli operators,Î is the identity operator, andẑ is the position operator for the z-component of the COM motion. Similarly, for odd-parity transitions it iŝ

B. Momentum representation
Due to the periodicity of potentials and couplings, the TDSE is most conveniently solved in the momentum representation [40,41]. Using we obtain the equations for the momentum-space spinor wavefunctions for the even-parity case, and similarly for the odd-parity case From these equations it is seen that in momentum representation the Hilbert space breaks up into subspaces . Subspaces with different p 0 do not couple to each other. Using the initial condition Ψ 1 (p, t = 0) = δ(p − p 0 − 2n 0 k) and Ψ 2 (p, t = 0) = 0, the evolution is restricted to the subspace for p 0 . We can then set and write Eq. (7) in matrix form as where For even parity we have . . . 0 1 0 while for odd parity we have Defining p 0 = β k and the two-photon recoil frequency ω k = 2 k 2 /m, the lattice-free energies in Eq. (11) can be written as where n = 0, ±1, ±2, · · · . For 85 Rb atoms in an optical lattice formed by a pair of counter-propagating 1064-nm laser beams, the case studied below, the two-photon recoil frequency is ω k = 2π × 8.300 kHz = 52.15 krad/s. From the above analysis we observe that states with momentum p couple only to the neighboring states p ± 2 k for both Ψ 1 (p) and Ψ 2 (p). The internal-state-conserving momentum coupling strengths are V 1 /2 and V 2 /2 for states |1 and |2 , respectively. The POL modulation induces internal-state-changing couplings between initial-and target-level momentum states that also differ by 2 k; those couplings have a strength determined by Ω.

C. Averaging over the thermal momentum distribution
For reference, we state the one-dimensional Maxwell distribution in velocity, momentum and our dimensionless scaled momentum β, The scaled-momentum distribution, f (β), can be written in terms of the width of the energy distribution in units of the thermal energy. Defining and ω T /ω k = 2.5. Hence, lattices and atom distributions that are about T = 1 µK deep or wide, respectively, exhibit quantum features such as tunneling and photonrecoil effects. Conversely, for T = 100 µK, ω T = 2π × 2.08 MHz and ω T /ω k = 250; for this temperature regime and lattice depth we expect convergence of classical and quantum treatments of the COM. We obtain the excited-state coefficients b n by solving the above equations for given lattice parameters over a range of detunings δ and over a range of initial scaled atom momenta β. For a given β, the initial condition for solving the TDSE, Eqs. (10), (11), and (12) or (13), is b n = 0, a n = δ n,n0 with n 0 = NINT(β/2) and p 0 = k(β − 2 n 0 ), where the function NINT(x) is the integer closest to x. The range of β is adapted to temperature; we use ranges β ∈ [−10, 10] for T = 1 µK and β ∈ [−80, 80] for T = 100 µK. For the time-dependence of the drive, Ω(t), we use Gaussian pulses with durations of tens to hundreds of µs. The drive pulse Ω(t) peaks at time t p /2, and the TDSE is integrated from t = 0 to t p . The target-state population is evaluated at t = t p and averaged over the Maxwell momentum distribution, yielding In the calculation, the duration t p is set sufficiently large such that the population P b (δ) becomes independent of t p . The population P b (δ) is equivalent to experimental spectra in which the Rydberg target state population is measured as a function of lattice modulation frequency, and is therefore a main result of this paper.

IV. BAND STRUCTURE MODEL
The above described method of finding the spectrum P b (δ) has the advantage that drive pulses Ω(t) of any time dependence can be modeled. Further, the spectrum usually exhibits some Fourier broadening, which occurs when the drive pulse is substantially shorter than the inverse Rabi frequeny (t p π/Ω), or a certain amount of saturation broadening, which occurs when the pulse area exceeds π (that is, if t p Ω π). Solving the TDSE is a good way to account for these effects.
For an approximate and fast solution, we use perturbation theory (Fermi's golden rule, FGR) to obtain transition probabilities between the Bloch states of the initial and target Rydberg levels in their respective opti-cal lattices (which have identical periods and Brillouin zones). The transition probabilities are averaged over the Maxwell momentum distribution along the same lines as described in Sec. III C. These results are a good approximation if the transition is not saturated (t p Ω π) and if COM coherences are not important, as discussed further below. In the following, the FGR method is briefly outlined.
Diagonalizing the Hamiltonian in Eq. (11) for both Rydberg levels |1 and |2 separately, we find the initial and target Bloch states where A and B stand for the initial and target Rydberg levels, respectively, i and i are band indices, and the quasimomentum p 0 ∈ [− k, k]. We typically obtain the Bloch states and their respective band energies, E a (i, p 0 ) and E b (i , p 0 ), on an equidistant grid of 200 p 0 -values. The squares of the transition matrix elements for a timeindependent Rabi frequency Ω in Eqs. (12) and (13) are for even-and odd-parity transitions, respectively, and with i denoting a target-and i an initial-level band index. The p 0 -dependence of the a-and b-coefficients is suppressed for brevity. Following FGR for the case of a harmonic drive, the transition rate from the initial Bloch state i into the target Bloch state i then is with energy detuning ∆E of the lattice-modulation drive from the transition energy between the initial and target Bloch states, For the energy density of states ρ(∆E) we use a Gaussian, ). There, the spectral width of the drive is chosen in the range 1 kHz σ E /h 100 kHz, in accordance with the Fourier widths for our typical drive-pulse durations.
If the drive pulse was a square pulse, the (non-saturated) FGR transition probability would be P (i , i, p 0 , δ) = R(i , i, p 0 , δ) t p , with pulse duration t p . The transition rate depends on time if the Rabi frequency Ω in Eq. (18) depends on time.
In the examples discussed in Sec. VI, the Rabi frequency Ω(t) has a Gaussian time dependence. In that case, P (i , i, p 0 , δ) is given by the integral of R(i , i, p 0 , δ, t) over the duration of the pulse.
It is then assumed that the initial COM states have a normalized thermal probability distribution for temperature T , W (n 0 , p 0 ) ∝ exp(−(p 0 + n 0 k) 2 /(2mk B T )), with p 0 ∈ [− k, k] and integer n 0 . For a sudden lattice turn-on, the thermally populated momentum base states |p 0 + 2n 0 k are projected into the basis of Bloch states, where the expansion coefficients a i,n0,p0 are known from the diagonalization of the Hamiltonian in Eq. (11). The FGR spectrum is then obtained from The quantity Q b (δ) is the FGR transition probability per atom, averaged over the initial thermal distribution of atoms over free-particle momentum states, for the case that the amplitude-modulated POL is suddenly turned on. This corresponds with the analysis performed in Sec. III.

V. SEMI-CLASSICAL MODEL
In the limit of temperatures and potential depths equivalent to energies much larger than the lattice recoil energy, ω k /4, semi-classical and quantum results should converge in certain aspects. To explore classical-quantum correspondence in the system, we use a model [26,30] in which the COM dynamics is treated classically by solving Newton's equations with a 6th-order Runge-Kutta routine. As the atoms move, the internal-state dynamics in the Rydberg-state space {|1 , |2 } is propagated quantum-mechanically, taking the explicit and implicit time dependence of the Rabi frequency into account. For even-parity transitions, for instance, the Rabi frequency in the atomic frame follows Ω(t) cos(2kz(t)), where Ω(t) is the real-valued, positive Rabi frequency from Eqs. (2) and (7), which is a Gaussian pulse in our examples, and z(t) is the classical COM position of the atom. Doppler effects arise from the sign flips of the cosine function that occur when the atoms pass through the inflection points of the cosine function. The sign flips are equivalent to phase jumps of π, which, when occurring at semi-regular time intervals in the moving frame of a hot atom moving through the lattice, cause a Doppler effect similar to the usual Doppler effect (see discussion in Sec. VI B). The internal evolution also depends on the detuning between the atomic transition and the POL AM modulation frequency, δ. In cases of non-magic lattices (V 1 = V 2 ) in the frame of reference of a moving atom with classical COM position coordinate z(t), the detuning further depends on position; in that case, the detuning for the internal quantum evolution is δ c = δ + (V 2 − V 1 ) cos(2kz(t)).
As the states |1 and |2 have generally different POL amplitudes V 1 and V 2 , a procedure is needed to compute a classical force. Here, we average the classical force over the two lattice potentials for the states |1 and |2 , using the time-dependent quantum-mechanical probabilities of the atom being in |1 or |2 as weighting factors. This method does not account for COM quantum effects such as tunneling, dispersion, and statedependent wave-packet splitting. This deficiency ties into the overall failure of the semi-classical description at low energies, and in certain other cases. Nevertheless, over wide swaths of parameter space the semi-classical model is quite successful. It is noted that in magic lattices (V 1 = V 2 ) the two coupled internal states have identical potentials, in which case the classical force simplifies to F (z(t)) = 2kV 1 sin(2kz(t)).
In accordance with the quantum models, in the semiclassical description initial atom velocities v(t = 0) are drawn from the Maxwell distribution given in Sec. III C, and initial positions z(t = 0) are random. The semiclassical spectrum, K b (δ), is given by the transition probability at the end of the interaction time t p , averaged over a classical thermal ensemble of initial atoms. In this work, we average over 10,000 to 100,000 trajectories.

VI. RESULTS
In general, it is desired to strike a balance between observing a transition with high signal-to-noise ratio and avoiding saturation broadening. Therefore, in our calculations we use Rabi frequencies that, for a given drive pulse shape and duration, yield a π pulse area for an atom at a COM location z where the Rabi frequency is maximal, Ω(t)dt = π. For our Gaussian pulses, Ω(t) = Ω 0 e −(t−tp/2) 2 /τ 2 0 . In the cases discussed below we quote the full width at half maximum of the pulse, which is τ FWHM = 2τ 0 ln(2), and the corresponding utilized value of Ω 0 that leads to the π pulse.

A. Structure of POL amplitude-modulation spectra
The range of pulse durations of interest is between 1 µs and 1 ms, as this range is accessible given the lifetimes of typical Rydberg states [20]. Within this range, we consider a pulse short if its Fourier bandwidth suffices to resolve the vibrational sidebands ∆ν = 0 in the spectrum from the central ∆ν = 0 peak, but the bandwidth is too large to resolve the recoil energy and anharmonic effects of the band structure of the COM (see Fig. 1). Conversely, a pulse is considered long if it resolves the anharmonicity-induced quantum structures in the vibrational sidebands. For two selected cases of short and long pulses, we will discuss spectra for several representative lattice depths and COM atom temperatures. Our initial discussion is focused on the TDSE and semi-classical models.
In Fig. 2 we show POL modulation spectra P b (δ), Q b (δ) and K b (δ), obtained by solving the TDSE, by computing FGR transition probabilities between Bloch states, and by performing semi-classical simulations, respectively. The transitions are even-parity and are driven by Gaussian field pulses with τ FWHM = 20 µs, corresponding to τ 0 = 12.0 µs. Our calculations are for ω k = 8.300 kHz, corresponding to 85 Rb atoms in a POL formed from two counter-propagating beams of 1064 nm wavelength, and lattice amplitudes V 1 = V 2 = 2π × 1.25 MHz (a case of a magic lattice). These conditions are similar to experimental work in [26,30]. The spectra exhibit the lowest-order allowed vibrational sidebands, ∆ν = ±2, as well as a Doppler-free central band ∆ν = 0. Since the COM oscillation frequency of the atoms at the bottoms of the wells is 2k V 1 /m = 2π × 144 kHz, the frequency separation between the ∆ν = ±2 vibrational sidebands is about 69 ω k /(2π). Inspecting Fig. 2, this corresponds to the separation between the outer fringes of the sidebands. The average separation is somewhat less, because atoms that are thermally excited into vibrational states of the COM motion above the vibrational ground state have a smaller oscillation frequency, leading to less separation. Anharmonicity-induced substructures are not resolved due to the FWHM width of the spectral density of the pulse, which is 3.77ω k = 2π × 31.3 kHz (in the case of low saturation, for the given pulse length). Also, higher-order sidebands, ∆ν = ±4, ±6, ... are too broad and weak to be observed.
For measurement and metrology purposes, the Doppler-free peak for ∆ν = 0 is of particular interest, which in Fig. 2 has a FWHM of 4ω k . This is in good agreement with the low-saturation Fourier width of the Gaussian pulse, 3.77ω k . Note that the presented case is for a π-pulse, and therefore the central peak is slightly saturation-broadened.

B. Doppler-free spectroscopy in modulated POL
The origin of the Doppler-free nature of the central peak becomes clear from the following semi-classical analysis. The classical COM of an atom oscillating within the center portion of a lattice well is approximately harmonic, z(t) = z 1 cos(2πf 1 t), with COM oscillation frequency f 1 = k V 1 /m/π and COM oscillation amplitude z 1 λ/8. Using the Jacobi-Anger expansion, the time dependence of the Rabi frequency in the frame of the atom is given by Ω a (t) = Ω(t) cos(2kz(t)) = Ω(t)[J 0 (2kz 1 ) cos(2kz 1 ) Here, the explicit time dependence is contained in the Gaussian envelope function Ω(t) = Ω 0 e −(t−tp/2) 2 /τ 2 0 with constant Ω 0 . The implicit time dependence in Eq. (20) arises from the atomic motion, z(t), and is contained in the Fourier series in the square brackets. The atomic motion generates Fourier components of the drive at even multiples of f 1 , which correspond with the vibrational signals for ∆ν = 0 and ±2 in Fig. 2. For the case of Fig. 2, Ω(t) is a slowly-varying real-valued envelope function with a duration of a few COM oscillation periods of the atoms in the wells. The ∆ν = 0 and ±2 components are resolved in Fig. 2 because the Fourier width of the envelope Ω(t) is less than the frequency separation 2f 1 between the components. This basic interpretation applies if a substantial fraction of atoms is trapped within the approximately harmonic regions of the POL wells, and if the spectrum of the pulse envelope Ω(t) is sufficiently narrow to resolve the vibrational sidebands. The treatment based on the band structure and FGR, visualized in Fig. 1, yields equivalent conclusions. Importantly, in Eq. (20) it is evident that the atom velocity plays no direct role in the spectrum. In particular, in magic lattices as in Fig. 2 the cental peak ∆ν = 0 is Dopplereffect-free and is unaffected by the anharmonicity of the lattice, making it ideal for high-precision spectroscopy of Rydberg-atom transitions. It is noted that Eq. (20) resembles the spectrum of frequency-modulated fields and the spectrum seen by optically driven ions oscillating in ion traps. These similarities can be born out more clearly in a photon picture of POL modulation spectroscopy, in which the atoms scatter a photon from one lattice-field mode into a counterpropagating field mode, where the mode frequencies differ by the POL modulation frequency (which equals the atomic-transition frequency and is orders of magnitude larger than f 1 and ω k ). The re-scattering is a stimulated process effected by the A 2 -term of the atom-field interaction. A detailed analysis of this picture is not of interest in the present paper.
Spectral broadening akin to the usual Doppler effect arises from atoms traversing over many lattice wells. Again, a semi-classical picture is well-suited to explain this effect. The phase of the drive term in Eq. (2), Ω(t) cos(2kz), undergoes a jump of value π at every inflection point of the lattice potential. For a hot atom moving at constant velocity v through the lattice, in the reference frame of the moving atom the phase of the drive field follows a step-function that is centered around the linear function φ(t) ≈ 2πz(t)/(λ/2) = 2kvt, equivalent to a Doppler shift of 2kv. This resembles the Doppler shift of stimulated Raman scattering between counterpropagating beams. It is noted, however, that POL modulation spectroscopy is fundamentally different from Raman spectroscopy, because it employs a first-order A 2process and not a second-order A · p-process.

C. Temperature insensitivity of spectroscopy in modulated POL
If the lattice depth 2 V 1 k B T , the initial temperature has only a minor effect on the spectrum. This is seen clearly in Fig. 2, where the spectra for T = 1 µK and 100 µK are quite similar. The heights and widths of the central peaks are near-identical in both cases. The temperature insensitivity results from the fact that atoms initially located away from a lattice minimum gain considerable potential energy at time t = 0, when the lattice is suddenly turned on. If 2 V 1 k B T , the initial potential energy dominates the initial kinetic energy, and the distributions of oscillation amplitudes z 1 are not very different (in the semi-classical model). Hence, for 2 V 1 k B T the signal strengths in the vibrational sidebands, ∆ν = ±2, are not very temperature-dependent. In the TDSE and band-structure models, the temperature insensitivity follows from the fact that the projection of a thermal ensemble of plane waves of the COM motion into the Bloch-state basis yields similar distributions as long as the potential energy is larger than the thermal energy of the ensemble prior to projection, leading to the same condition, 2 V 1 k B T .

D. Quantum-classical correspondence
Comparing the exact TDSE and the semi-classical spectra in Fig. 2, it is seen that the semi-classical model does quite well, even at a quantitative level. A leading deviation between TDSE and semi-classical results is that the vibrational sideband structure (for magic transitions) is symmetric in the semi-classical calculation, whereas it is asymmetric in the TDSE result. The symmetry in the semi-classical case directly follows from Eq. (20), where positive-and negative-frequency components of the drive have identical amplitude. After dropping the harmonic-COM approximation made in Eq. (20), this symmetry still holds. It is also noted that in the semi-classical model the red-and blue-detuned transitions between the Rydberg states |1 and |2 have no recoil effect on the classical COM motion, in accordance with the perfect symmetry between blue-and red-detuned vibrational sidebands in the semi-classical results.
In contrast, in both quantum treatments (TDSE and FGR) the spectra are non-symmetric. Using, for the sake of clarity, the notion of a harmonic COM motion, theẑ 2term in the expansion of the cosine in the drive term, Ω(t) cos(2kẑ), causes (most of) the ∆ν = ±2 vibrational sidebands in the spectrum. While the relevant matrix elements between COM states, ν −2|â 2 |ν and ν|â †2 |ν −2 , are symmetric, the population difference between the initial vibrational COM states |ν − 2 and |ν reduces the red-shifted sideband relative to the blue-shifted one by a Boltzmann factor of about exp(−2hf 1 /(k B T eff )), where T eff is an effective COM temperature after the sudden transfer of the atoms into the lattice. Generally, the vibrational spectral sidebands are the more asymmetric the lower the temperature and the shallower the lattices are.

E. Exact TDSE solution versus perturbative model
In Fig. 2, the perturbative FGR-model based on the band structure reproduces the TDSE results fairly well. As our drive pulses are near saturation of the transition, the FGR model generally overestimates the transition probabilities. This is qualitatively adjusted by scaling the FGR probabilities by a factor < 1. The FGR result exhibits asymmetries similar to those in the exact TDSE solution, as expected (see Sec. VI D). There are, however, deviations between the shapes of the FGR and the TDSE results. This is attributed to the fact that the FGR model does not account for any coherent transients, which can cause strong effects in the case of short drive pulses. This includes COM transients, which are generated by the sudden transfer of the atoms from free space into the lattice, which initiates a COM wave-packet motion. Due to the position dependence of the drive term, Ω(t) cos(2kz), the COM wave-packet maps onto a transient signal in the internal-state dynamics (which is our observable). The transients are the most pronounced if the pulse duration is on the order of the harmonic period of the COM motion, 1/f 1 . The cases discussed in the present and even more so in the next subsection are in this regime. While our FGR model cannot describe the transients, the transients are reproduced in large parts in the semiclassical model. This is not unexpected, because the semi-classical model incorporates a classical approximation of COM dynamics, the origin of the transients. Generally, the transients-related deviations between the FGR and the TDSE and semi-classical results are most visible in shallow lattices and at high temperature. In Fig. 2 the POL is deep enough that the deviations between the FGR lineshapes and the TDSE/semi-classical lineshapes are still fairly minor.

F. POL amplitude-modulation spectra in shallow lattices
Spectroscopic data often yield better results in less deep lattices due to a reduction in residual AC shifts, reduced photo-ionization losses, etc. In Fig. 3 we show results for the same parameters as in Fig. 2, with the exception that the lattice is only 1/100-th as deep. In this case, the POL has only two tightly-bound bands and an f 1 -value of only 14.4 kHz, corresponding to a round-trip time of 70 µs.
At T = 1 µK, quantum and semi-classical results still agree fairly well. The asymmetry of the TDSE result is much larger in the shallow lattice (Fig. 3) than it is in the deep one (Fig. 2). The asymmetry is somewhat reproduced in the FGR result, whereas the semi-classical spectra are perfectly symmetric and show a hint of resolving the ∆ν = ±2 sidebands. These trends in overall agreement and asymmetry behavior are expected from the discussions in the previous sections. Deviations between TDSE and FGR results are again attributed to the effects of COM transients on the internal-state dynamics.
At T = 100 µK, the TDSE and semi-classical results ( Fig. 3 (b)) agree very well, as the atoms are comparatively hot, leading to essentially classical COM dynamics. It is noteworthy that both results exhibit a central dip with a width near the Fourier limit. The dip is attributed to a rotary-echo behavior, which is most significant when τ 0 ∼ 1/(4f 1 ), as is the case in Fig. 3 (b). To explain the rotary-echo effect, we consider a drive-pulse duration on the order of 1/(4f 1 ). In that case, the classical motion of a large fraction of atoms covers about one half of the lattice period, over the duration of the drive pulse. Due to the cosine-dependence of the Rabi frequency, the atoms spend similar amounts of time in spatial regions with positive and with negative-valued Rabi frequencies.
The FGR model entirely fails to produce the echoinduced dip seen in the exact TDSE solution and the semi-classical result. This is due to the fact that our FGR model does not account for COM wave-packet dynamics and COM coherences at all. Any rotary-echo effects induced by COM motion and position-dependent sign flips of the Rabi frequency can therefore not be described with the FGR model. The absence of a dip in the FGR-spectrum in Fig. 3 (b) signifies that in certain cases the rotary-echo effect is critical in understanding POL modulation spectra. We notice that width of the peak in the lowertemperature case (Fig. 3 (a)) is substantially larger than the echo-induced dip in the high-temperature spectrum ( Fig. 3 (b)). Thus, suitable combinations of lattice depth and drive-pulse duration yield high-temperature spectra that are more suitable for high-precision spectroscopy than their low-temperature counterparts.

G. POL amplitude-modulation spectra at long interaction times
A prominent prospect of using POL modulation spectroscopy is to perform high-precision measurements with cold, long-lived circular-state Rydberg atoms. To that end, we now consider several cases of moderately-deep magic and non-magic POLs, and drive pulses with a FWHM field-pulse length of 400 µs. While this is still about an order of magnitude short of typical circularstate lifetimes in a 300 K thermal-radiation environment, it suffices for a discussion of quantum structures of COM vibrational side-bands, and of vibrationallyresolved Doppler-free spectroscopy in lattices that are non-magic. Here we choose a 10 µK temperature, a value attainable in rapid-cycle optical-molasses cooling [46].
In Fig. 4 we consider an even-parity magic transition in a moderately-deep POL with V 1 = V 2 = 250 kHz. As before, the Rabi frequency is chosen such that atoms pinned at a spatial Rabi-frequency maximum experience a π-pulse. It is seen that all models produce a Fourier- limited central peak, with virtually no signal background at detunings δ 5ω k .
Since in Fig. 4 the drive pulse is much longer than in Figs. 2 and 3, the ratio between pulse length and vibrational period is much larger (τ FWHM f 1 = 26), and the Rabi frequency is smaller. Therefore, the case in Fig. 4 is deeper within the validity range of FGR perturbation theory. As a result, in Fig. 4 the agreement between the exact solution of the TDSE and the FGR approximation is quite good. Both quantum models show a substructure of the ∆ν = ±2 vibrational sidebands that arises from level shifts of the excited COM levels. The reduced Fourier width afforded by the 400-µs long pulse length allows for the observation of the ν−resolved substructure, with each sub-peak of the sidebands characterized by a single ν−value. The vibrational splitting of the ∆ν = ±2 sidebands arises from the anharmonicity of the wells and is, unsurprisingly, well-represented by both the TDSE and the FGR-band-structure models.
Due to the softening of the POL potential near its maxima, the ∆ν = ±2 transitions for larger ν-values occur at smaller absolute values of the detuning δ. Defining ν min as the minimum of the coupled states ν and ν , it is seen in Fig. 4 that for ν min 4 the change in transition frequency is proportional to ν min . This trend is easily confirmed by considering the effect of the lowest-order non-linear correction of the trapping potential, which is ∝ −ẑ 4 , on the transition frequencies within the ∆ν = ±2 sidebands. We also see that the lines for larger ν min are broadened according to the widths of the corresponding lattice bands, as visualized by the square boxes in Fig. 4 for the ν = 7 → ν = 9 and ν = 8 → ν = 10 transitions and their conjugates. Due to the ensemble average that is being taken, the extrema of the band-energy differences, E b (p 0 ) − E a (p 0 ), which occur at p = 0 and p = ± k, produce enhanced signals near the edges of the corresponding spectral features.
The behavior of the strengths of the vibrationallyresolved lines in the quantum results for the ∆ν = ±2 sidebands is given by thermal populations and the squares of the COM matrix elements of ν , p 0 | cos(2kẑ)|ν, p 0 . For conditions as in Fig. 4, the thermal populations in the lowest few tightly bound states drop off slowly. Over the spatial range of their vibrational COM wave-functions one may approximate cos(2kẑ) ≈ 1 −ẑ 2 /2. Considering the usual expansion of z in raising and lowering operators, it is seen that the line strength scales as (ν min + 1)(ν min + 2) ≈ ν 2 av , where ν min is the smaller of ν and ν and ν av is the average vibrational quantum number, (ν + ν )/2. Classically, according to Eq. (20) the line strength of the ∆ν = ±2 transitions is ∝ (J 2 (2kz 1 )) 2 , which for small oscillation amplitudes z 1 scales as z 4 1 , which in turn scales as the square of the quantum-mechanical vibrational quantum number ν. We see that quantum and classical analysis give the same line-strength scaling ∝ ν 2 av for the vibrationally-resolved sub-lines within the ∆ν = ±2 sidebands. The height-above-base of the lowest few of the sub-lines in Fig. 4 clearly follows this trend; the scaling also explains why the ν min = 0 sub-lines are so weak.
The semi-classical spectrum K b (δ) for the ∆ν = ±2 sidebands in Fig. 4 shows excellent qualitative agreement with the δ-averaged quantum spectra, averaged over the ν-quantization in the sidebands. Since in the classical treatment the vibrational COM energy is not quantized, this type of agreement between quantum and classical analyses accords well with our expectations. H. POL amplitude-modulation spectra in non-magic lattices As seen in the previous section, long drive-pulse durations combined with moderately deep POL potentials allow us to resolve vibrational quantization in the ∆ν = ±2 sidebands. In non-magic lattices, this feature also extends to the ∆ν = 0 central band. This makes non-magic lattices with a homogeneous lattice-laser intensity distribution suitable for high-precision spectroscopy. As an example, in Fig. 5 we consider a case with parameters identical to those of Fig. 4, except that V 2 = 0.8V 1 . Is is seen that the central peak, ∆ν = 0, splits up into lines at frequencies δ/(2π) ≈ (ν min + 1/2)(f 2 − f 1 ), with ν min = 0, 1, 2... and using the harmonic approximation. To find the unshifted atomic resonance in a high-precision measurement, one may plot δ/(2π) against x = ν min + 1/2. Extrapolation to x = 0 yields a yintercept that marks the un-shifted atomic resonance. Also, the measured line spacing f 2 − f 1 allows for a calibration of lattice depths V 1 and V 2 .
Figs. 4 and 5 show that magic and non-magic POL modulation spectra can, in principle, be modeled in great detail. Thereby, the only narrow spectral feature that is suitable for Doppler-free high-precision spectroscopy and that is, at the same time, insensitive to lattice-depth inhomogeneities is the central (∆ν = 0) peak in magic lattices. All other Fourier-limited features in Figs. 4 and 5 are Doppler-free, but exhibit shifts that scale with linear combinations of the COM vibration frequencies, which in turn scale with the square root of lattice power. Assuming that the lattice intensity in the atom-field interaction region can be stabilized to within ∼ 1%, one may expect to be able to resolve the sub-lines within the various vibrational sidebands. In that case, the only free parame- ters to fit an entire experimental spectrum with a set of calculated spectra will be the lattice depths V 1 and V 2 , and a detuning offset of δ. As the ratio V 1 /V 2 is known (it only depends on the Rydberg levels, the lattice-laser wavelength and the beam angles [21]), there are only two independent fit parameters. Hence, a two-parameter fit should simultaneously yield a detailed match of all sub-Doppler lines. The two-parameter fit yields a latticedepth calibration and a result for the δ-offset. The latter amounts to an accurate measurement of the lattice-free atomic transition frequency that is corrected for latticeinduced transition shifts.
I. POL amplitude-modulation spectra for odd-parity transitions In Fig. 6 we finally consider a case of a magic lattice with odd-parity drive. In this case, the vibrational selection rules are ∆ν = ±1, ±3, ..., with corresponding vibrational shifts. Most importantly, the central Doppler-free feature (which in magic lattices is insensitive to latticedepth variations) is forbidden, and the lowest-order allowed vibrational bands are ∆ν = ±1. These are separated in frequency by 2f 1 . Given a drive pulse with a narrow enough spectrum, the bands are split into vibrational sub-lines due to the anharmonicity of the lattice. For a sufficiently well-controlled lattice intensity, POL modulation spectroscopy of odd-parity transitions could also be employed in high-precision spectroscopy work.

VII. CONCLUSION
We have presented three models of modulation spectroscopy of Rydberg atoms in ponderomotive optical lattices (POL), an application that harnesses the oftenignored A 2 -term in the atom-field interaction. Foremost, it has been stressed that the modulation-induced drive term generates a spatially periodic Rabi frequency with a step-function phase dependence on position. This peculiar phase behavior of the drive enables a novel type of sub-Doppler spectroscopy suitable for applications in high-precision measurement in Rydberg-atom systems. The vibrationally-resolved lines rely on the quantum entanglement between electronic and COM motion, as well as the peculiar couplings afforded by POL modulation.
In our case studies we have used modulation functions with a Gaussian time dependence, because the spectra are free of Fourier sidebands (that would arise in squarepulse drives, for instance). It is, in principle, fairly straightforward to realize modulation drive pulses with a time-dependent Rabi frequency. This can be done by implementing time-dependent POL amplitude modulation with electro-optic fiber modulators.
We have developed a picture of quantum-classical correspondence between the quantum and semi-classical models. We have found a rotary-echo phenomenon that arises from the interplay between the center-of-mass motion of the atoms and the spatial dependence of the modulation-induced Rabi frequency. Our perturbative quantum model fails to reproduce the rotary-echo effect, as expected. The rotary-echo effect can improve spectral resolution in some cases.
A variety of Fourier-limited, Doppler-free vibrational transitions is seen in fairly deep lattices, even at temperatures in the range of 100 µK. Magic Rydberg-atom optical lattices lead to the most robust spectroscopic structure. However, non-magic lattices with differences on the order of 20% between the lower-and upper-state potential depths should yield equivalent spectroscopic accuracy and precision, when using suitable fitting methods. Ponderomotive transitions can also be driven by spatial "shaking" of the lattice; this is possible by phase modulation of a lattice beam (as opposed to amplitude modulation of the entire lattice, the case studied in the present paper). In any case, POL modulation spectroscopy is expected to yield line widths in the kHz-range, opening venues for Rydberg-atom-based high-precision spectroscopy [36] and quantum simulators [37].