Modeling of the Near-Earth Low-Energy Antiproton Fluxes

The local interstellar antiproton spectrum is simulated taking into account antineutron decay, (He,p) interaction, secondary and tertiary antiproton production, and the solar modulation in the “force field” approximation. Inclusive invariant cross-sections were obtained through a Monte Carlo procedure using the Multistage Dynamical Model code simulating various processes of the particle production. The results of the simulations provided flux values of 4 ·10−3 to 10−2 and 10−2 to 1.7 ·10−2 antiprotons/(m2 s sr GeV) at energies of 0.2 and 1 GeV, respectively, for the solar maximum and minimum epochs. Simulated flux of the trapped antiprotons in the inner magnetosphere due to galactic cosmic ray (GCR) interactions with the atmospheric constituents exceeds the galactic antiproton flux up to several orders. These simulation results considering the assumptions with the attendant limitations are in comprehensive agreement with the experimental data including the PAMELA ones.


Introduction
The interest of antimatter component in the cosmic radiation ranges from the basics of cosmology [1] to a possible utility of antimatter fuel in future interplanetary missions [2].The antiprotons are perhaps the most intensively studied constituent of the antimatter accessible for observation.The experimental results of their measurements from [3] to [4] agree well with the hypothesis of their secondary origin [5].However, it does not exclude primary constituents produced, for example, in evaporation of primordial black holes or dark matter annihilations.
Before considering the existence of such exotic sources, it is imperative to evaluate precisely the antiproton production cross sections, uncertainties due to the propagation in the galaxy and heliosphere, and additional sources of antiproton production in the Earth's environment.The recent experimental data from different experiments at different solar activity phases provided a good opportunity for understanding the modulation process as well as the similarity in the proton and antiproton transport in the galaxy, interplanetary space, and the Earth's magnetosphere in spite of fluctuations in data and the model approximations.

Interstellar Antiproton Fluxes
In our simulation the antiproton local interstellar spectrum (LIS) F p (E p ) was obtained through the leaky-box model as a solution of the following system of integro-differential equations [6]: This considers the production of the secondary antiprotons Q 2 p by galactic cosmic ray (GCR) proton flux F p (E p ) and Advances in Astronomy subsequent tertiary antiprotons Q 3 p energy losses dE/dx in the interstellar matter, the flux decreases due to escape (λ esc ) and inelastic interaction (λ inel ).

Antiproton Production Spectrum.
The antiproton production spectrum, that is, the source function ) is a sum of the contributions from interactions of the GCR protons, He nuclei, and antiprotons with the interstellar H, He, and O nuclei in the interstellar matter.The corresponding densities n j are 1, 0.1, 8•10 −4 cm −3 [7].The antineutron production was also taken into account because they totally decay into antiprotons during confinement in the galaxy.The decaying antineutrons pass practically all their energy to the antiprotons with identical energy spectra.The production of antiprotons and antineutrons was simulated with a Multistage Dynamical Model (MSDM) Monte Carlo code [8].The code produces energy spectra and angular distributions of the reaction products (x, A) of an incident particle x with a target nucleus A together with total and inelastic cross sections and multiplicities.The projectile can be a hadron (n, n, p, p) or a meson (π + , π − , π 0 , K + , K − , K 0 ) with kinetic energies from 10 MeV up to 1 TeV.The target can be any nucleus with the atomic mass A ≥ 1.All the models included in the code were comprehensively tested by adjusting to the experimental data.
The MSDM code simulates all the stages of hadronnucleus and nucleus-nucleus interactions inside the target using an exclusive approach based on the models described in [9].It considers the cascade and precompound stages of the reaction as well as evaporation/fission, multifragmentation, and Fermi breakup of residual nuclei.At the cascade stage for projectile energies below 1 GeV, when consideration of only nucleons, pions, and D-resonances is sufficient, an original intranuclear cascade model [10] is used fitted to dynamics of formation and absorption of pions [11].For energies exceeding ≈10 GeV the independent quark-gluon string model is applied.
Hadron-hadron interactions for energies E h < 1 GeV are simulated on the base of parameterization of experimental data and πΔ-dynamics, for energies 1 GeV ≤ E h ≤ 10 GeV the improved version of the quark gluon string (QGS) model is used, and for E h > 10 GeV the model of independent QGS is used.
Figure 1 shows differential cross sections of the antiprotons dσ H (E p , E p )/dE p and antineutrons dσ H (E p , E n )/dE p for proton projectiles of various energies E p in (p, H) reaction (scattering of the data in Figures 1-4 is determined by only the statistic accumulated in Monte Carlo simulation.The overall statistics of several hundred millions events provides a necessary convergence of the final result).
One can see that the MSDM code provides some excess of the antiprotons over the approximation [12] widely used for simulations in GCR energy range.The spectra produced in (p, He) and (p, O) reactions are identical.The excess of the antineutrons over the antiprotons clearly seen in the figure is explained by a specific mechanism of their generation Antiproton multiplicity x100, [12] Antiproton multiplicity x100, MSDM Antineutron multiplicity x100, MSDM through the evaporation process favoring neutral particle production in the absence of the electromagnetic interaction.Inelastic, total (p, H) cross sections, and multiplicities are shown in Figure 2.  uniform distributions ∝ 1/projectile energy from [12].In contrast with (p, H) reaction, which produces practically identical spectra of p and n at all the projectile energies those spectra from ( p, H) reactions are similar only at higher projectile energies and are strongly different at low ones.As a parent LIS spectrum for the simulation a power low approximation F p (E p ) = 1.3710 4 E p −2.732 • (m 2 sr s GeV) −1 for protons and F He (E He ) = 7.0610 3 E He −2.699 • (m 2 sr s GeV) −1 for He [13] are used.The production spectra Q 2 p obtained from F p (E p ) are depicted in Figure 4.The cross section invariance for interaction of the GCR He with interstellar H was treated as a convolution of He LIS spectrum with dσ He (E p , E p )/dE p , where for He the E p is expressed as energy per nucleon.

