# Rejection of wavefront aberrations in an atomic gradiometer

Louis Pagot, Sébastien Merlet, Leonid A Sidorenkov, Franck Pereira dos Santos

► To cite this version:

Louis Pagot, Sébastien Merlet, Leonid A Sidorenkov, Franck Pereira dos Santos. Rejection of wavefront aberrations in an atomic gradiometer. 2025. hal-05354528

HAL Id: hal-05354528

https://hal.science/hal-05354528v1

Preprint submitted on 7 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Rejection of wavefront aberrations in an atomic gradiometer

LOUIS PAGOT, SÉBASTIEN MERLET, LEONID A. SIDORENKOV, AND FRANCK PEREIRA DOS SANTOS*

LTE, Observatoire de Paris, Université PSL, Sorbonne Université, Université de Lille, LNE, CNRS
61 Avenue de l'Observatoire, 75014 Paris, France

*franck.pereira@obspm.fr

**Abstract:** One of the main residual limitations of inertial sensors based on atom interferometry stems from laser beam distortions, which cause parasitic phase shifts and non-homogeneous matter-light couplings. Here we present numerical simulations, accompanied by analytical calculations, which quantify the impact of these effects in a cold atom gradiometer. We demonstrate that the propagation of interferometric laser beam aberrations, combined with initial asymmetry and significant time-of-flight expansion of the two atomic sources, limit the common-mode rejection of phase noise in a differential configuration. The resulting deviations in gravitational acceleration and its gradient are within reach of current experimental devices. Our study allows us to evaluate the surface quality requirements for retroreflective optics in cold-atom gradiometers of various baselines, and can be extended to other sensors based on different interferometer geometries.

## 1. Introduction

Since atom interferometry has been used to determine inertial quantities [1,2], this measurement principle has been applied to various areas of physics, as it allows for accurate and sensitive measurements [3-6]. Many systematic errors and noise sources can be eliminated in atomic interferometers by performing multiple measurements [7]. For example, the *k*-reversal technique allows to eliminate all systematics that are independent of the direction of the effective wave vector. Moreover, since atomic accelerometers, and gravimeters in particular, are limited in sensitivity by the vibrations of the retroreflection mirror, sensors based on a differential configuration, such as gradiometers, are marginally affected by vibrations that appear as a common-mode noise source [8-11].

The differential configuration also makes it possible to suppress systematic errors common to both interferometers, such as those related to delays between pulses, coupling inhomogeneities due to the finite width of the laser beam, or even Coriolis accelerations and wavefront aberrations, provided that the two atomic clouds have the same kinematic parameters. Thanks to the rejection inherent in this configuration, differential atomic interferometry is well suited to cutting-edge experiments. For example, beyond their direct application in geosciences, gradiometers are used to determine Newton's gravitational constant *G* [12, 13]. Other examples include long baseline atomic interferometers currently being developed for dark matter searches, as gravitational wave detectors [14], and to test the universality of free fall [15].

However, rejection based on this differential configuration is not necessarily perfect. The asymmetry of atomic sources and the evolution of the laser field during propagation will lead to residual contributions from the Coriolis effect on the one hand, and from wavefront aberrations on the other, which are the main sources of uncertainty in current experiments [7, 10, 16-18]. Due to their predominant contribution to the accuracy budget of the sensors, these wavefront aberrations have been the subject of numerous studies. In parallel with the development of theoretical frameworks [19, 20], experimental techniques are being elaborated and implemented to measure the effect of these aberrations in situ [21-23]. In some cases, the introduction of controlled distortions and the quantification of their impact with an atomic interferometer itself,
have made it possible to efficiently correct the associated systematic errors [19, 24–27].

In this article, we examine the impact of laser beam wavefront aberrations in a differential gradiometric configuration, based on numerical simulations and a theoretical framework previously developed to study the influence of optical aberrations on the accuracy of an atomic gravimeter [28]. We evaluate the phase residuals of these aberrations due to their finite rejection in a regime where the atomic cloud expands up to a size of the order of the laser beam width. Within this limit, we derive analytical formulas for the curvature of the Gaussian beam and for the distortions introduced by mirror surface defects, described by Zernike polynomials, which allows us to study the dependence of the bias on various experimental parameters. In addition, typical bias ranges are estimated in the presence of fluctuations in the initial transverse positions and velocities of the two atomic clouds, and in the specific cases of a 1-meter-long baseline gradiometer [9] and of a state-of-the-art 10-meter-high experiment [15, 29].

## 2. Gradiometer apparatus

The experiment considered here is a 2-meter-high cold atom gradiometer, illustrated in Figure 1. The output of the laser collimator is placed at height $z = 0$ m and the retro-reflecting mirror at the top, at position $z = 2$ m. Two $^{87}$Rb atomic clouds are prepared in magneto-optical traps at heights $z^{(l)} = 0.50$ m and $z^{(u)} = 1.50$ m, then further cooled in optical molasses down to a temperature of $\Theta = 2$ µK. At time $t_0 = 0$ s, they are launched upward at a speed of $1.29 \text{ m s}^{-1}$. After a Raman velocity selection pulse along the vertical axis, a Mach-Zehnder light pulse atom interferometer, starting at $t_1 = 15.8$ ms, is performed with a $\pi/2 - \pi - \pi/2$ sequence of Raman pulses separated by a time interval $T$ of typically 200 ms. In this case, the atomic clouds reach their apogee between the first and second pulses. The stimulated Raman transitions are performed using two laser beams with a wavelength $\lambda \approx 780.2$ nm. They arrive through a common collimator with two orthogonal circular polarizations and a radius at $1/e^2$ intensity of $w_0 = 5$ mm at position $z = 0$ m. A combination of a quarter-wave plate and the mirror located at the top of the experiment allows these polarizations to be exchanged for the reflected beams. During the interferometric sequence, the atoms fall freely in the Earth's gravitational potential modeled by the gravity acceleration $g = -9.81 \text{ m s}^{-2}$ and the vertical gravity gradient $\partial_z g = T_{zz} = 3.08 \times 10^{-6} \text{ s}^{-2}$ (or 3080 Eötvös, or 3080 E) at $z = 0$ m. To compensate for the Doppler effect and scan the interferometric fringes, the frequency difference between the two lasers is chirped at a rate $\alpha \approx 2\pi \times 25 \text{ MHz s}^{-1}$. Up to second order in the vertical gravity gradient [30, 31], the interferometric phase shift $\Delta\phi^{(j)}$ for the lower ($j = l$) and the upper ($j = u$) cloud in the ideal case of a plane wave laser beam is

$$ \Delta\phi^{(j)} = (k_{\text{eff}}g - \alpha)T^2 - k_{\text{eff}}T^2 T_{zz} \left( z_0^{(j)} + \bar{v}_0^{(j)}T - \frac{7}{12}gT^2 \right) \\ - k_{\text{eff}}T^4 T_{zz}^2 \left( \frac{7}{12}z_0^{(j)} + \frac{1}{4}\bar{v}_0^{(j)}T + \frac{1}{2}v_{\text{rec}}T - \frac{31}{360}gT^2 \right). \quad (1) $$

$k_{\text{eff}} \approx 2k = \frac{4\pi}{\lambda}$ is the effective wave vector of the stimulated Raman transition. $z_0^{(j)}$ and $v_0^{(j)}$ are respectively the vertical position and velocity of the atom at the beginning of the interferometer, and $\bar{v}_0^{(j)} = v_0^{(j)} + v_{\text{rec}}/2$ is the average velocity of the two paths just after the first interferometric pulse, with $v_{\text{rec}} = \frac{\hbar k_{\text{eff}}}{m}$ the two-photon recoil velocity of $^{87}$Rb atoms of mass $m$. The atomic distributions $w$ of each cloud are considered to be independent normal variables with a standard deviation $\sigma_{\rho}(0) = 0.3$ mm in position and $\sigma_{v} = \sqrt{\frac{k_B\Theta}{m}} = 13.8 \text{ mm s}^{-1}$ in velocity where $k_B$ is the Boltzmann's constant. After Raman selection along the vertical axis, the width of the vertical velocity distributions of the two clouds is reduced to approximately 6 mm s$^{-1}$. In addition, for the following analytical developments, we can consider that the Raman selection only modifies the longitudinal velocity distribution, since the cloud is relatively small compared to the size of
the laser beam during this pulse. The atomic cloud distributions after the Raman pulse are of the form

