UHECR anisotropy and extragalactic magnetic fields with the Telescope Array

. We study the energy-dependent distribution of ultra-high energy cosmic ray arrival directions with respect to luminous matter in the local Universe. We use a specially designed test statistic (TS) that is robust to uncertainties of the galactic magnetic field. We generate realistic mock UHECR sets assuming various injected compositions, and di ff erent strengths of the extragalactic magnetic field (EGMF). Applying the TS to both mock sets and the Telescope Array Surface Detector data we constrain, for a given EGMF strength, the UHECR injected mass composition at energies above 10 EeV. At highest energies, higher than 100 EeV, only heavy composition — iron or at least silicon — is compatible with the data, irrespective of the EGMF strength. We then compare the obtained results with the direct Telescope Array fluorescence measurements of the UHECR mass composition. Requiring that the TA composition measurements are compatible with the arrival direction distribution allows us to constrain the parameters of the EGMF. It appears that light composition, measured by TA at around 10 EeV, is compatible with the arrival directions if the EGMF has strength of order 1 nG.


Introduction
Ultra-high energy cosmic rays (UHECR) are particles of cosmic origin with extreme energies in excess of 10 18 eV that bombard the Earth atmosphere producing extensive air showers observed either by large arrays of particle detectors on the ground, or by the telescopes that detect the fluorescence light emitted by the shower, or by both, the latter being referred to as the hybrid mode. Despite several decades of observations, the nature and origin of these particles still remains unknown, undoubtedly due to the indirect character of the existing detection techniques.
The key questions that remain unsolved are the mass composition of the incident cosmic rays, and their sources. Both observational techniques -the ground array and the fluorescent telescopes -allow for a direct reconstruction of the arrival direction and, to a lesser extent, of energy of the incident particle. However, in both cases the mass of the incident particle has to be inferred from indirect observables through complicated modeling of the air shower development. This only works so far with large uncertainties, and on a statistical rather than event-by-event basis.
The arrival directions of UHECR are measured with good accuracy and are basically free from systematic errors. However, they do not necessarily point back to the UHECR sources because of the deflections in the Galactic and extragalactic magnetic fields. These magnetic fields are not known well enough to reconstruct true directions to the sources. The unknown charges of particles make the problem even more complicated, so the two problems are entangled.
Although the UHECR accelerators are not yet identified, their global distribution must follow the large-scale structure (LSS) of the ordinary matter. The latter is inhomogeneous at scales ≲ 100 Mpc which is a typical size of the UHECR collection region at the highest energies. Depending on the magnitude of deflections, these inhomogeneities may have an imprint on the distribution of UHECR arrival directions over the sky, which may be possible to identify in the data. This may give information on typical deflections, from which both the cosmic ray composition and the magnetic fields may potentially be constrained.
A way to quantify these arguments has been recently proposed in Ref. [1] which is based on measuring a typical UHECR deflection angle with respect to the LSS and comparing it to the predictions of various composition models. In this contribution we apply this method to the most recent TA data. This allows us to estimate or to constrain from below, depending on the energy, the typical UHECR deflections. Using these results, even after accounting for possible uncertainties, some light composition models can be ruled out. We will also see that, if the direct TA measurements of the composition are included, the conclusions about the strength of the extragalactic magnetic field can be made.

Experiment and data
Telescope Array [2,3] is a hybrid detector of UHECR located at 39.3 • N, 112.9 • W in Utah, USA. Its surface detector (SD) consists of an array of 507 scintillator stations placed in a square grid with the 1.2 km spacing covering the area of ∼ 700 km 2 . The atmosphere above the surface array is overlooked by 38 fluorescence telescopes arranged in three stations. The detector is taking data in a hybrid mode since May 2008. The duty cycle of the SD is about 95% [4].
In this analysis we use the TA SD data collected during 14 years of operation from 2008-05-11 to 2022-05-10. We use the standard TA data set with the reconstruction optimized for anisotropy analysis; the corresponding quality cuts are described in Refs. [4,5]. It contains events with energy E ≥ 10 EeV and zenith angle z ≤ 55 • . The angular resolution of these events is estimated as 1.5 • [6]. The energy resolution is 18% [4,5], while the systematic uncertainty of the energy scale is estimated as 21% [7].
In addition to standard cuts, we also eliminate the events that may have been induced by lightnings [8,9]. The latter are taken from the Vaisala lightning database compiled by the U.S. National Lightning Detection Network [10]. We correlated the list of the lightning events detected within 15 miles from the Central Laser Facility of the TA during the full time of TA operation with the list of TA events and removed all the TA events that occur within 10 minutes before or after the lightnings. This cut was shown to reduce the total exposure by less than 1% [11]. The total number of events in the resulting set is 5966; of these 19 events have energies exceeding 100 EeV. We implement the analysis method in three steps. First, for each injected mass composition model considered we generate a large mock set of realistic UHECR events, represented by their directions and energies. Then, we quantify the mean deflection of a given event set with respect to the Large Scale Structure of the Universe (LSS) using the specially designed test-statistics (TS) which is robust to the uncertainties of the magnetic fields. Finally, we compute the TS for each mock event set and for the experimental data set and quantify the compatibility of each mass composition model with the data.