Antiproton LIS.
In the simulation, the rigidity (R) dependent escape path length for antiprotons in the galaxy [14] is used: λ esc = 11.8β for R < 4.9 GV, λ esc = 11.8β(R/4.9) 0.54 for R ≥ 4.9 GV, the interaction length λ inel of antiprotons including annihilation is also simulated with the MSDM and the stopping power dE/dx is calculated utilizing standard procedure (e.g., [15]).Equations (1a), (1b), and (1c) is solved through an iteration procedure using the "Mathematica" package.The solution readily converges in the third iteration.As shown in Figure 4 the MSDM cross section provides about two-times larger tertiary output Q 3 p in the range of 0.3-3 GeV as compared to those from the uniform distribution used in [12].In the energy range of 0.04-2 GeV the LIS obtained with MSDM slightly exceeds that obtained with the approximation [12] with the maximum deviation of ≤40% at E p = 0.2 Gev.

Solar Modulation of the LIS.
The LIS modulation in the heliosphere is considered on the basis of a transport equation for the spherically symmetric case [16].The "force field" approximation inherently neglects the heliospheric gradient and curvature drifts but considers the diffusion, convection, and adiabatic deceleration: ( Here F HB , E HB , P HB , and F 1AU , E 1AU , and P 1AU are the proton flux, proton total energy in GeV, momentum in GV at heliospheric boundary (HB) and at the Earth's orbit (1 AU), respectively, m 0 is the proton rest mass in GeV.V is the average solar wind speed in 10 3 km/hr.Physical sense of the solution implies a conservation of the distribution function F/P 2 for particle energy decreases from E HB down to E 1AU in travel from heliosphere R HB to the Earth at 1 AU.
The heliospheric conditions are described by the "force field" parameter Φ determined by the solar wind speed V and the heliospheric boundary distance R HB .In our simulation we used A = 17, P c = 1.015GV [17], where it was shown that for the energies ≥0.02 GeV (2) approximates the exact solution of the equation in [16] for the power proton spectrum ∝ βE HB −γ (where β = P HB /E HB is the proton speed) and also consistent with the solar flare proton observations.BESS 95 + 97 [19]; BESS 00 [20], BESS 02 [21], BESS-Polar II [22], and PAMELA [4]."Force field" parameter Φ in GeV.

Advances in Astronomy
The Φ magnitude is determined from the best fit approximation with (2) of the observed proton spectrum F 1AU assuming the interstellar spectrum as F HB (E HB ) = 16470βE HB −2.76 protons/m s sr GeV.The fits results are Φ max = 0.964 GeV and Φ min = 0.368 GeV corresponding to solar maximum and minimum epochs.
The results of our simulated spectrum are shown in Figure 5 together with the antiproton fluxes obtained in different experiments performed at different solar minimum and maximum periods.The increases in the low-energy fluxes are provided by the higher fluxes of more energetic particles enriching the <1 GeV region due to adiabatic energy losses.The steeper the low-energy branch of the LIS spectrum the more pronounced the above-mentioned increases and [18].The results of the simulations provided flux values of ≈4•10 −3 to 10 −2 and 10 −2 to 1.7•10 −2 antiprotons/m 2 s sr GeV at energies of 0.2 and 1 GeV, respectively, corresponding to the solar maximum and minimum epochs.The curve for Φ = 1.5 is the lower limit for all the experimental data.It may correspond, for example, to V = 10 3 km/hour and R HB = 70 AU.

Geomagnetically Trapped Antiproton Fluxes
Another source of low-energy antiprotons in the magnetosphere is interaction of GCR with the Earth's residual atmosphere within geomagnetic trap.The source is quite the same as for the GCR antiprotons, but specific characteristics of the magnetic trap result in a significantly softer spectrum [23,24].
The flux captured within the plane of the geomagnetic equator can be simulated in frames of the radial diffusion model with an equation [25]: where μ is the antiproton magnetic moment, L is the McIlwain L-shell parameter in R Earth , f (L, μ) = F(L, E p )/ P(E p ) 2 represents the phase space distribution function (as in ( 2)), dμ/dt describes antiproton energy losses due to Coulomb scattering, λ inel describes inelastic and annihilation interactions as in (1a), the radial diffusion coefficient D LL is a sum of those of the magnetic D M LL and electric D E LL diffusion: The source functions Q 2 p (L, E p )and Q 3 p (L, E p ) are the same as in (1b) and (1c) but depend on L because of changing composition of the upper atmosphere [26] and on the cutoff rigidity modulating the parent CR proton and antiproton spectra.During the diffusion antiprotons born in situ at L-shell ≈1.2 sink in the underlying atmosphere or propagate to higher L-shells up to those where its rigidity is equal to that of the antiprotons and escape from the magnetosphere.
The Finite Difference Method (FDM) is utilized for the solution of the partial differential equation (4).The results of the simulations for the direct reaction are presented in Figure 6 for different antiproton energies.In the low-energy region the trapped population exceeds the CR antiproton one by a few orders of magnitude.
Equation ( 4) accounts for only the antiprotons produced in reaction in situ and trapped in the same place.One more source of the antiproton belt stems from decay of antineutrons born in another channel of the pair production reaction: (p + p, n + n + p + p) of GCR with the Earth's atmosphere constituents (this process was first considered in detail by [2] under a project supported with NASA's grant for search for effective space rocket fuel).In contrast to antiprotons, antineutrons owing to their zero charge travel far from their birthplace.They follow primarily the parent proton direction, but partially are backscattered in the atmosphere.Antiprotons resulted from decay of these albedo antineutrons acts as a supply for antiparticle radiation belts surrounding the Earth.This is the same CRAND (Cosmic Ray Albedo Neutron Decay) mechanism forming innermost proton radiation belt but with the antineutrons instead of the neutrons CRA ÑD Due to that the CRA ÑD flux can be calculated through normalization of the empirical CRAND model data by the antiproton-to-proton ratio for the antiproton source function: 1 antineutron for 10 5 up to 10 9 protons in dependence on the energy.This provides a few orders excess of this trapped population over the GCR antiproton flux in the lowenergy region.The belt is replenished every few years.L dependence of the antiproton flux from CRA ÑD source is shown in Figure 7. Maxima of the fluxes are located at the L = 1.3 − 1.6 [27] that is slightly higher being compared with that from the direct reaction (L = 1.25 in Figure 6); due to that this population was named as the "external antiproton belt" [23] due to its spatial position relative to the "inner belt" from the "in situ source" described above.
The flux magnitudes simulated for the PAMELA orbit in the South Atlantic Anomaly (SAA) region are compared with the experimental PAMELA results [28] and simulation [29] for L = 1.25 in Figure 8.The SAA region where the PAMELA trapped antiproton flux was measured (L = 1.15 − 1.4) is characterized by strong L and pitch-angle dependence of the trapped flux.Due to that the result of simulation there is very sensitive to small variation of the model parameters used.This may explain the difference between our simulation and that of [29].

