1 Introduction

Gravitational properties of antimatter have never been measured directly. A promising experimental method to do so consists in producing sufficiently cold antihydrogen atoms (\(\overline{\mathrm {H}}\)) and timing their free fall in the Earth’s gravity field. This approach is being pursued by AEGIS [1], ATHENA-ALPHA [2], ATRAP [3] and GBAR [4] collaborations.

In order to get the highest accuracy for measuring the free-fall acceleration \(\overline{g}\) of \(\overline{\mathrm {H}}\), one has to cool atoms down to low temperatures and to measure, or at least to deduce from design and calculations, the initial velocity distribution. We discuss here the method proposed by Walz and Hänsch [5], which is used in the GBAR project to reach very low temperatures: \({\overline{\mathrm {H}}}^+\) ions are trapped and cooled down to the lowest quantum state in a Paul trap and \(\overline{\mathrm {H}}\) is then produced by photo-detaching the excess positron. The photo-detachment pulse is the START signal for the free-fall timing measurement, while the STOP signal is provided by the annihilation of \(\overline{\mathrm {H}}\) atoms on a detection plate placed at a height \(H\) below the center of the ion trap.

The precision of this measurement depends crucially on the dispersion of vertical velocities before the free fall, which corresponds to the residual kinetic energy of the atoms after the cooling process. The aim of the present paper is to propose a new filtering method to further reduce the initial distribution of the vertical velocities and thus improve the accuracy in the measurement of \(\overline{g}\).

In Sect. 2 we justify our choice of characteristic values for the spatial localization of the initial atomic cloud by considering the spreading of the freely falling wavepacket of \(\overline{\mathrm {H}}\) in the gravitational field. We describe in Sect. 3 the new method for shaping vertical velocities of \(\overline{\mathrm {H}}\) in the quasi-classical approximation, and we show in Sect. 4 that the improvement of accuracy due to the velocity selection overcomes the degradation associated with the decrease of the statistics. We then present in Sect. 5 a quantum-mechanical description of the experiment in order to validate the quasi-classical estimations of the preceding sections. In Sect. 6 we list possible systematic effects and show that they scale down compared to those in the case of unrestricted free fall of \(\overline{\mathrm {H}}\). We then deduce the accuracy which could be reached on the measurement of \(\overline{g}\) under realistic experimental conditions.

We neglect throughout this paper systematic effects related to the energy-dependent probability of quantum reflection of \(\overline{\mathrm {H}}\) from the detection plate [6]. The atomic recoil in the photo-detachment process induces an additional velocity dispersion which is discussed in the last section on systematic effects.

2 Spreading of a freely falling wavepacket

In the simplified description presented in the introduction, the initial distribution at time \(t=0\) is the lowest quantum state in the Paul trap. This corresponds to a Gaussian wavepacket with vertical velocity dispersion \(\upsilon \) and vertical position dispersion \(\zeta \) reaching the minimum in the Heisenberg uncertainty relation:

$$\begin{aligned} m\upsilon \zeta =\frac{\hbar }{2} \end{aligned}$$
(2.1)

where \(\hbar \) is the reduced Planck constant and \(m\) the inertial mass of \(\overline{\mathrm {H}}\).

After their release from the trap at time \(t=0\), atoms start falling freely in the Earth’s gravity field until they reach the detection plate placed at a height \(H\) below the center of the trap. The time of fall is measured as the delay \(t\) from their release to their annihilation on the detection plate. The acceleration of gravity \(\overline{g}\) for antihydrogen is then deduced from the distribution of these fall times. This acceleration \(\overline{g}\) for antihydrogen is related to the analog quantity \(g\) defined for hydrogen by \(\overline{g}=Mg/m\), where \(M\) is the gravitational mass of \(\overline{\mathrm {H}}\).

We now discuss the distribution of free-fall times, assuming for simplicity that this distribution is determined by initial dispersions of the vertical velocity and position (other sources of uncertainty negligible). If the initial quantum state is poorly localized (large values of \(\zeta \)), then the spread of the fall times is too large, because of the initial dispersion of height. In the opposite case where the wavepacket is too localized (small values of \(\zeta \)), the spread of the fall times is too large, because of the initial dispersion of the vertical velocity. An optimum localization of the initial quantum state should be found as a compromise between these two limiting cases.

As the variations of position and velocity are uncorrelated in the initial wavepacket, a classical calculation gives the relative spread \(\Delta t\) of the free-fall times arising from both effects:

$$\begin{aligned} \frac{\Delta t}{t_H}&= \sqrt{{\left( \frac{\zeta }{2H}\right) }^{2}+{\left( \frac{\upsilon }{v_H}\right) }^2} \end{aligned}$$
(2.2)
$$\begin{aligned}&= \sqrt{{\left( \frac{\zeta }{2H}\right) }^{2}+{\left( \frac{\hbar }{2m v_H\zeta }\right) }^{2}}. \end{aligned}$$
(2.3)

The second of these relations uses (2.1) while the first one is valid even when \(\upsilon \) and \(\zeta \) do not reach the minimum in Heisenberg uncertainty relation. We have defined \(t_H=\sqrt{2H/\overline{g}}\) and \(v_H=\sqrt{2\overline{g}H}\) as the free-fall time and velocity for a free-fall height \(H\) with zero initial velocity. The optimum size of the initial state, which minimizes \(\Delta t\) in (2.3), is

