text:
# Penning collisions of laser-cooled metastable helium atoms

## I Introduction

Helium atoms in the metastable triplet state $`2^3S_1`$ (He\*) appear to be a good candidate for Bose-Einstein Condensation (BEC) according to theoretical predictions . The cross section for elastic collisions between spin-polarised metastable helium atoms is expected to be large, allowing efficient thermalization and evaporation in a magnetostatic trap, which is the standard technique to reach BEC . On the other hand, very high autoionization rates (Penning collisions) prevent reaching high densities of metastable helium atoms, both in the presence and in the absence of light, unless the sample is spin polarized.
If a metastable helium atom collide either with an other metastable atom, or with an helium atom excited in the $`2^3P_2`$ state, the quasi molecule formed can autoionize according to the following reactions:

$`\mathrm{He}(2^3S_1)+\mathrm{He}(2^3S_1)`$ $``$ $`\{\begin{array}{c}\text{He }(1^1S_0)\text{ + He}^+\text{ + }e^{}\hfill \\ \text{He}_2^+\text{ + }e^{}\hfill \end{array}`$ (3)
$`\mathrm{He}(2^3P_2)+\mathrm{He}(2^3S_1)`$ $``$ $`\{\begin{array}{c}\text{He }(1^1S_0)\text{ + He}^+\text{ + }e^{}\hfill \\ \text{He}_2^+\text{ + }e^{}\hfill \end{array}`$ (6)
A first experiment at subthermal energy ($`E=1.6`$ meV) with the metastable helium system was performed by M$`\ddot{\mathrm{u}}`$ller et al. , allowing the determination of the interaction potentials. Using those potentials the rate $`\beta _{SS}`$ for the reactions (3) has been calculated  to be a few $`10^{10}`$ cm3/s, which agrees with measurements performed in Magneto-Optical Traps (MOT) . According to theoretical predictions , the ionization rate corresponding to the reactions (3) should be suppressed by four orders of magnitude in a magnetostatic trap. Spin polarization of the atoms and spin conservation in the collisional process are the causes of this suppression, which makes the quest of BEC reasonable. Actually, a reduction of more than a factor of 20 in the two-body loss rate in an optically polarized sample was observed experimentally .

In presence of light exciting the transition $`2^3S_12^3P_2`$, the reaction (6) is dominant. “Optical collisions” with metastable helium atoms were measured to have surprisingly large cross sections when compared with alkali systems . The study of optical collisions is of fundamental importance in order to optimize the first step towards BEC, consisting in pre-cooling and trapping the atoms in a MOT. The goal is to transfer a cloud as dense as possible in a magnetic trap, in order to increase the elastic collision rate and start evaporation. The experimental study of optical collisions is the subject of this paper.

Several groups reported measurements of optical collisions rates, by studying losses in the MOT at small detunings  around $`5`$ MHz and at large detunings  at $`35`$ MHz and $`45`$ MHz. Measurements over a broad range of detunings, from $`5`$ MHz to $`20`$ MHz, were reported in  and the dependence of the loss rate on the intensity of the MOT laser beams was investigated. In reference  a theoretical model for optical collisions is also proposed predicting rates in good agreement with the measurements, but differing by more than one order of magnitude with all the other measurements previously quoted.

Our measurements are performed in a MOT loaded with $`10^9`$ atoms, at a peak density of $`10^{10}`$ atoms/cm3. With respect to previous works, we extend the measurements of the two-body loss rate to a wider range of detunings and intensities with a good precision, by measuring the number of atoms and the size of the trap using absorption techniques. Also, by measuring accurately the excited state population in each trapping condition, we are able to interpret our data with a simple model, expressing the two-body loss rate in terms of the excited state population and of a rate constant $`K_{sp}`$, found to be independent of the laser detuning and intensity.

Our experimental setup is described in section II, while in section III we explain our detection system and we give the working conditions and performance of our magneto-optical trap. In section IV we describe in detail the experimental procedure used to measure the two-body loss rate and the excited state population for different trapping conditions. The results are given in section V, and the conclusions in section VI.

## II Experimental set-up