Analysis
Let us describe the analysis steps in more detail. To generate realistic mock UHECR event sets, we first build a state-of-the-art model of the UHECR flux. Basically, we set all the model parameters but nuclei injected fractions to some reasonable values. The effect of the uncertainties in these parameters can be estimated in the very end of the analysis, by varying them for typical composition models. We left this issue for further studies.
For the flux model we assume the UHECR sources to trace the luminous matter distribution in the local Universe. To achieve this we attribute an equal intrinsic UHECR flux for each galaxy in a complete volume-limited sample. Technically, we use a highly-complete fluxlimited galaxy sample. We derive it from the 2MRS galaxy catalog [12] in which we are cutting out galaxies with mag > 12.5 and those located closer than 5 Mpc and farther than 250 Mpc from the Earth. To compensate the observational selection effects of the flux-limited sample we assign a larger flux to more distant galaxies. We assume that the sources beyond 250 Mpc are distributed uniformly and have same mean density as the sources within this distance. The procedure of the catalog construction is described in Ref. [13]. This source model successfully represent all the scenarios with sufficiently numerous sources (source number density n ≫ 10 −5 Mpc −3 ). The sensitivity of our method for the models with source densities of order 10 −5 Mpc −3 , that are not yet ruled out [14], will be given in further studies.
We specify the flux model further by building a separate flux map for each injected primary. We use 5 benchmark pure injected compositions: protons (p), helium (He), oxygen (O), silicon (Si) and iron (Fe). We set the injection spectrum for each primary by making a separate fit to the TA and Auger observed spectra [15]. Namely, we set power law spectra with indices −2.55, −2.20, −2.10 and without cut-off for p, He and O, respectively; power law with index −1.50 and with sharp cut-off at 280 EeV for Si; power law with index −1.95 and with sharp cut-off at 560 EeV for Fe. For each primary we take into account the effect of absorption in interstellar medium. The flux of secondary particles produced during propagation of injected primary through the intergalactic medium is taken into account for He and O but neglected for Si and Fe, see Ref. [15] details.
We also simulate the deflections in magnetic fields with the account of particle charge Z and its energy E. We take into account the particles deflection in regular [16] and random [17] galactic magnetic field (GMF) components. The deflections in the extra-galactic magnetic field (EGMF) are neglected for the basic flux model, this corresponds to magnetic field strength B EGMF ≪ 1 nG for the correlation length λ ∼ 1 Mpc. We also build a separate flux model with non-negligible deflections in EGMF to test its impact on results.
In order to make mock event sets similar to data sets we also take into account the experimental effects: the energies of the events are generated according to the observed TA spectrum [18] and we add effects of TA experiment energy resolution and exposure to the simulations. The examples of two flux maps used for mock event set simulations are shown in Fig. 1.
The test-statistics we use is based on the UHECR flux maps similar to those used for the mock event sets generation, but with smaller number of free parameters. Namely we use the same source catalog and the flux is attenuated as protons with ∼ E −2.5 injection spectrum without cutoff. The only free parameter of the TS is the magnitude of the uniform smearing of sources defined at the energy 100 EeV: θ 100 .
For each given value of θ 100 we build a set of maps Φ k (θ 100 , n) where n is the direction in the sky and k is the number of energy bin. The smearing of each map scales as 100 EeV/E k . The test statistics TS(θ 100 ) for a given event set with directions n i is defined in a following way: where the sum run over the events in the set i and energy bins k. For the particular analysis we use logarithmic energy bins of 0.05 decade width. The TS is normalized with a factor Φ iso (n i ) = Φ(∞, n i ) that corresponds to an isotropic distribution of events. This TS was first introduced in our study [1], where it was shown that for a large number of events it is distributed around its minimum according to χ 2 -distribution with one d.o.f. The position of the TS minimum θ min 100 for each event set has a meaning of the energy-rescaled average event deflection with respect to the sources in LSS. Therefore, for a large enough mock event set of a given composition the TS should have a deep and narrow minimum and the respective value θ min 100 characterizes this composition model. These values could then be compared with the TS(θ 100 ) distribution for the data sets. The test statistics Eq. 1 calculated for the data in five energy bins is shown in Fig. 2 as a function of θ 100 . Note, that the plotted bins are logarithmic and have 0.25 decade width, i.e. each of them contain 5 technical energy bins used for computation of Eq. 1. We limit the range of θ 100 by 200 • as set by the size of the TA field of view; beyond this value we consider configurations to be essentially isotropic. It can be seen that in all energy bins the small values of θ 100 are excluded, the strongest exclusion arising from the highest energy bin E > 10 20 eV. The low energy bins show indications of minima at 10 − 30 • which are not significant (recall that 1-σ interval corresponds to the change of the TS by 1). The minimum is more pronounced in the energy bin 10 19.75 − 10 20 eV where the TA data exhibit the hot spot of the size ∼ 25 • matching the position of the minimum in this bin. Surprisingly, the highest-energy bin E > 10 20 eV shows no minimum and prefers instead a nearly isotropic distribution with deflections θ 100 ≳ 50 • at 2σ CL. The positions of the minima of the TS (where they exist) and the 1-and 2-σ confidence intervals are shown in Fig. 3. In order to interpret this result in terms of composition models we generate, for each given composition, a large mock set of events as we describe in previous Section. For each of these sets we calculate the TS in exactly the same way as for the data. Since the mock sets are very large, the statistical uncertainties in the positions of the TS minima are negligible. The comparison of the model predictions with the data is shown in Fig. 4. The cases where the minimum for the model does not exist (like Si and Fe in low energy bins) are shown by upward arrows.