$$\begin{aligned} {\zeta }_{\mathrm{opt}}=\sqrt{\frac{{\hbar H}}{mv_H}}. \end{aligned}$$
(2.4)

It leads to an optimum resolution for the free-fall measurement:

$$\begin{aligned} {\left( \frac{\Delta t}{t_H}\right) }_{\mathrm{opt}}=\sqrt{\frac{\hbar }{2mv_HH}}. \end{aligned}$$
(2.5)

The larger the product \(mv_HH\) with respect to \(\hbar /2\), the better this optimal resolution is.

Better precisions are also obtained by increasing the fall height with the characteristics of the trap kept fixed. However, the setup size is limited by practical arguments involving price and space considerations. Note that Eq. (2.5) is translated in an uncertainty twice larger on the acceleration of gravity,

$$\begin{aligned} \frac{\Delta \overline{g}}{\overline{g}}=2\frac{\Delta t}{t} \end{aligned}$$
(2.6)

in the simple derivation presented here (a detailed analysis based on Monte-Carlo simulations is given in [4]).

With the typical numbers used for the design [4] of the GBAR experiment (\(H=0.3\) m so that \(v_H\approx 2.4\) m/s if \(\overline{g}\approx g\)), one obtains \(\zeta _\mathrm{opt}\approx 88\) \(\upmu \)m and \(\left( \Delta t/t_H\right) _\mathrm{opt} \approx 2.1{\times }{10}^{-4}\). If this optimum operation could be experimentally realized, the accuracy would reach \(\left( \Delta \overline{g}/\overline{g}\right) _\mathrm{opt} \approx 4.2{\times }{10}^{-4}\) for each detection of an annihilation event. With a total number of events \(N_{\mathrm {tot}} \approx 2.6{\times }{10}^{4}\), calculated for a typical measuring time of 1 month and an average production rate of 1 ultracold \(\overline{\mathrm {H}}\) atom per period of 100 s, this would lead to the resolution after one month:

$$\begin{aligned} \left( \frac{\Delta \overline{g}}{\overline{g}\sqrt{N_{\mathrm {tot}}}} \right) _\mathrm{opt}= \sqrt{\frac{2\hbar }{m v_H H N_{\mathrm {tot}}}} \approx 2.6 \times 10^{-6}. \end{aligned}$$
(2.7)

We have assumed there were no large systematic effects.

However, the size of the initial cloud used in the design of the GBAR experiment is far from this optimum. The Paul trap is characterized by its oscillation frequency \(\omega \), which fixes the velocity and position dispersions in the ground state:

$$\begin{aligned} \zeta =\sqrt{\hbar /2m\omega },\quad \upsilon =\sqrt{\hbar \omega /2m}. \end{aligned}$$
(2.8)

The mean kinetic energy in the ground state is then \(m \upsilon ^2/2=\hbar \omega /4\). Therefore the range of trap frequencies that can be used is limited by the residual kinetic energy of the atoms after cooling.

In GBAR, the considered frequency range is 0.1 MHz \(< \omega /2\pi <\) 1 MHz, so that one gets 0.22 \(\upmu \)m \(>\zeta >\) 0.07 \(\upmu \)m and 0.14 m/s \(<\upsilon <\) 0.44 m/s. This means that the initial cloud is smaller than the optimum by about three orders of magnitude. The resolution is thus limited by the dispersion of initial velocity:

$$\begin{aligned} \frac{\Delta \overline{g}}{\overline{g}\sqrt{N_{\mathrm {tot}}}}\approx \frac{2\upsilon }{v_H\sqrt{N_{\mathrm {tot}}}}. \end{aligned}$$
(2.9)

As it is not experimentally feasible to further cool down the ions to reach the optimum size of the initial cloud, we propose in this paper to select the initial vertical velocity of the atoms. This will improve the resolution after each annihilation event by a factor scaling as the reduced velocity range \(\Delta v /\upsilon \). The statistics is reduced by a factor scaling as \(\sqrt{N/N_{\mathrm {tot}}}{\propto } \sqrt{\Delta v /\upsilon }\) [see Eq. (3.1)] so that an overall improvement is expected. Also systematic uncertainties will decrease dramatically. The description of the shaping device and the evaluation of its performance are discussed in more detail in the next sections.

3 Shaping the distribution of vertical velocities of \(\overline{\mathrm {H}}\) in GBAR

The current design for GBAR is a classical free-fall experiment which aims at an accuracy of the order of 1 % [4]. With a quantum detection technique, one could get a significantly higher precision, in analogy to spectroscopy [7] or interferometry [8] of near-surface quantum states [9] of ultracold neutrons (UCNs) [10, 11]. However, these techniques require high energy resolution and sufficient statistics [12, 13]. The method that we propose in this paper is an intermediate step in this direction which is less precise than the full quantum detection technique but allows for better statistics and simpler design.

This method is analogous to the one used in the experiment on the observation of gravitational quantum states of ultracold neutrons [7, 1416]. The distribution of initial vertical velocities is shaped by selecting the atoms passing through a shaping device consisting of two disks. A scheme of principle of the shaping device where all useful quantities are defined is shown in Fig. 1. In the sequel of this section, a simple analysis of the problem is presented in terms of quasi-classical arguments, to be confirmed in the next sections. A more complete quantum-mechanical description is also available in papers devoted to ultracold neutrons [15, 1720].

Fig. 1
figure 1