A beam of metastable helium atoms is generated by a continuous high voltage discharge in helium gas, cooled to liquid nitrogen temperature. Radiation pressure on the metastable beam allows one to increase its brightness, and to deflect it from the ground state helium beam . The metastable atoms are then decelerated by the Zeeman slowing technique and loaded in a magneto-optical trap (MOT) in a quartz cell at a background pressure of $`5\times 10^{10}`$ torr. More details on the experimental setup will be given in a forthcoming paper . MOT parameters for optimal loading of the trap are listed in table I. For the laser manipulation of the atoms, we use the line at 1083 nm, connecting the metastable triplet state $`2^3S_1`$ to the radiative state $`2^3P_2`$. The saturation intensity $`I_{sat}`$ for this transition is 0,16 mW/cm2 and the linewidth $`\mathrm{\Gamma }/2\pi `$ is 1.6 MHz. Our laser system consists of a DBR laser diode (SDL-6702-H1) in an extended cavity configuration, injecting a commercial Ytterbium doped fiber amplifier (IRE-POLUS Group). The diode is stabilized by saturation spectroscopy at -240 MHz from resonance. At the fiber output we obtain 600 mW of power, in a TEM00 mode at the same frequency. The estimated linewidth is around 300 kHz. All the frequencies required for collimation, deflection, trapping and probing are generated by acousto-optical modulators in a double pass configuration, while we use directly part of the fiber output beam for slowing the atoms.

## III Detection system and characterization of the MOT

In order to fully characterize the cloud, we use a probe laser beam on resonance, whose diameter is about 1 cm, which is turned on 100 $`\mu `$s after the MOT field and light beams have been turned off. Our detection setup (see fig. 1) allows different measurements. With the combination of $`\lambda `$/2 plates and polarization beam splitter cubes, we can create either (i) a progressive wave, circularly polarized, passing through the atomic cloud towards a photodiode (PD1 in figure 1), giving the total absorption by the atoms, or (ii) a stationary wave, also circularly polarized, one arm of which is sent to a CCD camera, allowing spatially resolved absorption pictures of the cloud. A second photodiode (PD2 in figure 1) is used to collect the cloud fluorescence. We use the absorption photodiode PD1 to measure N, the number of atoms in the steady state of the MOT. The probe beam saturates the transition when the incident power exceeds 10 mW (see fig. 2). The maximum absorbed power is then $`P=Nh\nu \frac{\mathrm{\Gamma }}{2}`$. Our Watt-meter (Coherent lab-master) is calibrated to 3% accuracy and allows a rectilinear calibration fit of the photodiode voltage. We measure a maximum total absorption of 1 mW, corresponding to $`(1\pm 0.1)\times 10^9`$ atoms. We estimate the accuracy for the measurement of N to be about 10%.

The typical parameters of our magneto-optical trap with the operating conditions of table I are listed in table II.

We stress the fact that the case of He\* differs of that of alkalis, for which the imaging method gives a direct measurement of both the two-dimensional column density and the rms sizes of the MOT, by absorption of a brief and low intensity probe pulse ($`II_{sat}`$). In the case of He\*, the quantum efficiency of the CCD camera ($`10^3`$ at 1.083 $`\mu m`$) is too low to provide images with a sufficient signal to noise ratio. We need instead to illuminate the atoms with a 200 $`\mu s`$ pulse whose intensity is about 0.1 mW/cm2 ($`II_{sat}`$), and use a moderate magnification of 1/5. Another difficulty with He\* occurs from the large recoil momentum $`\mathrm{}k/m`$ (9.2 cm/s) due to the light mass of the atoms : the atoms are pushed out of resonance during the 200 $`\mu s`$ pulse if a traveling wave pulse is used. The solution we adopted is to illuminate the atoms in a standing wave with the set-up shown in figure 1. Though this scheme allows us to obtain pictures with a good contrast, the drawback is that the images obtained in the standing wave configuration for $`II_{sat}`$ are more difficult to analyze than in the low intensity case. In order to interpret the absorption pictures in the standing wave configuration, and for any saturation parameter, we developed a handy theoretical model (see appendix A) giving the column density of the atoms for each pixel of the CCD camera. The resulting density is then fitted by a Gaussian curve to extract the size of the cloud.

## IV Measuring the trap decay by fluorescence

Once the loading of the MOT is interrupted, the evolution of the number of trapped atoms $`N`$ is given by the following equation:

