Optimization of cw sodium laser guide star efficiency
R. Holzlöhner^{1}  S. M. Rochester^{2}  D. Bonaccini Calia^{1}  D. Budker^{2}  J. M. Higbie^{3}  W. Hackenberg^{1}
1  Laser Systems Department, European Southern Observatory
(ESO), KarlSchwarzschildStr. 2, 85748 Garching b. München, Germany,
http://www.eso.org/sci/facilities/develop/lgsf/
2 
University of California Berkeley, Department of Physics,
Berkeley, CA 947207300, USA
3 
Bucknell University, Department of Physics, 701 Moore Avenue,
Lewisburg, PA 17837, USA
Received 12 August 2009 / Accepted 1 October 2009
Abstract
Context. Sodium laser guide stars (LGS) are about to enter a
new range of laser powers. Previous theoretical and numerical methods
are inadequate for accurate computations of the return flux, hence for
the design of the nextgeneration LGS systems.
Aims. We numerically optimize the cw (continuous wave) laser format, in particular, the light polarization and spectrum.
Methods. Using Bloch equations, we simulate the mesospheric
sodium atoms, including Doppler broadening, saturation, collisional
relaxation, Larmor precession, and recoil, taking all 24 sodium
hyperfine states into account and 100300 velocity groups.
Results. LGS return flux is limited by ``three evils'': Larmor
precession due to the geomagnetic field, atomic recoil due to radiation
pressure, and transition saturation. We study their impact and show
that the return flux can be boosted by repumping (simultaneous
excitation of the sodium D_{2}a and D_{2}b lines with 1020% of the laser power in the latter).
Conclusions. We strongly recommend the use of circularly
polarized lasers and repumping. As a rule of thumb, the bandwidth of
laser radiation in MHz (at each line) should approximately equal the
launched laser power in Watts divided by six, assuming a
diffractionlimited spot size.
Key words: instrumentation: adaptive optics  methods: numerical  atmospheric effects  telescopes
1 Introduction
Laser guide stars (LGS) are becoming essential in providing artificial beacons for adaptive optics (AO) in large telescopes (Hubin 2009). The generation of 810 m class telescopes, such as the Very Large Telescope (VLT) on Cerro Paranal, Chile, or the Keck telescopes on Mauna Kea, has been retrofitted with LGS and is still operated in many observing programs with natural guide stars or no AO at all. The upcoming 30+ m telescopes, such as the 42m European Extremely Large Telescope (EELT) or the Thirty Meter Telescope (TMT), by contrast, are being designed from the start as adaptive telescopes and will require LGS in nearly all of their science operations. At the same time, sensing subaperture sizes of 2050 cm and nearkHz AO frame rates, designed to provide high Strehl ratios in the infrared or even in the visible, require unprecedented LGS brightness and hence laser power.
Sodium (Na) LGS at 589 nm are most commonly used because of the large fluorescence cross sectionabundance product of ca. 10^{11} cm m^{2}, a fluorescence wavelength in the visible, and the high altitude of the sodium layer around 80100 km compared e.g. to Rayleigh LGS, (Happer 1994, Table 1), which allows one to sense a large fraction of the turbulent atmosphere column above the telescope with a small number of guide stars. Unfortunately, powerful diffractionlimited laser beams at 589 nm are quite expensive to produce because of the lack of solidstate materials that amplify 589 nm or 1178 nm light. Recently, ESO has demonstrated a frequencydoubled narrowband Raman fiber laser emitting 25 W at 589 nm (Taylor 2009), and we expect this technology to significantly improve the experimental situation.
Careful dimensioning of the required laser power and optimization of the laser format (spectrum, polarization, spot size) are needed. Numerical simulations are necessary since the experimental situation has so far been unsatisfactory: the sodium layer and atmospheric parameters often fluctuate rapidly (Thomas 2008; Pfrommer 2009), few reliable powerful lasers at 589 nm have been available up to now, and LGS sky experiments lack commonly agreed on measurement standards. Sodium cells cannot easily simulate mesospheric conditions. This work provides optimization rules for the case of continuous wave (cw) lasers over a wide range of laser powers.
Sodium LGS take advantage of the dipole transition, known as the D_{2} line. The sodium ground state consists of two hyperfine multiplets with 8 magnetic substates combined, separated by 1.772 GHz, splitting the D_{2} line into the D_{2}a and D_{2}b transition groups, corresponding respectively to the F = 1 and F=2 Na ground states, where F is the total atomic angular momentum quantum number. The four multiplets ( ) are separated by only 16, 34, and 60 MHz (16 magnetic substates, see Ungar 1989). At mesospheric temperatures near 190 K, the D_{2}a and D_{2}b lines are Doppler broadened to about 1 GHz each, giving rise to the characteristic doublehump absorption profile (see for instance Bradley 1992, Fig. 11, or Fig. 1 of this work). When the D_{2}a transition is excited by circularly polarized light at high irradiance (i.e. optical power density, unit W/m^{2}), a large fraction of atoms is pumped into the substate, and the atoms cycle on the transition to the substate, with m the magnetic quantum number, so that sodium effectively becomes a twolevel system. This situation is twofold desirable because of the large transition cross section and the directional return light peak towards the laser emitter and the telescope. Numerous studies of optical pumping have been conducted, e.g. by Happer (1972) and McClelland (1985), and the process is well understood.
The effective sodium return flux depends on the environmental conditions such as collisions with constituent gases, temperature, and the geomagnetic field. At higher irradiance, the atoms are driven away from thermal equilibrium, and it is necessary to take effects into account such as saturation, optical pumping, and the recoil caused by radiation pressure. Understanding the complicated interplay of these effects and obtaining quantitative values of the fluorescence efficiency requires numerical simulations. A commonly used method is the solution of the density matrix evolution (Bloch equation) of a multilevel atom.
Milonni & Thode (1992) simplified the D_{2} scheme to a 2level Bloch model which they solve in time domain. Bradley (1992) simulated the full 24state density matrix, exciting the sodium by a train of short (nanosecond range) laser pulses like Milonni, using RungeKutta integration for one pulse period and exploiting the periodicity. Linear and circular light polarizations were treated. However, both works neglect the geomagnetic field. Morris (1994) studied frequencymodulated pulses over a wide range of pulse durations and linewidths up to 3 GHz, hence spanning the entire Doppler broadened D_{2} line, employing time domain integration, for both linearly and circularly polarized light. In the end of his paper, he briefly estimates the impact of the geomagnetic field, but he does not include Larmor precession terms into his Bloch equations.
To our knowledge, Milonni (1998,1999) has published the most advanced and detailed Blochequation simulation of sodium LGS to date, later generalized by Telle (2006,2008). In his 1998 article, Milonni treats the cases of laser pulses that are short, comparable, and long compared to the lifetime of ns, using numerical solution methods similar to Bradley's. His 1999 publication deals with cw excitation only and introduces spin relaxation and Larmor precession into the Bloch equations for the first time. He shows that the Larmor terms due to the geomagnetic field tend to redistribute the magnetic sublevel populations and hence impede optical pumping. All of the works cited in this and the previous paragraph solve the Bloch equations separately for a number of different sodium velocity groups (100400) and then perform a weighted average over the results, and all neglect recoil.
An alternative method of simulating atomic fluorescence is to use rate equations, either implemented as a set of differential equations, (Pique 2006), (Hillman 2008), or by employing Monte Carlo rate equation techniques^{}. In the latter, one follows a more intuitive approach and tracks the time evolution of a single atom, executing photon absorption and emission events according to steadystate cross sections and branching ratios using a random number generator, while accounting for elapsed physical time and recording the atomic state history. After a sufficient time span has been simulated, one divides the total number of spontaneously emitted photons by to obtain the fluorescent flux. The advantage of this method is that it is easy to understand, and complex physical events such as collisions with different kinds of particles with or without spin exchange can be modeled in a transparent way. Also, the transit of optically pumped atoms across the Doppler spectrum due to velocitychanging collisions and recoil can be easily modeled, hence different velocity groups can be properly coupled, including their atomic polarization exchange (using Bloch equations, such coupling requires the simultaneous solution of sets of 24^{2}dimensional equations per velocity group, which we describe in this work for the first time in application to sodium LGS, to our knowledge).
A serious disadvantage of rate equations is that they rely on steadystate atomic cross sections, limiting the scope to the simulation of events slow compared to the sodium transition lifetime, (Milonni 1992), and the combination of the excitationemission (Rabi) cycling with Larmor precession is hard to implement correctly. Furthermore, when compared to atomic Bloch equations, rate equations neglect the coherences (offdiagonal terms of the density matrix that describe transverse atomic polarization), which is problematic when modeling Larmor precession or linearly polarized light (the latter even for B=0). Finally, Monte Carlo rateequations need considerable CPU time in order to converge.
There are also approaches that mix aspects of Bloch and rateequation codes, such as BEACON, which has been adapted to model twostep sodium excitation for polychromatic LGS, (Bellanger 2004). BEACON neglects atomic collisions and recoil. Guillet de Chatellus (2008) reports that it requires on the order of 24 h per run on a 2.6 GHz processor, and that the agreement with a pure rateequation model can be good in certain cases.
Throughout this paper, we highlight what may be called the ``three evils'' of sodium LGS, ordered by decreasing importance
 1.
 Larmor precession,
 2.
 Recoil (radiation pressure),
 3.
 transition saturation (stimulated emission),
