Scattering of wave packets with phases

A general problem of $2\rightarrow N_f$ scattering is addressed with all the states being wave packets with arbitrary phases. Depending on these phases, one deals with coherent states in $(3+1)$ D, vortex particles with orbital angular momentum, the Airy beams, and their generalizations. A method is developed in which a number of events represents a functional of the Wigner functions of such states. Using width of a packet $\sigma_p/\langle p\rangle$ as a small parameter, the Wigner functions, the number of events, and a cross section are represented as power series in this parameter, the first non-vanishing corrections to their plane-wave expressions are derived, and generalizations for beams are made. Although in this regime the Wigner functions turn out to be everywhere positive, the cross section develops new specifically quantum features, inaccessible in the plane-wave approximation. Among them is dependence on an impact parameter between the beams, on phases of the incoming states, and on a phase of the scattering amplitude. A model-independent analysis of these effects is made. Two ways of measuring how a Coulomb phase and a hadronic one change with a transferred momentum $t$ are discussed.


Non-plane-wave states
In a quantum theory of scattering, the in/out-states are commonly chosen as delocalized plane waves. This model allows one to tremendously simplify the calculations; however, its limits of applicability are not explicitly articulated in the overwhelming majority of textbooks. Under standard conditions, finite sizes of the wave packets, their spreading during a collision, and finiteness of an interaction region do not play any essential role, especially for ultrarelativistic energies. There are, however, important exceptions.
The first example in which this model fails to work is collision of beams with the large impact parameters -the so-called MD-effect, observed at the collider VEPP-4 in Novosibirsk [1]. A somewhat similar effect is collision processes with a t-channel singularity when initial particles are unstable [2,3]. It is beams' finite sizes that provide natural regularization of the divergence. The next illustration is neutrino oscillations -a phenomenon that is intrinsically spatially and temporarily localized [4][5][6].
If colliding particles' wave fronts are neither plane nor gaussian, even approximately, then the plane-wave approximation is no longer valid either. The so-called vortex (or twisted) particles with orbital angular momentum (OAM) relative to a propagation axis and the Airy beams represent the simplest examples of such non-plane-wave states. They were shown to be solutions of the wave equations [7][8][9][10][11][12][13][14], and the corresponding beams of photons, electrons, and neutrons were generated in recent years [15][16][17][18][19][20]. Vortex electrons with the kinetic energy of 200 − 300 keV can be focused to a spot of anÅngström size [21], their OAM can be as high as ℓ = 200 [22], their magnetic moment increases proportionally to ℓ [13], and this brings about new effects in the electromagnetic radiation [23]. Such photons and electrons were also proved useful for optical manipulation [14], for probing phase of a transition amplitude, for creating pairs entangled in their OAM [24][25][26][27][28], etc.
Besides that, several groups have recently managed to create even more sophisticated photonic quantum states, including those combining features of the vortex-and Airy beams [29][30][31]; and one can await generation of the corresponding states of massive leptons and hadrons in near future. Although these novel beams differ from the coherent states, as they possess a distinct set of quantum numbers, the difference between them, mathematically speaking, lies only in the phases of their wave functions ψ(p).
A consistent relativistic scattering theory beyond the plane-wave approximation is absent by now, even though a number of non-plane-wave solutions for relativistic wave equations have long been known [32] and specific calculations were made [1][2][3][4][5][6]. For vortex beams, a corresponding ad hoc formalism has been recently developed by Ivanov and Serbo [25,26]. Generalization of their procedure for other quantum states (e.g. for coherent states, the Airy beams, etc.) may nevertheless represent a challenge. It is highly desirable therefore to have at hand an approach that would enable us to study scattering of the wave packets whatever their wavefront is, that is, for arbitrary phases of their wave functions. This work aims at developing such a method by generalizing the customary (plane-wave) S-matrix formalism. In doing so, we follow the work [1] in which the incoming particles are described by their Wigner functions.

Why Wigner functions?
The reader may well ask why we should deal with the Wigner formalism, especially when we are already going to treat sophisticated non-plane-wave effects. Indeed, although this approach was successfully applied both for non-relativistic-and relativistic scattering problems (see [33] and [1], respectively), it does not seem to have drawn much attention. The answer is that this formalism turns out to be the most convenient and elegant tool in the paraxial regime when σ p ≪ p , complementary and alternative to the wave-function approach (cf. [5,6,34]).
To put it in more detail, consider a generic matrix element which is a functional of the plane-wave one S (pw) (p ′ ) integrated over all pairs of momenta, p and p ′ , with some weights. Such products of amplitudes cannot, as a rule, be reduced to the standard traces and the terms such as a pair of spinors u(p)ū(p ′ ) have to be expanded over a complete set of 16 Dirac matrices. This makes the customary calculation procedure rather cumbersome and technically challenging (although not impossible per se).
Conversely, both the momenta coincide in the plane-wave limit when S (pw) f i (p)| 2 . Therefore when the packets are narrow, σ p ≪ p , the following inequality holds true: A density matrix written in these variables, (p + p ′ )/2 and p − p ′ , is called the Wigner function, introduced by Wigner in 1932 [35]. It is this inequality that allows one to develop a perturbation theory in which the ratio σ p / p serves as a small parameter, the in/out-states are described by their Wigner functions, and the number of events, a luminosity, and the cross section are expanded into σ p / p -series. When σ p → 0, the particle's wave front turns flat and the phase of its wave function does not contribute to the observables. Conversely, for any finite σ p the cross section gets corrections that depend explicitly on phases of the incoming states.

What other non-plane-wave effects can one expect?
Naively one would think that non-plane-wave corrections to the (plane-wave) cross section are attenuated as λ 2 c /σ 2 x ≪ 1 where λ c = /(mc) is a Compton wave length of a particle and σ x ∼ 1/σ p . As we shall demonstrate, at least some of them indeed are. There are, however, corrections with a more complex model-dependent behavior. For instance, for any finite σ p the observables grow dependent not only on the absolute value of the amplitude |M (pw) f i | but also on its phase ζ f i . The problem of determining the phase of scattering amplitudes is of high importance for hadronic physics and has a long history (see, for example, [36]). The ratio ρ(s, t) = Re M is calculated in different models, including Regge approaches, and it is extracted from elastic pp− and pp-collisions. It is known that at the small transferred momenta |t| this function is small, |ρ(s, t)| ≪ 1, on the energy scale from several GeV to several TeV [37][38][39][40], thus manifesting the high value of the phase ζ f i . A proper analysis of the elastic proton-proton collisions at the Large Hadron Collider (LHC) at √ s = 7 TeV by the TOTEM collaboration [41] has also shown that the real part of the amplitude can dominate for the large transferred momenta [42]. Finally, just a few months ago the TOTEM collaboration managed to estimate the amplitude's phase at √ s = 8 TeV in the interval of |t| from 6 × 10 −4 GeV 2 to 0.2 GeV 2 with a much better accuracy than before [43].
As we have recently shown [44] and demonstrate in more detail in this paper, scattering beyond the plane-wave approximation does not allow one to extract the Coulomb-or hadronic phase itself, but it does allow us to estimate how these phases change with t or with the scattering angle θ sc . We make an analysis that does not depend on phases of the incoming states. First steps towards this direction have been taken by Ivanov for vortex beams in [27].
The paper is organized as follows. In Sec.2 we give general formulas for the probability and the cross section within the Wigner formalism. We calculate several different Wigner functions in Sec.3, including the ones of a vortex particle and of an Airy state. In the latter case we compare an exact Wigner function with the corresponding approximate expression derived when σ p / p ≪ 1. In Sec.3.7 we give a general approximate formula for such a function of a wave packet with an arbitrary complex phase. The Lorentz invariant generalizations are given in Sec.3.8. Sec.4 is devoted to derivation of a general probability formula when the instates represent the wave packets with phases. The corresponding generalization for beams is presented in Sec.4.5. The first non-vanishing corrections to the plane-wave results are given in Secs.4.6, 4.7. A QED example is given in Sec.4.8. The effects of the scattering amplitude's phase are discussed in Sec.5. In Sec.6 we solve a similar problem but when the out-states are described as the Bessel ones with some OAM. We summarize in Sec. 7.
In what follows the terms "a wave function" and "a phase" will refer mostly to the momentum representation, and we mark σ ≡ σ p everywhere. As we do not put any limitations on the scattering amplitude, we shall use the term "scattering", for the sake of conciseness, for all the elastic and inelastic 2 → N f processes. The units = c = 1 are used.