$$\frac{dN}{dt}=\alpha N\beta n^2(𝐫,t)d^3𝐫$$
(7)
where $`n(𝐫,t)`$ is the atomic density at position r and time t, $`\alpha `$ is the decay rate due to collisions between trapped atoms and the residual gas, and $`\beta `$ is the two body intra-MOT loss factor. Assuming that the spatial distribution is independent of the time evolution of the number of atoms, which is valid at low enough densities, one can write the density as

$$n(𝐫,t)=\frac{N(t)}{(2\pi )^{\frac{3}{2}}\sigma _x\sigma _y\sigma _z}e^{\frac{x^2}{2\sigma _x^2}\frac{y^2}{2\sigma _y^2}\frac{z^2}{2\sigma _z^2}}$$
(8)
At low enough pressure and high enough density, losses due to background gas are negligible, so that the equation reduces to

$$\frac{dN}{dt}=\beta \frac{N^2(t)}{(4\pi )^{\frac{3}{2}}\sigma _x\sigma _y\sigma _z}$$
(9)
whose solution is

$$N(t)=\frac{N(t_0)}{1+\frac{\beta }{2\sqrt{2}}n(\mathrm{𝟎},t_0)(tt_0)}$$
(10)
where $`t_0`$ is the initial time. In order to follow the evolution of the number of trapped atoms, we monitor the fluorescence decay of the MOT with a photodiode (PD2 in figure 1). As the fluorescence signal is proportional to the number of atoms, we obtain a fluorescence decay curve reproducing equation (10), which we fit to get the parameter $`\beta n(\mathrm{𝟎},t_0)`$. In order to determine $`\beta `$, one still has to measure $`n(\mathrm{𝟎},t_0)`$, which means that one has to measure the rms size of the cloud along the three directions and the initial number of atoms $`N(t_0)`$.
Our goal is to measure the loss rate for a wide range of detunings and intensities. The experimental procedure, divided in three successive steps, is the following.
(1) First, we load the trap for 1 s at $`\delta =45`$ MHz and at the highest intensity in the trapping beams (I/Isat=50 per laser arm). Then, we stop the loading by blocking the slowing beam with a mechanical shutter. 20 ms later, we “compress” the MOT by suddenly changing its detuning and intensity using acousto-optical modulators. We record the fluorescence signal during this procedure. A typical fluorescence curve is shown in figure 3.

The loading is stopped at t=-20 ms and the photodiode signal drops by a factor of 2 at t=-10 ms because the background light from the slowing beam is blocked. The fluorescence is greatly enhanced in the beginning of the compression phase at t=0 ms, as expected when the detuning is set closer to resonance (the detuning is set here to -20 MHz), but decays to almost zero in about 100 ms because of the two-body losses. Figure 4 shows the time evolution of the size of the cloud during this phase of compression, showing that 10 ms are enough to reach the new equilibrium size. Thus, we extract the parameter $`\beta n(\mathrm{𝟎},t_0)`$ from a fit of the fluorescence decay starting from $`t=t_0=10`$ ms. At this very time we measure the sizes of the MOT along x and y and the number of atoms in order to calculate $`n(\mathrm{𝟎},t_0)`$.

(2) Then, the sizes along the weak axis of the magnetic field gradient are measured by absorption on the CCD camera as explained in section III. Figure 5 shows the rms size along x for various laser detunings and intensities. The size along z (strong axis of the quadrupole field) is inferred from measurements of the sizes along x and y with a magnetic field gradient $`b`$ twice as large. We find a typical size along z 20% smaller than along the weak axes of the quadrupole. We did not correct the sizes for the expansion of the cloud during the pulse lasting 200 $`\mu s`$, as this would have required the measurement of the temperature for all the detunings and intensities. Nevertheless, we performed some time of flight measurements, giving temperatures ranging from 0.3 mK at -10 MHz to 1 mK at -40 MHz, from which we estimate that the sizes are overestimated at most by 5 % at -25 MHz and by 15 % at -5 MHz. In addition, we measured the statistical error on the sizes to be relatively small at large detunings, 2 to 3%, but larger at small detunings (about 10% at -5 MHz). This is due (i) to the poor spatial resolution of our imaging system (pixel dimension $`80\mu m\times 130\mu m`$), and (ii) to a low signal to noise ratio for small detunings where the loss rate is larger, as most of the atoms are lost during the compression phase.