There is some ambiguity in the sodium LGS community about the term ``saturation'', for it is sometimes used to describe only transition saturation, and sometimes the depopulation of the F=2 upper ground state towards F=1, known as downpumping. Once an atom is in the F=1 state, it can only be excited if the D_{2}b line is pumped as well, which is known as repumping (either by allocating 1020% of the laser power to the D_{2}b line, or by widening a single laser line to 2 GHz). Downpumping in the absence of repumping becomes more severe with increasing laser irradiance and reduces the return flux long before the onset of stimulated emission, hence true transition saturation. Spinexchange collisions can bidirectionally exchange populations between the ground states and within them (transitions between the ( ) and ( ) ground substates are particularly strong). Repumping has already been experimentally demonstrated to be able to boost the LGS return flux by a factor of 1.6, (Telle 2008 finds a factor 2.4 through Bloch simulations), using two separate laser beams at different frequencies), and we show in this paper that more than a factor of 3 can be achieved.
To date, the most powerful sodium LGS system is installed at the Starfire Optical Range (SOR) at the Kirtland Air Force Base near Albuquerque, New Mexico, (Denman 2006a), fed by a singlefrequency (10kHz linewidth) cw laser that emits 50 W with circular or linear polarization. Due to its location, the median atmospheric seeing at SOR is significantly worse than at most sites of large astronomical telescopes, causing enlarged LGS spots in the mesosphere, and consequently limiting the irradiance. The upcoming generation of civilian 20Wclass LGS in astronomy, such as for the Adaptive Optics Facility of UT4 of the VLT^{}, is therefore expected to venture into unprecedented mesospheric irradiance regimes (much higher laser irradiances have of course been applied to gas cells, albeit at higher gas pressures). Nextgeneration lasers will emit circularly polarized cw radiation, use repumping, have a linewidth of a few MHz, and are expected to return on the order of photons/s/m^{2} on the ground. Quantifying and optimizing their return flux vs. the laser parameters is one of the principal purposes of this paper.
So far, we have only focused on photon return flux. What is really desired when designing AO systems, though, is a bright guide star, hence high luminosity concentrated in a small spot size. Compared to the uplink laser irradiance in the sodium layer, the above mentioned saturation effects spatially broaden the LGS return fluorescence distribution by emphasizing the low irradiance regions and dimming the peaks. We can show using physicaloptics simulations, Holzlöhner (2008), that this effect increases the instantaneous spot sizes on a wavefront sensor by about 0.1'', hence it is not negligible. At large angles between the laser beam and the geomagnetic field, downpumping can be mitigated by repumping.
In this work, we present a Blochequation method that models any alkali atom taking into account spontaneous and stimulated emission, Larmor precession due to the geomagnetic field, arbitrary elliptical light polarization, recoil, on the order of 100 coupled velocity groups with velocitychanging collisions and spin exchange, finite atomic dwell time in the beam (atom replacement), arbitrary laser bandwidth, and repumping. We neglect nonlinear Zeeman shifts and hyperfine coherences since we found them to have a small effect on the result. In contrast to the above cited Bloch simulation publications, we directly compute the steady state solution, which is more efficient than time domain solutions (a single run takes about 2 s on a modern PC). The program is written in Mathematica and based on the publicdomain Atomic Density Matrix package^{}.
In order to validate our results, we have also implemented a Monte Carlo rateequation simulation called Exciter in Matlab, whose results we compare with the Blochequation method.
Section 2 describes the Blochequation method, Sect. 3 gives details about the simulation parameters, Sect. 4 presents the results, suggesting optimal cw sodium LGS formats, and we conclude in Sect. 5.
2 Bloch equations
In order to calculate the observed fluorescence from mesospheric sodium atoms, the evolution of the atoms is modeled using the optical Bloch equations for the atomic density matrix. The density matrix describes the statistical state of an ensemble of atoms in the state space of the Na D_{2} transition. In order to account for atoms with different Doppler shifts, the density matrix is also considered to be a function of atomic velocity along the laser beam propagation direction. (An additional degree of freedom is included to account for laser line broadening as discussed below.) The calculation is semiclassical in the sense that while the atoms are treated quantum mechanically, the light fields are treated classically (the effect of spontaneous decay must be included phenomenologically since we do not apply field quantization). Because the density matrix describes all populations of, and coherences between, the 24 Zeeman sublevels making up the ground and excited states, the calculation describes, in principle, all saturation and mixing effects for essentially arbitrarily large optical and magnetic fields. (In practice, certain coherences in the system are negligible under our experimental conditions and can be neglected in order to increase the computational efficiency.)
In order to perform numerical calculations, the velocity dependence of the density matrix is discretized to describe an appropriate number of velocity groups, each with a fixed longitudinal velocity. Because coherences between atoms with different velocities can be neglected, the complete density matrix can be thought of as a collection of separate but coupled density matrices, each of dimension .
The evolution of the density matrix is given by a generalization of the Schrödinger equation:
where is the total Hamiltonian, with H_{0} the Hamiltonian for the unperturbed energy structure of the atom, the Hamiltonian for the interaction of the electric dipole of the atom with the electric field of the light, the Hamiltonian for the interaction of the magnetic moment of the atom with the local magnetic field , where Js is Planck's constant, and the square brackets denote the commutator. The term in Eq. (1) represents phenomenological terms added to account for relaxation processes not described by the Hamiltonian (Budker 2002). In our case these relaxation processes include spontaneous decay (omitted from the Hamiltonian due to the semiclassical approximation), collisional spin relaxation (``Sdamping'') proportional to (Happer 1987), and the exit of atoms from the light beam due to motion of the atoms and the beam. In addition, there are terms included in to describe changes in atomic velocity due to collisions and lightinduced recoil, as well as an effective relaxation rate that describes dithering of the laser phase in order to simulate a finite bandwidth. These terms are described in more detail below. Each relaxation process described by includes a corresponding repopulation process, so that the trace over the density matrix for all velocity groups is conserved, corresponding to conservation of the total number of atoms. The repopulation process describing the entrance of atoms into the beam is independent of and so is written as a separate term .
Velocitychanging collisions (vcc) are treated as hard collisions in which the velocity of the colliding atom is rethermalized in a Maxwellian distribution (no speed memory). The internal state of the atom is assumed to be unchanged.
Lightinduced recoil is described phenomenologically by causing a fraction of the excitedstate atoms in each velocity group to be transferred upon decay into the next higher velocity group. Here is the recoil velocity and is the width of the particular velocity group. This model relies on the fact that cm/s (equivalent to a Doppler shift of 50.004 kHz) is much smaller than the typical value of .
In order to simulate a finite bandwidth laser, a form of phase dithering is used (frequency or amplitude dithering can also be employed). To avoid resorting to a timedomain calculation, the dithering is implemented in the spatial domain: density matrices are written for two ``regions'' with light fields that are out of phase with each other, and relaxation terms are included that transfer the atoms between the regions (this doubles the size of the system of equations). The model is that of a laser beam with very fine ``speckles'' of different phases. The result is an effective laser spectrum of Lorentzian shape with a width proportional to the transfer rate between the regions. This method has been verified by comparison to a timedomain model (implemented for a nuclearspinless system) in which the light frequency randomly changes with a Lorentzian distribution. Identical results from the two methods are obtained for the case in which the rate that the light frequency changes is faster than the natural decay rate.
Equation (1) supplies a linear system of differential equations for the density matrix elements, known as the optical Bloch equations. Thinking of as a column vector of density matrix elements, the Bloch equations can be written as , where A and b are a matrix and vector, respectively, that are independent of . The vector b corresponds to and A to the rest of the righthand side of Eq. (1).
The laser light field has a frequency component tuned near the D_{2} transition group (D_{2}a), and may have an additional ``repump'' component tuned near the transition group (D_{2}b). Thus the matrix A has components that oscillate at each of these frequencies. Under the rotating wave approximation (Corney 1977), the overall optical frequency is removed from A. However, the beat frequency between the two lightfield components remains. This beat frequency can also be removed from the Bloch equations in our case: each frequency component interacts strongly with one transition group and very weakly with the other, so the weak coupling can be neglected for each transition. If, in addition, the small magneticfieldinduced mixing between the two hyperfine ground states is neglected, the beat frequency can be entirely removed from the evolution equations. This makes A timeindependent for cw light. To find the steadystate density matrix, we can set and solve the linear system . The vectors and b have 322 elements per velocity group (576 if hyperfine states are not neglected), so that the sparse linear equation system has dimension 32 50065 000 in practice.
To solve the Bloch equations for a particular set of experimental parameters, we first choose an appropriate set of velocity groups. Since the signal is strongly peaked for atoms whose Dopplershifted resonance frequency is near the light frequency, we can obtain more accurate results for a given number of velocity groups if narrower bins are used for resonant atoms, and wider for offresonant. We have two methods for doing this.
The first method is to choose two fixed bin sizes, one narrow and one wide, and the number of narrow bins to cluster near each resonance. The wide bins are then used to take up the rest of the Doppler distribution. This method is useful when we don't know beforehand what the spectrum of the signal in velocity space is.
If we have an estimate of the spectrum (obtained using the first method), we can refine it using the second method, which takes advantage of this knowledge. We create a weighting function consisting of three terms: a constant term, which tends to make equalsized bins, a term proportional to the spectrum, which makes more bins where the signal is large, and a term proportional to the magnitude of the second derivative of the signal, which makes more bins where the signal changes rapidly as a function of velocity. The bin sizes are then found by dividing the integral of the weighting function evenly into the chosen number of bins.
The linear system is solved using the implementation of the iterative BiCGSTAB method (stabilized biconjugate gradient, van der Vorst 1992) builtin to Mathematica. This is a Krylov subspace method in which an initial guess is improved by minimizing the residual over a subspace with dimension much smaller than that of the full system. The rate of convergence of the method is increased by premultiplication with a blockdiagonal preconditioner (approximate inverse of A), obtained by setting all terms that connect density matrix elements from different velocity groups to zero, and then inverting the block for each velocity group.
The fluorescent photon flux per solid angle emitted in a given direction can be found from the steadystate solution for as the expectation value of a fluorescence operator (Corney 1977).
3 Simulation parameters
3.1 Determination of parameters
Many physical constants and the atomic level diagram of Na have been summarized by Steck ^{} and will not be repeated here. Table 1 lists further simulation parameters, and we walk through it to explain some quantities. We work in MKS units, except that we show magnetic field strength in Gauss ( T) and atomic cross sections in squared centimeters.
The launched laser power P equals the laser device output beam power, diminished by optical losses in the beam train and launch telescope (LT). The value of is the ``50% power in the bucket'' irradiance I_{P/2}, as explained in the next subsection. We mostly deal with circularly polarized light (the ellipticity angle denotes respectively LH/RH circular and linear polarization; is the polar angle in the Poincaré sphere, and is the ellipticity, hence the majorminor axes ratio of the polarization ellipse). The repumping fraction q is the fraction of the total laser power allocated to the repump beam: the D_{2}a beam power is (1q)P and the power in the D_{2}b beam (tuned above the D_{2}a frequency) is qP.
Table 1: Simulation parameters and their standard nominal values.
The geomagneticfield strength B has a strong impact on the return flux. Its value varies considerably over the world and can be computed for different mesospheric altitudes using the International Geomagnetic Reference Field model^{}. Cerro Paranal in northern Chile, the location of the VLT and the reference site for this work (24.6S, 70.4W), has G at 92 km altitude, Mauna Kea (Hawaii) has G, and the Starfire Optical Range (Albuquerque, New Mexico) has G, about twice the field strength at Paranal.
The atmospheric transmission at 589 nm was measured at Paranal to be 0.89 in photometric nights, (Patat 2004). We use a value that is 5% lower to account for higher aerosol levels. The mesospheric temperature, as well as partial gas densities, can be derived using the MSISE90 model^{}.
Our values for the sodium layer centroid altitude and column abundance are obtained from studies taking place for over 30 years in São Paulo, (Simonich 1979; Moussaoui 2010), a site whose latitude differs only by one degree from the ESO Paranal Observatory. We believe that the sodium layer parameter statistics on the seasonal and daily variations are valid for Paranal.
Atomic collisions have a significant effect on the sodium
states and hence on the LGS return flux. Since mesospheric sodium is
rarefied (the total mass of global mesospheric sodium is about
600 kg), NaNa collisions are negligible compared to NaN_{2} and NaO_{2}
collisions. Most of these collisions are binary (collision of two
molecules). One important effect of collisions is to change the
velocity of the atoms, causing diffusion of optically pumped atoms in
velocity space, sometimes called vdamping (Happer 1987). Since the masses of N_{2} and O_{2}
molecules are comparable to that of Na atoms, one can assume that every
collision completely randomizes their velocity. The collision rate of a
gas of particle mass M_{1} with another gas type of particle mass M_{2} and number density n_{2} is given by (Wright 2004)
where J/K is Boltzmann's constant and is the collisional cross section with the effective particle radii r_{1}, r_{2}. Measuring these radii for velocitychanging collisions is difficult; here we just assume the VanderWaals radii of pm, pm, and pm^{}, (Fishbane 2005). The effect of other gas species is negligible. With these numbers, we have, for example, and s) at the sodium centroid ( ).
The sodium layer has a median FWHM thickness of 11.1 km and its median centroid lies at km (Moussaoui 2010). The gas pressure decreases exponentially with altitude and the collision rate varies across the layer by about one order of magnitude, as shown by the MSISE90 atmospheric model. We compute the mean collision rate of Na with N_{2} and O_{2}, based on a table of , , and as functions of altitude, weighted by the assumed Gaussian sodium density distribution. The result, which we will use throughout this work, is s), which is three times higher than Milonni's (1999) assumption of s). Since the sodium abundance and layer thickness and altitude are highly variable and its profile often deviates significantly from Gaussian, one cannot expect a high accuracy in this parameter.
The other important relaxation mechanism beside vdamping is spinexchange relaxation or Sdamping (Happer 1987, Sect. 13), in particular for NaO_{2} collisions. Spin relaxation time constants have been measured between rubidium and metastable triplet helium (He, Dmitriev 2008), as well as between rubidium and H_{2}, O_{2}, and N_{2}, (Nagengast 1998) and sodium and various gases (Ramsey 1964; Kartoshkin 1998). A major difficulty with such measurements is that the overwhelming contribution to Sdamping of Na in the mesosphere is due to collisional spin exchange with O_{2}, however, in gas cells O_{2} oxidizes Na quickly and hence this particular cross section is hard to determine experimentally. Theoretical calculations of the cross section involve BornOppenheimer molecular potential curves of doublet/quartet surfaces for NaO_{2}, analogous to the singlet/triplet curves for NaNa and have not yet been carried out to our knowledge.
We estimate cm^{2} at 185 K, based on spinexchange cross section measurements of NaHe and O_{2}He. However, only 1/2 of this cross section is effective in our case (Dmitriev 2008, Eq. (3)), although the details of the crosssection calculation are under investigation at the moment. Setting in Eq. (2), we obtain s) at 92 km altitude, which is close to the value s) that Milonni (1999) finds through fitting to experiment. Note, however, that Milonni's initial guess was cm^{2}, and he does not apply the scaling factor of 1/2. Performing the same sodium density weighted averaging over altitude as above, we obtain s), which will be used throughout this work. We will discuss the sensitivity of the Na return flux to variations in and in the following section.
The rms lateral velocity describes the sodium atom exchange into and out of the beam, where cm is the median FWHM mesospheric speckle diameter (see the following section), and is the atom exchange rate. We assume that consists of four components that we sum in quadrature since in general they have uncorrelated directions:
 (a)
 gas diffusion orthogonal to the beam;
 (b)
 mesospheric wind orthogonal to the beam;
 (c)
 beam wander caused by atmospheric turbulence; and
 (d)
 LGS beam slewing due to star tracking.
 (a)
 the diffusion coefficient of Na in air (mostly N_{2}) can be calculated using the ChapmanEnskog formula (Cussler 1997), yielding cm^{2}/s, close to Milonni's assumption of cm^{2}/s. The effective lateral diffusion velocity across the laser speckle is hence m/s, where we apply the factor since we consider only the diffusion orthogonal to the beam. Note that scales like . At an altitude of about 104 km the mean free path of sodium exceeds d and the atomic motion within the beam is no longer diffusive (the mean lateral ballistic velocity is then 370 m/s);
 (b)
 we use the horizontal wind model HWM07 which computes zonal and meridional winds at selected mesospheric altitudes and times (Drob 2008). We find a median horizontal wind speed magnitude at Paranal of 20 m/s. Nocturnal fluctuations by a factor of 2 are common;
 (c)
 physicaloptics simulations show an rms beam wander in the mesosphere of about 20 cm ( ), (Holzlöhner 2008), a value supported by analytical approximations (Andrews 2005, Ch. 12, Eq. (51)). The beam wanders on the time scale of (atmospheric wind speed) 23 cm/10 m/s = 23 ms, leading to a beam wander velocity of 8.7 m/s;
 (d)
 the tracking speed of the laser beam in the mesosphere at is about m/s, where is the lineofsight distance to the guide star centroid with the geometrical length extension factor (airmass) , equaling the secant of the zenith angle .