Probability and cross section
Let us consider a generic (not necessarily elastic) scattering of two wave packets with N f final plane waves with the momenta p 3 , p 4 , ..., p N f +2 . As demonstrated by Kotkin et al. [1], the scattering probability can be represented as a functional of the (generalized) cross section dσ(k, p 1,2 ) and a function that we shall denote L(k, p 1,2 ) and call the particle correlator : where the amplitudes no longer depend on the normalization volume V , and n(r, p, t) is a (bosonic part of a) particle's Wigner function with the following properties: Note that in this approach we do not need fermionic Wigner functions (see, for example, [45,46]), because the spin parts of all wave functions are factorized and enter into the amplitude T (pw) f i . That is why the Wigner functions that we shall use in this paper are Lorentz scalars. The function dσ(0, p 1,2 ) coincides with a conventional definition of the plane-wave cross section. The dependence on k appears because of translational non-invariance of the Wigner function and in the plane-wave case this invariance is recovered (see below). For the wellnormalized Wigner functions, the probability (3) represents an unambiguous quantity in a sense that it does not depend on the auxiliary normalization variables, such as V and T . When making a comparison with the plane-wave case, it is convenient, however, to define also the effective cross section by dividing the probability by a luminosity factor L, Note that the quantities dW, L, and dσ are Lorentz-invariant, whereas the correlator is not. Probability formula when all the particles, including the final ones, are not plane waves but rather some generic wave packets with the phases can be derived following the very same procedure as described in Ref. [1]. For a special case with only two final states, the result reads: where L represents a 4-particle correlator, dn f is an integration measure for the final states, p and k denote the sets of all vectors: p 1 , p 2 , p 3 , p 4 and k 1 , k 2 , k 3 , k 4 , respectively. This expression can be useful when the final detected state is characterized not with the 3-momenta and probably with spins, but with a different set of quantum numbers. In this case the measure dn f may also include discrete variables such as the OAM. A need for spatial-and temporal localization of the detected states appears, for instance, in the theory of neutrino oscillations [4][5][6] or in QED calculations with the twisted photons [24].

Generalities
Before we turn to scattering, let us take a closer look at the Wigner functions and their properties. If the system is in a pure state with a wave function ψ(p), then its Wigner function can be found as follows (see, for example, a good pedagogical introduction by Case [47] or the textbook [48]) where In what follows we shall derive and analyse several examples of such functions. Since the integral in (8) cannot be always evaluated exactly, we shall develop a method for obtaining a suitable approximate expression, applicable when σ/ p ≪ 1. The literature on the Wigner functions of optical vortices, Airy beams and their generalizations is extensive -see, for example, Refs. [49][50][51][52]. Here we do not intend to make a comprehensive comparison of our results with those obtained in optics.
In this paper we employ only wave packets with a Gaussian envelope in momentum representation, bearing in mind that a wide class of packets can actually be approximated by this form (see discussion in Ref. [5]). Conversely, the corresponding wave functions in configuration space may appear to be non-Gaussian, because of the phases. Parameters characterizing the packets are: a mean momentum p , a momentum uncertainty σ, and a phase ϕ(p). While the latter is Lorentz invariant, a dispersion σ 2 generally transforms under the Lorentz transformation not as a scalar, but as a 3-vector: We also introduce a 3-tensor σ ij = diag{σ x , σ y , σ z }, which can be non-diagonal but still symmetric in an arbitrary frame of reference. For a boost along, say, z axis we have with σ ′ i in a frame where the packet is at rest on average. It is much more illustrative, however, to start with the customary "non-relativistic" expressions for packets with the only one σ. That is to say, we first work in the frame of reference in which Lorentz invariant generalizations for packets with σ → σ ij are straightforward, and it will be done at the very last stage in Sec.3.8.

Wigner function of a Gaussian wave packet
Let us consider a Gaussian wave packet with the following wave function: where r 0 denotes initial conditions. Lorentz invariance of the normalization requires that σ 3 (≡ detσ) transform as an inverse volume: in accordance with Eq.(10). Substituting (12) into Eq. (8), we see that the main contribution to the integral comes from the small values of k: |k| σ (recall the Ineq.(2)): Therefore, we can make the following expansion Within this accuracy, we arrive at the simple expression (compare with that for coherent states [33,48,53]) where r = r 0 + u(p)t.
As a cautionary remark, we note that this is not yet the true mean path of the system, as it still depends on p, not p . Although the function (16) itself does not spread, it is consistent with the momentum-coordinate uncertainty relations. It also satisfies the continuity equation, and represents a relativistic (3 + 1) D generalization (that is why n(r, p, t) ≤ 2 3 ) of the corresponding function of a quantum oscillator, and when → 0 it describes the so-called quasi-classical trajectory-coherent state of a boson (see, for example, [32,34,48]). The higher-order terms that we have neglected in (15), give corrections to (16) of the order of σ 4 . Indeed, taking into account the next term in the expansion (15), we obtain the corresponding correction to Eq.(16), where the second term in parentheses is supposed to be small compared to unity, and that is why this function stays everywhere positive in all the orders in σ. In the plane-wave limit with σ → 0, we find from where we get the following useful rule (recall Eq. (13))

Wigner function of a Gaussian beam
Having obtained the everywhere positive Wigner function of a Gaussian wave packet, we can now derive a similar expression for a beam of N b identical non-interacting particles. Such a function is a result of the statistical averaging of the wave-packets with some distribution over the packets' centers, f (r 0 ): and we imply that the distance between the particles, ∼ σ b /N b , does not exceed the coherence length of one wave-packet 1/σ, that is, For a Gaussian beam with we find: where Note that for a beam, the momentum uncertainty σ and the spatial width σ b are two independent parameters, and in the overwhelming majority of practical cases Say, for the LHC proton beam with σ b ∼ 10µm and σ/ p 1% [54], we have σσ b > 10 11 and, therefore, Σ 2 ≈ 1/σ 2 b . The case with σσ b ∼ 1 is realized, for instance, for 300-keV electrons focused in a spot of 1Å [21]; and now Σ 2 ≈ 1/(2σ 2 b ).