$$w(\vec{r}, \vec{v}) \, dx \, dy \, dz \, dv_x \, dv_y \, dv_z = \frac{e^{-\frac{x^2+y^2}{2\sigma_{\rho}^2(0)}} e^{-\frac{v_x^2+v_y^2}{2\sigma_{\rho}^2}}}{2\pi\sigma_{\rho}^2(0) 2\pi\sigma_{\rho}^2} \, dx \, dy \, dv_x \, dv_y \, \eta(z, v_z) \, dz \, dv_z \quad (2)$$

<figure><img src="image_2.png" /><figcaption>Fig. 1. Gradiometer configuration. The upper and lower atomic clouds are prepared in 3D MOTs separated by 1 m, launched at a velocity $v_z = 1.29 \, \text{m s}^{-1}$, and then an atomic interferometer is performed by successively applying three Raman pulses $\pi/2 - \pi - \pi/2$ during free fall. During the time of flight, the clouds expand ballistically and eventually reach the size of the laser beam in the transverse direction.</figcaption></figure>

Numerically, Monte-Carlo simulations are performed with 25 000 to 50 000 atoms in each atomic cloud and 11 to 21 frequency chirp rates $\alpha$ to scan the interferometer fringes. The parameters obtained by fitting the fringes are eventually averaged over 50 runs. The transition probabilities of the two interferometers are adjusted to the functions

$$P^{(j)} = P_m^{(j)} - \frac{C^{(j)}}{2} \cos \left( \Delta\phi^{(j)} + \delta\phi^{(j)} \right). \quad (3)$$

$P_m^{(j)}$ is the mean transition probability, $C^{(j)}$ is the fringe contrast, and $\delta\phi^{(j)}$ is a deviation from the ideal value $\Delta\phi^{(j)}$ given in (1). The bias on the vertical gravity gradient and gravity acceleration for each interferometer is calculated as

$$\delta T_{zz} = - \frac{\delta\phi^{(l)} - \delta\phi^{(u)}}{k_{\text{eff}} T^2 (z_0^{(l)} - z_0^{(u)})} \quad \text{and} \quad \delta g^{(j)} = \frac{\delta\phi^{(j)}}{k_{\text{eff}} T^2}. \quad (4)$$

Analytically, the general expression of the transition probability for an atom with initial coordinates $(\vec{r}, \vec{v})$ is

$$p^{(j)}(\vec{r}, \vec{v}) = p_m^{(j)}(\vec{r}, \vec{v}) - \frac{c^{(j)}(\vec{r}, \vec{v})}{2} \cos \left( \Delta\phi^{(j)} + \delta\phi^{(j)}(\vec{r}, \vec{v}) \right). \quad (5)$$
The total interferometer phase corresponds to a linear combination of the laser phase acquired by the atom at each pulse: $\phi_1^{(j)}(\vec{r}, \vec{v}) - 2\phi_2^{(j)}(\vec{r}, \vec{v}) + \phi_3^{(j)}(\vec{r}, \vec{v})$. It is here written as the sum of two terms. $\Delta\phi^{(j)}$: the phase contribution for a plane wave laser beam at the average coordinates of the atomic cloud, which is independent of the atom. $\delta\phi^{(j)}(\vec{r}, \vec{v})$: is related to the deviation of the wavefront with respect to the plane, and depends on the kinematic parameters of the atom. $c^{(j)}(\vec{r}, \vec{v})$ is the contrast for this atom. The transition probability for the entire cloud is obtained by integrating expression (5) over the atomic cloud distribution (2). Assuming small phase biases $|\delta\phi^{(j)}(\vec{r}, \vec{v})| \ll 1$ rad, the approximate expression of contrast [32] is

$$C^{(j)} = \int w(\vec{r}, \vec{v}) c^{(j)}(\vec{r}, \vec{v}) dx dy dz dv_x dv_y dv_z, \qquad (6)$$

and that of phase shift is

$$\delta\phi^{(j)} = \frac{1}{C^{(j)}} \int w(\vec{r}, \vec{v}) c^{(j)}(\vec{r}, \vec{v}) \delta\phi^{(j)}(\vec{r}, \vec{v}) dx dy dz dv_x dv_y dv_z. \qquad (7)$$

Calculating the expression (7) requires knowing the expression of the contrast $c^{(j)}(\vec{r}, \vec{v})$ up to a multiplicative constant. Furthermore, in the simple case where all atoms have a similar interferometric contrast $c^{(j)}(\vec{r}, \vec{v}) \approx \text{cte}$, expression (7) simplifies to the phase bias averaged over the atomic cloud distribution: $\delta\phi^{(j)} = \int w(\vec{r}, \vec{v}) \delta\phi^{(j)}(\vec{r}, \vec{v}) dx dy dz dv_x dv_y dv_z$. In this limit of uniform contrast, this relationship can be used to derive the contribution of a Gaussian beam curvature (8) [20] or that of aberrations described by a Zernike polynomial $Z_n^m$ [28]. In the latter case, for a Zernike polynomial defined on a disk of radius $R$, the phase shift is a function of: the amplitude of the aberration, the ratio between the atomic transverse size and the typical transverse length of the aberration $\sigma_\rho/l_{xy}$ with $l_{xy} = R/(n+1)$, and the ratio between the distance along the optical axis and the typical propagation length of the aberration $\Delta z/l_z$ where $l_z = 2kl_{xy}^2$.

### 3. Finite size gaussian beam

In the limit of uniform contrast, the phase bias introduced by the curved wavefront of a Gaussian beam [20], for centered clouds, is

$$\delta\phi_{\text{curv}} = 2(\beta_1 \sigma_{\rho}^2(t_1) - 2\beta_2 \sigma_{\rho}^2(t_2) + \beta_3 \sigma_{\rho}^2(t_3)). \qquad (8)$$

$\beta_i = k [R^{-1}(z_i) - R^{-1}(z_{i,\text{ref}})]/2$ is proportional to the difference between the inverse of the radius of curvature of the direct and reflected beams. Together with the contribution of the Gouy phase, the gravity biases for the two interferometers of the gradiometer are represented by dashed lines in Figure 2. These analytical calculations of the biases for the two interferometers diverge from the Monte-Carlo simulations (data points in Figure 2) as the time interval $T$ between pulses increases, since the assumption of uniform contrast, which requires an atomic cloud smaller than the waist of the laser beam, breaks down. As atomic clouds have normal distributions in transverse position and velocity, their size over time is $\sigma_\rho(t) = \sqrt{\sigma_\rho^2(0) + \sigma_v^2 t^2}$, which, for a free fall time of $2T \approx 400$ ms, corresponds to $\sigma_\rho = 5.5$ mm, which is similar to the waist $w_0$ of the Gaussian beam. We therefore expect the phase bias to also depend on the ratio between the size of the atomic cloud and the size of the beam: $\sigma_\rho/w_0$ and to reproduce previous results when this ratio is negligible.

To model the contrast loss due to a non-negligible ratio $\sigma_\rho/w_0$ in a simple way, we neglect the impact of inhomogeneities at the first and second pulses, and only take them into account at the final pulse, where the cloud is largest. In addition, since the initial longitudinal distributions of the atomic clouds are independent of the transverse distributions, other effects related to
<figure><img src="image_3.png" /><figcaption>Fig. 2. Gravity bias for the two interferometers of the gradiometer. Each data point corresponds to 50 Monte-Carlo simulations with 400 000 atoms to limit the error bars for T close to 25 ms. Solid (dashed) lines correspond to the analytical expression that (do not) take into account the contrast weights in (7).</figcaption></figure>

longitudinal coordinates, such as Doppler shift, are assumed to contribute only to a multiplicative constant in the expression of the contrast. Thus, using the analytical expression from [32], the contrast is approximated by

$$c(x, y, v_x, v_y) \propto \sin \left( \mathcal{A} \exp \left( -2 \frac{(x + v_x t_3)^2 + (y + v_y t_3)^2}{w_0^2} \right) \right) \quad (9)$$