3.2 Mesospheric spot size
Given a certain launched laser power, the spot size and shape
determine the mesospheric irradiance and thus have a strong impact on
the return flux, particularly in the presence of strong optical pumping
and for small magnetic field polar angles .
The instantaneous mesospheric laser spot size has been simulated using physical optics (Holzlöhner 2008).
Atmospheric turbulence produces a fast changing (timescale a few
milliseconds) speckle pattern on the sky due to diffraction, where the
speckles have a FWHM diameter of about
and w is the 1/e^{2} irradiance radius of the laser beam at the projector. The number of speckles and their beam wander are governed by the seeing, and the pattern changes on the time scale of a few milliseconds. Note that the medium/longterm spot size as observed on the ground is significantly larger than d, and it does, in contrast to d, depend directly on the Fried length r_{0} (Fried 1966), which is a size scale over which atmospheric phase variations remain below . For a typical r_{0}=525 cm, assuming an outer turbulence scale of L_{0}=25 m, there is a single strong central speckle and a couple of much weaker satellite speckles.
The LGS return flux depends in general nonlinearly on the mesospheric laser irradiance I,
and the return light has to pass the turbulent atmosphere again on the
downlink to reach the telescope; thus the irradiance distribution I(x,y),
as it would appear on a screen at the mesosphere, cannot be directly
observed. To account for the effect of turbulence on the return light
and characterize speckle pattern statistically, we collect statistics
over many simulated realizations of I(x,y) for different turbulence phase screens. We rewrite the usual optical power integral
(4) 
as an integral over irradiance up to the maximum value , where A is area orthogonal to the beam and is the atmospheric transmission along the line of sight (zenith angle ). The quantity dA/dI can be thought of as a histogram of mesospheric laser irradiance. We define I_{P/2} by
(5) 
representing an irradiance of ``50% power in the bucket'', which will be convenient to compute laser efficiencies later (note that usually ). In order to model a realistic LGS beam scenario, we assume a D=40 cm LT, a launched Gaussian beam of 1/e^{2} (in irradiance) radius w=0.72D/2=14.4 cm, two different seeing values of 0.6'' and 1.0'' (site monitor seeing at 500 nm and zenith), and total wavefront errors (WFE) of the launched beam of 70 nm rms () and 100 nm rms (). We simulate the actual mesospheric irradiance distribution and I_{P/2} using physical optics, and record the equivalent beam diameter FWHM that a single Gaussian spot of the same mesospheric power and same I_{P/2} would have, as well as the speckle size d, as in Holzlöhner (2008). Table 2 summarizes the results. The mean of I_{P/2} on lines 2 and 4 equals 46 W/m^{2}, which we will use as the reference irradiance in this work. The median values of d in the last column are somewhat above the analytical prediction of Eq. (3) of 21.2 cm, mostly due to the finite WFE.
Table 2: Mesospheric spot parameters for a 40 cm launch telescope under standard conditions (P=20 W, )
3.3 Flux quantities
In the following, we define five flux quantities used in this article: our simulations compute the return flux per solid angle as the number of photons per atom and per unit time spontaneously emitted in the direction of the launch telescope, as observed in the mesosphere, with units of photons/s/sr/atom;
 in order to compare different laser formats at similar irradiances, one can divide by the laser irradiance (in W/m^{2}) in the mesosphere to derive the specific atomic return flux in units of photons/s/sr/atom/(W/m^{2}), where area is measured in the mesosphere orthogonal to the beam. Most of our results are expressed in this quantity;
 to compute the expected return flux on a receiver at a distance L from the fluorescing atoms, we begin by integrating over the receiver area A orthogonal to the beam, subtending the solid angle A/L^{2} to obtain the quantity ;
 the photon flux on the detector