A scheme of principle of the proposed shaping device: an \(\overline{\mathrm {H}}\) atom is released from the Paul trap (central spot) and it bounces a few times on the mirror surface of the bottom disk (arrows); if it scatters on the rough top surface, it annihilates (lightnings); otherwise, it escapes from the aperture between the two disks, and it falls to the detection plate where it annihilates (lightning on the detection plate). \(R\) is the radius of the bottom and top disks, \(r\) is the radius of central openings in the disks, \(h\) is the distance between the top surface of the bottom disk and the bottom surface of the top disk, \(H\) is the distance between the top of the detection plate and the top of the bottom disk, \(L\) is the horizontal distance between the initial spot and the detection point

In the zone between the two disks, atoms with sufficiently small vertical velocities bounce on the bottom mirror disk due to the high efficiency of quantum reflection in the Casimir–Polder potential [6]. If the top surface of the mirror disk is flat, smooth and horizontal, the horizontal velocity component as well as the total energy of the vertical motion do not change and atoms thus pass through the shaping device with high probability. This last statement would be precisely valid for ideal quantum reflection from the mirror surface; otherwise corresponding corrections have to be taken into account (more discussion below). On the other hand, atoms with large vertical velocities rise in the Earth’s gravity field to the height of the rough surface of the top disk and scatter non-specularly on this surface. As this scattering mixes horizontal and vertical velocity components, it leads to rapid loss of scattered \(\overline{\mathrm {H}}\) through annihilation on the top or bottom disk.

A few remarks are useful at this point: (1) the shaping device has to be coupled with the Paul trap (not shown on the figure); this point is not discussed in this paper except for the role of the openings left in the center of the disks for operating the Paul trap; note that the disks may consist of several sectors not covering the complete \(2\pi \) horizontal angle in order to include the Paul trap in the overall design; (2) annihilation events are assumed to be detected with position-sensitive and time-resolving detectors; this will allow one to account for the time spent in the shaping device (see below); (3) due to the cylindrical symmetry of the device, all atoms with small enough vertical velocity components and any value and direction of the horizontal velocity component can pass through it with high probability.

In order to describe the operation of the shaping device, we follow the possible classical trajectories of atoms from the initial point where they are released to the points where they annihilate. As the size of the initial spot (discussed in Sect. 2) is much smaller than any other characteristic size of the shaping device, it plays no role in the following. We assume the initial spot of atoms to be placed at the height \(H\) of the top surface of the bottom disk (origin for altitude placed at the detection plate). In a first step, we let the radius \(r\) of the central opening tend to zero and the radius \(R\) of the disk tend to infinity. Disregarding the losses due to imperfect quantum reflection on the mirror disk, we obtain the fraction of atoms going through the angular acceptance of the shaping device as

$$\begin{aligned} \frac{N}{N_{\mathrm {tot}}}\approx \frac{\Delta v}{\upsilon }\sqrt{\frac{1}{2\pi }} \end{aligned}$$
(3.1)

where \(\upsilon \) is the standard deviation of the Gaussian distribution of the vertical velocities and \(\Delta v\) the range of vertical velocities fitting the aperture of the shaping device. With the geometry sketched in Fig. 1, the latter corresponds to atoms with vertical velocities \(0<v<\Delta v\) with \(\Delta v\) deduced from the energy needed to rise the height from \(H\) to \(H+h\) in the gravity field:

$$\begin{aligned} \Delta v=\sqrt{2\overline{g}h}. \end{aligned}$$
(3.2)

Note that the fraction of atoms going through the angular acceptance of the shaping device changes as a function of the height of the initial spot above the mirror disk and as a function of the radius \(r\); therefore Eq. (3.1) has to be modified for other positions of the spot. Also Eq. (3.1) has been written in the limit of a good velocity selection \(\Delta v < \upsilon \), which entails through (3.2) that \(h\) has a maximum value \(h_{\max }\):

$$\begin{aligned} h<h_{\max }=\frac{{\upsilon }^{2}}{2\overline{g}}. \end{aligned}$$
(3.3)

With the GBAR numbers considered above, the maximum value \(h_{\max }\) lies in the interval 1–10 mm. If this condition is not obeyed, Eq. (3.1) has to be replaced by the appropriate integral.

We now take into account the finite values of the radii of the central openings \(r\) and of the disks \(R\). In order to do it properly, we have to consider the shape of the angular distribution of initial velocities. The operation of the Paul trap may indeed require anisotropy to be introduced between horizontal and vertical directions. This can be described by the ratio \(\varepsilon \) between frequencies of operation in horizontal and vertical directions \({\omega }_{\mathrm {hor}}={\varepsilon \omega }\) (\(\omega \) is the frequency already introduced for operation of the vertical trap). This ratio should be in the interval \(2<\varepsilon <4\) for a proper operation of the Paul trap [21]. Using the same reasoning as in the preceding section, we deduce that the horizontal dispersions are

$$\begin{aligned} \upsilon _{\mathrm {hor}}=\upsilon \sqrt{\varepsilon },\quad \zeta _{\mathrm {hor}}=\zeta /\sqrt{\varepsilon }, \end{aligned}$$
(3.4)

where \(\upsilon \) and \(\zeta \) are the dispersions already introduced for operation of the vertical trap.

We can now discuss the role of the finite radius \(r\) of the central openings. We want to avoid extra loss of statistics at the entrance of the device, and thus we choose \(r\) small enough so that the angular divergence there fits the angular acceptance of the shaping device:

$$\begin{aligned} \frac{h}{r}>\frac{\Delta v}{\upsilon \sqrt{\varepsilon }}, \quad r<r_{\max } =\frac{\upsilon \sqrt{\varepsilon h}}{\sqrt{2\overline{g}}} =\sqrt{\varepsilon h h_{\max }}. \end{aligned}$$
(3.5)

To write these relations, we have neglected the effect of gravity on the short distance \(r\) and used the value in (3.4) of the root-mean-square (rms) dispersion of the horizontal velocity.

We then consider the role of the finite radius \(R\) of the disk, using the following classical arguments. We want to produce an efficient loss of atoms having too large velocities with respect to the designed angular acceptance of the shaping velocity, and thus we choose \(R\) large enough so those atoms efficiently touch the top disk. Saying that they touch it at least once, this implies that the time \(T\) they spend in the zone between the disks is about two times larger than the time \(t_h=\sqrt{2h/\overline{g}}\) corresponding to a free fall on a height \(h\):

$$\begin{aligned}&T=\frac{R}{\upsilon \sqrt{\varepsilon }}>2t_h = 2\sqrt{\frac{2h}{\overline{g}}}\nonumber \\&\quad \Rightarrow \quad R>{R}_{\min }= \frac{4\upsilon \sqrt{{\varepsilon h}}}{\sqrt{2\overline{g}}}=4r_{\max }. \end{aligned}$$
(3.6)

Again, we have used the dispersion (3.4) of the horizontal velocity to calculate the time \(T\) spent in the shaping device for an atom with the rms velocity. Of course, the time \(T\) depends on the actual horizontal velocity (not its rms value) so that a value larger than that calculated in (3.6) is required to produce an effective shaping for the whole distribution.

We want also to stress that the time \(T\) appears as a systematical delay in the free-fall timing experiment so that its knowledge is crucial for accuracy. Here the fact that annihilation event detectors are position sensitive is important. Measuring the horizontal distance \(L\) between the initial spot and the detection point indeed gives the actual horizontal velocity of the atom \(L/T_{\mathrm {tot}}\) with \(T_{\mathrm {tot}}\) the time between escape from the trap and annihilation on the detector and allows one to correct the timing measurement for the time spent in the shaping device \(T=R T_{\mathrm {tot}}/L\).

At the exit of the shaping device, the height lies in the interval \(\left[ H,H+h\right] \), while the vertical velocity lies in the interval \(\left[ -\Delta v,+\Delta v\right] \). As discussed in the next section, this affects the resolution of the timing measurement in the same manner as the dispersion of velocities did affect the free-fall measurement discussed in Sect. 2. In order to optimize the various parameters, in particular the value of the radius \(R\), we have to simulate the whole experiment, that is, the photo-detachment, the passage through the shaping device, the free fall from its output slit to the detection plate, the timing of annihilation events, and the correction from the time spent in the device. In the present paper, we use simpler arguments to estimate the resulting accuracy of the measurement.

4 Estimation of statistical uncertainty

At this point, we have all information needed to give a simple estimation of the statistical accuracy in this experiment. To this aim, we use the analogy with the free-fall timing measurement to write the relative spread of the free-fall times as [compare with (2.2)]:

$$\begin{aligned} \frac{\Delta t}{t_H}=\sqrt{\alpha {\left( \frac{h}{2H}\right) }^{2} +\beta {\left( \frac{\Delta v}{v_H}\right) }^{2}}; \end{aligned}$$
(4.1)

\(\alpha \) and \(\beta \) are dimensionless numbers smaller than unity describing the shapes of position and velocity distributions at the output slit of the shaping device. For simplicity, we have assumed that these distributions are uncorrelated and we have considered that the correction for the time \(T\) spent in the shaper has been done. As \(\Delta v=\sqrt{2\overline{g}h}\) and \(v_H=\sqrt{2\overline{g}H}\) with \(h{\ll }H\), it follows that the relative spread \(\left( \Delta t/t_H\right) \) is dominated by the effect of velocity dispersion and can be written as:

$$\begin{aligned} \frac{\Delta t}{t_H}\approx \sqrt{\frac{{\beta h}}{H}}. \end{aligned}$$
(4.2)

This corresponds to an accuracy \(\left( \Delta \overline{g}/\overline{g}\right) \approx 2\sqrt{{\beta h}/H}\) for each detection of an annihilation event. We then obtain the resolution after one month of measurement, taking into account that the number of events is reduced by the velocity selection [compare with (3.1)]:

$$\begin{aligned} \frac{\Delta \overline{g}}{\overline{g}\sqrt{N}} =2\sqrt{\frac{\beta h}{H}}\sqrt{\frac{\upsilon \sqrt{2\pi }}{N_\mathrm {tot}\Delta v}} =2 \left( \frac{{\pi h}{\beta }^{2}{\upsilon }^{2}}{\overline{g}{H}^{2}N_{\mathrm {tot}}^{2}} \right) ^{1/4}. \end{aligned}$$
(4.3)

It is instructive to compare this resolution with the analogous result obtained without the velocity selection mechanism. The improvement is described by the ratio of (4.3) to (2.9):

$$\begin{aligned} 2 \left( \frac{{\pi h}{\beta }^{2}{\upsilon }^{2}}{\overline{g}{H}^{2}N_{\mathrm {tot}}^{2}} \right) ^{1/4} \left( \frac{2\upsilon }{v_H\sqrt{N_{\mathrm {tot}}}} \right) ^{-1} =\left( \frac{2{\pi h}{\beta }^{2}}{h_{\max }} \right) ^{1/4}. \end{aligned}$$
(4.4)