Wigner function of a Bessel state
A Bessel state is characterized with the following quantum numbers: the longitudinal momentum p , an absolute value of the transverse momentum κ, the energy ε(κ, p ), and a projection of the OAM onto the propagation axis, L z ≡ ℓ (see, for example, [12,13,25,55]). The wave function is where φ p is the azimuthal angle, and one can use the following rules The corresponding Wigner function can be found exactly applying Eq. (8): where and the Heaviside function is extended so that Θ(0) = 1. This Wigner function can be negative, it is invariant under the Lorentz boosts along the z axis, and it coincides up to a common factor with the so-called Wolf function in optics (see Eq. (22) in Ref. [52]). As can be easily shown, as should be for an azimuthally symmetric Bessel state. Hence, the maximum value of ℓ is The singularity ξ → 0 is integrable, and all the properties of a generic Wigner function (5) hold true. The lack of time dependence means the absence of spreading, the well-known feature of the Bessel beams. Unlike the wave function (28), the Wigner function (30) can possess any transverse momentum up to κ, and it is just maximized when p ⊥ reaches κ.

Wigner function of a packet with OAM
If an OAM eigenstate has a distribution over the transverse momentum, then such a state, unlike the pure Bessel one, has a finite OAM dispersion (or the OAM spectrum) when κ = 0 -see, for example, Refs. [56][57][58]. When calculating the corresponding Wigner function with the use of Eq.(8), we need to deal with the factor where φ ± are the azimuthal angles of the vectors p ± k/2. By analogy with Eq.(15), we make an expansion of this over the small k: This leads to the following result with r ⊥ being a transverse part of (17). When the longitudinal momentum has also the same uncertainty, we arrive at the generalization of the coherent state (16): As can be seen, the wave packet with the OAM implies, similar to the pure Bessel state, a finite transverse momentum and in the plane-wave limit, σ → 0, the OAM vanishes. We would like to emphasize that even for a state with p ⊥ = 0, the mean absolute value of the transverse momentum appears to be non-vanishing, as can be readily checked (see Eq. (40) in [59]). That is why it is sometimes helpful to think of σ as of the transverse momentum. The maximum value of the OAM is Such a state with p ⊥ = 0 has an azimuthally symmetric distribution of the intensity with a central minimum (see the left panel in Fig.1), thus representing a well-normalized generalization of the pure Bessel state. Unlike the Wigner function of the latter, however, Eq.(37) is positive for all the values of ℓ. The right panel in Fig.1 describes a similar wave packet with a non-vanishing transverse momentum. Such states can be useful for quantum entanglement in the OAM [60] and for probing the phase of the scattering amplitude in a collision experiment with vortex particles (see Sec.5).
Performing the Weyl transformation forL z operator, as explained for instance in Refs. [47,48], we find the OAM expectation value calculated with the use of these functions: as should be. And of course the OAM-dependent term in (37) does not change the mean trajectory r calculated with this Wigner function. If the OAM is quantized not relative to the z-axis, but along a unit vectorn, then one should make the following substitution in the Eqs. (36), (37): This allows one to use packets with the so-called orbital helicity for which n = p / p 2 [26].

Wigner function of an Airy particle
The normalized wave function of an Airy particle is parameterized with a 2D vector 1 ξ = {ξ x , ξ y } [7,11], which transforms as coordinates under the Lorentz boosts: The corresponding exact Wigner function is found with the use of Eq. (8): This function may become negative together with the arguments of the Airy functions. It is instructive to obtain the corresponding approximate expression by employing the small-k expansion, As a result, we arrive at a simple everywhere-positive function: This demonstrates explicitly that the negative values of the Wigner function are connected with the non-Gaussian O(σ 4 )-terms that we have neglected. An important distinction between the Airy phase and that of the vortex particle is that the fourth derivative of the former vanishes. That is why the neglected O(σ 4 )-term is the only correction to this Wigner function. From Eq.(44) we infer: In Figs.2,3 spatial distributions of the Airy wave packet are depicted for both the expressions: the exact-and approximate one. As can be seen, a transit from the not-everywherepositive Wigner function (42) to the everywhere-positive one (44) implies smoothing out the fast Airy oscillations for negative values of x, y, and the approximate Wigner function works better when ξ ξ max ∼ 1/σ. Since the total areas under the both surfaces coincide, there are also the regions where the black curve in Fig.3 may exceed the blue one (namely, when x = 0, y = 0).
In collision of an Airy particle with another wave packet, these fast oscillations in Fig.3 may play role when the latter is focused in a spot comparable to the oscillation period. Study of these effects obviously lies beyond the current experimental possibilities. In other words, although the differences in Figs.2,3 can be seen with the naked eye, the function (44) nevertheless represents a very good approximation in the scattering problems.

Wigner function of a packet with an arbitrary phase
Comparing Eq.(44) with the one for the vortex state (37), one can notice that they can be easily generalized for the wave function with an arbitrary, but regular complex phase ϕ(p), Then employing the same expansion, we find the following Wigner function: The similar expression for a Gaussian beam reads as: with the vector r b pointing to the center of the beam at t = 0. It is instructive to write down explicitly the condition of smallness of the higher-order terms that we neglected in Eq. (47): For Airy beams, this yields just σ 2 ≪ p 2 ∼ p 2 . This inequality is obviously violated for the vortex particle with p ⊥ = 0 being a singularity in Eq. (35). Consequently the expansion (47) is not applicable when p ⊥ → 0 or rather p 2 ⊥ ≪ σ 2 . However the Wigner function (37) itself is exponentially suppressed in this case as exp{−ℓ 2 σ 2 /p 2 ⊥ }, and that is why the formulas (37), (48) can still be used in the entire p-domain.
As before, the O(σ 4 )-corrections to these functions are responsible for the possible negativity. However in the paraxial regime with σ ≪ p , they are small, and a series like (20) could never make the Wigner function not-everywhere-positive. To put it simply, in scattering problems the possible negativity of the Wigner functions can reveal itself beyond the perturbative regime, that is, when σ ∼ p . This dictates focusing of the beam to a spot of which seems to be feasible only for cold systems with p ≪ m or σ b ≫ λ c where λ c = 1/m ≡ /(mc) is a Compton wave length of a particle. In the overwhelming majority of practical cases, therefore, the everywhere-positive Wigner functions can be used with a good accuracy.

Lorentz invariant generalizations
Lorentz invariance of the Wigner functions can be restored by making the following substitution in the wave packet (12): Then we have instead of Eq. (12): Here and in what follows The terms like ipr 0 can be made explicitly invariant by substituting −ipr 0 → i(εt 0 − pr 0 ) ≡ i(pr 0 ) with t 0 = 0 in our case. The corresponding generalization for coherent states (16) and (48) is and for a beam characterized with a symmetric matrix σ bij we arrive at the following formula instead of Eqs. (26), (49): when the matrices σ, σ b are diagonal. An analogous expression for an Airy particle can be readily guessed from Eq.(42) and we encourage the reader to make this generalization.