(unit photons/s/m^{2}) at the sodium column abundance
(unit atoms/m^{2}) equals
(6)
The airmass X appears in the numerator because the sodium column along the laser beam scales like the airmass (note that the distance L grows by the same factor). We do not account for laser power depletion with increasing propagation distance in the sodium layer because only about m m % of the laser photons interact with any sodium atoms under the standard conditions of Table 1, the other 96% travel unused into space;  in order to compare LGS systems at similar powers but
different laser formats, beam spot sizes in the mesosphere, atmospheric
transmission, and zenith angle, we can divide
by the laser power in the mesosphere
to arrive at the figure of merit quantity
in units of photons/s/W/(atoms/m^{2})
where appears quadratically in the denominator of the last term, accounting for both up and downlink. The quantity is hence the photon flux on the ground, corrected for its dependence on sodium centroid height and abundance, airmass, atmospheric transmission, and, to first order, launched laser power (since depends strongly on , we will only compute it for in this work). We will compare LGS laser beam efficiencies based on . Note that in our definition is not a slope efficiency (which would be proportional to ); see also d'Orgeville (2000).
4 Results
We compute using the method and parameters described above. Using the standard conditions of Table 1, we obtain ph/s/sr/atom/(W/m^{2}). In the following, we first provide some numerical validations including scans of , followed by optimization studies of when varying light polarization ellipticity angle , repumping frequency offset , repumping power fraction q, and laser linewidth . Unless otherwise noted, the standard conditions of Table 1 apply.
4.1 Some numerical validations
Figure 1: Absorption cross section vs. detuning from the Na D_{2}a line center. Solid blue: simulation at I=0.001 W/m^{2} (q=0), dashed red: analytical cross section, green: simulation of effective cross section at standard conditions (I=46 W/m^{2}, q=0.12). Inset: reference coordinate system. 

Open with DEXTER 
Figure 1 shows a comparison between the analytical Na absorption cross sections given by Eqs. (10), (11) in (Milonni 1998) at T=185 K (dashed red line, the sum of two Gaussians) and the simulated effective cross section
(8) 
where is the center frequency of the D_{2}a line, with c the vacuum speed of light, and W is the actual rate of spontaneous emissions at irradiance I (solid blue curve; the inverse of W is known as the cycle time). The abscissa shows detuning of the (main) laser line from the D_{2}a line center. The agreement is excellent. We add the caveat that this agreement is necessary, but not sufficient to prove the validity of our algorithm since effects like Larmor precession, recoil, and stimulated emission have no effect at low irradiance as the Na atom is in thermal equilibrium.
The steadystate spontaneous emission rate W cannot exceed ns), and hence tends to zero in the limit of infinite I. At I=46 W/m^{2} and q=0.12 (green curve), the reduction in is still modest, however. We render the tail of the green curve towards the D_{2}b line center at 1.772 GHz dashed, where our computational simplification of letting the ``D_{2}a laser line'' only excite the D_{2}a Na transitions, and analogously allowing the ``repumping'' line to only excite the D_{2}b transitions, whenever repumping is used (q > 0), breaks down. Once the detuning approaches the D_{2}b line center, this assumption obviously becomes invalid. By contrast, the blue curve was computed for q=0 and without using this simplification, and it is valid for any frequency offset.
The inset in Fig. 1 sketches the definition of the spherical angles and of the vector in a coordinate system where the laser beam is projected along the zaxis. The major axis of the polarization ellipse for noncircular polarization is parallel to x.
Figure 2 shows the simulated atomic velocity distribution under standard conditions, except that Larmor precession and repumping are absent ( ). The abscissa shows relative atomic velocity away from the receiver in frequency units (proportionality constant ). Recoil is usually associated with a frequency redshift, whereas we display the shift towards the blue. We choose this convention since a recoiled atom becomes resonant with a blueshifted laser line.
The blue curve represents the occupation histogram from the Monte Carlo rateequation simulation Exciter after sampling 10^{6} different atomic velocities within the Doppler distribution ( Monte Carlo steps simulating 26.6 s of physical time). The green curve is a Gaussian with the theoretical Doppler FWHM width of 1.033 GHz at K, representing the thermal equilibrium, and the magenta curve shows the occupation probability obtained from the Blochequation simulation. Finally, the cyan curve depicts for comparison the simulated return flux spectrum in the atomic frame, which in the present case of singlefrequency excitation is close to the sodium natural line shape of a Lorentzian with a FWHM of MHz (plotted at arbitrary vertical scale).
The effect of spectral hole burning is quite striking; in fact the occupation at the D_{2}a line center is depleted to 64% below the green Doppler curve. The atomic population, as seen from the telescope, is blueshifted within about one velocity group. We observe that spectral hole burning due to recoil is in general exacerbated if repumping is applied, presumably due to the larger number of excitations per time. On the other hand, hole burning is mitigated when the laser bandwidth is extended at constant laser power because of the reduced spectral irradiance.
In Exciter, we increment the Doppler shift of the simulated atom by 50 kHz each time a spontaneous emission occurs, which is correct on average. The agreement between Exciter and the Bloch equations is excellent, giving us confidence that the Gaussian velocity distribution and recoil are properly modeled.
At I=46 W/m^{2}, Exciter simulates a spontaneous emission every 9.5 s on average. Of all emissions, 10.2% are stimulated and 89.8% are spontaneous, and 70.6% of all atomic transitions occur along , indicating efficient optical pumping. We conclude that besides suffering from recoil, we are also starting to saturate this transition. At I=100 W/m^{2}, 17.9%, and at I=1000 W/m^{2}, 55.2% of all emissions are stimulated, respectively. On average, one spontaneous emission occurs every 5.1 s and 1.1 s for these respective irradiances.
Figure 2: Atomic velocity distribution under standard conditions (but ), blue: Exciter, green: Doppler profile, magenta: Bloch equations, cyan: return light spectrum in the atomic frame. Inset: zoom on line center. 