The best accuracy is therefore achieved for smaller slit sizes. We take the value \(H=0.3\) m chosen in the current design of GBAR, the worst case of \(\beta =1\) and a velocity dispersion \(\upsilon =\) 0.44 m/s and discuss three cases corresponding to decreasing values of \(h\):

  1. 1.

    Equation (4.4) shows that \(h\) should be smaller than \(\approx h_{\max }/2\pi \) for the shaping device to improve the resolution of the experiment. We choose as an example \(h=1\) mm, so that the statistics is \(N\approx 3.3\times 10^3\). The opening radius has to be smaller than \(r_{\max }\simeq 3.2 \sqrt{\varepsilon }\) mm and the disk radius should be larger than \(R_{\min }\simeq 13 \sqrt{\varepsilon }\) mm. The statistical accuracy is then \(\Delta \overline{g}/\overline{g}\sqrt{N} \approx 2.0\times 10^{-3}\). Note that, for a conducting mirror and a maximal vertical velocity \(\sqrt{2 g h}\approx 0.14\) m/s, the reflection probability for an atom is \(78\,\%\) [6]. To simultaneously improve the resolution and reduce losses from annihilation on the bottom mirror, we move to smaller values of \(h\).

  2. 2.

    For \(h<50\) \(\upmu \)m, the atom flux through the slit can no longer be evaluated from classical arguments and the quantum behavior of \(\overline{\mathrm {H}}\) in the slit between the disks has to be taken into account [7, 14, 15]. At the boundary \(h=50\) \(\upmu \)m, the statistics is \(N \approx 7.3\times 10^2\) and the statistical accuracy is \(\Delta \overline{g}/\overline{g}\sqrt{N} \approx 1.0{\times }{10}^{-3}\). The opening radius has to be smaller than \(r_{\max }\simeq 0.7 \sqrt{\varepsilon }\) mm. Note that the reflection probability for an atom with the maximal velocity \(\sqrt{2 g h}=3.1\times 10^{-2}\) m/s is 94 % for a conducting mirror.

  3. 3.

    For \(h< 20\,\upmu \)m, only atoms in the lowest quantum state can pass through the slit. The reflection probability approaches unity in this case which also corresponds to the highest accuracy for the free-fall timing measurement. This quantum limit is analyzed in Sects. 5.2 and 5.3.

The first and second cases provide more comfortable conditions for merging the proposed shaping device and the Paul trap, as well as better statistics. In this discussion, we have disregarded several factors which may decrease statistics (annihilation of \(\overline{\mathrm {H}}\) in the bottom disk, non-perfect merging of the angular acceptance of the optical device and the incoming beam of \(\overline{\mathrm {H}}\), quantum reflection of \(\overline{\mathrm {H}}\) from the reference plate, etc.). These factors have to be evaluated at a later stage.

5 Quantum-mechanical description

We now perform a quantum-mechanical description of the experiment, which will turn out to reproduce the main features and estimations of the quasi-classical treatment given above.

5.1 Free fall of a wavepacket

We consider the free fall of a pre-formed quantum wavepacket of \(\overline{\mathrm {H}}\) in the Earth’s gravity field, and we estimate the accuracy of the corresponding time-of-fall measurement. We know that the initial state \({\Psi }_{0}\left( z\right) \) of the wavepacket is a Gaussian function centered in the vertical direction \(z\) around the height \(H\) of the center of the trap, with the vertical position dispersion given by (2.8):

$$\begin{aligned} \Psi _0(z)= \left( \frac{m\omega }{\hbar \pi }\right) ^{1/4} \exp \left( -\frac{m \omega }{2\hbar }\left( z-H\right) ^2\right) . \end{aligned}$$
(5.1)

This wavefunction is calculated prior to the release, at a time where the gravity is compensated by the trap. After the photo-detachment event, the atom is suddenly released and its state is modified by the free fall in the gravity field.

This evolution is given by the propagation equation:

$$\begin{aligned} \Psi \left( z,t\right) = \int _{-\infty }^\infty G\left( z,z',t\right) \Psi _0\left( z'\right) \mathrm {d}z' \end{aligned}$$
(5.2)

where \(t\) is the free-fall time and \(G\) the propagator:

$$\begin{aligned} G\left( z,z',t\right)&= \sqrt{\frac{m}{2 i \pi \hbar t}} \exp \left[ \frac{im}{2\hbar t}\left( z-z' +\frac{\overline{g}t^2}{2}\right) ^2\right] \nonumber \\&\times \exp \left[ \frac{m\overline{g}zt+\frac{1}{6} m \overline{g}^2 t^3}{i\hbar } \right] . \end{aligned}$$
(5.3)

Integrating (5.2) for the initially Gaussian wavepacket (5.1), one gets

$$\begin{aligned} \Psi \left( z,t\right)&= \left( \frac{m\omega }{\hbar \pi (1+i\omega t)^2}\right) ^{1/4} \exp \left[ \frac{m\overline{g}zt+\frac{1}{6} m \overline{g}^2 t^3}{i\hbar } \right] \nonumber \\&\times \exp \left[ -\frac{m \omega }{2\hbar (1+i\omega t)} \left( z-H+\frac{\overline{g}t^2}{2}\right) ^2\right] .\nonumber \\ \end{aligned}$$
(5.4)