.1 Generalities
Having studied properties of the Wigner functions, we intend now to substitute these formulas into the general expression for the number of events from Sec.2 and then to expand it into series with a small parameter of σ/ p ≪ 1. As we shall see, the first term in this expansion represents the conventional plane-wave result, and the corrections to it embrace effects of finite monochromaticity of the incoming beams and of their spreading with time, of a finite impact-parameter, of phases of the incoming states, as well as of the general phase of the scattering amplitude. The three latter effects appear thanks to finite overlap of the incoming wave packets, which is proportional to This constant serves as a small parameter for such interference phenomena. As before, we start in the frame of reference in which σ ij = σδ ij . Generalizations for arbitrary frames will be made in Sec.4.7. In this case the overlap is described by the following matrix: which is symmetric but can be non-diagonal. Expansion of the probability and the cross section into series with a small α −1 turns out to be not a Lorentz-covariant procedure, but it becomes so in the relativistic case. Figure 4: Scattering with the in-states characterized by the mean momenta, momentum uncertainties, phases (will be added in Sect.4.2), and probably by the spins. Depending on the phases ϕ 1,2 (p 1,2 ), the wave front may become non-gaussian.

Benchmark case: 2 wave packets → plane waves
Let us start with collision of two wave packets with N f final plane waves (see Fig.4). Our current goal is to derive formulas for the observables that are reduced to the customary plane-wave expressions in the corresponding limit, σ 1 , σ 2 → 0, and can also be generalized when the in-states carry phases, as well as for collisions of beams. To this end, we take first the simplest (3 + 1) D coherent states 2 (16). Making use of the equality we arrive at the following formula for the correlator in Eq. (3): is the relative impact parameter of two particles at t = 0.
The small values of k give the main contribution to the integral in (3) and, consequently, one can make the following expansions 3 , where . We imply that the amplitude be a smooth-and analytical function of its arguments.
In what follows we shall need the integrals As these states are approximate, we neglect the terms O(σ 4 ) from the very beginning; see also Sec.4.4. 3 We shall not deal with the non-plane-wave matrix element S f i any longer, that is why below we shall omit where in our case we find: with α from Eq. (56). Note that by virtue of the equality ku 2 = ku 1 , one could have written instead u 2,i u 2,j in B ij . Using integral representations of the delta-functions, we obtain the following expression for the integral over k in Eq. (3): (63) and the integral over t ′ is gaussian and can be easily evaluated. After some algebra, we arrive at the following probability formula: Note that the first correction to |T f i | 2 vanishes when b = 0. An attentive reader might have already noticed that in the plane-wave regime, when σ 1 = σ 2 ≡ σ → 0, the following limits hold true: that is why we wrote down in (64)  If the wave packets do not spread much during the collision (say, for ultra-relativistic neutrinos [5]), that is, then one can neglect |β| compared to α, and the "standard" energy delta-function, δ(ε 1 (p 1 )+ ε 2 (p 2 ) − ε f ), is recovered. This is so, in particular, in the plane-wave regime: where |∆u| = (u 1 ( p 1 ) − u 2 ( p 2 )) 2 . If Eq. (22) is then applied, we have σ 2 /(2π) → 1/(2V 2/3 ), and the plane-wave cross-section is obtained dividing (66) by 1/(2V 2/3 |∆u|) instead of T /V (see Ref. [61]). This implies the so-called "time-to-space conversion" 4 [4], T = V 1/3 /(2|∆u|). One can avoid these (purely technical) subtleties by dealing with the cross section (6) instead, which is obtained by dividing (64) by the luminosity. Within the same accuracy, the latter is where the first terms in the exponent describes overlap of the incoming states. In the planewave limit, we have exactly the factor L → σ 2 υ( p 1 , p 2 )/(2π|∆u|) that cancels the extra term in (66) and leaves us with the conventional plane-wave cross section.
We would like to emphasize that the probability formula (64) includes plane-wave processes as a special case. Not only does this expression describe all the well-known specifically quantum phenomena, such as recoil-and spin-flip effects for instance, it also describes quantum effects that have no classical counterpart and vanish in the plane-wave approximation, although we have used the everywhere-positive Wigner functions in the derivation.

2 packets with phases → plane waves
Now let both the in-states possess phases, ϕ 1 (p 1 ) and ϕ 2 (p 2 ). Calculating the particle correlator in Eq.(3) with the Wigner functions from (48), we note that the phase terms do not depend on the integration variables and, consequently, the final results for the correlator, the luminosity, and for the probability are given by the Eqs. (59), (67), and (64), respectively, simply with the following substitution: We shall call this vector b ϕ the effective impact parameter. For the vortex-and Airy particles we have respectively, see Eqs. (37), (44). An envelope determining dependence of the probability (64) For vortex particles with a very small σ, this is proportional to with ℓ max from Eq. (38). If only one of the in-states is a plane wave with σ → 0, then and the phase-dependent terms (69) vanish anyway. The probability in this case does not depend on a phase of the second (non-plane-wave) state at all. In other words, the observables become sensitive to the in-states' phases if and only if the normalized wave-packets are used for both of the in-states. This observation generalizes the analogous conclusion for vortex particles [25] and can be easily understood: dependence upon the phases appears thanks to an overlap of the in-states, which are not orthogonal. When at least one of them is a plane-wave or even a non-normalized pure Bessel-or Airy state, then this overlap vanishes together with α −1 from Eq.(56).
We return now to the correction to |T f i | 2 in (64) that is linear in b ϕ and note that when the packets possess phases this correction survives even when b = 0. It depends on phases of the in-states as well as on the overall phase of the scattering amplitude, ζ f i (p 1 , p 2 ). Indeed, if one represents the amplitude as follows then one can obtain the following simple expression: where the function is odd in b ϕ . The second correction to |T f i | 2 in (64), on the contrary, does not depend on the effective impact parameter: where an imaginary part on the left-hand side vanishes when convoluted with a (m, n)symmetric expression in (64). In order to quantify an effect of the phase ζ f i , we define the following asymmetry: which vanishes in the plane-wave limit. Explicit formulas for the first correction dσ (1) to the plane-wave cross section dσ (pw) and for the asymmetry will be given hereafter.
There are two ways how one can change the sign of b ϕ : • If the incoming states are just wave packets with no phases whatsoever, the sign of b ϕ = b (the latter is to be small, b 1/σ, but non-vanishing) can be changed by replacing the initial wave packets: r 0,1 ↔ r 0,2 • Conversely, when b = 0 and the in-states possess phases, the change b ϕ → −b ϕ can be achieved by inverting a sign of parameters, such as the OAM ℓ or the vector ξ for Airy beams. The phases must contain only odd degrees of the parameters, which is the case for both types of the states.
It is clear that in both these scenarios we need to deal with the realistic beams of N b ≫ 1 particles and of a width σ b instead of single wave packets.