Open with DEXTER 
Figure 3: Specific return flux . Lines: Bloch equations, symbols: Monte Carlo rateequation simulation Exciter. a) Green, crosses: ; red, diamonds: ; blue, squares: (standard conditions); magenta, circles: (all with circular polarization and ). b) Same as a), but for linear polarization. c) Red and magenta lines, diamonds and circles: as in a); black, crosses and cyan, squares: same, respectively, without recoil. 

Open with DEXTER 
Figure 3 shows three semilogarithmic plots of for the standard conditions of Table 1 unless noted otherwise, but neglecting exchange with thermal Na atoms outside the beam ( , which is not implemented in Exciter, but causes only a small difference). We hence simulate a singlefrequency laser ( ) tuned to the peak of the D_{2}a line. Solid curves show obtained from the Bloch method, and the symbols indicate the results of Exciter. The green and red curves and symbols represent the case of laser beam parallel to the geomagnetic field () with and without repumping, hence in the absence of Larmor precession, respectively. The blue and magenta curves and symbols are the same as green and red, respectively, but with the laser orthogonal to the field ( ). Figure 3a shows the case of circular, and Fig. 3b that of linear polarization, the blue and magenta curves for (note the difference in vertical scale).
We can make a number of interesting observations. First, the
impact of the magnetic field is profound and reduces the return flux
strongly, at some irradiances by several times. At very low irradiance (I=10^{2} W/m^{2}),
the atom is in thermal equilibrium and all magnetic sublevels are
nearly equally populated. Optical pumping sets in with increasing I if using circular polarization, but Larmor precession is powerful enough to completely suppress it at
,
as evident from the monotonically falling magenta curve. Conversely, in the absence of Larmor precession
strongly peaks near I=2 W/m^{2}. The Larmor precession period equals
(9) 
where is the Bohr magneton and g_{F} is the Landé factor that depends on F ( for the F=1, 2 sodium ground states, respectively). On Paranal, s (159 kHz) and at SOR we find s (333 kHz, B=0.48 G).
The green and blue curves portray the ``healing'' effect of repumping (q=0.12, i.e., 12% of the laser power shifted to D_{2}b). In the absence of Larmor precession, the peak is more pronounced and shifted to higher I (green). The highest impact occurs when the influence of the magnetic field is highest ( ), where remains shallow until 80 W/m^{2} and then decays gently. This finding refutes the often heard notion in the LGS community that ``one should always stay away from the circular polarization saturation irradiation ( W/m^{2} to avoid saturation''.
At very low I, repumping is ineffective and even slightly decreases the return flux due to the smaller cross section of the D_{2}b transition. At high I, repumping can readily compensate the degrading effect of Larmor precession; in its absence, decays to zero. For linear polarization, this decay even occurs for (red curve) due to depopulation of the upper ground state. Conversely, circular polarization ``rescues'' many atoms into the transition, in which they are safe from downpumping.
Note that besides increasing the return flux by several times, a flatter function also leads to smaller observed spot sizes since spatial powerbroadening is reduced, a very welcome bonus.
The agreement between the Bloch calculation and Exciter is very good in almost all cases, given the completely different nature of the two methods, with the Bloch code usually yielding the higher values of . Exciter finds the magnetic resonance bump at smaller I due to difficulties with the proper simulation of two concurrent effects on similar, but not equal, time scales in our Monte Carlo rateequation. Also the agreement between the green curve and the green crosses is somewhat poor around I=0.220 W/m^{2}, presumably due to the simplified way in which Exciter models Sdamping and/or the absence of coherences.
Figure 3c shows the impact of recoil. The red and magenta curves and symbols are as in subplot a), while black and cyan crosses and squares denote the same, respectively, but neglecting recoil. The magnetic resonance is better visible in the cyan than in the magenta curve. Recoil leads to a significant reduction in above 2 W/m^{2}.
Figure 4a shows as a function of the Larmor frequency (blue curve) at I=1 W/m^{2}. The red dashed curve, overlapping the solid blue curve, is a fit function composed of the sum of two Lorentzians of different widths, centered at the origin, plus a constant term. The (half)width of the narrower Lorentzian ( ) is determined by the Sdamping resonance, and that of the broader Lorentzian ( ) by the velocitychanging collision rate, and the widths change proportionally when varying or . However, both resonances are somewhat powerbroadened. The geomagnetic field is strong enough at I=1 W/m^{2} to place us on the lowest terrace of the blue curve, as indicated by the vertical lines. The blue diamonds show the result of Exciter for comparison, also exhibiting the terraces.
Figure 4b shows the corresponding relative populations of the five S, F=2 upper ground state sublevels ( ), where the solid lines indicate Bloch equations and dotted lines Exciter. For Hz, the sublevel populations collapse to the same value due to Larmorinduced sublevel mixing. Conversely, for Hz, the populations diverge and the laser pumps the m=2 sublevel (blue) most strongly. When increasing the irradiance to I=46 W/m^{2} (not shown in the plot), the curves in b) would look similar, but the divergence point shifts upward to 1 kHz. Furthermore, when plotting the sublevel populations at I=46 W/m^{2} as a function of , one observes that the magnetic field at Paranal is just strong enough to collapse the populations for . In other words, by scanning across the sky with a narrowband 20Wclass laser, we can observe sodium excitation all the way from full optical pumping to no pumping at all, leading to a severe return flux penalty, as shown in the following figures.
Figure 4: Magnetic field impact. a) Blue curve: as a function of the Larmor frequency at I=1 W/m^{2}. Red dashed curve (almost obscuring the blue curve): fit (sum of two Lorentzians). Blue diamonds: Exciter. Solid green (dashed magenta) lines: at Paranal (SOR). b) F=2 fractional sublevel populations, solid (dashed) curves: Bloch (Exciter) simulation. 

Open with DEXTER 
Table 3: Sensitivity of to a 1% perturbation in some simulation parameters under standard conditions.
Table 3 lists the sensitivity of w.r.t. a 1% perturbation in some simulation parameter p, more precisely , where all other parameters are those of Table 1. A value of b% in the last table column thus indicates that in some range around the chosen value of p. We have selected only those parameters that we will not study in greater detail in the following subsections (except I). In addition, we have excluded those parameters around which is stationary ( ), and also those parameters whose influence on the observed return flux is obviously linear, such as .
From Table 3, we notice that with the parameters of Table 1 the sensitivity to changes in I is small. A comparison with Fig. 3b (dashed blue curve) shows that is very shallow near I=46 W/m^{2} (the slope is actually weakly positive due to the magnetic resonance bump), meaning that , i.e. the absolute photon return grows linearly with the irradiance. The dependence of on B, however, is strong: Switching from B=0.23 (Paranal) to B=0.48 (SOR, factor 2.11) decreases by a factor of 0.76! The influence of temperature is significant as well since the width of the velocity distribution scales like , and conversely its normalization (the number of Na atoms per velocity group) scales like . However, increasing from 185 K to 200 K, which may be a typical seasonal variation, reduces by only 3.7%, and hence common temperature variations do not directly influence much (we note that temperature influences the sodium abundance). Both and have a small influence on , mainly due to the presence of repumping, as Milonni (1999) has also noted (for q=0, the sensitivities are 0.049% and 0.56%, respectively).
We can draw two lessons here: There is little point in spending much effort trying to improve the spinexchange cross section estimate, and, secondly, repumping makes more robust towards some parameter variations, leading to higher stability in numerical simulations and possibly also in experiment. The variation of with is much larger than with , and at first sight surprisingly, the derivative is positive. We explain the positivity by the fact that collisions mitigate spectral hole burning due to recoil (neglecting recoil, the derivative is 0.089%). Finally, we have constants and (more precisely kHz, with the atomic mass of Na kg), in order to demonstrate the importance of proper recoil modeling at I > 1 W/m^{2}.
In the following subsection, we will show that the laser parameters of Table 1 are close to optimal, given the other conditions of the table. We focus on the case where the laser beam is directed orthogonal to the geomagnetic field (the worst case, for which LGS lasers must be designed to achieve a given return flux requirement).
4.2 Optimization of
4.2.1 Polarization ellipticity angle
Figure 5: as a function of laser polarization ellipticity angle , where (0) marks circular (linear) states. The curves indicate different combinations of the magnetic field polar angle and azimuth angle . Note that is highest for circular polarization in all cases. 