$\mathcal{A}$ is the area of the third pulse at the center of the laser beam, which we take as $\pi/2$ in the following. With this expression for the loss of contrast, the calculations of the integrated contrast (16) and the contribution of the Gaussian beam curvature (17) are presented in **Appendix: Calculations with contrast weights**. The analytical series, evaluated up to $n = 5$ [33], are represented by solid lines in Figure 2 and exhibit a behaviour that is in good agreement with the Monte-Carlo simulations. We attribute the difference between this simple analytical model and the exact numerical simulation mainly to the non-uniformity of the coupling during the second pulse. Compared to the results corresponding to the uniform coupling case, displayed as dashed lines, we observe effects of the order of tens of nm s⁻² on the gravity bias and of tens of Eötvös on the gravity gradient, which are within reach of current experiments [10, 11, 17].

## 4. Mirror with Zernike polynomial surface

In Figure 2, the rejection of wavefront aberrations of the Gaussian beam is not ideal, even though the Rayleigh length $z_R \approx 100$ mm, which is the typical propagation length of the gaussian beam for the considered $w_0 = 5$ mm, is much larger than the typical longitudinal lengths of the system. We will now study the effect of an additional wavefront aberrations caused by optical surface defects of the retro-reflective mirror. As a simple example, we first assume that the Gaussian beam, the mirror surface, and the two clouds are coaxially centered and that the atomic clouds have no initial transverse velocity. The mirror surface is represented by a Zernike polynomial defined on a reference disk of diameter $2R = 27$ mm, of the order of the diameter of our mirror surfaces [28]. Only Zernike polynomials with rotation invariance $Z_n^0$, with $n$ even, should be considered, as the others have at least one axis of antisymmetry, which results in a null contribution to the
interferometer phase. In this case, an analytical approximation is developed, see **Appendix: Calculations with contrast weights**, similar to that presented in [28], which also takes into account the non-uniformity of contrast due to the last pulse, as in the previous case of the Gaussian beam curvature. This analytical approach can be adapted to more complex cases, such as an initially off-centered cloud. With a first-order Taylor expansion of the wavefront aberration and using the series S defined in expression (18), the phase bias induced by a mirror with surface $S_{\text{mir}} = \gamma Z_n^0$ can be approximated by