Assuming that all atoms annihilate instantaneously when they touch the detection plate at \(z=0\), we deduce that the distribution for annihilation times is given by the flux \(\mathcal {F}(t)\) of atoms passing through the plane at height \(z=0\), which is also the opposite of the current (downward velocities have negative values):

$$\begin{aligned} \mathcal {F}(t)&= -j(0,t)=-\frac{\hbar }{m} \text {Im}\left( \Psi (0,t)^* \frac{\partial }{\partial z} \Psi (0,t) \right) , \nonumber \\&= \sqrt{\frac{m\omega ^5 t^2}{\hbar \pi (1+\omega ^2 t^2)^3}} \left( H + \frac{\overline{g}t^2}{2} + \frac{\overline{g}}{\omega ^2} \right) \nonumber \\&\times \exp \left[ -\frac{m\omega }{\hbar (1+\omega ^2 t^2)} \left( \frac{\overline{g}t^2}{2} - H \right) ^2 \right] . \end{aligned}$$
(5.5)

This probability distribution is shown in Fig. 2 for an initially Gaussian wavepacket dropped from 30 cm, in the two cases of an initial size typically expected for the GBAR expected (upper plot) and the optimal size discussed above (lower plot).

Fig. 2
figure 2

The distribution of annihilation times of \(\overline{\mathrm {H}}\) falling from the height \(H=30\) cm; in the upper plot, the initial state is a Gaussian of width \(\zeta =70\) nm, a typical value expected in the GBAR experiment; in the lower plot, it is a Gaussian with the optimal width \(\zeta _\mathrm{opt}=88\) \(\upmu \)m. For comparison, the time scale is the same on both graphs. A zoom on the peak is shown in the inset for the lower plot

The optimal case (lower plot) leads to an extremely narrow time distribution, with a peak having the Gaussian shape deduced by expanding at lowest order in \((t-t_H)\) the distribution (5.5) :

$$\begin{aligned} \mathcal {F}(t) \underset{t\approx t_H}{\simeq } C \exp \left[ -\frac{(t-t_H)^2}{2 \Delta t^2} \right] . \end{aligned}$$
(5.6)

The width \(\Delta t\) of the distribution agrees with the classical result (2.3):

$$\begin{aligned} \Delta t= \sqrt{ \frac{\hbar (1+\omega ^2 t_H^2)}{2 m\omega \overline{g}^2 t_H^2} } = t_H \sqrt{{\left( \frac{\zeta }{2H}\right) }^{2}+{\left( \frac{\hbar }{2m v_H\zeta }\right) }^{2}}.\nonumber \\ \end{aligned}$$
(5.7)

The upper plot in Fig. 2, which corresponds to the typical numbers of the GBAR design, leads to a much broader distribution and shows a deformed shape with respect to a Gaussian distribution. As already discussed, this is a consequence of the large dispersion of initial vertical velocities.

5.2 Gravitational quantum states in the shaping device

We come now to the discussion of the shaping device in the regime where quantum gravitational states play an important role. The wavefunction of the atoms can thus be developed over the basis of eigenstates \(\Psi _n\) with energies \(E_n\) in the gravity field, here calculated above a perfectly reflecting mirror [17],