Alternative representation of the probability formula
When deriving the general probability formula (64), we used the approximate expressions for the Wigner functions and for the correlator, neglecting O(σ 4 )-corrections. Here we show how to obtain a formula, equivalent to (64) within this accuracy, but written in a more compact fashion, when the assumption of small σ 1,2 is made only once and at the very last stage 5 . Such a procedure allows us to come to the final result quicker, however, all the properties of the Wigner functions remain hidden here. We start with the exact (i.e. not expanded over the small σ) Wigner function, and the exact correlator: This yields the following (exact) expression for the probability: Now we expand all the functions in series over the small k up to the 2nd order inclusive, thus keeping terms not higher than O(σ 2 ), and then integrate over k, similar to the procedure in Sec.4.2. The result is: where we have used the following representation: and we have also denoted: The matrix B is still defined by the Eq. (62). In contrast to Eq.(64), in Eq.(79) the phase ζ f i also enters the exponent, and there are no terms linear in b in the pre-factor. As we have already neglected the O(σ 4 ) corrections when deriving both of these expressions, they can be easily shown to be equivalent within this accuracy.

Generalization for beams
The probability formula derived above describes scattering of two wave packets with the spatial widths ∼ 1/σ 1,2 . However the real beams with N b particles are many orders of magnitude wider in the majority of cases: say, for the LHC beam 1/σ is less than a femtometer and σ b ∼ 10µm. Taking as an example beams with the Gaussian distributions from Eq.(25), we perform statistical averaging of the particle correlator (77). The result is where the relative impact parameter of the two beams, is introduced and The corresponding number of events is That is to say, one just needs to replace σ 1,2 in α from Eq.(56) with Σ 1,2 in the probability formula (Eq.(64) or (79)): When σ b ≫ 1/σ, the effective values of k are |k| 1/σ b , and it is much smaller than those for p: |p| σ. As a result, the exponential envelope in the probability formula looks as follows (see Eq.(69)): The maximum values of the parameters entering the phases ϕ follow from the inequality: In particular, for beams with a phase vortex we come to the very same estimate (38), ℓ max ∼ p ⊥ σ b ∼ σσ b . Similar considerations for Airy beams, yield (cf. with Eq. (45)) Unlike the OAM's effective value for vortex beams, ξ max gets smaller when σσ b ≫ 1.

First correction to the plane-wave cross section
The integrals in Eqs. (64), (79) can be evaluated numerically for a specific model of T f i . It would be much more illustrative, however, to have at hand a purely analytical modelindependent expression for corrections to the conventional plane-wave results. To this end, we expand Eq.(79), or rather its generalization for beams, into Σ 1,2 -series when spreading of the packets is small (but not vanishing) according to Ineq. (65). Then we expand all the functions under the integral into p 1,2 − p 1,2 series, also up to the 2nd order inclusive, and finally integrate over p 1,2 and t. In doing so, one can notice that and by virtue of this we obtain the following intermediate result: where dN kin (σ 1,2 → 0) = dN (pw) and the second term, dN int , describes interference of the incoming packets. Neglecting the higher-order terms, we can simply take the integrand in dN int in the points p 1,2 . In doing so, however, we need to be cautious since: • The functions ∂ϕ 1,2 /∂p 1,2 | p 1,2 = p 1,2 may not be analytical everywhere or, in other words, in the expansion may not be small in the entire p domain. Say, for vortex beams it is big when p ⊥ ≪ ℓ/σ b . The probability itself, however, is exponentially suppressed in this case and that is why one can still use the expansion but keep in mind that p ⊥ > ℓ/σ b . The neglected small momenta contribute to the higher powers of ℓ/(σσ b ) < 1 (or to dN (n) , n ≥ 2) only.
• Parameters of the functions ϕ 1,2 can also depend on σ and σ b , which is the case, say, for Airy beams with ξ max = ξ max (σ, σ b ). As a result, for vanishing transverse momentum the higher-order σ-corrections to ϕ 1,2 may give contributions to the leading order: ξ 3 x,y p 2 x,y = ξ 3 x,y p x,y 2 + ξ 3 x,y σ 2 /2. When p x,y → 0, the seconds term survives.
Thus when making the small-∆p expansion, one should either suppose that the azimuthal asymmetry is broken from the very beginning, i.e. u ⊥,1,2 = 0, or keep the higher σ-terms in ϕ 1,2 so that we can always return to the special case of a vanishing transverse momentum. It is the latter scenario in which the mean value b ϕ , i.e. b ϕ weighted with the Gaussians ( √ πσ 1,2 ) −3 exp{−(p 1,2 − p 1,2 ) 2 /σ 2 1,2 }, comes into play. After the integration, we arrive at the following result for the first correction: (90) with C ij andb ϕ from Eq.(81) and p 1,2 = p 1,2 is implied everywhere. As a last step, we can also write down the first corrections to the luminosity and to the cross section, which are supposed to be small: The former is and everywhere p 1,2 = p 1,2 is implied. This expression is even in b ϕ , it also contains "kinetic" terms that are due to finite sizes of the packets and an interference term, proportional to α −1 and depending on the particles' phases. The analogous correction to the cross section is Here, as distinct from the number of events, dependence on the effective impact-parameter survives in the linear terms only, and when ζ f i = 0 (in some models on a tree-level, for instance) this correction depends neither on the impact-parameter nor on the particles' phases. Note that the term linear in b ϕ , ∂ ∆p ζ f i can be either positive or negative and this brings about the non-vanishing scattering asymmetry defined in Eq. (75). As we shall see, this term also breaks an up-down symmetry in angular distributions of the final particles.

Relativistic generalizations
Now let us return to an arbitrary frame of reference with σ → σ ij (recal Eq. (10)). Formulas for the correlator, Eqs. (59) and (77), and for the probability, Eqs. (64) and (79), stay valid with the following substitutions: and with the following matrix α instead of the scalar from Eq. (56): , or for a beam: where Σ ij is from Eq.(55). Then for a new matrix B we find: with β from Eq. (62). As a result, is denoted.
As neither the small "parameter" of this expansion, α −1 , nor the condition of small spreading (65) is Lorentz-invariant (actually the latter inequality is automatically fulfilled for ultrarelativistic particles), the resultant first correction to the number of events, unlike its plane-wave counterpart, turns out to be non-invariant either: where p 1,2 = p 1,2 is to be put. The terms neglected in Eq.(97) are O(γ −4 ) and that is why Lorentz invariance is restored in the relativistic case. Indeed, let in a frame where the packets are at rest on average σ x ∼ σ y ∼ σ z ≡ σ. Then in a collider frame of reference with we have: which is invariant under boosts along the collision axis and where spreading of the packets is neglected. The (diagonal) matrix α −1 ⊥ in this frame is Note that the corrections due to the finite sizes of the wave packets, σ 2 /m 2 , are positive. The corresponding expressions for the luminosity are and for the cross section we obtain the following formula: or for relativistic energies in the collider frame (100): which is also z-invariant and dσ (pw) is from Eq.(91). The first rows in (101) and (105) contain "kinematic" terms due to the finite width of the packets, while the second ones depend on derivatives of the amplitude in C ij and on the phase ζ f i . Note that the term b ϕ α ⊥ −1 b ϕ is absent in dσ (1) /dσ (pw) and that is why when ζ f i = 0 the first correction to the cross section, unlike the one to the number of events (101), depends neither on the impact-parameter nor on the phases of the incoming particles. In other words, the effective cross section turns out to be less sensitive (than the number of events) to the possible spatial inhomogeneity of the colliding wave fronts.