(3) Finally, to determine the number of atoms that were still trapped at $`t_0=10`$ ms, we simultaneously switch off the magnetic field and set the trapping beams on resonance at $`t_0`$, instead of letting the trap decay as in figure 3. The laser intensity is set to a high enough value to strongly saturate the transition. We get a peak of fluorescence, whose amplitude is proportional to the number of atoms. We compare it with the peak obtained with the same procedure but for the MOT in the best loading conditions of figure 2, for which we measured the number of atoms precisely. From this comparison, we infer the number of atoms at $`t=t_0`$ in the compressed MOT, and thus determine $`n(\mathrm{𝟎},t_0)`$.
This measurement also gives access to the value of the average population of the excited state $`\pi _p`$. Indeed, $`\pi _p`$ is given by

$$\frac{F}{F_{max}}=\frac{\pi _p}{1/2}=2\times \pi _p$$
(11)
where $`F`$ is the fluorescence signal we measure in the compressed MOT at $`t_0`$, and $`F_{max}`$ the fluorescence signal at resonance, when the transition is saturated, and $`\pi _p`$ expected to be 1/2.
Figure 6 shows the results of the fluorescence measurements, giving $`F_{max}/F`$ as a function of the inverse of the laser intensity I for various detunings. It is interesting to note that the inverse of $`F`$ is found to vary linearly with the inverse of $`I`$.

Following , the fluorescence of N atoms in the compressed MOT can be modeled by the following equation:

$$F=\eta Nh\nu \frac{\mathrm{\Gamma }}{2}\frac{C_1\frac{I}{I_{sat}}}{1+C_2\frac{I}{I_{sat}}+4\frac{\delta ^2}{\mathrm{\Gamma }^2}}$$
(12)
where $`\eta `$ is the detection efficiency, $`I`$ is the total intensity of the six MOT beams, and $`C_1`$ and $`C_2`$ phenomenological factors. $`C_1`$ and $`C_2`$ would be 1 for a two-level atom, but they are expected to be smaller for an atom placed at the intersection of 6 differently polarized laser beams, as happens in a MOT. In reference , $`C_1`$ and $`C_2`$ are found to be equal, and slighly larger than the average of the squares of the Clebsch-Gordan coefficients over all possible transitions. For a $`J=1J=2`$ transition, this average is 0.56. We can rewrite equation (12) as

$$\frac{F_{max}}{F}=\frac{C_2}{C_1}+\frac{1+4\frac{\delta ^2}{\mathrm{\Gamma }^2}}{C_1}\frac{I_{sat}}{I}$$
(13)
where $`F_{max}=\eta Nh\nu \frac{\mathrm{\Gamma }}{2}`$.
The results of figure 6 show a good agreement with (13). But, $`C_2`$ and $`C_1`$ are not found equal, and both depend on the detuning. For example, $`C_1`$ is found to be 0.58, 0.48, 0.46, 0.44, 0.22 for $`\delta =`$-25, -20, -15, -10, -5 MHz respectively. We stress the fact that, for the fluorescence at resonance, and for full saturation, C1 and C2 are expected to be equal.

## V Results

The results of the Penning collisions rate $`\beta `$ are shown in figures 7 and 8.

Figure 7 presents the loss parameter $`\beta `$ as a function of the laser intensity for different detunings $`\delta `$, from -30 to -5 MHz. The uncertainty of the measurements varies from 25 $`\%`$ for large detunings to 60 $`\%`$ for small detunings. For all detunings, $`\beta `$ increases with power, which shows that S-P collisions are dominant.