$$ \delta\phi^{(j)} \approx (-1)^{\frac{n}{2}} k_{\text{eff}} \gamma \left[ \cos\left(\frac{\Delta z_{1,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_1) - 2\cos\left(\frac{\Delta z_{2,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_2) \right. \\ \left. + \cos\left(\frac{\Delta z_{3,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_3) \right], \quad (10) $$

with $l_z^{(n)} = 2k(l_{xy}^{(n)})^2$, $l_{xy}^{(n)} = R/(n+1)$ and $\Delta z_{i,\text{mir}}^{(j)}$ the distance between the mirror surface and the position of the atomic cloud during the i-th pulse. From (10), the resulting bias on the vertical gravity gradient is

$$ \delta T_{zz} \approx (-1)^{\frac{n}{2}} \frac{\gamma}{T^2} \sum_{i=1}^{3} \xi_i \frac{\cos\left(\frac{\Delta z_{i,\text{mir}}^{(l)}}{l_z^{(n)}}\right) - \cos\left(\frac{\Delta z_{i,\text{mir}}^{(u)}}{l_z^{(n)}}\right)}{z_0^{(l)} - z_0^{(u)}} S(t_i), \quad (11) $$

with $\xi_1 = \xi_3 = 1$ and $\xi_2 = -2$.

<figure><img src="image_4.png" /><figcaption>Fig. 3. Biases in gravity acceleration (a) and its vertical gradient (b) for the $2T = 400$ ms gradiometer configuration, with a gaussian beam of waist $w_0 = 5$ mm and a mirror surface proportional to a Zernike polynomial with rotation invariance $S_{\text{mir}} = \frac{\lambda}{200}Z_n^0$ defined on a disk of diameter $2R = 27$ mm. Each data point is generated from 50 simulations with 50 000 atoms. The crosses correspond to the evaluation of equation (10), to which the bias of the Gaussian beam calculated with equation (17) is added. The series are expanded up to the 5th order.</figcaption></figure>

The results of the numerical simulations for a free fall time $2T = 400$ ms are represented as data points in Figure 3. The calculations of analytical expressions (10) and (11) up to $n = 5$, in addition to the analytical calculations of the Gaussian beam effect, are represented by crosses. For comparison purposes, the amplitude of the mirror aberration is the same for all polynomials, $\gamma = \lambda/200$, which corresponds to a typical upper limit for a real mirror [28]. We note that as long as the amplitude is sufficiently small, i.e. $k_{\text{eff}} \gamma \ll 1$ rad, linearization with respect to
the amplitude of the aberrations is valid for the full range of Zernike polynomials considered here, since the simulated biases closely correspond to analytical expressions (10) and (11). The cosines in expression (11) result from the propagation of aberrations, and factorisation with their difference is specific to ideal centered case considered here. When the propagation is neglected or lower-order aberrations are considered, i.e., when the gradiometer baseline is negligible compared to the propagation length: $|z_0^{(l)} - z_0^{(u)}| \ll l_z^{(n)}$, the rejection can be considered ideal. In this case, the bias on the vertical gravity gradient in Figure 3 (b) is of -5.6(1.3) E, which corresponds to the contribution of the Gaussian beam reflected on a perfectly flat mirror of -5.5(5) E in Figure 2. On the other hand, rejection deteriorates when these two lengths are comparable, which, with the parameters used, corresponds to $n \approx 37$, of the order of the threshold observed in the simulations. For higher-order aberrations, the rejection depends on both the cosines difference and the envelope defined by the series $S$ in expression (11). For $(\frac{\sigma_\rho(t_2)}{R/(n+1)})^2 \ge 10$ equivalent to $n \ge 16$, the contribution of wavefront aberrations of the last two pulses of the interferometer is strongly averaged, so that the effect of the first pulse dominates. This is clearly visible in Figure 3 (a), as the zero of the lower interferometer for $n \approx 54$ equivalent to $l_z^{(n)} = 0.97$ m $\approx 2\Delta z_{1,\text{mir}}/\pi$, corresponds to the first cancellation of the cosine in (10). By replacing the series $S(t_1)$ with $\exp(-\frac{1}{2} \left( \frac{\sigma_\rho(t_1)}{l_{xy}^{(n)}} \right)^2)$, its expression in the uniform contrast limit, the maximum bias on the vertical gravity gradient is expected for $n \approx 65$ before decreasing to zero, as in the simulations shown in Figure 3 (b).

## 5. Sensitivity to initial transverse conditions

Until now, atomic clouds were considered to be coaxially centered with the Gaussian beam and the mirror surface. In practice, the initial positions of the clouds could be shifted in the transverse plane. They could also have non-zero initial transverse velocities, which would result in variable position offsets over throughout the interferometer. These position offsets must be compared to the typical transverse length of the aberrations under consideration. For instance, in the extreme case where the typical transverse length $l_{xy}$ of the aberrations is much larger than the typical offsets under consideration, the sensitivity to these initial conditions can be neglected.

To evaluate the sensitivity of the experiment to initial transverse conditions, Monte-Carlo simulations are repeated with similar atomic clouds (same widths for position and velocity distribution) whose average transverse conditions $(\langle \vec{r} \rangle(0), \langle \vec{v} \rangle(0))$ are shifted relative to the previous co-centered configuration where we had $(\langle \vec{r} \rangle(0), \langle \vec{v} \rangle(0)) = (0, 0)$. Assuming that there is no preferred direction for the fluctuations in the transverse coordinates, and that the two atomic clouds behave independently, the average initial positions $(\langle x \rangle(0), \langle y \rangle(0))$ and average initial velocities $(\langle v_x \rangle(0), \langle v_y \rangle(0))$ of each cloud are drawn from independent centered normal distributions. These distributions are characterized for positions by a standard deviation $\sigma_\rho(0) = 0.2$ mm and for velocities by a standard deviation $\sigma_{v_{x,0}, v_{y,0}} = 1$ mm s$^{-1}$. These fluctuation amplitudes are realistic compared to those measured in experiments [18, 34].

Since the clouds are no longer coaxially centered, Zernike polynomials with anti-symmetric axes will have an effect on measurement bias. An example is shown in Figure 4, where the mirror surface is $S_{\text{mir}} = \frac{d}{200} Z_{25}^{-1}$. The resulting vertical gravity gradient distribution is characterized by a mean of -6.2 (7) E and a standard deviation of 16.3 (5) E. The uncertainty on the standard deviation is defined as the squared root of the statistical contribution [35] and the uncertainty resulting from the fringe fitting. The mean value is consistent with -5.5(5) E, the bias calculated for the interferometer of $2T = 400$ ms with the atomic clouds being coaxially centered with the Gaussian beam reflected on a perfectly flat mirror in Figure 2. Since the main contribution to the bias is imprinted during the first Raman pulse of the interferometer and $\sigma_{v_{x,0}, v_{y,0}} t_1$ is an order
<figure><img src="image_5.png" /><figcaption>Fig. 4. Simulation of a gradiometer with $2T = 400$ ms, $w_0 = 5$ mm, and with (a) a mirror surface $\frac{d}{200}Z^{-1}_{25}$. (b) 500 simulations are performed with different initial mean transverse positions chosen from $N(0, \sigma_{x_0,y_0} = 0.2 \text{ mm})$ and initial mean transverse velocities chosen from $N(0, \sigma_{v_{x_0,y_0}} = 1 \text{ mm s}^{-1})$. Each simulation consists of 50 evaluations of the interferometers phase with 25 000 atoms. The biases are plotted as a function of the difference between the initial mean positions of the two clouds.</figcaption></figure>

of magnitude smaller than $\sigma_{x_0,y_0}$, the bias depends mainly on the difference between the initial average positions of the two clouds, as represented in Figure 4 (b). The axis of antisymmetry of the mirror surface is reflected in the sign of the bias of the vertical gravity gradient.

These simulations are repeated for different Zernike polynomials [36], and the standard deviations obtained are reported in Figure 5. For lower-order aberrations, the impact is lower since the typical lengths of the clouds $\sigma_\rho, \sigma_{x_0,y_0}$ are small compared to the characteristic size $l_{xy}^n$ of the aberrations. With these initial transverse fluctuations, the bias on the vertical gravity gradient for Zernike polynomials with rotation invariance $m = 0$ and indices $n \le 20$ can vary by a few Eötvös. This effect is significant compared to the ideal case illustrated in Figure 3 (b), where only the effect of the gaussian beam is initially visible. For high-order aberrations, with an index $n$ larger than $\approx 100$, the sensitivity returns to zero due to the averaging of aberrations with smaller characteristic transverse lengths $l_{xy} = R/(n+1)$.

The dependence of the Zernike polynomials near the center of the disk is

$$Z_n^{|m|} \propto \rho^{|m|} \cos(|m|\theta) + O(\rho^{|m|+2} \cos(|m|\theta)). \quad (12)$$

This local behaviour explains why $Z_n^{|m|=1}$, which is locally linear in one direction, is more sensitive to fluctuations in the initial transverse conditions than $Z_n^{|m|=0$ and $Z_n^{|m|=2$, which are locally quadratic. Moreover, near the center, the coefficient of the quadratic term of $Z_n^{|m|=2$ is actually smaller, by a factor of 1/2, than that of $Z_n^{|m|=0$, which corresponds approximately to the ratio of the calculations in Figure 5. For Zernike polynomials with higher $m$ indices, which are locally flatter, the biases caused by wavefront aberrations will thus be even less pronounced. This definition of the sensitivity of the measurement to mirror surface aberrations described by Zernike polynomials in the linear amplitude regime could be used to define the mirror specifications.

To calculate analytical estimates of the standard deviation for clouds that initially have non-zero transverse conditions, analytical formulas (10) and (18) must be extended. As indicated in [28], within the limit of uniform contrast, for a given Zernike polynomial $Z_n^m$ and for an atomic cloud close to the center of the mirror surface, i.e. $\sqrt{\langle x \rangle^2 + \langle y \rangle^2} \ll R$, the exponential term describing the averaging of aberrations is identical to that of a centered cloud and for the Zernike polynomial $Z_n^0$, so that the bias is simply multiplied by the value of the Zernike polynomial defined at the average position of the cloud $(\langle x \rangle, \langle y \rangle)$. Assuming that this relation remains valid in the present
<figure><img src="image_6.png" /><figcaption>Fig. 5. Standard deviation of the vertical gravity gradient $\sigma (\delta T_{zz})$ for a mirror surface $S_{mir} = \frac{d}{200} Z_n^m$ and fluctuations in initial transverse positions: $N (0, \sigma_{x_0,y_0} = 0.2 \text{ mm})$ and velocities: $N (0, \sigma_{v_{x,0},v_{y,0}} = 1 \text{ mm s}^{-1})$. Each data point is calculated over 500 runs, as in Figure 4 (b), with error bars corresponding to $\pm 3$ times the uncertainty [35]. Analytical estimates are calculated from 20 000 random initial conditions.</figcaption></figure>

case, since the displacements considered are small compared to the mirror radius and the beam waist: $\sigma_{x_0,y_0}, \sigma_{v_{x,0},v_{y,0}} t_3 \ll R, w_0$, the formula for the bias (10) is extended to

$$ \delta\phi^{(j)}[Z_n^m] = k_{\text{eff}}\gamma \left[ \begin{aligned} & Z_n^m(\langle x(t_1) \rangle, \langle y(t_1) \rangle) \cos\left(\frac{\Delta z_{1,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_1) \\ & -2Z_n^m(\langle x(t_2) \rangle, \langle y(t_2) \rangle) \cos\left(\frac{\Delta z_{2,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_2) \\ & +Z_n^m(\langle x(t_3) \rangle, \langle y(t_3) \rangle) \cos\left(\frac{\Delta z_{3,\text{mir}}^{(j)}}{l_z^{(n)}}\right) S(t_3) \end{aligned} \right]. \quad (13) $$

($\langle x(t_i) \rangle, \langle y(t_i) \rangle$) are the average transverse positions of the cloud at the time of the $i^{th}$ pulse. The standard deviation of the vertical gravity gradient values estimated using this formula is displayed as solid lines in Figure 5 on 20 000 random initial conditions. The results agree well with the numerical simulations. The largest deviations are observed for low $n$ indices, which can be explained by the replacement of Zernike radial polynomials by Bessel functions, which is more accurate for high $n$ indices.

## 6. Scaling

In this section, we discuss how biases on the measured gravity acceleration and its vertical gradient vary as a function of experimental parameters.

First, scaling with the interferometer geometry parameters, the free-fall time and the baseline of the gradiometer, can be deduced either from expression (7), by majoring the phase of each atom by the peak-to-valley of the phase fluctuations or from equation (13) as Zernike polynomials and cosines are bounded by 1 in absolute value. Denoting $k_{\text{eff}}\Delta\gamma$ the maximum amplitude of the
phase fluctuations at the three pulses, the biases in gravity accelerations and vertical gradient are bounded by

$$|\delta g^{(j)}| \le \frac{4\Delta\gamma}{T^2} \quad \text{and} \quad |\delta T_{zz}| \le \frac{4\Delta\gamma}{T^2(z_0^{(u)} - z_0^{(l)})} \qquad (14)$$

The factor 4 results from the sum of the laser phases at the three pulses, the phase of the laser at the mirror pulse being imprinted twice. Thus, increasing the free-fall time and the baseline of the gradiometer allows to reduce the biases in this worst case scenario. Though, even with the chosen parameters, $2T = 400$ ms and $z_0^{(u)} - z_0^{(l)} = 1$ m necessitating a 2-meter-high experiment, expression (14) with $\Delta\gamma = \lambda/200$ corresponds to a bias on the vertical gravity gradient of 390 E. Hopefully, experimentally the bias is lower as correlation between the phases fluctuations with large transverse aberrations length and efficient averaging of those with typical transverse length smaller than the atomic cloud size allows for a decrease of the wavefront aberrations effect as shown in Figure 3.

To study the scaling of the wavefront aberrations effect with the kinematic parameters of the cloud, we keep the geometry of the gradiometer defined above and we chose to examine the vertical gravity gradient bias for mirror surfaces described by Zernike polynomials $Z_{34}^0$ and $Z_{54}^0$ with the same amplitude $\lambda/200$, which cause biases of the order of tens of Eötvös, as shown in Figure 3 (b). For atomic clouds axially co-centered with the mirror surface and the laser beam, and with no average transverse velocity, the bias on the gravity gradient is calculated as a function of the initial size $\sigma_\rho(0)$ of the cloud for different temperatures in the range [2 nK, 2 $\mu$K]. In order to focus on the contribution of the mirror surface, the contribution of the Gaussian beam is eliminated by considering the difference with simulations performed with a perfectly flat mirror $S_{\text{mir}}^{\text{flat}} \propto Z_0^0$. The absolute value of the difference is represented in solid symbols in Figure 6. Finally, for clarity, for each temperature, the error bars defining the uncertainties in the simulation results have been connected and the resulting areas filled with the corresponding color, as shown in Figure 6. These uncertainties represent the standard deviation of the results of 50 experiments with 50 000 atoms in each clouds.

<figure><img src="image_7.png" /><figcaption>Fig. 6. Absolute value of the vertical gravity gradient bias for a mirror surface (a) $S_{\text{mir}} = \frac{\lambda}{200}Z_{34}^0$ and (b) $S_{\text{mir}} = \frac{\lambda}{200}Z_{54}^0$ without fluctuations in initial conditions (positions and velocities) as a function of the square of the initial cloud size $\sigma_\rho(0)$ and for different temperatures $\Theta \in \{2, 20, 60, 2000\}$ nK. The shaded areas represent uncertainties. Each value is determined from 50 simulations of the experiment with 50 000 atoms in each cloud.</figcaption></figure>

Regardless of temperature, all simulations in Figure 6 show an exponential decrease proportional to the square of the initial size $\sigma_\rho(0)$ of the clouds, and the slope in Figure 6 (b) for $Z_{54}^0$ is higher
than in Figure 6 (a) for $Z_{34}^0$. This agrees with the analytical expression (18), as this decrease is faster for aberrations with a smaller typical transverse size (with $l_{xy} \propto 1/(n+1)$), due to averaging over the size of the cloud. Thus, when the clouds initial size $\sigma_\rho(0)$ is large enough so that the effects of the mirrors surface are averaged to zero, only the effect due to the curvature of the gaussian beam remains and, as it has been subtracted in Figure 6, the difference is compatible with zero within the error bars.

The temperature dependence is more complex. While the temperature-dependent curves are ordered in Figure 6 (b) for the mirror surface proportional to $Z_{54}^0$, this is not the case in Figure 6 (a) for the mirror surface proportional to $Z_{34}^0$. At lower temperatures, the cloud expands less during its free fall, so that the averaging during the second and third pulses is similar to that of the first pulse. Depending on the cosine terms related to propagation in equation (11), these contributions can either increase or decrease the total bias. As illustrated in Figure 6, the uncertainties represented by the shaded areas increase as the clouds get larger and experience more inhomogeneities from the laser beam.

## 7. Simulations with real mirror surface

So far, we have examined a simplified model based on a mirror surface described by a single Zernike polynomial with a fixed amplitude. In this section, the contribution to the gravitational acceleration and gradient bias is calculated from an actual mirror surface measured using a Fizeau interferometer [37]. The surfaces shown in Figures 7 (a) and 7 (b) are the same as in [28], and their power spectral density [38] (PSD) is illustrated in Figure 7 (d).

As in section **Sensitivity to initial transverse conditions**, 500 numerical simulations are used to evaluate the sensitivity of the measurements to typical variations in the initial transverse positions $\sigma_{x_0,y_0} = 0.2$ mm and transverse velocities $\sigma_{v_{x,0},v_{y,0}} = 1$ mm s$^{-1}$. In addition, since the mirror surfaces have been decomposed into Zernike polynomials [28], estimates of the mean and standard deviation of the bias distributions [40] can also be calculated from the Zernike coefficients using formula (13). The results of the simulations and analytical estimates with the

<table><caption>Table 1. Mean and standard deviation of the bias on the gravity acceleration and its vertical gradient with mirror surfaces SN1 and SN2.$^{a,b}$</caption><thead><tr><th></th><th>SN1 num.</th><th>SN1 Zer.</th><th>SN2 num.</th><th>SN2 Zer.</th></tr></thead><tbody><tr><td>⟨δg<sup>(l)</sup>⟩ (nm s<sup>-2</sup>)</td><td>7.6(1)</td><td>10.7(7)</td><td>6.1(2)</td><td>5.6(2)</td></tr><tr><td>σ(δg<sup>(l)</sup>) (nm s<sup>-2</sup>)</td><td>1.7(1)</td><td>3.3(7)</td><td>3.2(1)</td><td>5.5(2)</td></tr><tr><td>⟨δg<sup>(u)</sup>⟩ (nm s<sup>-2</sup>)</td><td>18.1(2)</td><td>22.2(7)</td><td>8.7(2)</td><td>9.8(2)</td></tr><tr><td>σ(δg<sup>(u)</sup>) (nm s<sup>-2</sup>)</td><td>5.1(1)</td><td>6.4(7)</td><td>3.7(1)</td><td>5.6(2)</td></tr><tr><td>⟨δT<sub>zz</sub>⟩ (E)</td><td>-10.5(3)</td><td>-11.4(9)</td><td>-2.7(2)</td><td>-4.2(3)</td></tr><tr><td>σ(δT<sub>zz</sub>) (E)</td><td>5.4(1)</td><td>7.2(9)</td><td>5.0(2)</td><td>7.8(3)</td></tr></tbody></table>

$^a$The results of numerical simulations (num.) are calculated over 500 runs. The results obtained with Zernike decomposition of the surfaces (for polynomials with $|m| \in \{0, 1, 2\}$ up to order 2500 in ISO-14999) and the extrapolated equation (13) (Zer.) are calculated over 20000 runs.

$^b$Uncertainties in (Zer.) result from the calculation of the coefficients $(c_i \pm \delta c_i)_{1 \le i \le 2500}$ of the Zernike polynomials. As changing the coefficients leads to correlated errors for every atoms, the uncertainty is calculated by performing the quadratic sum of (13) with amplitudes $(\delta c_i)_{1 \le i \le 2500}$.

surfaces of the two mirrors are summarized in Table 1. The average biases obtained from the
<figure><img src="image_8.png" /><figcaption>Fig. 7. Map of the real mirror surfaces SN1 (a) and SN2 (b) measured using a Fizeau interferometer, with their PSD (d). (c) Example of a mirror surface generated numerically using the PSD model, similar to [39], represented by the black dashed line in (d).</figcaption></figure>

numerical and analytical calculations are marginally consistent, given their uncertainties. As for the standard deviations, those estimated analytically are generally larger than those simulated, which we attribute to the overestimation in the formula (13) of the effect of lower-order Zernike polynomials, as illustrated in Figure 5. Nevertheless, when the Zernike decomposition of the mirror surface is available, equation (13) can be used to quickly estimate the bias induced by surface defects as well as an upper limit of the expected standard deviation for given initial transverse distributions. The average bias and associated standard deviation of the order of ten Eötvös on gravity gradient measurements with mirror surfaces 7 (a) and 7 (b), correspond to the typical errors encountered in existing experiments [10].

## 8. Simulations for a 10 m gradiometer

To conclude this study on wavefront aberrations, we perform simulations in a state-of-the-art 10-meter-long gradiometer [15, 29]. The beam has a radius at $1/e^2$ in intensity of $w_0 = 2.3$ cm, while the initial transverse size is $\sigma_\rho(0) = 300$ µm and the expansion temperature is 50 nK for
both clouds. The mirror is placed at a height $z = 10.10$ m, its diameter is $2R = 16.4$ cm and its surface, assumed to be isotropic, is generated according to the PSD theoretical model in Figure 7 (d). This theoretical model, similar to that of [39], is significantly better than the real mirror surfaces SN1 and SN2 shown in Figures 7 (a) and 7 (b), although it corresponds to existing mirrors [39, 41]. As in previous simulations, the two atomic clouds are prepared at two positions separated by the baseline length. They are then released, with zero initial average longitudinal velocity, and fall into the gravitational field. At the end of the interferometer sequence, which lasts $2T = 1$ s and starts at $t_1 = 15.8$ ms, they have traveled an approximate distance of 5 m. In the following simulations, where the distance between the two clouds will be modified, the upper cloud is always assumed to be prepared at a height $z = 10$ m.

<figure><img src="image_9.png" /><figcaption>Fig. 8. Standard deviation of the vertical gravity gradient for different baseline lengths. (blue dots) The standard deviation is calculated over 100 mirror surfaces for coaxially centered clouds. Standard deviation obtained for clouds with initial transverse fluctuations (orange squares) {$\sigma_{x_0,y_0} = 0.2$ mm, $\sigma_{v_{x,0},v_{y,0}} = 1$ mm s$^{-1}$} and (green triangles) {$\sigma_{x_0,y_0} = 0.1$ mm, $\sigma_{v_{x,0},v_{y,0}} = 0.5$ mm s$^{-1}$} and averaged over 20 mirror surfaces.</figcaption></figure>

First, in the ideal case where both atomic clouds are co-centered with the mirror and the laser beams, i.e. $\langle \vec{r} \rangle (0) = \vec{0}$, and have no average transverse velocity, $\langle \vec{v} \rangle (0) = \vec{0}$, the bias on the vertical gravity gradient is calculated for different baseline lengths in the range [0.5, 5] m. The typical magnitude of the bias for a random mirror surface with the same PSD is evaluated by calculating the standard deviation over 100 mirror surfaces. It is displayed as blue dots in Figure 8. Remarkably, in the range of baselines explored here, its scaling is not proportional to the inverse of the baseline length. This can be explained by the bias expression (11), since the difference between the cosine terms that account for the propagation of the aberrations also depends on the baseline. For the PSD model considered here, even in this ideal case, the bias on the gravity gradient is typically of the order of 0.1 E, which corresponds to the order of magnitude of the accuracy required for the most demanding scenario in [15].

On top of the bias in the ideal case due to the specific defects of the mirror surface, additional contributions arise from fluctuations in the initial transverse conditions of the atomic clouds. As before, typical fluctuations in position $\langle \vec{r} \rangle (0) \sim N(0, \sigma_{x_0,y_0} = 0.2$ mm) and in velocity $\langle \vec{v} \rangle (0) \sim N(0, \sigma_{v_{x,0},v_{y,0}} = 1$ mm s$^{-1})$ are assumed. Since typical defect fluctuations are similar for mirror surfaces with the same PSD, the phase bias standard deviation caused by fluctuations
of initial conditions should also be similar for different mirror surfaces of this type. We therefore calculate the standard deviations over 100 runs with fluctuating positions $\langle \vec{r} \rangle (0)$ and velocities $\langle \vec{v} \rangle (0)$ for 20 different synthetic mirror surfaces. These 20 values are eventually averaged and represented by the orange squares in Figure 8. The contributions of the ideal case of a gaussian beam reflected on a perfectly flat mirror are not subtracted, as they are two orders of magnitude smaller. Notably, the standard deviation related to these initial transverse conditions fluctuations decreases linearly with the baseline length $L$. Adjusting these points to the function $\sigma_{\text{eff}}[\sigma_{x_0,y_0}, \sigma_{v_{x,0},v_{y,0}}]/(LT^2)$, represented by the orange dashed line in Figure 8, gives an effective height fluctuation $\sigma_{\text{eff}}[200 \ \mu\text{m}, 1 \ \text{mm s}^{-1}] = 0.1700(3) \ \text{nm}$. Similarly, simulations with fluctuations in position $\sigma_{x_0,y_0} = 0.1 \ \text{mm}$ and in velocity $\sigma_{v_{x,0},v_{y,0}} = 0.5 \ \text{mm s}^{-1}$ were performed, although in this case, 20 different mirror surfaces were used for each baseline length. The corresponding results, represented by green triangles in Figure 8, still show behaviour inversely proportional to the baseline length, and the fit gives an effective height fluctuation $\sigma_{\text{eff}}[100 \ \mu\text{m}, 0.5 \ \text{mm s}^{-1}] = 0.091(2) \ \text{nm}$. This decrease in the mean standard deviation is mainly due to the decrease in the initial velocity distribution, as simulations with the parameters $\sigma_{x_0,y_0} = 0.1 \ \text{mm}$ and $\sigma_{v_{x,0},v_{y,0}} = 1 \ \text{mm s}^{-1}$, not shown in Figure 8, give an effective height fluctuation $\sigma_{\text{eff}}[100 \ \mu\text{m}, 1 \ \text{mm s}^{-1}] = 0.161(3) \ \text{nm}$.

Thus, in this 10-meter-long configuration and with a mirror whose PSD is defined in Figure 7 (d), it is possible to achieve a gravity gradient accuracy of a few 0.1 E, which corresponds to the most stringent requirement in [15], given that the fluctuations in the initial transverse conditions are less than the realistic values considered here [18, 34].

## 9. Conclusion

We have demonstrated that the rejection of wavefront aberrations in an atomic gradiometer is finite, and that the residual bias on the measurement of the vertical gravity gradient is within the reach of current experiments. Besides, in an interferometric configuration where an atomic cloud expands up to the size of the laser beam, contrast non-uniformity plays a crucial role in the phase bias of the interferometer. This applies to both the gaussian beam curvature and the aberrations caused by mirror surface defects. Within this limitation, the analytical framework developed in [28] has been extended to allow rapid estimation of the induced biases and can be adapted to atomic interferometers of different geometries. In addition, simulations with fluctuations in the initial transverse conditions of the atomic clouds allowed us to estimate the sensitivity of the measurement to the initial kinematic parameters of the atomic sources. Eventually, in a state-of-the-art configuration, with a high quality mirror such as that designed for optical gravitational wave detectors [39, 41], we show that wavefront aberrations could cause a bias of the order of 0.1 E on the measurement of the vertical gravity gradient, which corresponds to the accuracy requirement targeted in some experiments [15].

Since we focused on the effects of wavefront aberrations, we did not take into account the Coriolis effect, which also contributes to the interferometer phase [10, 16, 31], particularly when considering fluctuations in the initial transverse velocity. Moreover, for simplicity, the impact of detection on the interferometer phase was omitted. In particular, in the expression (1), the term proportional to $\propto T_{zz}^2 T_{\text{det}}$, where $T_{\text{det}}$ is the delay between the last pulse and detection, has been neglected [30]. More importantly, the finite size of the detection area may also play a role, when certain classes of atoms are not detected or have lower weights [17, 32]. Furthermore, our simulations are based on calculating laser beam distortions at the average longitudinal position of the atoms at each pulse. Since the typical separation distance along the optical axis is $\Delta z = v_{\text{rec}}T$, aberrations with a comparable typical propagation length correspond to a typical transverse length $\sqrt{\Delta z/k_{\text{eff}}}$, which corresponds to 20 $\mu$m in the case of the 10-meter-long gradiometer. This value is smaller than the size of the pixels used to represent the mirror surfaces and is an order
of magnitude smaller than the transverse size of the atomic clouds considered, so their effect, once averaged, should be negligible. However, for initially smaller clouds, e.g. those prepared in a dipolar trap, the corrections could be significant. Moreover, although Zernike polynomials offer the advantage of analytical approximations and enable us to understand the role of key experimental parameters, characterizing optics surfaces at length scales well below the radius on which these polynomials are defined is a tedious task. An alternative is to use the PSD to characterize the mirror surface. In this case, and within the limits of small defects, an analytical approach based on Sinusoidal-Gaussian beams [42] for the propagation of aberrations might be developed. This would be particularly useful for defining specifications of the optics used in experiments, especially since the PSD can be reconstructed from different measurements at different length scales [39, 41, 43].

Although the primary objective of this work is to improve the characterization of our gradiometry experiment [9], it can be adapted to other experiments with different geometries, such as a four-pulse cold atom gyroscope [44]. In this case, even if the atomic cloud is ideally located at the same mean position in the laser beams during the ascending and descending pulses, and if wavefront aberrations are rejected at the first order, residual effects are expected due to the expansion of the cloud during the interferometer sequence. Additionally, other atomic interferometers may have similar configuration to the one considered here, i.e., where due to expansion, the size of the atomic cloud reaches the size of the laser beam. This is all the more true, on the one hand, the power of the laser is limited, particularly in the context of space missions [45], and its size must therefore be restricted in order to maintain sufficient intensity and, on the other hand, in order to improve the sensitivity of the experiment, a larger interferometer surface area and therefore a longer evolution time are sought. Finally, since laser beam distortions are today one of the main factors limiting the accuracy of existing experiments [10, 17, 23], it is important to be able to model and measure their impact in order to improve the error budget of current and future experiments [14, 45].

## Appendix: Calculations with contrast weights

The contrast of the interferometer for a single atom is assumed to be of the form

$$c = f \cdot \sin \left( Ae^{-\frac{2 \rho^2(t_3)}{w_0^2}} \right) \qquad (15)$$

$A$ is the area of the third and last pulse of the interferometer, typically $\pi/2$, and $f$ is a function that does not depend on the transverse coordinates of the atom. The average contrast [46] of the interferometer is then obtained by integrating this expression over the atomic cloud distribution

$$C = \frac{1}{4\pi^2\sigma_{\rho}^2(0)\sigma_v^2} \int e^{-\frac{x^2+y^2}{2\sigma_{\rho}^2(0)}} e^{-\frac{v_x^2+v_y^2}{2\sigma_v^2}} \cdot c \cdot dx dy dv_x dv_y \cdot \eta(z, v_z) dz dv_z \\ = C_z \cdot \sum_{n=0}^{+\infty} \frac{(-1)^n A^{2n+1}}{(2n+1)!} \cdot \frac{1}{1 + (2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2} \equiv C_z \cdot C_{\perp}. \qquad (16)$$

$C_z$ is the integration of the function $f$ over the longitudinal distributions. The phase biases are calculated using expression (7), assuming that the wavefront aberrations depend only on the transverse positions of the atoms, which is valid as long as the size of the cloud along the optical axis is much smaller than the typical propagation length $l_z$ of the aberrations [28].
The contribution of the curvature of a Gaussian beam wavefront $\delta\phi[\rho(t)] = \beta\rho^2(t)$ is

$$ 
\begin{aligned}
\delta\phi &= \frac{\beta C_z}{4\pi^2\sigma_{\rho}^2(0)\sigma_v^2 C} \int e^{-\frac{x^2+y^2}{2\sigma_{\rho}^2(0)}} e^{-\frac{y^2+v^2}{2\sigma_v^2}} \cdot \rho^2(t) \sin\left(\mathcal{A}e^{-\frac{2\rho^2(t_3)}{w_0^2}}\right) \cdot dx dy dv_x dv_y \\
&= \frac{2\beta\sigma_{\rho}^2(t)}{C_{\perp}} \sum_{n=0}^{+\infty} \frac{(-1)^n \mathcal{A}^{2n+1}}{(2n+1)!} \cdot \frac{\left(1+(2n+1)\left(\frac{2\sigma_{\rho}(0)\sigma_v(t_3-t)}{w_0\sigma_{\rho}(t)}\right)^2\right)}{\left[1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2\right]^2}.
\end{aligned}
\quad (17) $$

To the extent that the transverse selection is negligible $\rho(t_3) \ll w_0$, the expression converges toward $2\beta\sigma_{\rho}^2(t)$, which is obtained with uniform contrast [20].

The contribution of a wavefront aberration described by a Zernike polynomial, as in [28], assuming that the Zernike polynomial defined on a disk of radius $R$ can be replaced by a Bessel function: $\delta\phi[\rho(t)] = \beta Z_n^0(\frac{\rho(t)}{R}) \approx \beta J_0(\frac{\rho(t)}{l_{xy}})$ with $l_{xy} = R/(n+1)$, is

$$ 
\begin{aligned}
\delta\phi &= \frac{C_z}{4\pi^2\sigma_{\rho}^2(0)\sigma_v^2 C} \int e^{-\frac{x^2+y^2}{2\sigma_{\rho}^2(0)}} e^{-\frac{y^2+v^2}{2\sigma_v^2}} \cdot \beta J_0\left(\frac{\rho(t)}{l_{xy}}\right) \sin\left(\mathcal{A}e^{-\frac{2\rho^2(t_3)}{w_0^2}}\right) \cdot dx dy dv_x dv_y \\
&= \frac{\beta}{C_{\perp}} \sum_{n=0}^{+\infty} \frac{(-1)^n \mathcal{A}^{2n+1}}{(2n+1)!} \cdot \frac{\exp\left(-\frac{1}{2}\left(\frac{\sigma_{\rho}(t)}{l_{xy}}\right)^2 \cdot \frac{1+(2n+1)\left(\frac{2\sigma_{\rho}(0)\sigma_v(t_3-t)}{w_0\sigma_{\rho}(t)}\right)^2}{1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2}\right)}{1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2} \equiv \beta S(t).
\end{aligned}
\quad (18) $$

Within the uniform contrast limit $\rho(t_3) \ll w_0$, the expression converges to $\beta e^{-\frac{1}{2}\left(\frac{\sigma_{\rho}(t)}{l_{xy}}\right)^2}$, which is the one obtained in [28]. $\beta$ contains the dependence on the amplitude of the aberration and on the ratio between the position along the optical axis and the typical propagation length of the Zernike aberration $l_z = 2kl_{xy}^2$.

In the two previous expressions (17) and (18), the square of the cloud size is multiplied by the same correcting term

$$ 
\begin{aligned}
\sigma_{\rho}^2(t_i) \cdot \frac{1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2 \cdot \left(\frac{\sigma_{\rho}(0)\sigma_v(t_3-t_i)}{\sigma_{\rho}(t_i)\sigma_{\rho}(t_3)}\right)^2}{1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2} \\
= \sigma_{\rho}^2(t_i) \cdot \left(1 - \frac{(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2 \cdot \left(\frac{\sigma_{\rho}^2(0)+\sigma_v^2 t_i t_3}{\sigma_{\rho}(t_i)\sigma_{\rho}(t_3)}\right)^2}{1+(2n+1)\left(\frac{2\sigma_{\rho}(t_3)}{w_0}\right)^2}\right).
\end{aligned}
\quad (19) $$

This effective reduction of the cloud transverse size corresponds to the fact that atoms with higher transverse velocities participate less to the interference process since they are on the edge of the laser beam during the final recombination pulse, where the atom-light coupling is lower. In particular, if the initial size of the cloud can be neglected with respect to the thermal expansion of the cloud, expression (19) can be simplified as

$$ \frac{(\sigma_v t_i)^2}{1 + (2n + 1) \left( \frac{2\sigma_v t_3}{w_0} \right)^2}. \quad (20) $$
For example, considering the contribution of defocus $Z_2^0$ corresponding to $l_{xy} = R/3$, the exponential in expression (18) can be linearized supposing $\sigma_\rho \ll R$. For low temperature ($2\sigma_v t_3 \ll w_0$) the contribution is linear with the temperature, and for higher temperature ($2\sigma_v t_3 \approx w_0$) the effect is reduced as seen in Figure 3 of [17]. The analytical formulas are derived supposing that contrast losses due to the cloud transverse size at the first and second pulses are negligible.

**Funding.** The authors acknowledge the support from a government grant managed by the Agence Nationale de la Recherche under the Plan France 2030 with the reference "ANR-22-PETQ-0005" (project QAFCA).

**Acknowledgment.** The authors thank Arnaud Landragin for fruitful discussions. Optical elements in Figure 1 were taken from ComponentLibrary by Alexander Franzen, licensed under a Creative Commons Attribution-NonCommercial 3.0 Unported License.

**Disclosures.** The authors declare no conflicts of interest.

**Data availability.** Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

## References

1. F. Riehle, T. Kisters, A. Witte, et al., "Optical ramsey spectroscopy in a rotating frame: Sagnac effect in a matter-wave interferometer," Phys. Rev. Lett. 67, 177-180 (1991).

2. M. Kasevich and S. Chu, "Atomic interferometry using stimulated raman transitions," Phys. Rev. Lett. 67, 181-184 (1991).

3. J. J. Hudson, D. M. Kara, I. J. Smallman, et al., "Improved measurement of the shape of the electron," Nature 473, 493-496 (2011).

4. P. Cladé, F. Nez, F. Biraben, and S. Guellati-Khelifa, "State of the art in the determination of the fine-structure constant and the ratio h/mu," Comptes Rendus Physique 20, 77-91 (2019). The new International System of Units / Le nouveau Système international d'unités.

5. L. Antoni-Micollier, D. Carbone, V. Ménoret, et al., "Detecting volcano-related underground mass changes with a quantum gravimeter," Geophys. Res. Lett. 49, e2022GL097814 (2022).

6. Q. d'Armagnac de Castanet, C. Des Cognets, R. Arguel, et al., "Atom interferometry at arbitrary orientations and rotation rates," Nat. Commun. 15, 6406 (2024).

7. A. Louchet-Chauvet, T. Farah, Q. Bodart, et al., "The influence of transverse motion within an atomic gravimeter," New J. Phys. 13, 065025 (2011).

8. M. J. Snadden, J. M. McGuirk, P. Bouyer, et al., "Measurement of the earth's gravity gradient with an atom interferometer-based gravity gradiometer," Phys. Rev. Lett. 81, 971-974 (1998).

9. R. Caldani, K. X. Weng, S. Merlet, and F. Pereira Dos Santos, "Simultaneous accurate determination of both gravity and its vertical gradient," Phys. Rev. A 99, 033601 (2019).

10. W. Lyu, J.-Q. Zhong, X.-W. Zhang, et al., "Compact high-resolution absolute-gravity gradiometer based on atom interferometers," Phys. Rev. Appl. 18, 054091 (2022).

11. C. Janvier, V. Ménoret, B. Desruelle, et al., "Compact differential gravimeter at the quantum projection-noise limit," Phys. Rev. A 105, 022801 (2022).

12. J. B. Fixler, G. T. Foster, J. M. McGuirk, and M. A. Kasevich, "Atom interferometer measurement of the newtonian constant of gravity," Science 315, 74-77 (2007).

13. G. Rosi, F. Sorrentino, L. Cacciapuoti, et al., "Precision measurement of the newtonian gravitational constant using cold atoms," Nature 510, 518-521 (2014).

14. A. Balaz, D. Blas, O. Buchmueller, et al., "Long-baseline atom interferometry," (2025).

15. J. Hartwig, S. Abend, C. Schubert, et al., "Testing the universality of free fall with rubidium and ytterbium in a very large baseline atom interferometer," New J. Phys. 17, 035011 (2015).

16. S.-Y. Lan, P.-C. Kuan, B. Estey, et al., "Influence of the coriolis force in atom interferometry," Phys. Rev. Lett. 108, 090402 (2012).

17. R. Karcher, A. Imanaliev, S. Merlet, and F. P. D. Santos, "Improving the accuracy of atom interferometers with ultracold sources," New J. Phys. 20, 113041 (2018).

18. Q. Luo, L. Li, H. Zhou, et al., "Evaluating the effect of wavefront aberrations in an atom gravimeter*," Metrologia 62, 055006 (2025).

19. S. Bade, L. Djadaojee, M. Andia, et al., "Observation of extra photon recoil in a distorted optical field," Phys. Rev. Lett. 121, 073603 (2018).

20. J. M. Cervantes and E. Gomez, "Effect of an aperture in atomic gravimetry," J. Opt. Soc. Am. A 41, 881-891 (2024).

21. S. Seckmeyer, H. Ahlers, J.-N. Kirsten-Siemß, et al., "Spatially resolved phase reconstruction for atom interferometry," EPJ Quantum Technol. 12, 34 (2025).

22. J. Junca, J. Kitching, and W. McGehee, "Wavefront mapping for absolute atom interferometry," (2025).
23. S. Gaudout, R. Si-Ahmed, C. Debavelaere, et al., "Probing the spatial distribution of k-vectors in situ with bose-einstein condensates," arXiv (2025). Version arXiv: v3.

24. V. Schkolnik, B. Leykauf, M. Hauth, et al., "The effect of wavefront aberrations in atom interferometry," Appl. Phys. B **120**, 311–316 (2015).

25. M.-K. Zhou, Q. Luo, L.-l. Chen, et al., "Observing the effect of wave-front aberrations in an atom interferometer by modulating the diameter of raman beams," Phys. Rev. A **93**, 043610 (2016).

26. A. Trimeche, M. Langlois, S. Merlet, and F. Pereira Dos Santos, "Active Control of Laser Wavefronts in Atom Interferometers," Phys. Rev. Appl. **7**, 034016 (2017).

27. W.-J. Xu, J. Liu, Q. Luo, et al., "In situ measurement of the wavefront phase shift in an atom interferometer," Phys. Rev. Appl. **22**, 054014 (2024).

28. L. Pagot, S. Merlet, and F. P. D. Santos, "Influence of optical aberrations on the accuracy of an atomic gravimeter," Opt. Express **33**, 18843–18854 (2025).

29. P. Asenbaum, C. Overstreet, T. Kovachy, et al., "Phase shift in an atom interferometer due to spacetime curvature across its wave function," Phys. Rev. Lett. **118**, 183602 (2017).

30. A. Peters, K. Y. Chung, and S. Chu, "High-precision gravity measurements using atom interferometry," Metrologia **38**, 25 (2001).

31. J. M. Hogan, D. M. S. Johnson, and M. A. Kasevich, "Light-pulse atom interferometry," (2008).

32. P. Gillot, B. Cheng, S. Merlet, and F. Pereira Dos Santos, "Limits to the symmetry of a Mach-Zehnder-type atom interferometer," Phys. Rev. A **93**, 013609 (2016).

33. The series comes from the development of the sinus in the contrast expression. The truncation of the sine series at $n = 5$ is an 11th-order polynomial, and its difference with the sine function on $[0, \pi/2]$ is less than $10^{-7}$.

34. A. Gauguet, B. Canuel, T. Lévèque, et al., "Characterization and limits of a cold-atom sagnac interferometer," Phys. Rev. A **80**, 063604 (2009).

35. The variance estimator of $(X_i)_{1\le i\le n}$ is chosen to be $\overline{V}_n = \frac{1}{n}\sum_i (X_i - \overline{X}_n)^2$, with $\overline{X}_n = \frac{1}{n}\sum_i X_i$ the mean estimator. Statistical uncertainty is defined as the square root of the mean squared error (MSE) of the estimator. For the mean estimator: MSE$(\overline{X}_n) = \sigma^2(X)/n$. For the variance estimator: MSE$(\overline{V}_n) = (\mu^{(4)}(X) - \sigma(X)^4)/n + O(n^{-2})$, with $\mu^{(4)}(X) = \mathbb{E}[(X - \mathbb{E}[X])^4]$ is the 4th centered moment [47]. The standard deviation estimator is defined as $\overline{\sigma}_n = \sqrt{\overline{V}_n}$, with the associated uncertainty $\epsilon(\overline{\sigma}_n) = \epsilon(\overline{V}_n)/(2\overline{\sigma}_n)$.

36. The Zernike polynomials $Z_n^{|m|}$ and $Z_n^{-|m|}$ differ only by a rotation. Thus, simulations are performed for only one of the two polynomials.

37. P. F. Forman, "The Zygo Interferometer System," in *Interferometry*, vol. 0192 G. W. Hopkins, ed., International Society for Optics and Photonics (SPIE, 1979), pp. 41–49.

38. The 1D PSD is defined as $\text{PSD}_{1D} = f \int_0^{2\pi} d\theta \text{PSD}_{2D}(f, \theta)$ with $f$ the spatial frequency modulus, such as $\sigma_{\text{RMS}}^2 = \int_0^{+\infty} df \text{PSD}_{1D}$. A radial Hann window is used for calculations [43].

39. E. Hirose, D. Bajuk, G. Billingsley, et al., "Sapphire mirror for the kagra gravitational wave detector," Phys. Rev. D **89**, 062003 (2014).

40. The mirror surfaces were decomposed over the first 2500 Zernike polynomials in the ISO-14999 indexing scheme. Analytical estimates of the biases were calculated using Zernike polynomials $Z_n^m$ with azimuthal indices $|m| \in \{0, 1, 2\}$, as adding those with $|m| = 3$ results in changes smaller than the error indicated in parentheses in Table 1.

41. Y. Drori, J. Eichholz, T. Edo, et al., "Scattering loss in precision metrology due to mirror roughness," J. Opt. Soc. Am. A **39**, 969–978 (2022).

42. L. W. Casperson, D. G. Hall, and A. A. Tovar, "Sinusoidal-gaussian beams in complex optical systems," J. Opt. Soc. Am. A **14**, 3341–3348 (1997).

43. T. D. B. Jacobs, T. Junge, and L. Pastewka, "Quantitative characterization of surface topography using spectral analysis," Surf. Topogr. Metrol. Prop. **5**, 013001 (2017).

44. R. Gautier, M. Guessoum, L. A. Sidorenkov, et al., "Accurate measurement of the sagnac effect for matter waves," Sci. Adv. **8**, eabn8009 (2022).

45. T. Lévèque, C. Fallet, J. Lefebvre, et al., "Carioqa: Definition of a quantum pathfinder mission," (2022). ArXiv 2211.01215.

46. The series in equation (16) is a hypergeometric function $_1F_2$.

47. J. Garnier, S. Méléard, and N. Touzi, *Aléatoire* (Les Éditions de l'École polytechnique, Palaiseau, 2019).