QED example: ee → e ′ e ′
In order to illustrate the non-plane-wave effects, let us take a head-on elastic scattering of electrons in quantum electrodynamics. The tree-level amplitude [61], has only an arbitrary constant phase due to the electrons' spinors u(p), and that is why one can put ζ f i = 0. We shall work in the collider frame (100) with the identical incoming beams of the width σ b ≈ 1/σ. The "kinematic" term in (101) is where λ c = 1/m ≡ /(mc) is the electron's Compton wave length. Even for the electrons focused to a spot of an Angstrom size [21] this ratio is of the order of 10 −6 .
The main correction to the plane-wave result therefore is For the small t, we find where the leading term, must be small by definition, although on the tree-level −t < 4m 2 . Thanks to the latter inequality, this correction can be much larger than (107) for the following scattering angles: For the same electrons with the kinetic energy of 300 keV from the Ref. [21] (regardless of the OAM), we get Measuring angular distributions of the scattered electrons at such angles is definitely challenging although not impossible. When the incoming states possess phases the number of events also depends on the effective impact-parameter: When we collide two Bessel beams with ℓ 1 = ℓ 2 ≡ ℓ, p ⊥ = 0 and σ 1 = σ 2 = σ, this correction becomes which is mostly negative as ℓ < ℓ max ≈ σσ b and where γ 0 ≈ 0.577 is the Euler's constant. Calculating the mean value of b 2 ϕ from Eq.(68) with b = 0 one should recall that p ⊥ > ℓ/σ b (see the discussion before Eq.(90)).
Whereas the correction (108) increases the number of events, the second one, Eq.(114), diminishes it. As a result, these two contributions can nearly compensate each other. In Fig.5 we show angular distributions of the scattered electrons with these corrections taken into account. As can be seen on the left panel, in a region where the correction is big but the perturbation theory still works, its contribution can reach the values of 10 − 20%.
On the tree-level, neither of the corrections depend on the azimuthal angle; however the amplitude's finite phase, ζ f i , as we shall demonstrate hereafter, restores this dependence.

Scattering asymmetry
When the phase ζ f i is non-vanishing, which is true on the loop level in QED or in the more sophisticated theories like quantum chromodynamics, the cross section also depends on signs of the incoming particles' phases (say, on OAM of the vortex beams) and also on the azimuthal angle. In order to quantify this effect we substitute the first non-plane-wave correction (104) into the asymmetry formula, Eq. (75). We arrive at the following compact expression: or when Σ x ≈ Σ y ≈ Σ z ≡ Σ: with b ϕ from Eq.(68) and Σ 1,2 from Eq.(83). This formula could have actually been guessed from the symmetry considerations. Indeed, for our kinematics the asymmetry can depend only upon the following three vectors: ∆u, b ϕ , ∂ ∆p ζ f i , and, simultaneously, it must be a linear function of the two latter ones. The only true scalar that satisfies these criteria is Eq.(116). We would like to stress, however, that this formula was obtained in the lowest order of the perturbation theory in α −1 and that is why |A| ≪ 1 or, at the best, |A| 1. Otherwise these expressions are inapplicable. Consider a 2 → 2 process (not necessarily elastic: say, pp → X, ee → X, etc.) with m 1 = m 2 in the collider frame of reference (100). Using the standard invariant variables, one can write down the derivative ∂ ∆p in this frame as follows For azimuthally symmetric dispersion with Σ x ≈ Σ y ≡ Σ and α −1 ⊥,11 ≈ α −1 ⊥,22 ≡ α −1 in (102), we get We shall stick to this model in what follows. When the incoming packets' widths are the same, σ 1 ≈ σ 2 , σ b,1 ≈ σ b,2 ≡ σ b and Σ 2 ≈ 1/σ 2 b (recal Eq.(83)), then As b ϕ ≡ b ϕ ⊥ , the asymmetry (118) is odd with respect to Therefore, the amplitude's phase ζ f i violates an up-down symmetry in angular distributions of the final particles, if this symmetry takes place without the phase of course. Note that for the strictly forward scattering, p 3 → p, the asymmetry vanishes.
As we have already mentioned in Sec.4.6, the averaging of b ϕ has appeared because some phases ϕ 1,2 may not be analytical in the entire p-domain, but contain a finite number of removable singularities. Say, for vortex beams with ϕ = ℓφ the derivative is not analytical for a vanishing transverse momentum. This singularity is removable and the mean value of this, simply vanishes when p ⊥ → 0. Note that this asymmetry is a purely quantum effect that vanishes in the plane-wave limit and might seem to be counter-intuitive from a classical perspective. Indeed, for a pair of azimuthally symmetric wave packets their substitution clearly does not alter the (classical) cross section. It is violated when either the packets are not-azimuthally symmetric (the 2 nd scenario) or the particles themselves have some inner structure (atoms, ions, hadrons). It is the latter case in which the phase ζ f i comes into play.

1 st scenario: off-center collision of Gaussian beams
As discussed in Sec.4.3, there are two ways how one can measure the asymmetry in a collision experiment. In the first one with two phaseless Gaussian beams collided at a nonvanishing impact-parameter one can put b otherwise the number of events is exponentially suppressed. Moreover, as clear from Eq.(118), one should not necessarily swap the beams: it is enough to have a non-vanishing b and then to compare angular distributions of the scattered particles in the upper-and in the lower semi-spaces. Their difference reveals itself in the asymmetry, which is It is only linearly attenuated with σ b and its pre-factor has a simple sin θ sc cos φ sc dependence upon the scattering angles θ sc ≡ θ 3 , φ sc ≡ φ 3 . Any deviation of the measured asymmetry from this dependence would be an evidence of a non-trivial phase ζ f i (s, t).
Further simplifications are possible for elastic scattering in the relativistic case with and now (compare this with Eq.(110)) where λ c = 1/m is the Compton wavelength of the incoming particle. We see that the asymmetry is only linearly attenuated by λ c /σ b , and it gets bigger when the momentum uncertainty of the beams approaches m and |t| becomes greater than 4m 2 (unlike the correction (110)). As we know, it is exactly when loop diagrams become significant. Assuming that the phase is a fast function of the scattering angle θ sc , but a slow one of p, we get the formula which shows how the phase changes with the scattering angle. Since in our approximation p ≈ ε = γm, one can re-write this formula as follows: The factors λ c /σ b and γ −1 ∂ ζ f i /∂θ sc are Lorentz invariant separately, and for protons the former is of the order of 10 −10 for moderately relativistic beams focused in a spot of ∼ 1µm and of the order of 10 −8 for protons with p ≈ 2 MeV and focused to σ b 10 nm [62,63]. The estimate (125), however, is inapplicable for such non-relativistic particles. Conversely, in collision of electrons the ratio λ c /σ b becomes bigger than 10 −3 for 300-keV electrons focused in a spot of the order of 1Å [21] (regardless of the OAM), although for such intermediate energies the formula (125) can be used only for qualitative analysis. According to West and Yennie, for a Coulomb phase on a one-loop level [36] 1 γ and hence with α em ≈ 1/137. This estimate is in accordance with that of the recent paper [64]. In this scenario, we bring two sub-nm-sized electron beams into collision (note that in this case 1/σ ∼ σ b ), slightly off-center, and that is why one ought to be able to control their relative position with the accuracy better than 0.5Å. Then angular distributions of the scattered electrons are measured and compared in the upper-and in the lower semi-spaces. Their difference reveals itself in the asymmetry and its conservative estimate for the scattering angles of θ sc ∼ 10 −2 − 10 −1 is which is in principle measurable with high statistics. One could further increase it by performing measurements at yet smaller scattering angles or by making the impact parameter very large, b ≫ σ b . In the latter case, however, the price is a drop in the number of events.
Returning to the elastic scattering of protons, little can be said, unfortunately, in a modelindependent way about the factor in the left-hand-side of (126). The TOTEM collaboration has managed to perform measurements at the scattering angles lower than 10 −5 at √ s = 8 TeV [43], which yields γθ sc ∼ 10 −2 − 10 −1 , and the hadronic phase ζ f i itself, unlike the Coulomb one, is not attenuated by a small parameter α em → α s , as scattering within a diffraction cone is not described by the perturbation theory.
As an example let us take three following models for the hadronic phase employed in the recent experiment at the LHC [43]: Bailly et al. [65], where ρ = ReM f i /ImM f i ≡ ρ(t = 0). Taking the same parameters as in [43], that is, √ s = 8 TeV, ρ = 0.1, t 0 = −0.5 GeV 2 , τ = 0.1 GeV 2 , t d = −0.53 GeV 2 , ζ 1 = 800, κ = 2.311, ν = 8.161 GeV −2 , we can estimate the asymmetry in Eq.(123). For the proton beam's width of σ b ∼ 10 µm we get the results shown in Fig.6. The derivative of the phase itself is no longer small, but the asymmetry is suppressed by the following factor: Although at the small transferred momenta interference of the hadronic phase with the Coulomb one can become prominent [43], the asymmetry (unlike the one for electron scattering) stays too small due to the large σ b . Thus the effects of the amplitude's phase are governed by width of the colliding beams.