Fig. 8 shows the loss parameter as a function of detuning for a fixed intensity ($`\frac{I}{I_{sat}}=80`$). For the same reason, the rate increases when the detuning goes to zero, as the population in the P state increases. Our results for $`\beta `$ agree with previous results  within the given error bars, extending the measurements to a wider range of parameters. For example, at -5 MHz and in an intensity range for which $`\beta `$ is not expected to vary strongly ($`I=140`$ to $`200I_{sat}`$), Kumakura et al.  find $`\beta =(4.2\pm 1.2)\times 10^8`$ cm3/s, Browaeys et al.  $`\beta =2\times 10^8`$ cm3/s with an uncertainty of a factor 2 and Tol et al.  $`\beta =(1.3\pm 0.3)\times 10^8`$ cm3/s. Our measurement $`\beta =(3.5\pm 1.4)\times 10^8`$ cm3/s agrees best with . One should also note that we find neither a decrease of $`\beta `$ for high intensities at small detunings, nor a decrease of $`\beta `$ at small detunings for a given intensity : this differs from the results of . In fact, we find that $`\beta `$ increases with intensity at small detunings, and also increases with decreasing detunings at a given intensity. We also disagree with the results of  where much smaller rates are found.
Finally, we also measured the loss rate in the trapping conditions ($`\delta =45`$ MHz, $`I=310I_{sat}`$) : the decay rate of the number of atoms was found to be $`\beta n(\mathrm{𝟎})=30`$s-1 at a density of $`10^{10}`$ at/cm3, which gives $`\beta =3\times 10^9`$ cm3/s.
One can further analyze these data following the simple model of  which relates the decay constant $`\beta `$ to the constant rate coefficients $`K_{ss}`$, $`K_{sp}`$ and $`K_{pp}`$ and to the populations of the excited and ground state levels, $`\pi _p`$ and $`\pi _s`$ respectively:

$$\beta =K_{ss}\pi _s\pi _s+2K_{sp}\pi _s\pi _p+K_{pp}\pi _p\pi _p$$
(14)
Experiments  or theory  have shown that the contributions $`K_{ss}\pi _s^2`$ and $`K_{pp}\pi _p^2`$ to the total rate $`\beta `$ are smaller than the $`K_{sp}`$ term by approximately two orders of magnitude.
From the measurements of the fluorescence signal in figure 6, we derive $`\pi _p`$ for each experimental point, as $`F/F_{max}`$ in eq. (11) is equal to $`2\times \pi _p`$.

In figure 9, we then plot $`K_{sp}`$ for the ensemble of our data. We do not see clear evidence for a dependence of $`K_{sp}`$ with the detuning or the intensity within the dispersion of our data. To a good approximation, we estimate then that $`K_{sp}`$ is actually constant in the explored range of parameters: $`K_{sp}=(1.0\pm 0.4)\times 10^7`$ $`\text{cm}^3`$/s, with a dispersion that roughly agrees with the error bars we claim. This result agrees with the first measurement ever performed , but the precision is now much improved. It also agrees well with the measurements of  where the authors found $`K_{sp}=(8.3\pm 2.5)\times 10^8`$ $`\text{cm}^3`$/s, assuming that for their parameters ($`\delta =5`$ MHz and $`I=30`$ mW/cm2), $`\pi _s=\pi _p=0.5`$.
An important point is that, in contrast with the measurement of the fluorescence at resonance where the transition is assumed to be saturated, $`\pi _p`$ in the compressed MOT never reaches 0.5 in our measurements : even for the smallest detuning and the highest intensity, $`\pi _p`$ is only 0.2. This explains why the results for $`\beta `$ in figure 7 at $`\delta =5`$ MHz strongly increase for increasing intensity over the whole explored range.

## VI Conclusion

We measured the absolute two-body loss rate between metastable atoms in a magneto-optical trap as a function of detuning and intensity. We extended the range of these parameters and compared the results to those of previous measurements, mostly performed at small detunings. Using a new experimental approach, we obtained reliable values for the two-body loss rates with an improved accuracy as compared to most earlier results. In the region of overlap of parameters, we find a good agreement with previous measurements, within the quoted uncertainties. We find a loss rate monotonically increasing as a function of intensity and decreasing with detuning. Our measurements are interpreted with a simple model, giving a rather constant loss rate $`K_{sp}`$, with an average value of $`(1\pm 0.4)\times 10^7`$ cm3/s, as already found in the very first measurement of . We believe that the quality and the extended range of our measurements should motivate more theoretical work, in order to understand better the peculiar dynamics of Penning collisions between metastable helium atoms in the presence of light.

Acknowledgments: The authors wish to thank C. Cohen-Tannoudji for helpful discussions and careful reading of the manuscript.

a Permanent address: Laboratoire de Physique des Lasers, UMR 7538 du CNRS, Université Paris Nord, Avenue J.B. Clément, 93430 Villetaneuse, France.

b Permanent address: Institute of Opto-Electronics, Shanxi University, 36 Wucheng Road, Taiyuan, Shanxi 030006, China.