Open with DEXTER 
Figure 6: Contour plot of in ph/s/sr/(W/m^{2}) for magnetic azimuth at I=46 W/m^{2}. 

Open with DEXTER 
Figure 5 plots as a function of the laser polarization ellipticity angle . To our knowledge, this is the first publication to discuss for LGS. The axis limits of the quantity denote circular polarization (the poles of the Poincaré sphere; in our convention, pumps towards increasing m), and indicates linear polarization. The curves pertain to different combinations of and , the polar angle and azimuth of the vector in a system where the laser points along z (note that for , the azimuth becomes meaningless at any light polarization. Conversely, when using any laser that is not purely circular polarized, does influence the return flux, which is often forgotten. When using purely linear polarized light, , is equivalent to setting B=0 within our convention, while induces the strongest effect from Larmor precession.)
The cosinelike shape of presents an initially gradual decrease of the return flux from , meaning that high polarization purity is not required when pumping the sodium with circularly polarized light. If one is willing to accept a decrease in of 5% along the solid red curve, it is enough to maintain , i.e., down to an ellipticity angle of . The polarization extinction ratio (PER) in dB is given by and assumes the values of 0 and for circular and linear polarization, respectively. The condition of implies dB. This insensitivity is good news for the design of launch telescope optics.
Figure 6 shows a contour plot of in ph/s/sr/(W/m^{2}) for (worst case) under standard conditions. Current LGS systems, even SOR's FASOR laser, operate at lower ``P/2'' irradiances than I_{P/2}=46 W/m^{2} and/or at linear polarization, for which varies slowly. At higher I and circular polarization, is a steeper function (along the right vertical plot edge). For , also varies slowly near , indicating again that the circular polarization purity does not have to be very high.
4.2.2 Repumping frequency offset
Figure 7: as a function of repumping laser line frequency offset at I = 1 W/m^{2} (dashed green), I = 10 W/m^{2}(dashdotted blue), and I = 46 W/m^{2} (solid red). 

Open with DEXTER 
Figure 8: Contour plot of in ph/s/sr/(W/m^{2}). 

Open with DEXTER 
Figure 7 shows as a function of , the repumping frequency offset between D_{2}a and D_{2}b, relative to 1.7178 GHz, for two singlefrequency laser lines and q=0.12 at I=1 W/m^{2} (dashed green), I = 10 W/m^{2} (dashdotted blue), and I = 46 W/m^{2} (solid red). The repumping from requires excitation to the states (P,F'=1) or (P,F'=2) since no transitions are allowed.
Figure 8 is a contour plot of , zooming into the region of MHz. In order to harness the full improvement of repumping, one needs to adjust with a precision of a few MHz, particularly when at intermediate irradiance ( W/m^{2}). At larger irradiance, powerbroadening washes out the peak (a laser beam has irradiances in the mesosphere, located along a horizontal line in the diagram). The tolerance in is so narrow because repumping is most efficient if a given Na atom that has been downpumped to the F=1 lower ground state can be reexcited before its velocity changes due to a collision, for which the frequency offset of the repumping laser line should lie between the transition frequencies minus the transition frequency that the main laser line is tuned to, i.e., close to 1717.8 MHz.
4.2.3 Repumping power fraction
Figure 9: as a function of repumping power fraction q for different field angles . Dashdotted red curve: standard conditions, solid blue: (circular polarization, no Larmor precession), dashed black: standard conditions, but B=0.48 G (SOR), solid green: linear polarization ( ). 

Open with DEXTER 
Figure 10: a), b) Contour plots of in ph/s/sr/(W/m^{2}) at Paranal for linear and circular polarization, respectively; d) the ratio at Paranal; e) same as b), but for SOR geomagnetic field; c), f) repumping gain at Paranal and SOR, respectively, circular polarization. Blue crosses: I=46 W/m^{2}, q=0.12. 

Open with DEXTER 
Figure 9 portrays the calculated for different magnetic field polar angles and field strengths B. The dashdotted red curve is for standard conditions of Table 1, the solid blue curve for (i.e., circular polarization and no Larmor precession), the dashed black curve is for standard conditions, but with the magnetic field at SOR, and the solid green curve is for linear polarization (worst case ).
All curves peak around q=825%, and the maximum occurs at larger q the stronger the magnetic field effect is. Note that peaks at q > 0 even in the absence of Larmor precession (blue curve) since there is always a chance that atoms decay to the lower ground state. Moreover, is very steep at q=0, rewarding even the weakest repumping, which is the reason why broadline lasers (2 GHz) often do reasonably well, despite their poor spectral overlap with (see Fig. 1). Our simulations, however, indicate that a narrowband laser with 12% repumping beats the return flux of a singleline 2 GHzlinewidth laser by a factor of 3.7 under standard conditions (see next subsection).
For circular polarization we find at I=46 W/m^{2}, underlining again that repumping is extremely beneficial, in particular in future LGS systems involving high irradiance.
Figure 10 shows six semilogarithmic contour plots of (note the different color scales of the plots). The blue crosses indicate the point (I=46 W/m^{2}, q=0.12). Plot a) exhibits for linear polarization at Paranal for , which only varies weakly with B. At the blue cross, ph/s/sr/(W/m^{2}) (for q=0, it is only 44 ph/s/sr/(W/m^{2})!). Plot b) shows the same for circular polarization ( ph/s/sr/(W/m^{2}) at the blue cross, which is close to the maximum of in Plot b) and lies within the magnetic resonance ``island''). Plot d) shows the ratio . Circular polarization at only improves the return flux at Paranal if W/m^{2} (at smaller , there is a much higher gain, as indicated by Fig. 3). The improvement at the blue cross is 1.52 (with the Bfield at SOR, it would be only 1.2, but can reach 2.0 for smaller ). We have limited the plot range in the lower right corner of d) to 1.8 for better rendering; in reality, the value exceeds 3. Plot e) shows the same as b), but for the Bfield of SOR. As noted in Fig. 9, the optimal q lies closer to 19% under this condition. Note that even at q=0.5 (equal power in D_{2}a and D_{2}b), ph/s/sr/(W/m^{2}) is not much below its optimum of 201. Plot c) depicts the gain in when turning on repumping at standard conditions, i.e., . Plot d) shows the same for SOR. At the blue crosses equals 3.6 and 4, respectively.
4.2.4 Laser linewidth
Figure 11: a) Contour plot of in ph/s/sr/(W/m^{2}) for standard conditions (Paranal). b) Same, but for the geomagnetic field of SOR. The optimal spectral irradiance is near 10 and 17 W/m^{2}/MHz, respectively (dashed blue lines). 