2 nd scenario: colliding beams with phases
Within the second scenario, we start with a head-on collision of two vortex beams with the impact-parameter b = 0, the phases ϕ 1,2 = ℓ 1,2 φ 1,2 , ℓ ≡ ℓ z , and opposite signs of their orbital helicities [25]. The spatial distribution of such beams is no longer Gaussian but a doughnut-shaped one with a minimum on the collision axis (see Fig.1). As before, we need not necessarily to swap the beams or even change the signs of their OAM. It is enough to compare angular distributions of the final particles in the upper-and in the lower semi-spaces, that is, when φ sc → φ sc ± π.
Still working in the frame (100), we find with the help of Eq.(120): This vector vanishes, together with the asymmetry, when either the total OAM is zero, ℓ 1 + ℓ 2 = 0, or one collides the azimuthally symmetric Bessel beams with u ⊥ = 0 (which is implied in the frame (100)). Clearly, what is happening here is that in order to have a non-zero A the azimuthal symmetry of the problem must be broken already in the initial state, exactly as in the previous scenario. When the impact-parameter is vanishing, violation of the (initial) azimuthal symmetry can be achieved by shifting a phase vortex off the beam's symmetry axis. When dealing with the holograms (as in Refs. [16,17]), a shift of a fork dislocation off the beam center provides a (small) azimuthal asymmetry or, in other words, a non-vanishing transverse momentum (see details, for example, in Ref. [60]). The probability density for such a state is depicted in the right panel of Fig.1. Such a shift is to be small, δρ σ b , δp ⊥ 1/σ b , δθ ∼ 1/(pσ b ) and it is made to opposite directions for both beams. To put it simply, a non-vanishing transverse momentum plays in this scenario the same role as does a finite impact parameter in the previous section.
By using Eq.(118), we arrive at the following estimate for the asymmetry: The major difference between this expression and Eq.(121) is appearance of the factor ℓ 1 + ℓ 2 , which can be very large. It might seem therefore that the second scenario with ℓ 1,2 ≫ 1 provides a much higher value of the asymmetry. However the price for this increase is again a drop in statistics due to the factor exp{−ℓ 2 1,2 /(2σ 2 b,1,2 p 2 1,2,⊥ )} in the number of events. By analogy with Eq.(124), the factor determines sensitivity to the asymmetry in relativistic case. The maximum value of OAM for which the number of events is not suppressed is ℓ max ∼ p ⊥ σ b ∼ σσ b (see Eq. (38)), and that is why exactly as in the previous scenario.
Since the production of twisted hadrons with azimuthally non-symmetric profiles seems to be more technologically challenging than it is for electrons, we turn to elastic scattering of the latter particles. In order to maximize the effect, one can take again 300-keV twisted electrons focused to σ b ∼ 1Å, with the monochromaticity 6 σ/p 1%, and σσ b ∼ ℓ max ∼ 1. For measuring the asymmetry, angular distributions of the scattered electrons are to be compared in the two experiments with ℓ 1,2 = 1 and ℓ 1,2 = −1, respectively. Alternatively, one can carry out only one experiment with ℓ 1,2 = 1 when comparing angular distributions in the upper-and in the lower semi-spaces. The numerical estimate (128) stays valid. Since for such a study we need vortex electrons with the azimuthally asymmetric profiles, we would also like to find such states for which the requirement of a non-vanishing transverse momentum can be relaxed.
As can be readily seen, it is the case for Airy states as their azimuthal distribution itself is highly asymmetric. For collision of two such beams with u ⊥ = 0, the phases ϕ = (ξ 3 x p 3 x + ξ 3 y p 3 y )/3 (see Eq.(41)), and the opposite signs of their parameters ξ 1 = −ξ 2 ≡ ξ = {ξ, 0, 0}, we find: where we have used p 2 x = p x 2 + σ 2 /2. The typical values of ξ follow from the factor exp{−Σ 2 (σ 2 ξ 3 /2) 2 /2} in the probability formula: see Eq.(87). In any case, this yields the same p 3 /σ b factor in the asymmetry as in Eq.(121) and λ c /σ b for relativistic elastic scattering. Therefore the use of Airy beams leads to the very same predictions for the asymmetry as in the previous examples.
Moreover, one could think of such a phase ϕ(p) that maximizes the asymmetry. Within the paraxial case with Σ ≪ p , however, the phases are limited by Ineq.(86). That is why the asymmetry stays O(λ c /σ b ) for all the other types of non-plane-wave states as well.
The idea of using vortex states for probing the amplitude's phase was put forward by Ivanov [27]. By analogy with his work, let us consider now scattering of a light particle by a heave one (say, ep → X, γp → X) with σ 1 /σ 2 ≪ 1. Working in a frame in which the longitudinal momentum of the heavy particle is zero, we assume the light one to be in the pure Bessel state with u ⊥,1 = 0. We obtain that the asymmetry, does not depend on the OAM ℓ 1 of the light particle and, compared to Eq.(131), it has an additional small factor σ 1 /σ 2 , which is less than 10 −3 for available beams. This factor also appears for the Airy beams when p 1 = −p 2 but σ 1 ≪ σ 2 . That is why the higher values of the asymmetry favor the case with σ 1 ∼ σ 2 , in accordance with the Ref. [27]. The difference between the two methods described above can be elucidated by comparing two ways of colliding two rubber balls. If the balls are pumped up well, they are azimuthally symmetric and in order to violate this symmetry in scattering we need to collide them slightly off-center. Conversely, when the balls are deflated they are most likely no longer azimuthally symmetric and that is why they can collide even at a zero impact parameter. One simply needs to imagine a wave packet with a non-trivial wave front instead of such a deflated ball.
Concluding, scattering experiments probing the Coulomb phase, albeit being on the frontiers of technology, can be carried out at the modern electron microscopes, both with the Gaussian beams and with the vortex-and/or Airy ones if they are focused to a spot of the order of or less than 1Å in diameter. Predictions for the hadronic (or relative) phase are less encouraging, due to the small ratio λ c /σ b , and inevitably model-dependent.