c Permanent address : Dept. of Physics, Univ. of Perugia, Via Pascoli, Perugia, Italy; Lens and INFM, L.go E. Fermi 2, Firenze, Italy

d Present address : Universit$`\ddot{\mathrm{a}}`$t Hannover, Welfengarten 1, D-30167 Hannover, Germany.

e Permanent address : TIFR, Homi Bhabha Road, Mumbai 400005, India.

Unité de Recherche de l’Ecole Normale Supérieure et de l’Université Pierre et Marie Curie, associée au CNRS (UMR 8552).

## A Model of the absorption

In this appendix we describe the method we used to quantitatively interpret the absorption images of the atomic cloud when a standing wave configuration of the probe beam is used, and for an arbitrary saturation parameter. We describe the atoms as two-level atoms characterized by a non linear susceptibility:

$$\chi =n(x,y,z)\left[\frac{d^2}{\mathrm{}ϵ_0}\frac{\delta i(\mathrm{\Gamma }/2)}{(\mathrm{\Gamma }/2)^2+\delta ^2+|\mathrm{\Omega }|^2/2}\right]$$
(A1)
where $`n(x,y,z)`$ is the atomic density, $`d`$ the atomic dipole, $`\delta `$ the detuning, $`\mathrm{\Gamma }`$ the inverse lifetime of the excited state and $`\mathrm{\Omega }`$ is the Rabi frequency given by

$$\frac{\mathrm{}\mathrm{\Omega }}{2}=d^{(+)}E=^{(+)}e^{i\omega t}+c.c.$$
(A2)
where E is the electric field. The direction of propagation of the beam is $`z`$ and the field is supposed to be uniform in the plane (x,y). The propagation of the field is then described by the Maxwell equations:

$$\left[\mathrm{\Delta }+k_0^2(1+\chi )\right]\mathrm{\Omega }(z)=0$$
(A3)
where $`k_0`$ is the wavevector of the light.
The principle of the model is to use the slowly varying envelope approximation generalized to the case of a standing wave. We then decompose the probe beam field as:

$$\mathrm{\Omega }(z)=A_+(z)e^{ik_0z}+A_{}(z)e^{ik_0z}$$
(A4)
where $`A_+`$, $`A_{}`$ are the slowly varying envelopes of the wave going towards positive $`z`$ and negative $`z`$ respectively. A similar decomposition holds for the nonlinar susceptibility of the atoms:

$$\chi (z)=\chi _0(z)+\chi _+(z)e^{2ik_0z}+\chi _{}(z)e^{2ik_0z}+\mathrm{}$$
(A5)
where $`\chi _0`$, $`\chi _+`$ and $`\chi _{}`$ are slowly varying envelopes, and where we neglect terms in the expansion describing generation of frequencies others than the probe frequency via the non linear interaction.
If we insert the expansions (A4) and (A5) into the propagation equation (A3) and use the rotating wave appoximation, we obtain a set of two coupled differential equations for the slowly varying field amplitudes $`A_+`$, $`A_{}`$. By splitting the complex amplitudes into modulus and phase:

$$A_+=|A_+|e^{i\varphi _+}A_{}=|A_{}|e^{i\varphi _{}}$$
(A6)
and by introducing the phase difference $`(\varphi _+\varphi _{})`$ in the definition of the slowly varying susceptibilities $`\chi _+`$ and $`\chi _{}`$:

$$\chi _+=\stackrel{~}{\chi }_+e^{i(\varphi _+\varphi _{})}\chi _{}=\stackrel{~}{\chi }_{}e^{i(\varphi _+\varphi _{})},$$
(A7)
one can write :

$`{\displaystyle \frac{d|A_+|}{dz}}`$ $`=`$ $`{\displaystyle \frac{k_0}{2}}\left(\mathrm{Im}\stackrel{~}{\chi }_+|A_{}|+\mathrm{Im}\stackrel{~}{\chi }_0|A_+|\right)`$ (A8)
$`{\displaystyle \frac{d|A_{}|}{dz}}`$ $`=`$ $`{\displaystyle \frac{k_0}{2}}\left(\mathrm{Im}\stackrel{~}{\chi }_{}|A_+|+\mathrm{Im}\stackrel{~}{\chi }_0|A_{}|\right).`$ (A9)
By using expressions (A1) and (A4), the quantities $`k_0\mathrm{Im}\stackrel{~}{\chi }_+`$, $`k_0\mathrm{Im}\stackrel{~}{\chi }_{}`$ and $`k_0\mathrm{Im}\stackrel{~}{\chi }_0`$ are readily calculated:

$`k_0\mathrm{Im}\stackrel{~}{\chi }_0`$ $`=`$ $`{\displaystyle \frac{3\lambda ^2}{2\pi }}n(x,y,z)\alpha f_0`$ (A10)
$`k_0\mathrm{Im}\stackrel{~}{\chi }_+`$ $`=`$ $`k_0\mathrm{Im}\stackrel{~}{\chi }_{}={\displaystyle \frac{3\lambda ^2}{2\pi }}n(x,y,z)\alpha f_1`$ (A11)
where

$`\alpha `$ $`=`$ $`{\displaystyle \frac{(\mathrm{\Gamma }/2)^2}{(\mathrm{\Gamma }/2)^2+\delta ^2+\frac{1}{2}\left(|A_+|^2+|A_{}|^2\right)}}`$ (A12)
$`f_0`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{1ϵ^2}}};f_1={\displaystyle \frac{1f_0}{ϵ}}`$ (A13)
$`ϵ`$ $`=`$ $`{\displaystyle \frac{|A_+||A_{}|}{(\mathrm{\Gamma }/2)^2+\delta ^2+\frac{1}{2}\left(|A_+|^2+|A_{}|^2\right)}}.`$ (A14)
As a last step we eliminate the atomic density $`n(x,y,z)`$ from the equations by changing variable:

$$Z(z)=_{\mathrm{}}^zn(x,y,z^{})𝑑z^{}$$
(A15)
and we obtain the final coupled equations:

$`{\displaystyle \frac{d|\stackrel{~}{A}_+|}{dZ}}`$ $`=`$ $`{\displaystyle \frac{3\lambda ^2}{4\pi }}\alpha \left(f_1|\stackrel{~}{A}_{}|+f_0|\stackrel{~}{A}_+|\right)`$ (A16)
$`{\displaystyle \frac{d|\stackrel{~}{A}_{}|}{dZ}}`$ $`=`$ $`{\displaystyle \frac{3\lambda ^2}{4\pi }}\alpha \left(f_1|\stackrel{~}{A}_+|+f_0|\stackrel{~}{A}_{}|\right),`$ (A17)
where:

$$\stackrel{~}{A}_{}=A_{}/(\mathrm{\Gamma }/2)\stackrel{~}{A}_+=A_+/(\mathrm{\Gamma }/2).$$
(A18)
For $`\delta =0`$ and in the limit of small saturation parameters, one has $`\alpha =1`$, $`f_01`$, $`f_10`$ and one recovers the usual decoupled equation for low saturation absorption. We have now to solve equations (A16) and (A17). More precisely we wish to calculate the column density

$$Z^{\mathrm{}}=_{\mathrm{}}^+\mathrm{}n(x,y,z^{})𝑑z^{}$$
(A19)
for each effective pixel (x,y) of our image of the cloud. For each effective pixel, we can measure the initial conditions:

$`|\stackrel{~}{A}_+|^2(Z(\mathrm{})=0)`$ $`=`$ $`I_i`$ (A20)
$`|\stackrel{~}{A}_{}|^2(Z(\mathrm{})=0)`$ $`=`$ $`I_f`$ (A21)
corresponding respectively to the intensity of the probe beam before passing through the cloud, or equivalently without the atoms, and to the intensity of the probe beam that passed through the atomic cloud. For symmetry reasons, the column density (A19) is given by $`2Z_0=Z(0)`$, where $`Z_0`$ verifies

$$|\stackrel{~}{A}_+(Z_0)|^2=|\stackrel{~}{A}_{}(Z_0)|^2.$$
(A22)
For each pixel (x,y), we then integrate equations (A16)-(A17) numerically using the initial conditions (A20)-(A21) until $`|\stackrel{~}{A}_+(Z)|^2=|\stackrel{~}{A}_{}(Z)|^2`$. The corresponding $`Z`$ multiplied by 2 gives the column density. Note that, contrarily to what happens in the low saturation regime, we here need the values $`I_i`$ and $`I_f`$ separately (and not only their ratio), which implies a calibration of our CCD camera.