Open with DEXTER 
So far, we have only simulated singlefrequency line lasers. We now turn to lasers with Lorentzian lines of finite bandwidth. Figure 11 shows two contour plots of for conditions at Paranal and SOR. The specific return is highest near the dashed blue lines, hence near 10 W/m^{2}/MHz (100 W/m^{2} per natural linewidth of sodium of 9.8 MHz) for Paranal and near 17 W/m^{2}/MHz for SOR. For irradiances below about 20 W/m^{2}, this linearity does not hold anymore and will yield higher return flux. Conversely, at W/m^{2}, is no longer negligible compared to the Doppler width, and hence the absorption cross section falls off in the wings of the laser spectrum (our simulations for q > 0 lose accuracy for wide bandwidths, see the discussion of Fig. 1). In the upper left region of the plot, decreases due to low pumping that is too weak to overcome the population equalization from Larmor precession, and in the lower right we suffer from saturation. We note that d'Orgeville (2000) has conducted a related study (d'Orgeville's Fig. 8 differs from ours in that its horizontal axis shows laser power and not irradiance, and it is rendered in double logarithmic scale).
Based on the dashed blue lines in Fig. 11, the optimum FWHM laser bandwidth can be estimated by , where is in MHz, I_{P/2} in W/m^{2} with W/m^{2}, and B in Gauss. We expect this optimal spectral power density to only depend on the relaxation time scales. This scaling relationship is good news for future LGS systems, since it means that we can achieve high return flux even at very large irradiances, as long as we keep increasing the laser bandwidth. Lasers with microsecond pulses designed for mesospheric spot tracking or lidar with duty cycles of 1:201:100, and hence high peak powers, may take advantage of this possibility.
4.3 Optimization of laser beams
In this subsection, we extend the analysis to entire laser beams by carrying out weighted integrals of for (Gaussian) spot profiles, from which we compute (Eq. (7)).
To begin, we conduct a comparison with the beam return flux computed by Milonni (1999). We apply our Bloch code to Milonni's Na conditions of B=0.5 G, q=0, , , P=1 W, W/m^{2}, circular polarization ( ), no recoil, , K, , , s)) and obtain ph/s/W/(atoms/m^{2}), in contrast to Milonni's value ph/s/W/(atoms/m^{2}) from Eq. (30) of that paper. If we switch to s) (longdashed curve in Milonni's Fig. 6a), we obtain ph/s/W/(atoms/m^{2}). For the combination s), s), we find ph/s/W/(atoms/m^{2}), underlining the importance of properly modeling diffusion in atom velocity space. We suspect that the discrepancy in the flux results is partly due to Milonni's (1999) greatly simplified Sdamping formula Eq. (18), which tends to overestimate spin relaxation.
Figure 12: Beam efficiency as a function of FWHM laser bandwidth for a P=20 W laser at Paranal (solid red curve, standard conditions), Mauna Kea (dashed green), and SOR (dashdotted blue), depending on the local geomagneticfield strengths (P=20 W, ). Dotted magenta (solid black): Laser at Paranal for P=217 W (P=435 W). 

Open with DEXTER 
Figure 12 shows a plot of the beam efficiency for different geomagnetic field strengths at (solid red: Paranal, dashed green: Mauna Kea, dashdotted blue: SOR). The dotted magenta and solid black lines show for Paranal at P=217 W ( I_{P/2} = 500 W/m^{2}) and P=435 W ( I_{P/2}=1000 W/m^{2}), respectively. At these irradiances, the peaks at 38 MHz and 75 MHz, respectively, without losing much of its magnitude compared to P=20 W. We expect that pulsed lasers with pulse durations much longer than 16 ns can enjoy this high specific return. Note that the numerical similarity of ph/s/W/(atoms/m^{2}) with typical values of in ph/s/sr/(W/m^{2}) is a mere coincidence, and the two quantities should not be confused.
We add the caveat that the laser bandwidth is defined here as the shortterm bandwidth, as opposed to the apparent bandwidth obtained from a longterm spectral measurement of a laser whose central frequency is slowly drifting, or being chirped.
5 Conclusions
In this article, we describe a modeling method to accurately determine the photon return flux from sodium LGS using Bloch (density matrix) equations. Numerical simulations are crucial, since it is not easy to achieve the desired optimization experimentally.
We summarize some optimization conclusions for cw sodium LGS lasers as follows:
 1.
 The next generation of 20 Wclass (launched) cw laser sodium LGS will achieve unprecedented mesospheric irradiances near I_{P/2}=50 W/m^{2} (I_{P/2} = ``50% power in the bucket irradiance'') in good seeing conditions, using launch telescopes with 3050 cm clear aperture diameter.
 2.
 Atomic Larmor precession due to the geomagnetic field B can completely suppress optical pumping and thus strongly reduce the LGS return flux, depending on the orientation between the laser beam and the field lines. LGS systems must be dimensioned for the beam pointing that yields the lowest return, i.e., laser beam orthogonal to .
 3.
 In the future LGS irradiance regime, circular polarization will significantly elevate the return flux, with a rather relaxed polarization extinction ratio of dB. The polarization chirality (sense) does not matter.
 4.
 Repumping (excitation of the D_{2}b line) with 1020% of the laser power in combination with circular polarization increases the return flux further by a factor of 34. Repumping can moreover alleviate spatial powerbroadening and hence reduce the observed LGS spot size.
 5.
 The optimum FWHM laser bandwidth can be estimated for a Lorentzian line shape by , where is in MHz, I_{P/2} in W/m^{2}( W/m^{2}), and B in Gauss. As a rule of thumb, the laser bandwidth in MHz should approximately equal the launched laser power in Watt divided by six, assuming a diffractionlimited spot size in good seeing. Employing lasers with a single line of 23 GHz bandwidth is usually several times less efficient.
 6.
 A laser with the above properties can achieve a specific return of photons/s/W/(atoms/m^{2}) at a lowB location such as in northern Chile and photons/s/W/(atoms/m^{2}) at a highB location such as New Mexico (all values for laser beam orthogonal to ).
We also compare the Bloch simulation results with our Monte Carlo rateequation simulation Exciter and find good agreement. However, such Monte Carlo rateequation methods are based on numerous assumptions, and they always need to be validated against more rigorous methods such as a Bloch equation calculation. Moreover, our Bloch code is orders of magnitude faster than Exciter, with the runtime advantage quickly growing with the desired level of accuracy. Our Bloch simulation code can simulate any alkali and is publicly available on our websites^{}.
In this work, we treat only steadystate sodium excitation with cw lasers. Pulsed lasers, with pulse durations close to or shorter than the sodium lifetime of 16 ns, will be dealt with in a forthcoming publication. However, we surmise that our present results are valid for microsecond laser pulses that are useful for mesospheric spot tracking and lidar. Onsky experiments to validate our results with a mobile 20 Wclass laser unit are being prepared at ESO at the moment, and are expected to yield results in 2010.
In conclusion, our simulations indicate that we have only begun to realize the full capability of sodium LGS, to be harvested by upcoming generations of laserassisted AO, both with cw and pulsed lasers. The future of sodium LGS, so to speak, looks very bright.
AcknowledgementsWe thank Craig Denman, Paul Hillman, Peter Hoffmann, Victor Kartoshkin, Ed Kibblewhite, Peter Milonni, Steve Morgan, Brian Patton, John Telle, Yan Feng, Nourredine Moussaoui, and Luke Taylor for their valuableinput. We would like to give special credit to Will Happer, the father of sodium LGS, for his generous advice. This work was supported in part by the NGA NURI program.
References
 Andrews, L. C., & Phillips, R. L. 2005, Laser Beam Propagation through Random Media (SPIE Press, Bellingham), 2nd edn. Bellanger, V., Courcelle, A., & Petit, A. 2004, Comp. Phys. Comm., 162, 143 [CrossRef] [Google Scholar]
 Bradley, L. C. 1992, J. Opt. Soc. Am. B, 9, 1931 [NASA ADS] [CrossRef] [Google Scholar]
 Budker, D., Gawlik, W., Kimball, D. F., et al. 2002, Rev. Mod. Phys., 74, 1153 [NASA ADS] [CrossRef] [Google Scholar]
 McClelland, J. J., & Kelley, M. H. 1985, Phys. Rev. A, 31, 3704 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Corney, A. 1977, Atomic and Laser Spectroscopy (Oxford: Oxford University Press) [Google Scholar]
 Cussler, E. L. 1997, Diffusion. Mass Transfer in Fluid Systems 2nd edn. (Cambridge: Cambridge University Press) [Google Scholar]
 Denman, C. A., Drummond, J. D., & Eickhoff, M. L. 2006, Proc. SPIE, 6272, 62721L [CrossRef] [Google Scholar]
 Dmitriev, S. P., Dovator, N. A., & Kartoshkin, V. A. 2008, Tech. Phys. Lett., 34, 693 [NASA ADS] [Google Scholar]
 Drob, D. P., et al. 2008, J. Geophys. Res., 113, A12304 [NASA ADS] [CrossRef] [Google Scholar]
 Drummond, J., Denman, C., & Hillman, P. 2007, ed. S. Ryan, 1215, E67 [Google Scholar]
 Fishbane, P. M., Gasiorowicz, S. G., & Thornton, S. T. 2005, Physics for Scientists and Engineers (PrenticeHall, Upper Saddle River) [Google Scholar]
 Fried, D. L. 1966, J. Opt. Soc. Am., 56, 1372 [NASA ADS] [CrossRef] [Google Scholar]
 Guillet de Chatellus, H., Pique, J.P., & Moldovan, I. C. 2008, J. Opt. Soc. Am. A, 25, 400 [NASA ADS] [CrossRef] [Google Scholar]
 Happer, W. 1972, Rev. Mod. Phys., 44, 169 [NASA ADS] [CrossRef] [Google Scholar]
 Happer, W., & Van Wijngaarden, W. A. 1987, Hyperfine Interactions, 38, 435 [NASA ADS] [CrossRef] [Google Scholar]
 Happer, W., MacDonald, G. J., Max, C. E., & Dyson, F. J. 1994, J. Opt. Soc. Am. A, 11, 263 [NASA ADS] [CrossRef] [Google Scholar]
 Hillman, P. D., Drummond, J. D., Denman, C. A., & Fugate, R. Q. 2008, Proc. of the SPIE, 7015, 70150L [CrossRef] [Google Scholar]
 Holzlöhner, R., Bonaccini Calia, D., & Hackenberg, W. 2008, Adaptive Optics Systems, Proc. of the SPIE, 7015, 701521 [Google Scholar]
 Holzlöhner, R., Bonaccini Calia, D., & Hackenberg, W. 2009, ESO Technical Report VLTTREESO118753940, Issue 1 (available upon request) [Google Scholar]
 Hubin, N. 2009, Proc. of AO for ELT conference, abstract 118 [Google Scholar]
 Kartoshkin, V. A. 1998, Opt. and Spectr., 85, 177 [NASA ADS] [Google Scholar]
 Milonni, P. W., & Thode, L. E. 1992, Appl. Opt., 31, 785 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Milonni, P. W., Fugate, R. Q., & Telle, J. M. 1998, J. Opt. Soc. Am. A, 15, 217 [NASA ADS] [CrossRef] [Google Scholar]
 Milonni, P. W., Fearn, H., Telle, J. M., & Fugate, R. Q. 1999, J. Opt. Soc. Am. A, 16, 2555 [NASA ADS] [CrossRef] [Google Scholar]
 Morris, J. R. 1994, J. Opt. Soc. Am. A, 11, 832 [NASA ADS] [CrossRef] [Google Scholar]
 Moussaoui, N., Holzlöhner, R., Hackenberg, W., & Bonaccini Calia, D. 2009, A&A, 501, 793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Moussaoui, N., Clemesha, B., et al. 2010, A&A, in press [Google Scholar]
 Nagengast, W., Nass, A., Grosshauser, C., Rith, K., & Schmidt, F. 1998, J. Appl. Phys., 83, 5626 [NASA ADS] [CrossRef] [Google Scholar]
 d'Orgeville, C., Rigaut, F., & Ellerbroek, B. L. 2000, Proc. SPIE, 4007, 131 [NASA ADS] [CrossRef] [Google Scholar]
 Patat, F. 2004, The Messenger, ESO Publication, 118, 11 [Google Scholar]
 Pique, J.P., Moldovan, I. C., & Fesquet, V. 2006, J. Opt. Soc. Am. A, 23, 2817 [NASA ADS] [CrossRef] [Google Scholar]
 Pfrommer, T., Hickson, P., & She, C.Y. 2009, Geophys. Res. Lett., 36, L15831 [NASA ADS] [CrossRef] [Google Scholar]
 Ramsey, A. T., & Anderson, L. W. 1964, Il Nuovo Cimento, 32, 1151 [Google Scholar]
 Ressler, N. W., Sands, R. H., & Stark, T. E. 1969, Phys. Rev., 184, 102 [NASA ADS] [CrossRef] [Google Scholar]
 Simonich, D. M., Clemesha, B. R., & Kirchhoff, V. W. 1979, J. Geophys. Res., 84, 1543 [NASA ADS] [CrossRef] [Google Scholar]
 Taylor, L., Friedenauer, A., et al. 2009, Proc. CLEO 2009, postdeadline PDA.7 [Google Scholar]
 Telle, J., Drummond, J., et al. 2006, Proc. of the SPIE, 6215, 62150K1 [CrossRef] [Google Scholar]
 Telle, J., Drummond, J., Hillman, P., & Denman, C. 2008, Proc. of the SPIE, 6878, 68780G1 [CrossRef] [Google Scholar]
 Thomas, S. J., Gavel, D., Adkins, S., & Kibrick, B. 2008, Proc. SPIE, 70155L [CrossRef] [Google Scholar]
 Ungar, P. J., Weiss, D. S., Riis, E., & Chu, S. 1989, J. Opt. Soc. Am. B, 6, 2058 [NASA ADS] [CrossRef] [Google Scholar]
 van der Vorst, H. A. 1992, SIAM J. Sci. Comp., 13, 631 [CrossRef] [MathSciNet] [Google Scholar]
 Wright, M. R. 2004, An Introduction to Chemical Kinetics (Wiley) [Google Scholar]
Footnotes
 ... techniques^{}
 ``Sodium Laser Guide Star Return Flux Study for The European Southern Observatory'', contract study for ESO, ETREKIB2970001, 26 August 2008 (available upon request).
 ... VLT^{}
 Adaptive Optics Facility (AOF) for the VLT, European Southern Observatory (ESO), http://www.eso.org/sci/facilities/develop/ao/sys/.
 ... package^{}
 Rochester, S. M., Atomic Density Matrix package for Mathematica, version 09.08.07 or later, available at http://budker.berkeley.edu/ADM/.
 ... Steck ^{}
 Steck, D. A., ``Sodium D Line Data'', version 2.1.1, 30 April 2009, http://steck.us/alkalidata/.
 ... model^{}
 IGRF Release 2005 geomagnetic model online at NOAA website http://www.ngdc.noaa.gov/IAGA/vmod/.
 ... model^{}
 MSISE90 atmospheric model online at NASA website http://ccmc.gsfc.nasa.gov/modelweb/atmos/msise.html.
 ... pm^{}
 Atomic beam scattering experiments yield significantly larger cross sections. We will provide details in a forthcoming publication. Online atomic radius list http://www.ccdc.cam.ac.uk/products/csd/radii/.
 ... websites^{}
 ESO Laser Systems Department website http://www.eso.org/sci/facilities/develop/lgsf/; LGS extension to Atomic Density Matrix package, http://budker.berkeley.edu/ADM/ (requires Mathematica v. 6 or higher).
All Tables
Table 1: Simulation parameters and their standard nominal values.
Table 2: Mesospheric spot parameters for a 40 cm launch telescope under standard conditions (P=20 W, )
Table 3: Sensitivity of to a 1% perturbation in some simulation parameters under standard conditions.
All Figures
Figure 1: Absorption cross section vs. detuning from the Na D_{2}a line center. Solid blue: simulation at I=0.001 W/m^{2} (q=0), dashed red: analytical cross section, green: simulation of effective cross section at standard conditions (I=46 W/m^{2}, q=0.12). Inset: reference coordinate system. 

Open with DEXTER  
In the text 
Figure 2: Atomic velocity distribution under standard conditions (but ), blue: Exciter, green: Doppler profile, magenta: Bloch equations, cyan: return light spectrum in the atomic frame. Inset: zoom on line center. 

Open with DEXTER  
In the text 
Figure 3: Specific return flux . Lines: Bloch equations, symbols: Monte Carlo rateequation simulation Exciter. a) Green, crosses: ; red, diamonds: ; blue, squares: (standard conditions); magenta, circles: (all with circular polarization and ). b) Same as a), but for linear polarization. c) Red and magenta lines, diamonds and circles: as in a); black, crosses and cyan, squares: same, respectively, without recoil. 