$$\begin{aligned} \Psi _n(z)=\frac{1}{\sqrt{l}}\frac{\text {Ai}(z/l-\lambda _n)}{\text {Ai}'(-\lambda _n)}, \quad E_n=mgl\lambda _n. \end{aligned}$$
(5.8)

The typical scale \(l\) of the gravitational quantum states is

$$\begin{aligned} l=\left( \frac{\hbar ^2}{2 m^2 \overline{g}} \right) ^{1/3} \approx 5.9\,\upmu \mathrm {m}, \end{aligned}$$
(5.9)

and the quantized energy levels are determined by the zeros of the Airy function \(\text {Ai}\):

$$\begin{aligned} \begin{array}{ll} &{}\!\!\!\text {Ai}(-\lambda _n)=0 \\ &{}\!\!\!\lambda _1\approx 2.34,\;\lambda _2\approx 4.09,\; \lambda _3\approx 5.52, \ldots . \end{array} \end{aligned}$$
(5.10)

The high-\(n\) states are given by the asymptotic law

$$\begin{aligned} \lambda _n\underset{n\rightarrow \infty }{\approx }\left( \frac{3\pi }{2}\left( n-\frac{1}{4} \right) \right) ^{2/3}. \end{aligned}$$
(5.11)

The selectivity of the shaping device is based on the sharp dependence of the transmission of eigenstates \(\Psi _n\) versus the height \(h\) of the slit. The detailed formalism in [17] leads to propagation through the device described by the following propagator:

$$\begin{aligned} K(z,z',t)=\sum _n \Psi _n(z) \Psi _n(z') \exp \left[ \frac{(E_n-i\Gamma _n) t}{i\hbar } \right] .\nonumber \\ \end{aligned}$$
(5.12)

The width \(\Gamma _n\) of level \(n\) becomes large for high values of \(n\) [17], as explained by the following qualitative interpretation. When the spatial dispersion \(l\lambda _n\) of the state \(\Psi _n\) is smaller than the slit size \(h\), the overlap with the absorber is small and the atom has a high probability to pass through the device (\(\Gamma _n\) small). On the other hand, when \(l\lambda _n\) is larger than \(h\), the overlap of the wavefunction with the absorber is significant and atoms have a high probability to be absorbed (\(\Gamma _n\) large).

Fig. 3
figure 3

Transmission of first, \(n=1\), and second, \(n=2\), gravitational states through a shaping device with a length of 5 cm

As a quantitative illustration, Fig. 3 shows the probability of transmission for atoms in the two lowest gravitational states, \(\Psi _1\) and \(\Psi _2\), when the length of the shaping device is \(R-r=5\) cm and the roughness amplitude of the top absorber is 1 \(\upmu \)m. A slit size \(h=24\) \(\upmu \)m provides 72 % transmission probability for the first state but only 0.3 % for the second state. This implies that a nearly pure ground state or a superposition of a few lowest gravitational states can be prepared by a suitable choice of the parameters of the shaping device.

5.3 Free-fall experiment after the velocity shaping

The output of the velocity shaping device is a superposition of gravitational quantum states \(\Psi _n\), determined by the propagator (5.12) calculated for a time \(t=(R-r)/v_\mathrm {hor}\) for an atomic horizontal velocity \(v_\mathrm {hor}\). This shaped superposition then falls freely to the detection plate so that the time distribution of annihilation events depends on the properties of the shaped state. We stress again at this point that this presupposes that the time \(R/v_\mathrm {hor}\) spent in the shaping device, and before its entrance, is corrected in the data analysis, \(v_\mathrm {hor}\) being deduced from the position of the annihilation event.

The spatial and velocity dispersions of the state \(\Psi _n\) can be expressed in terms of \(\lambda _n\) [22]:

$$\begin{aligned} \Delta z_n=\frac{2l\lambda _n}{3 \sqrt{5}},\quad \Delta v_n=\frac{\hbar }{ml}\sqrt{\frac{\lambda _n}{3}}. \end{aligned}$$
(5.13)

In contrast with the case of Gaussian wavepackets discussed above, these dispersions do not reach the minimum in the Heisenberg inequality. Furthermore, \(\Delta v_n\) and \(\Delta z_n\) increase simultaneously as functions of \(n\). The dispersion of the annihilation time (after correction of the time spent in the device) is thus given for the state \(\Psi _n\) by (2.2) with \(\zeta ,\upsilon \) replaced by \(\Delta z_n,\Delta v_n\) :

$$\begin{aligned} \frac{\Delta t}{t_H} = \sqrt{\frac{l^2 \lambda _n^2}{45 H^2}+\frac{l \lambda _n}{3H}}{\approx } \sqrt{\frac{\lambda _n l}{3H}}. \end{aligned}$$
(5.14)

As \(l\lambda _n\sim h\ll H\), the initial velocity spread still dominates the uncertainty on the annihilation time. It follows that the dispersion of these times is determined by \(\Delta v_n\) and scales as \(\sqrt{\lambda _n}\).

In order to get an estimate of the dispersions, we assume that the state in the shaper is an incoherent superposition of the quantum states which fit in the slit. It follows from the arguments in the preceding section that the quantum states which fit in the slit correspond to

$$\begin{aligned} n\le n_{\max }, \quad l\lambda _{n_{\max }}\approx h. \end{aligned}$$
(5.15)

We then deduce the dispersion of the annihilation times as

$$\begin{aligned} \frac{\Delta t}{t_H} = \sqrt{ \sum _{n} \pi _n \frac{l\lambda _n }{3H}}, \end{aligned}$$
(5.16)

where \(\pi _n\) is the population in the state \(\Psi _n\). As the slit size is small compared with the incoming wavefunction size, we expect that the states are equally populated among the fitting gravitational quantum states, so that \(\pi _n\approx 1/{n_{\max }}\) for \(n\le n_{\max }\), \(\pi _n\approx 0\) otherwise. In the quasi-classical limit where \(n_{\max }\gg 1\), we can use the asymptotic expression (5.11) for \(\lambda _n\) and replace the sum by an integral to find

$$\begin{aligned} \frac{\Delta t}{t_H} \approx \sqrt{\frac{ h}{5H}}. \end{aligned}$$
(5.17)

This expression scales like the classical result (4.2) with \(\beta \) now specified to be \(1/5\).

The preceding argument disregards the coherence between the components \(\Psi _n\) in the superposition prepared by the shaping device. This approximation can be justified qualitatively by considering that the effects of coherence are washed out in the averaging associated with free-fall propagation as well as the horizontal velocity dispersion. However, it cannot be considered as exact, and it will have to be confirmed by more precise simulations, to be published in forthcoming papers.

Exact quantum calculations can be performed for the special case of an initial state for free fall prepared by the shaper as the ground gravitational state \(\Psi _1\). The initial velocity distribution, shown in Fig. 4, has a width \(\Delta v \approx 9.5\) mm/s. This is 30 times larger than the optimal velocity spread \(\upsilon _\mathrm{opt}\approx 0.36\) mm/s, but two orders of magnitude smaller than the initial velocity spread in the GBAR experiment. The exact quantum evolution of this initial wavepacket is then obtained by integrating the propagation Eqs. (5.25.3). The annihilation time distribution calculated in this manner is shown in Fig. 5. Its spread is in excellent agreement with the prediction \(\Delta t = t_H \sqrt{{l\lambda _1}/{3H}} \simeq 0.97\) ms, deduced from (5.14). As a comparison, this spread was of the order of 45 ms for the free-fall measurement performed without velocity shaping. The improvement reflects the velocity selection by the shaping device, which is only partly balanced by the degradation of the statistics (as discussed above).

Fig. 4
figure 4

The velocity distribution in the ground gravitational state \(\Psi _1\)

Fig. 5
figure 5

The distribution of arrival times of \(\overline{\mathrm {H}}\) falling from the height of 30 cm, assuming that the initial wavepacket has been shaped into the gravitational ground state \(\Psi _1\)

6 Estimation of systematic effects

For our proposal to be useful as an improved option of the GBAR measurement, one must ensure that there are no large systematic uncertainties which could contribute at a level comparable to the estimated statistical uncertainty of \({\ 10}^{-3}\).

We first examine the additional velocity dispersion caused by the photo-detachment recoil. As discussed in [23], the vertical velocity dispersion due to the absorption of the photon and the positron emission can be kept small (\(\sim 0.5\) m/s) by using a horizontal polarized laser beam with an energy tuned at around \(\Delta E\approx 10\) \(\upmu \)eV \(\approx 0.1\) cm\(^{-1}\) above threshold. The photo-detachment cross section near threshold follows the Wigner law and can be estimated by using the available information in the literature to be \(\sigma =6.8\times 10^{-26} (\Delta E/1 \text { cm}^{-1})^{3/2}\approx 2\times 10^{-27}\) m\(^2\) [2427]. With a \(P=1\) W laser beam tuned close the threshold energy \(E_T=6,\!083\) cm\(^{-1}=0.76\) eV focused on an area \(A=10\) \(\upmu \)m \(\times 10\) \(\upmu \)m covering the Paul trap center, the photo-detachment rate is \(R=\sigma P/A E_T=130\) s\(^{-1}\).

In GBAR, antihydrogen ions can be produced only every 110 s, the ejection period of the antiproton decelerator at CERN. This time is sufficient to photo-detach the excess positron with high efficiency. The method is to illuminate the ion during a short enough time so as to define the start time with high precision, at a low enough repetition rate so that in case of successful photo-detachment, the free fall is completed before the next laser shot. For example, since the free-fall time on 30 cm is only 250 ms, laser shots of 100  \(\upmu \)s duration at a repetition rate of 2 Hz during 100 s allows the start time to be known with enough precision (\(4\times 10^{-4}\)), it also avoids ambiguity on identifying the successful shot, and it leads to a photo-detachment efficiency larger than of 90 %.

Since the velocity dispersion induced by the atomic recoil is of the same order as that from the confinement in the Paul trap, one would not gain by trying to get closer to the optimal cloud size. Finally, this effect is equivalent to a slightly warmer antihydrogen cloud, which changes the effective value of the frequency \(\omega \) to be used in the calculations, without affecting the principle of the method.

A careful analysis of other systematic effects has to be performed in the future, in particular for the following list of possibilities:

  1. 1.

    Uncertainty of shaping/measuring the distribution of the vertical velocity components of \(\overline{\mathrm {H}}\) within the range of acceptance of the two-disk system.

  2. 2.

    Finite positioning and timing resolution for the detection of annihilation events.

  3. 3.

    Accuracy and reliability for the correction for the time spent in the shaping device.

  4. 4.

    Diffraction of atoms on the mirror edges.

  5. 5.

    Residual electromagnetic effects, and in particular the patch effect on mirror surfaces.

  6. 6.

    Defects of mechanical alignments, such as inclinations of the disks and detection plate.

  7. 7.

    Finite precision of production and adjustment of optical elements.

  8. 8.

    Vibrations able to cause parasitic transitions between gravitational quantum states.

Monte-Carlo simulations of the experiment are under way; they take into account photo-detachment, coupling of the shaping device with the Paul trap and detector vessel, as well as points 1) and 2). For most of these systematic effects, one may also rely on the experience accumulated in experiments with UCNs [7, 9, 14, 15]. We note that the main systematic uncertainties (in particular 1) are proportional to the ratio \(h/H\), and thus they decrease strongly when slit heights are decreased. We therefore think that the control of these systematic effects will be improved at small slit heights.