2 packets with phases → 2 Bessel states
Describing detected states as plane waves, we lose information about all the other quantum numbers the evolved state may possess per se. The simplest example here is a Compton backscattering of an optical twisted photon by an ultra-relativistic plane-wave electron [24]. In this case, one may use an orthonormal set of Bessel beams for describing the outgoing particles, which can be treated as entangled in their OAM [28].
Bearing this in mind, let us choose now two out-states as the pure Bessel ones with the OAM ℓ 3 and ℓ 4 , and the Wigner functions from Eq. (30). Both the incoming particles are still described with the Gaussian wave packets (16) which can be later generalized to possess complex phases. For the sake of simplicity, we shall quantize both the final OAM relative to the same z-axis, which means that the scattering angles should be smaller than the momentum's conical angles. Generalization for the case with the so-called orbital helicity [26] is straightforward when the final states are also described as wave packets rather than as the idealized Bessel states. This problem, however, requires tedious calculations and will be tackled elsewhere. We shall also stick to the model with σ ij = σδ ij throughout this section.
The integration measure for the final states in the probability formula (7) can be chosen as follows (see, for example, [24,25]): At this point an important remark is in order. Using the detected non-plane-wave states with a definite set of quantum numbers, we imply that there exists an appropriate detector, which is sensitive to this set. For twisted photons, an OAM-sensitive "detector" may be thought of as a combination of a computer-generated hologram projecting the twisted state back onto the fundamental mode with ℓ = 0, a mono-mode fiber, and a "usual" detector -a setup routinely used in quantum optics with the parametric down-converted twisted photons (see, for example, [60,67]). For electrons and other massive vortex particles, an analogous registration scheme may also include a pair of appropriate holograms and a CCD camera.
Let us now derive a probability formula analogous to Eqs.(64),(79) but with the final Bessel states. First we find the following relation for the Wigner function (30): (136) where sin ξ = 1 − (p ⊥ /κ) 2 . This yields the following formula for the correlator from Eq. (7): where we have used k 1 + k 2 − k 3 − k 4 = 0. Now we return to Eq. (7), integrate over k 3 and k 4 , and then, similar to the procedure employed in the previous sections, we make an expansion of dσ (k, p) over the small k 1 , k 2 . The integral over k 2 is then eliminated with the use of the delta-function, δ(k 1 +k 2 −k 3 −k 4 ), and for the k 1 -dependent part we get the integral similar to (63) with the same matrix B ij , but this time with • Even in the lowest σ-order, there is a contribution from the phase of the scattering amplitude.
It must be noted, however, that the last effect takes place solely because of the non-normalizable nature of pure Bessel beams. If we had used the well-normalized wave packets with the OAM (37) instead, we would have come to |T f i (p)| 2 in the leading order, as in Eq.(64).
Note that the final Bessel states are monochromatic: they posses definite energy but undetermined azimuthal component of the momentum p φ , according to the angular momentumangle uncertainty relations [69][70][71][72]. The OAM-entanglement appears because the azimuthal components of the wave functions do not factorize in |ℓ 3 , ℓ 4 , that is, as a result of p φinterference. We would like to emphasize, however, that this happens only when both the incoming states are described with the well-normalized wave packets. Indeed, even if the instates have no OAM whatsoever, their OAM spectra are finite (see, for example, [57,58]) and, as a result, there is some overlap between both the azimuthal distributions. It is precisely this overlap that makes the final state a coherent superposition of the one-particle ones and the entanglement measure (144) finite.

Summary and outlook
We have developed the relativistic scattering theory beyond the plane-wave approximation in the paraxial regime when the incoming packets are narrow in the momentum space. The Wigner formalism turns out, therefore, to be the rather powerful tool that allows one to study effects accessible neither in the plane-wave approximation nor in the quasi-classical regime. These non-plane-wave effects are brought about by a finite overlap of the incoming wave packets. Depending on the phases, the packets represent the coherent states, the vortex beams carrying orbital angular momentum, the Airy beams, as well as their various generalizations. We have derived the general model-independent expressions for the probability and for the cross section and, when the non-plane-wave effects are small, have also obtained the analytical formulas for the first corrections to the plane-wave results.
In the latter case, along with the "kinematic" terms, λ 2 c /σ 2 b , there are also corrections that depend on the amplitude's model, that is, dσ (1) ∝ f (s, t) λ 2 c /σ 2 b . In a region of parameters where the function f is large, this correction is no longer vanishing and it can reach the values up to about 10 − 20%. For scattering of electrons with intermediate energies in QED, this happens at the scattering angles of a few tenths of a degree and at yet smaller angles for relativistic case.
Perhaps the most compelling finding of this study appears to be the azimuthal up-down asymmetry brought about by the scattering amplitude's phase. It is only linearly attenuated by the small parameter λ c /σ b . We have discussed two methods for probing this effect in experiments either with the conventional Gaussian beams or with such novel states as the vortex particles and the Airy beams. For Coulomb phase, the asymmetry is bigger than 10 −3 − 10 −4 for beams of the modern electron microscopes with the energies of hundreds of keV or less. For hadronic phase in proton-proton collisions, the similar effects are much weaker due to the ratio λ c /σ b ∼ 10 −11 for beams at the LHC.
The Wigner formalism is thus alternative and complementary to such well-developed quasi-classical methods as the trajectory-coherent approach with being the small parameter [34] or the operator method in which the ultrarelativistic motion is also implied [61]. Neither of these approaches has an advantage of explicit Lorentz invariance, although covariant generalization of the Gaussian wave packets seems to be feasible [73,74].
As we have also demonstrated, in elastic scattering of two particles at least one of which carries OAM the final pair gets OAM-entangled. In addition to the standard optical technique of the parametric down conversion, such a scattering (or annihilation) could become another method for obtaining the OAM-entangled beams, not only of photons but of the massive particles (including hadrons) as well. The somewhat similar spin-entanglement of final electrons has been recently studied experimentally in [75]. The analogous procedure can also be applied to other non-plane-wave beams with other sets of quantum numbers. It is of general interest, therefore, to generalize these results when the final states also represent wave packets with phases. This implies that they are detected with the appropriate apparatus and localized spatially as well as temporarily.