Open with DEXTER  
In the text 
Figure 4: Magnetic field impact. a) Blue curve: as a function of the Larmor frequency at I=1 W/m^{2}. Red dashed curve (almost obscuring the blue curve): fit (sum of two Lorentzians). Blue diamonds: Exciter. Solid green (dashed magenta) lines: at Paranal (SOR). b) F=2 fractional sublevel populations, solid (dashed) curves: Bloch (Exciter) simulation. 

Open with DEXTER  
In the text 
Figure 5: as a function of laser polarization ellipticity angle , where (0) marks circular (linear) states. The curves indicate different combinations of the magnetic field polar angle and azimuth angle . Note that is highest for circular polarization in all cases. 

Open with DEXTER  
In the text 
Figure 6: Contour plot of in ph/s/sr/(W/m^{2}) for magnetic azimuth at I=46 W/m^{2}. 

Open with DEXTER  
In the text 
Figure 7: as a function of repumping laser line frequency offset at I = 1 W/m^{2} (dashed green), I = 10 W/m^{2}(dashdotted blue), and I = 46 W/m^{2} (solid red). 

Open with DEXTER  
In the text 
Figure 8: Contour plot of in ph/s/sr/(W/m^{2}). 

Open with DEXTER  
In the text 
Figure 9: as a function of repumping power fraction q for different field angles . Dashdotted red curve: standard conditions, solid blue: (circular polarization, no Larmor precession), dashed black: standard conditions, but B=0.48 G (SOR), solid green: linear polarization ( ). 

Open with DEXTER  
In the text 
Figure 10: a), b) Contour plots of in ph/s/sr/(W/m^{2}) at Paranal for linear and circular polarization, respectively; d) the ratio at Paranal; e) same as b), but for SOR geomagnetic field; c), f) repumping gain at Paranal and SOR, respectively, circular polarization. Blue crosses: I=46 W/m^{2}, q=0.12. 

Open with DEXTER  
In the text 
Figure 11: a) Contour plot of in ph/s/sr/(W/m^{2}) for standard conditions (Paranal). b) Same, but for the geomagnetic field of SOR. The optimal spectral irradiance is near 10 and 17 W/m^{2}/MHz, respectively (dashed blue lines). 

Open with DEXTER  
In the text 
Figure 12: Beam efficiency as a function of FWHM laser bandwidth for a P=20 W laser at Paranal (solid red curve, standard conditions), Mauna Kea (dashed green), and SOR (dashdotted blue), depending on the local geomagneticfield strengths (P=20 W, ). Dotted magenta (solid black): Laser at Paranal for P=217 W (P=435 W). 

Open with DEXTER  
In the text 
Copyright ESO 2010