Results
It is tempting to conclude from Fig. 4 that predominantly light compositions are excluded, but this is not so straightforward without examining mixed compositions in various combinations, which we do not attempt in this contribution. Let us instead focus on the highest-energy bin E > 10 20 eV where conclusions can be drawn from already existing results. As can be seen from Fig. 4, in this energy bin no pure composition fits the data, and even the heaviest -pure iron -is disfavored by the data at about 2σ level. The same is clearly true for an arbitrary mixture. We conclude therefore that the isotropy of the data at E > 10 20 eV indicates a very heavy composition at these  energies. It should be stressed that this result is obtained with the account of the deflections in the GMF. The latter is not well known; it can be checked however that the conclusion is robust: changing the GMF model to the one of Ref. [19] affects the model predictions only slightly, and thus the conclusion remains. We will see shortly that assuming a non-negligible extragalactic MF also does not change the picture. Let us now focus on the low energy bins. Here it is not possible to conclude that the composition should be heavy without the detailed study of the mixed composition models. It may appear that the data still exclude the light p-He mix. This conclusion, however, is not robust with respect to possible effects of the EGMF. While the EGMF is not known, there exist experimental upper bounds at the level of 10 −9 G assuming maximum correlation length of 1 Mpc [20]. If the EGMF is pushed to this extreme value, the model predictions change. This change is shown in Fig. 5 where the predictions for pure proton and helium compositions are compared to the data in two cases: zero EGMF, and maximum EGMF satisfying the current limits. one can see that inclusions of deflections in the EGMF can make the light composition models compatible with the TA data.
This argument can be inverted if extra information concerning the composition is included. In the energy bin 10 19 −10 19.25 eV there exist a TA measurement of the composition from the X max distribution of the fluorescence detector (FD) data [21]. This measurement indicates a light composition consisting of the proton-helium mix. If this measurement is combined with the above measurement of event deflections, together they imply that the EGMF should be large, at its current experimental upper-limit.
It remains to be checked whether the large EGMF, as large as indicated by the combined TA SD and FD data, changes the conclusion about a heavy composition at the highest energies. From Fig. 5 one can see that in the highest energy bin the light elements are still not deflected enough to explain the isotropy of the events. Fig. 6 shows how the predictions for the pure heavy composition (Si, Fe) change upon the inclusion of the maximum EGMF. In the highest energy bin, the predicted deflection for iron become compatible with the data, while those for silicone become marginally compatible. The conclusion about heavy composition at the highest energies remains, therefore, unchanged.

Conclusions
To summarize, we have introduced a quantitative measure of the overall cosmic ray deflections in a given data set with respect to the large-scale structure on the Universe, θ min 100 defined in terms of Eq. 1. We then calculated this quantity for the most recent TA SD data set in 5 energy bins at E > 10 19 eV. The result is shown in Fig. 3. This quantity can be compared to the predictions of different composition models. Despite the uncertainties inherent in these models, two conclusions can be robustly deduced: • In the highest energy bin E > 10 20 eV only the heaviest composition -iron if no EGMF is assumed, or possibly silicone in the case of a maximum EGMF -gives deflections large enough to be compatible with the measured value. • In the lowest energy bin 10 19 − 10 19.25 eV, the independently measured (from the TA FD data) light composition is only compatible with the observed large deflections if a maximum allowed EGMF is assumed.
We see therefore that the observed high isotropy of the UHECR and, consequently, their large deflections, indicate the heavy composition of UHECR at highest energies, and large, at the level of existing experimental limit, extragalactic magnetic fields.