7 Conclusion

In this paper, we have proposed a new method for shaping vertical velocities of antihydrogen atoms in the timing experiment to be performed by the GBAR collaboration [4]. We have given first estimations of the corresponding statistic uncertainties and listed possible systematic effects. The conclusion of these preliminary estimations, to be confirmed by further analysis, is that the accuracy in the measurement of the free-fall acceleration \(\overline{g}\) of \(\overline{\mathrm {H}}\) atoms could be pushed below 10\(^{-3}\) in realistic experimental conditions.

Statistical uncertainties in the experiment are improved for smaller slit heights, which lead to better defined vertical velocities of \(\overline{\mathrm {H}}\). This means that a better selection of the range of the vertical velocities overweighs the loss in statistics. Systematical uncertainties are expected to decrease even more dramatically for smaller heights of the slit between the two disks in the proposed experimental design. In the optimum experiment where the atomic wavepacket is shaped to the lowest quantum state, the effective temperature corresponding to the vertical motion of \(\overline{\mathrm {H}}\) is as low as 10 nK.

These preliminary estimations have to be confirmed by more complete simulations. We are currently working to develop a fully quantum treatment of the shaping device as well as a complete Monte-Carlo simulations.

Let us also mention that an even better accuracy could in principle be obtained by studying interference effects in the time-of-arrival distribution of a coherent superposition of a few lowest-lying gravitational quantum states [12, 13].