Conclusions
The present simulation of the fluxes of interstellar origin incorporating solar modulation is attempted to explain the recent measurements of antiprotons at solar maximum and minimum.Particularly for the possible excess of the <1 GeV interstellar antiproton observations, initially the simulation considered the tertiary and antineutron decay antiprotons of the LIS source.The interaction cross sections by the MSDM Monte Carlo code provided a slightly larger antiproton flux in the energy range of 0.1-1 GeV compared to the [12] approximation.Then the "force field" solution for the solar modulation with rigidity dependence in compliance with    the LIS and the 1 AU spectra showed satisfactory agreement between the simulations and the balloon results at the solar maximum and minimum periods.
The simulation of geomagnetically trapped antiproton fluxes generated by GCR in the atmosphere showed a few orders of magnitude higher flux at L = 1.2 regions being compared with the interstellar antiprotons at sub-GeV energies.This agrees well with the last PAMELA experiment [28].

Figure 1 :
Figure 1: Antiproton (aP) and antineutron (aN) differential cross sections for (p, H) reaction.Symbols mark the MSDM results; thick solid lines represent the Tan and Ng [12] approximation.

Figure 2 :
Figure 2: Inelastic, total (p, H) cross sections, and multiplicities simulated with the MSDM code compared with the approximation from [12].

Figure 3 Figure 3 :
Figure 3: Antiproton and antineutron differentia cross-section of p + H reaction simulated with the MSDM code.

Figure 4 :
Figure 4: The interstellar secondary antiproton production spectra simulated using the MSDM and [12] cross sections.

1 Figure 6 :
Figure 6: The L-dependence of the trapped antiproton flux spectra for various energies of antiprotons.The horizontal segments of the curves correspond to the direct penetration of the cosmic ray flux to the magnetosphere.

Figure 7 :
Figure 7: The L-dependence of the trapped antiproton fluxes spectra in the Earth's magnetosphere from the CRA ÑD source for various energies.

Figure 8 :
Figure 8: Comparison of results of simulations with data of PAMELA experiment: ε is the criterion of stable trapping, P is the particle moment and B is geomagnetic field strength module.