A personal decadal survey in ultrahighenergy cosmic-ray observatories

Toshihiro Fujii (fujii@cr.scphys.kyoto-u.ac.jp)
The Hakubi Center for Advanced Research, Kyoto University
ICRR seminar, 27th December 20181

Due to a significantly improved statistics of ultrahigh-energy cosmic rays
(UHECRs) by the Telescope Array Experiment and Pierre Auger Observatory,
we firmly confirmed a suppression of the energy spectrum at the highest energies
and observed intriguing large-scale anisotropies in arrival directions of UHECRs.
We also encountered a gradually transition to a heavier composition, a deficit of
the number of muons in simulations and a lack of desired small-scale anisotropies.
In this talk, I highlight recent results of the two observatories including on-going
updates and then address scientific goals and requirements for future UHECR
observatories in next decade. I introduce three ideas as a personal decadal survey:① a fine-pixel fluorescence telescope for low-energy extension, ② a layered
Water-Cherenkov detector array for sub-EeV anisotropy and ultrahigh-energy
photon search, and ③ a low-cost fluorescence telescope array suitable formeasuring the properties of UHECRs with an unprecedented aperture.
2the energy dependencies of anisotropies. We emphasize here that the energy spectrum, thenumber of cosmic ray particles per time in a unit area from a given direction in a given
energy range is, by definition, a function of the direction. The measurement of the full-sky
energy spectrum by the future Auger and TA will make a crucial contribution to identifying
the sources
2.5 of ultra- high energy cosmic rays.

Energy spectrum of cosmic rays

E

J(E)

Hillas&plot

1 particle/
km2/century

1 particle/m2/yr

Ultrahigh-energy cosmic
rays (UHECRs)
1 EeV= 1018 eV

V. Verzi et al., PTEP, 2A103 (2017)

Starburst
galaxies

γ : Lorentz factor of shock
Z : atomic number
B : magnetic field strength
R : sizeA. M. Hillas, Astron. Astrophys., 22, 425 (1984)

the energy dependencies of anisotropies. We emphasize here that the energy spectrum, thenumber of cosmic ray particles per time in a unit area from a given direction in a given
energy range is, by definition, a function of the direction. The measurement of the full-sky
energy spectrum by the future Auger and TA will make a crucial contribution to identifying
the sources
2.5 of ultra- high energy cosmic rays.Cosmic'microwave'

Energy spectrum of cosmic rays

E

J(E)

background'radia5on'(CMBR)'

Cosmic'Ray

1,2,3 : m=0, γg=2.7, 4 : m= 4, γg=2.45

Earth

Cosmic'Ray25-1log10J(E)E ,m sec sr eV

2

4

-1

1,2,3

E2× Flux (TeV cm -2 s-1)

supernova remnants,
galactic center?
W44

IC443

W44

IC443

10

Extragalactic
origin?

-10

17
68% CL= 2.9 PeV

90% CL
cut,p
95% CL
cut,p= 0.6 PeV

Model: Diffuse emission E
Model: Diffuse emission E

48

10
10-10

1047
10-11

0

Fitted π decay model

1046

Derived Proton spectrum

-12

10

VERITAS (Acciari et al. 2009)
MAGIC (Albert et al. 2008)

8

10

9

10

10

Proton Spectrum E2 dN/dE (erg)Gamma-ray flux E2 dN/dE (erg cm-2 s-1)

Model: Diffuse emission Ecut,pIC 443

V. Berezinsky etal., hep-ph/0107306 (2001)

3

1

H.E.S.S. collab., Nature
531, 476 (2016)

Model (best fit): Diffuse emission

1049

4

23

10-12

W44

3
E J(E)

10-11

Diffuse emission (× 10)

10-9

3

-2

Galactic origin

24

-13

= 0.4 PeVHESS J1745-290

1

10Energy (TeV)

Figure 3: VHE -ray spectra of the diffuse emission and HESS J1745-290. The Y axis shows fluxes multiplied bya factor E2 , where Eais the energy on the X axis, in units of TeVcm 2 s 1 . The vertical and horizontal error bars show
the 1 statistical error and bin size, respectively. Arrows represent 2 flux upper limits. The 1 confidence bands ofthe best-fit spectra of the diffuse and HESS J1745-290 are shown in red and blue shaded areas, respectively. Spectral
parameters are given in Methods. The red lines show the numerical computations assuming that -rays result fromthe decay of neutral pions produced by proton-proton interactions. The fluxes of the diffuse emission spectrum andmodels are multiplied by 10.

b

19

20

log10E, eV

FLUX MAP ABOVE 8 EeV

E > 8 EeV

Where is thements. The curves show the predicted diﬀerential spectra for the uniform distribution of sou
maximum energy
or without evolution. The case without evolution (m = 0, γ = 2.7) is given by curves (1of cosmic ray,for maximum generation energy E= 3 · 10 eV, 1 · 10 eV and ∞, respectively. Th
Emax?
curve 4 describes the evolutionary model with m = 4, γ = 2.45 and E= ∞.

FIG. 3. UHECR spectrum
as observed in Akeno (triangles) and AGASA (filled circles
Equatorial coordinates

Sgr A*

max

20

21

g

10

10

11

10

21

g

9

Sgr A*

18

2

max

12

10

Energy (eV)Figure 3: Proton and gamma-ray spectra determined for IC 443 and W44. Also shown are
the broadband spectral flux points derived in this study, along with TeV spectral data points for
IC 443 from MAGIC (29) and VERITAS (30). The curvature evident in the proton distributionat ∼ 2 GeV is a consequence of the display in energy space (rather than momentum space).
Gamma-ray spectra from the protons were computed using the energy-dependent cross section
parameterized by (32). We took into account accelerated nuclei (heavier than protons) as wellas nuclei in the target gas by applying an enhancement factor of 1.85 (33). Note that modelsof the gamma-ray production via pp interactions have some uncertainty. Relative to the model
adopted here, an alternative model of (6) predicts ∼ 30% less photon flux near 70 MeV; the twomodels agree with each other to better than 15% above 200 MeV. The proton spectra assume
average gas densities of n = 20 cm−3 (IC 443) and n = 100 cm−3 (W44) and distances of 1.5
kpc (IC 443) and 2.9 kpc (W44).

Figure 1: VHE -ray image of the Galactic Centre region. The colour scale indicates counts per 0.02 ⇥0.02 pixel.
Left panel: The black lines outline the regions used to calculate the CR energy density throughout the central molecularzone. A section of 66 is excluded from the annuli (see Methods). White contour lines indicate the density distributionof molecular gas, as traced by its CS line emission30 . The inset shows the simulation of a point-like source. Rightpanel: Zoomed view of the inner ⇠ 70 pc and the contour of the region used to extract the spectrum of the diffuse
emission.

Fermi-LAT collab. Science 339, 807 (2013)

V. Verzi et al., PTEP, 2A103 (2017)

Galactic centerWe can fit the Akeno-AGASA data in both cases, with and without evolutio

spectra without evolution, m = 0 can fit the data starting from relatively high18
Collab.
Science1266(2017)EPierre≥ 1 · 10Auger
eV. Thefit needs
γg = 357,
2.7. The
curves1, 2 and 3 in Fig.3 show the spec
9

diﬀerent Emax equal to 3 · 1020 eV, 1 · 1021 eV and ∞, respectively. The fit without eUHECR detections
Longitudinal Development
Fluorescence
detector (FD)~10% duty cycle

Xmax

Lateral Density Distribution
CORSIKA software

Surface detector
array (SD)~100% duty cycle

https://web.ikp.kit.edu/corsika/movies

5

UHECR observatories
Telescope Array
Experiment (TA)Utah, USA2
700 km2
TA×4 → 3,000 km

Pierre Auger
Observatory (Auger)
Malargue, Argentina2dŚĞƐƚƌĂƚĞŐǇ
3,000 km ŽĨƚŚĞhƉŐƌĂĚĞ

Google EarthAugerPrime:
additional scintillatorDĞĂƐƵƌĞǁŝƚŚƚŚĞWŝĞƌƌĞƵŐĞƌKďƐĞƌǀĂƚŽƌǇƵŶƚŝůϮϬϮϱ͘;DKhƐŚĂǀĞďĞĞŶƐŝŐŶĞĚŝŶEŽǀϮϬϭϱͿ

6

occurs at 1019.62 ± 0.02 eV, while the corresponding break in TA is seen at 1019.78 ± 0.06 eV, a factor
Calibration of TA FD with ELS
Bokkyun Shinof 1.4 ± 0.2 higher.

Energy spectrum at the highest energies
M. Ave et al. / Astroparticle Physics 42 (2013) 90–102(FY)
beam. The emittance (e) and Twiss parameters (a,b ,g) areFluorescence
related with particleyieldposition(X)
10% energy scale differenceofintheTA/Auger0 ) by e = b X 2 aXX 0 + gX 02 . The emittance and Twiss parameters were measured
Table
6and direction(XSystematic
uncertainties: 14%
in Auger,
21% in TA
Systematic uncertainties on the fluorescenceby
yield
measuredlaser calibration,scanning
the withquadrupole
magnets (QM) current and determining the beam size using Kakimotoa beam et al. [5]
TA
.
Y

101

laser
337

E J(E) / ( eV2 km-2 sr-1 yr-1 )E J(E) / ( eV2 km-2 sr-1 yr-1 )spot screen 15 cm downsteam of the second QM. The result of the beam spread measurement is
Calibration of TA FD with ELS
Nagano et al. [6]Data selection and background subtraction
1.0%
shown
in
Table
1.
38
38
r N2 10Integrating sphere efficiency
Geometry
Laser probe calibration
Calibration sphere transmission
Simulation of energy deposit
Monte Carlo statistics

10
1.0%
FLASH Coll. [7]Axis 0.9% ab [rad/m ] g [m/rad] e [m rad]PMTs in the two BRM telescopes that see the ELS beam. Q is defined as the beam charge [pC]
Vertical 0.3%-0.265
1.320
0.811
8.076e-06
AIRFLY Coll.
5.0%
measured by the core monitor.
this measurement
Horizontal0.8%-0.369
3.525
0.322
2.743e-06
The beam measurement by the FC is destructive, so we measure
the
Q [8]
of an ELS run by using
MACFLY
Coll.
2.0%
Table 1: Emittanceand Twiss parameterscharacterizethethenon-destructive
core thatmonitor(CM)
the beamcurrent
TAmeasuring
SDELS
Fullbeam.
Sky(ICRC2015)of each ELS shot. The CM is
1.0%
TA SD Full Sky (ICRC 2015)
thewithquenchingcross
sections(see
alsoe↵ectsLefeuvre
et al. [51]),[9]calibration
used as a relative monitor
its absolute(Esensitivity
calibrated
by dedicatedrunsthatusing are properly accounted for in Auger
rescaledby -5.2%)
5.8%
37
the FC before each ELS
run. One of theWe
results of the
FC-CMcalibration
run is shownin Figurethe
2. Kakimoto et al. and Airfly absolute
Auger SD Full Sky (ICRC 2017)
experiment.
that20%
di↵erence
between
10
Auger note
SD FullSkythe(ICRC
2017)
Waldenmaier et al. [10]AugerPrecision measurement offluorescence technique

The
Q BRM
for the FD
ELS simulation is obtainedby makingby
a separate
MC run with the FC installed
3. ELS Beam Detection by theTA+5.2%)
FYs is outside(Eofrescaled
the range
definedby the uncertainties stated by the two measurements,
3

3

Total
1037

Bokkyun Shinin the beam line (and destroying the beam at the FC).
10%
[38],
respectively.
3 [36]
4 and 53.9% 6
7 shots into
8 the air
9 at 0.5 10
ments are not statistically competitive. Nevertheless,a ELS
is composed
700 ELS
Hz. The intensity
The AF photonsthey
fromprovide
the ELSEachbeam
arerundetectedby twoofofapproximately
the BRM FD telescopes.
The total
Y337 from[photons/MeV]
useful cross-checkof the independenceofthefluorescence
yield
to
of
each
shot
is
slowlychangedduringarunusually
covering
the
range
20
pC. The
18
18.5
19
19.5
20
20.5
18
18.5
19
19.5
20pC to 200 20.5
field of view (FoV) of the two telescopes is 3 - 33 in elevation, and 18 in azimuth. A typical
the type of relativistic particle originating
the wasenergyInscatter
plotELSof SFADC
andQ
obtainedeach shot
is20
shown
in160
FigureThe value of S is obtained
M. beamAve
etal.forAstropart.Phys.
42(2013)
log ELS(E/eV)
log
data-settakendeposit.
with about
700
shots
with
chargevariedfrompC (E/eV)topC, 4. 90–102
Fig.
10.
Experimental
results
on
Y
.
For
some
experiments,
thefluorescence yield
337
10
10
deed, measurements performed withthese
secondary
beams
were
as theTheslope
ofofathelinearfit
tobandtheELS
scatterplot.map
about 20group),
minutes ofoperation.
Figure3 337(left)
isan
event
two telescopes,nm
is derived
from ofthetheintegrated
yield measured between &300 to
D. Ivanov et al. (Spectrum inworkingICRC
2017
consistent within their 6% statisticaluncertainty
the
400 nm(see text
for details).
andthe color axiswith
shows
thefluoresintegrated numberof FADCcounts.
3cence yield of Eq. (25).
× 10
Data
The different contributions to the systematic uncertainty are
Model-F
Electronlightsource(ELS)
ELS result
Model-A
experiments,
weconverttheintegrated
yield into a yield of thesummarized in Table 6. Taking these contributions as uncorrelated,
5000
Model-K30
12000
337 nm band by using the spectrum measured by AIRFLY [12]. Rea total systematic uncertainty of 5.8% on Y laser
isobtained.
337
25
10000
sults forY 337 are compared in Fig. 10. The AIRFLY measurement is
40 MeV, 109
4000
compatible with previous measurements and presents the smallest
8000
electrons,
160 pCyield measurement20
8. Combined
fluorescence
3000 It is also compatible with the detailed study by [27],
uncertainty.
6000
15
where it is argued that some of these measurements should be corThe fluorescence yield measurements presented
in
Sections
6
4000
10
2000
rectedfor a systematic bias in their calculation of the energy
ChDiff
ChMirrlaser
and 7 – Y 337 , Y 337 and Y 337 – are found to be in good agreement
2000
5deposit.
within their uncertainties. They are based on very different calibra1000
0
135
140
145
150
155a
24
25
26
27
28
29
30
31
32
tion light sources – photons from Cherenkov emission
or fromAzimuth [deg]Time [us]0nitrogen laser – whose corresponding systematic uncertainties
9. Conclusions
0 20 40 60 80 100 120 140 160 180 200
Figure
3:
Left:
Event
map
of
ELS
by
512
PMTs
which
is
made accumulated
over 700
shots in one
are largely uncorrelated. Thus, an appropriate weighted average
V. Verzi et al., PTEP. 12A103 (2017),
Beam Charge[pCELS/ Pulse]run. Right:
FADC
waveformof the (ICRC2015)
bottom center PMT.
IntegrationofE↵ect
the FADC
signalisperformed
over
7
Fig.
8
of
the
changeofthe fluorescence
yield(FY)etinreconstruction
BK
Shin
et al., PoS
640
of these measurements yields an improved
result.
T. Fujiial.,theProcof ICRC 2017 of the FD
We have performed a precise measurementof the absolute
1200

800

600

FADC countElevation [deg]1000

FADC [FADC count / Pulse]

Figure 1: Energy spectrum measurements by the Auger [8] and TA [9] surface detectors. Left:
Using energy scales of Auger and TA. Right: TA energy scale is reduced by 5.2% while Auger
energy scale is increased by 5.2%.
As pointed out in Section 1, although the TA and Auger techniques of reconstructing SD
event energies are very similar, there do exist differences in their respective instruments and themethods of how the final primary energies are assigned. The systematic uncertainty in the overall
energy scale is 14% for Auger and 21%Figure
for4:TA,
exposure
Scatter while
plot of beamthechargeuncertainties
vs. FADC for one ELS due
data-set.to
Thethedata anda fit to the data and
are
400

200

Entire Sky Spectra

Declination dependenceofδ toenergy•
Restrict[-15 ,24.8spectra] rangeo••oEntire sky of TA and Auger

Excludes TA hot spot

Inependence of exposure on
declination (aka “1/ω method”):
1
1
𝐽Common/ 𝐸 =Declination
Band
ΔΩΔ𝐸
ω(δ )
(UHECR 2016 proceedings)
Common declination band

10

Y. Tsunesada, D. Ivanov in ICRC 20178

Xmax moments

Mass composition using Xmax

(combining HeCo and FD-standard)

760

19

27

87

80

395
289

572

780

758

800

170
132

820

801

<Xmax> (g/cm2)

31

Proton dominated at 1018.4 eV,
Narrowing σ(Xmax)

740
720
700
data
QGSJet II-04 proton
QGSJet II-04 helium
QGSJet II-04 nitrogen
QGSJet II-04 iron

680
660
640
18.2

18.4

18.6

18.8

19

19.2

19.4

19.6

Xmax moments(combining HeCo and FD-standard)

19.8
log (E/eV)
10

Figure 14. Mean Xmax as a function of energy as observed by Telescope Array in BR/LR hybrid mode over
8.5 years of data collection. The numbers above the data points indicate the number of events observed.
The gray band is the systematic uncertainty of this analysis. Reconstructed Monte Carlo of four di↵erent
primary species generated using the QGSJet II-04 hadronic model are shown for comparison.

18.2 ≦ log(E) < 18.3, N = 801

18.8 ≦ log(E) < 18.9, N = 132

nitrogen, but statistics in the data there are very poor. Care must be taken in interpreting Figure 14,
since hXmax i by itself is not a robust enough measure to fully draw conclusions about UHECR
composition. When comparing hXmax i of data to Monte Carlo, in addition to detector resolution andsystematic uncertainties in the data which may hinder resolving the between di↵erent elements withrelatively similar masses, the issue of systematic uncertainties in the hadronic model used to generate
the Monte Carlo must also be recognized. This will be discussed in Section 5. Referring back to
Figures 12 and 13, we can see that though the hXmax i of the data in Figure 14, lies close to QGSJet II04 helium, the (Xmax ) of the data is larger than the helium model allows for energy bins with gooddata statistics. For this reason, we will test the agreement of data and Monte Carlo by comparing
not just hXmax i and (Xmax ), but by using the entire distributions. The elongation rate of the data
shown in Figure 14 found by performing a 2 fit to the data is found to be 56.8 ± 5.3 g/cm2 /decade.
The 2 /DOF of this fit is 10.67/9. Table 4 summarizes the observed first and second moments ofTA’s observed Xmax for all energy bins.

19.4 ≦ log(E) < 19.9, N = 19
14

5. STATISTICAL HYPOTHESIS TESTS

5.1. Method

TA collab. ApJ, 858, 76(2018)

M. Unger et al., ICRC 2017, J. Bellido et al., ICRC 2017 9
14

AugerTAs ( Xmax ) > s ( Xmax ) for lg( E/eV) = 18.5

19.0

TA’s Xmax analysis assuming Auger’s mixed composition
840

Auger

TA

(Auger2 + TA2)

sys. errors:

1/2

70

TA, ApJ (2018)

60

σ (X max) [g cm-2]

800〈 X max 〉 [g cm-2]

TA

80

AugerMix, PRD (2014) ⊗ TA

820

Auger

preliminarysys. errors:

780
760
740

50
40
30
20

720

EPOS- LHC

700
18.2

18.4

18.6

18.8
19.0
lg(E FD / [eV])

19.2

19.4

19.6

AugerMix, PRD (2014) ⊗ TA
TA, ApJ (2018)

10
0

18.2

18.4

18.6

18.8
19.0
lg(E FD / [eV])EPOS- LHC
19.2

19.4

19.6

A. Yushkov et al., (Mass Composition WG), UHECR 2018
A. Yushkov

Auger – TA working group on mass composition

11
10

own with red line, iron MC is shown with blue line.

max

780

650

max

650

760

Mass composition analysis using surface detector array
iron

600

740

600

720

18.2

550
17

17.5

18

18.5
850
log(E/eV)

iron

18.4

19

19.5

18.6

550
17

20

17.5

18

Auger FD ICRC17 (prel.) ± stat.
Auger SD ICRC17 (prel.) ± stat
Same as Fig. 15 including additional
X max measure± sys.
800

<ln A>4 Fe
3In Fig. 17, the data of Fig. 15 are compared withmeasurements made with 650
the fluorescence detectors[22].
The agreement is good; the results from the surface
detector alone are statistically 600
stronger and extend to higher
energies.

Si
N

2
1 He

1018

1017

0 p
4

19
log10 E, eV

19.5

20

Y. Zhezher et al. in ICRC 2017, arXiv: 1808.03680

21
Are cosmic rays accelerated to 10 eV (= 1 ZeV)?

8. Average atomic mass hln Ai in comparison with the
ZeVatronthe Universe?
scope Array hybrid
results in
[16].

3

log(E/eV)80

40

A. Interpretation of the measurements in termof average mass
30

A comparison with hadronic models allows the ex
20sion of P.
theSanchez-Lucasaverage depth ofet shower
maximain termal., in ICRC 2017,
the natural logarithmof the96,122003
atomic masshln 10Ai, follo
Phys.Rev.D(2017)
the procedure discussed in Sec. V. The evolution of has a function of energy is shown in Fig. 18. In0 the en

E [eV]

Fe

1019

1020

1017

Xmax ICRC 2015
Xmax 1500 m array
Xmax 750 m arrayFelines: air shower simulations using post-LHCh
4

QGSJetII-04

3

EPOS-LHC

N
2

N
2

He
1

He
1

p

0
1017

90
19.5

FIG. 17. Evolution of hX max i as a function of energy.
70 Thecompares the mean values of the X max distributions measur
the fluorescence and surface detectors of the60
Pierre A
Observatory. In most cases, the uncertainties are smaller
the size of the symbols.
50〈ln A〉

18.5〈ln A〉

18

Xmax ICRC 2015
Xmax 1500 m array
Xmax 750 m array

19

σ(Xmax) [g/cm ]5

TA SD, QGSJET II-03
TA hybrid, QGSJET II-03

18.5

2

Preliminary

FIG. 16.
ments from the surface detectors above and below 3 EeV. (Inset)
Three energy bins have been included below 3 EeV using the dataof the 1500 m array, and two measurements added above 3 EeV
750 is good agreement between
use the data of the 750 m array. Theremeasurements in the overlap region. The brackets correspond tosystematic uncertainties.
700〈Xmax〉 [g/cm2]

6

EPOS-LHCQGSJetII-04

1018

1019
E [eV]

1020

p

0
1017

1018

1019
E [eV]

1020

11

p+

!

CMB

p+⇡

0

p+

CMB

!

n+⇡

+

No GZK γ and ν at the highest energies+&

upper limits 95% CL&

1

Y 201010−1

10−5

Z-burstTD
SHDM ISHDM IIGZK proton I
GZK proton II

10−6

HP 2000
A 2002

Hy 2011

Hy 2016+syst.

Single flavour, 90% C.L.
νe : νµ : ντ = 1 : 1 : 1

e + 3⌫

Enrique Zas

ANITA-II (2010) (x 1/3)

Cosmogenic ν models
p, Fermi-LAT, E =1019 eV (Ahlers '10)min
p, Fermi-LAT, E =3 1017 eV (Ahlers '10)minp, FRII & SFR (Kampert '12)
p or mixed, SFR & GRB (Kotera '10)
Fe, FRII & SFR (Kampert '12)

Auger 1 Jan 04 - 31 Mar 17

Astrophysical sources

IceCube (2015) (x 1/3)

E2 dN/dE [ GeV cm-2 s-1 sr-1 ]upper limits 95% CL-2Integral photon flux Eγ > E0 [ km sr-1 yr -1]

UHE neutrinos at Auger

AGN ν (Murase '14)

10−7

TA 2015

10−2

10−8

10−3

10−9

SD 2015

18

10

19

10

20

10

Pierre Auger collab., JCAP04 (2017) 009

E0 [eV]

1017

1018

Eν [eV]

1019

1020

1021

Figure 2: Integral upper limit (at 90% C.L.) for a diffuse neutrino flux of UHE dN/dEn = kE 2 given as
E. Zas, Proc. of ICRC 2017a normalization, k, (straight red line), and differential upper limit (see text). Limits are quoted for a single
flavor assuming equal flavor ratios. Similar limits from ANITAII [8] and IceCube [9] are displayed alongwith prediction for several neutrino models (cosmogenic [10, 11, 12], astrophysical [13].)Top-down models are ruled out.

3.2 Limits to point-like sources of UHE neutrinos
Augerlimits become sensitive to GZK-ν and [16γ of 30]

The Earth-skimming channel is very effective at converting the tau neutrinos into exiting tau

12Energy Number(EeV) of events

Fouriercoefficient aa

Fouriercoefficient ba

AmplituderaPhaseϕa (°)

Probability
P (≥ ra)statistical fluctuations or indicate the presFLUX
MAP ABOVE
8 EeV
ence of additional
structuresat smaller
angular(a= (100°,require−24°). at
In least
galactic
coordinates,d, dd)would
scalestwice
as many Equatorial coordinates
thedirection of this dipole is (‘, b) = (233°,
events.
−13°). This dipolar pattern is clearly seen intheflux map infor
Fig.the2. Toorigin
establish
Implicationsof whether the
departuresfroma perfect
high-energy
cosmicrays dipole are merely
statisticalfluctuations
or indicate
theseenpresThe anisotropy
we have foundshould beinenceof additional
structures
angularthe context
of related
results at
atsmaller
lower energies.
scalesrequireat least oftwice
as many
Above awould
few PeV,
the steepening
the cosmic-ray
events.
energy spectrum has been interpreted as being

13

Large/intermediate scale anisotropies

þ0.0064 to 881,7010.001 ± 0.005
0.005 ± 0.005
0.005 %0.002
80 ± 60
0.60
Table
1.
First
harmonic
inright
ascension.
Data
are
from
the
Rayleighanalysis
of
the
first
.....................................................................................................................................................................................................................
≥8−0.008 for± 0.008
0.046 ±bins.
0.008
0.047 þ0.008
100 ± 10
2.6 × 10−8%0.007
harmonic
in32,187
right ascension
the two energy

Auger dipole: E > 8 EeV, 4.7% dipole with 5.2σ.....................................................................................................................................................................................................................

Energy Number
Fourier
Fourier
Amplitude
Phase
Probability(EeV) 2. Three-dimensionalof events
coefficient
coefficient bDirectionsrof
ϕa (°)(≥ ra) inTable
dipoleaareconstruction.
arePshownaa dipole components
equatorial coordinates.
4 to 881,7010.001 ± 0.005

0.005 ± 0.005

80 ± 60

0.60≥8

þ0.0060.005 %0.002

32,1870.046 ± 0.008

100 ± 10

2.6 × 10−8

Energy(EeV)

Dipolecomponent dz

0.047 þ0.008%0.007.....................................................................................................................................................................................................................

−0.008 ± 0.008

due to efficientfor
escapeofparticlesfrom the galImplications
theoriginof
Galactic centeraxy
and/orbecauseof theinability of the sources
high-energy cosmic
rays
to accelerate cosmic rays beyond a maximum
þ0.007þ17The
anisotropy
we
have
found
shouldbeseenin357,AllSky
Surveywith
TA&PAO
4 to 8−0.024 ± 0.009
0.006%0.003
0.025þ0.010−75%8
80 ± 60
value
of
E/Z.
The
originoftheparticles
remains%0.007
PierreAugercollab.
Science1266 (2017)
Table
2.
Three-dimensional
dipolereconstruction.
Directionsofdipolecomponents
are
shown
in
.....................................................................................................................................................................................................................
Oversamplingwith energies.
20°-radius circle
thecontextAlthoughof related
results
at lowerþ0.011
þ0.013
þ12
unknown.
supernova
remnants
are often≥8−0.026
±
0.015
0.060
0.065
−24
100
±
10
%0.010%0.009%13equatorial
coordinates.
.....................................................................................................................................................................................................................
Above
a few
the steepeningof thecosmic-ray
discussed
asPeV,
sources,
evidence has
beenreportedenergyspectrumhasgalactic
been interpreted
as being
for a sourcein thecenter capableof✦ TAmodulation
due to efficientescapetoofPeV
particlesfrom
thegalHotspot:
Eright> 57
EeV,Dipole[TA collab.
790:L21(2014)]
accelerating
particles
energies(31).
Diffuof the azimuthof ApJL,
the arrival
direction
ascensionof3.4σ100° ±anisotropy
10°. Dipole
DipoleEnergy is atDipole
Dipolerightof the
axy and/or
becauseofgalaxy
the inabilityof thesourcessive
escapefrom
theisexpected
tolead toshowerrather
than
in
terms
of
the
right
asThe(EeV)
maximumof the modulation
for the 4 dEeV<
componentdz
component
amplitudeddeclinationd(°)
ascension
a
(°)
⊥
d
d
todipolar
acceleratecosmic withrays abeyond
a maximumacomponent
maximumnear thecension.
The
results
are
b
E
<
8
EeV
bin,
at
80°
±
60°,
is
compatible
with=−0.013±0.005
inϕ✦ TA (7þ0.00757 EeV)+
þ0.010 Auger(10þ17 years, 157 events
years,
109
eventsabove
to 8−0.024 ± 0.009
0.006%0.003± 60
value of center
E/Z. The
origin oftheThisparticles
remains%0.007galactic
direction(32).
is compatible
the4one
determinedin the higher-energy
bin,0.025the4 EeV < E −75< 8%8
EeV bin and b80.....................................................................................................................................................................................................................
ϕ = −0.014 ±þ0.011
þ0.013
þ12
aboveoversampling15
18 are often
unknown.
Althoughsupernova
remnants≥8 57itEeV),−0.026±20°
0.015circle
0.060
0.0650.008
100probabilities± 10%0.010%0.009 in the −24
%13
to
10
eV rangewithresults
obtainedinthe10
althoughhas
highuncertainty
andtheampliE
≥
8
EeV
bin.
The
.....................................................................................................................................................................................................................
discussed
sources,
evidence
has been(15,
16, 23,as24,33), which
providevaluesreported
for thetude is not statistically significant. Table S1 showsthat larger or equal absolute values for bϕ arise
for a sourcein ascension
the galactic
centercapableofphase
in right
close
to thatof theresults
obtainedunder
the stricter
triggerfrom anisotropic distribution are 0.8% and
✦that
E
>
57
EeV,
no
excess
from
the
Virgo
cluster
Northern
TA :Diffu7 years 109 events (>57EeV)
accelerating
particles(31).
ofazimuth of the arrival direction of themodulation
is atascension eventsof 100°gained± 10°.
No
correction
for PeV energies
= to
266°.
galactic center,
aGC K.
condition andforright
the additional
8%,therespectively.
Southern
Auger: 10ofyears
157 events(>57EeV)
Kawata
et
al.,
Proc.
ICRC
2015
E 2018scale February
difference
The
AstrophysicalJournalLetters,
853:L29(10pp),
1
sive
escapefromthegalaxyis
expectedtoleadtoshowerthan in terms
thedominant
right asTheof themodulation
for the4 EeV <
Models proposing
a galactic
originupistotheaftermaximum
relaxing thetriggerare entirely
consistent
Underrather
the assumption
thatofthe
b/w TA and PAO!!
Southern
hotspotseen
at Cen A(Pre-trial ~3.6σ)adipolarcomponent
withamaximumnearthecension.
The
resultsarebE<pattern
8each
EeV other.
bin,correlation
at 80° ± 60°, iscompatible
with ApJL,± 0.005in
12
ϕ = −0.013
highest observed energies (34, 35) are in increasing
withcosmic-ray
based onprevious✦ Flux[PierreAuger collab.
853:L29anisotropy(2018)] is dipolar,
galactic
centerdirection (32).
This
is compatible
theFigure
one determinedin the higher-energy
bin,
the
4 EeV
E < 8 that
EeVthebineffects
and bofϕhigher-order= −0.014 ±tension withobservations.
If thegalactic
sources
1 shows the distributionof the normalstudies
that< found
15
18
10 aboveeV rangewithresultsto
obtainedin the10 toalthoughithas
highabove
uncertainty
andthe ampli0.008
in theE≥ significant
8 EeV bin.inThe
probabilities
postulated
acceleratecosmic
rays
EeV
ized
rate
of
events8EeV
asafunctionofmultipoles
are
not
this
energy
range
✦tude
With
a
flux
patternofstarburst
galaxies,
isotropyof(15,
16, 23,such
24, 33),
provide values
for thenot statisticallysignificant.function
Table S1 showsthatlargerequal
absolute
valuesandfor itsbϕ direcarise
energies,
as whichshort gamma-ray
burstsorrightisascension.
The sinusoidal
corre(25, 29,30),orthedipolecomponents
UHECR
is disfavored
with
4.0σ
above coordinates39
phase
in right
ascension
close
the
that
results
obtainedunder thetriggerfromisotropic
distribution
are
0.8% andhypernovae,
weredistributed
in tothethatdiskofof thetion inanequatorial(aspondsto thefirst
harmonic;
thestricter
distribution
isconfidenced, dd) can be= 266°.
galacticaGCcomponent
condition
foratheadditional
events gained8%,
respectively.
EeV andgalaxy, acenter,
dipolarof anisotropy isestimatedfromcompatible
withdipolarmodulation:
c2/n =
Models withproposinga galactic
upexisting
to theafter
relaxing
trigger are entirely
consistent
Under the assumption that the dominantpredicted
an amplitude
thatorigin
exceeds
10.5/10for thethe
first-harmonic
curve andc2/n =
highestenergies(24,(34,33).
35)Inarethisin increasing
withcosmic-ray
anisotropy
is dipolar, based on previous
bounds observed
at EeV energies
sense, the45/12
for aother.
constant
functionfraction(where n isandthe 12.9°✦ each9.7%
anisotropic
angularscale
tensionwith
observations.
thedipolegalactic
sourcesFigureof
1 showsof thenormalstudies that found that thereffects
of higher-orderaconstraint
obtainedhere onIftheamplitudenumberdegreestheof distribution
freedom, equal
to thenumd⊥ ≈
postulated
to4 acceleratecosmicrays above
EeVa
ized
rate
of events
above
8minus
EeV asa functionofmultipoles are not significant
in this energy range
hcos di(Table
2)
for
EeV
<
E
<
8
EeV
further
disfavorsberofpoints
in
the
plot
the
numberof✦right
The
other
three
flux patterns:
2.7σ–3.2σenergies,
such galactic
as short
gamma-ray
burstsor
ascension.
Thesinusoidal
function corre(25, 29, 30), the dipole components and its direcpredominantly
origin.
This tension
could
parametersof thefit).
hypernovae,
distributed
thetion in equatorial coordinates(ad, dd) can bespondsto the first harmonic;
theequatorial
distributionisbϕ
be alleviatedwere
if cosmic
rays atina thefewdisk
EeVofwere
The distributionof events incoor.....................................................................................................................................................................................................................

Dipolecomponent d⊥

Dipoleamplitude dDipoledeclination dd (°)

Dipole rightascension ad (°)

Significance(pre-trial)

Downloaded from http://science.sciencemag.org/
on September 22,
Downloaded from http://science.sciencemag.org

9

Aab e

Large scale anisotropy above 4 EeV

8−0.007±0.008−0.002±0.008

12−0.014±0.008
0.006±0.008

0.38
0.80

0.08
0.45

14

Figure 3. Evolution with energy of the amplitude (left panel) and direction (right panel) of the 3D dipole determined in different energy bins above 4EeV. In the sky
map in Galactic coordinates of the right panel the dots represent the direction toward the galaxies in the 2MRS catalog that lie within 100Mpc, and the cross indicatesthe direction toward the ﬂux-weighted dipole inferred from that catalog.

energies above 32EeV. Table 4 reports the results for thecorresponding azimuth analysis in these new energy bins.
3.2. Reconstruction of the CR Dipole
We now convert the harmonic coefﬁcients in R.A. and in
azimuth into anisotropy parameters on the sphere, assumingﬁrst that the dominant component of the anisotropy is the
dipole d . The ﬂux distribution can then be parameterized as afunction of the CR arrival direction û as' (uˆ)  '0 (1The AstrophysicalJournal, 868:4 (12pp), 2018 November 20
ovember
20

Aab et al.

Aab et al.

arrival directions,forofthefour energy
for 10 degrees
freedom)and ﬁrstbins
plusabove
second4EeV.
harmonics(χ /dof=1.02,
2(χ /dof=0.44,
for 8 degrees
of freedom),
respectively.
832,1871−0.008±0.008−0.008±0.008
0.046±0.008
0.047
100±10 0.046±0.008
2.6×10−8

arrival directions, for the four energy bins above 4EeV.

d · uˆ).

(4 )In this case, the amplitude of the dipole component along theTable 1
Table 1
rotation axis of Earth, dz, that in the equatorial plane, d⊥, andResults of the First- and Second-harmonic Analyses in R.A.
sults of the First- and Second-harmonic Analyses in R.A.
the
R.A. and decl.
of the dipole Bdirection, (αd, δd), are related to
B
B
B
B
B
B
B
B
B
Figure1.
Distribution
in
R.A.
of
the
normalized
ratesof
eventswithenergyK k (n)Energy (EeV)
Events bkk
K k (n)ak
bk P (.rk )rkP (.rk )akrk
the ﬁrst-harmonic amplitudes in R.A. and azimuth throughabove 8EeV. The black solid and blue dashed lines show the distributions
4–80.001±0.005
81,701
1 in Galactic
0.001±0.005
0.005±0.00580±60 events in windows
0.60
0.005±0.005
0.005
80±60
0.60
Mapscoordinatesof thebetweenthe0.005
numberof observedof 45° and those expected for an isotropic distribution ofobtained Figure
from the4.
weighted
Fourier analysis
corresponding
to a ratio
ﬁrst harmonic(The
PierreAuger
Collaboration
2015b)−0.001±0.005 70±80 0.001±0.005
0.002
70±80coordinates of the 0.94
−0.001±0.005 2
0.001±0.005 2
0.002
0.94
Figure 4. Maps in Galactic
ratio between the number of observed events in windows of 45° and those expected for an isotropic distribution of
0.013±0.008

0.012±0.008 2

0.013±0.008 21±12
0.018

0.012±0.008
0.065

0.047
0.018

100±1021±12dz −8
2.6×10
b1G
0.065

,

Table 5
cos
ℓ
sin
R

§
Table 5obs
Three-dimensional
Dipole Reconstruction for Energies above 4EeV
0.009. Regarding the second harmonic, none of the amplitudes
B
Three-dimensional Dipole
Reconstruction
for Energies
above 4EeVr
Energy
(EeV)
1
found
are
signiﬁcantly
different
from
zero.
d⊥ ,
dz
d
αd (deg)δd (deg)
Table 2
Table 2d?
Medianlts of the First- and Second-harmonic
Analyses thein Azimuth
Resultsof theFirstand Second-harmonic
Analysesin AzimuthInterval
Figure 2 displays
maps,
inequatorial
coordinates,ofthecos E§
Energy (EeV)
0.007
0.010
17
d
d
d
α
(deg)
δ
(deg)
4–8
5.0
−0.024±0.00980±60
0.006
0.025

75
⊥
z
d
d

0.003

0.007

8
B
G
G
insideEnergykakG (EeV) exposure-weighted
bkG average of the
P (.
∣akG∣) a top-hat window
akG ﬂuxbPkG(.∣bkG∣)
P
(
.
∣
a
∣)
P
(
.
∣
b
∣)
,
B
K

0.013
0.011
12
k
k
d
8
11.5
−0.026±0.015100±10
0.0650.009
2413
1 0.0600.010
Interval
Medianof radius 45°,
appreciate0.045
the large-scalefeatures,−0.010±0.005−0.013±0.005
0.009
4–8
1 so as to better−0.010±0.005−0.013±0.005
0.045
0.0090.013 ⎞
0.015
16
8–16
10.3
−0.008±0.017104±11
0.058⎛0.0590.008
816
d
0.007
0.010
17
0.011
z
0.002±0.005 for the 4–8−0.002±0.005
0.69 inenergyand E0.698 EeV. An
excess2 bins [4, 8]EeV
0.002±0.005−0.002±0.005
0.69
0.69 0.025 0.025
5.0−0.024±0.009
0.006
75
15
0.03
0.003
Ed  arctan(580±60) 0.100.02⎟ , 0.007−0.08±0.03
501416–32
20.2
82±20 80.065⎜
0.018
0.013
0.011
12
⎝
⎠
0.07
28
0.05
d
the ﬂux fromdirections
the predominantfeature?
32
39.5−0.08±0.07
115±35
0.110.03
46260.08
8 −0.014±0.008
11.5 is0.38−0.026±0.015100±10
0.065
0.060
2413−0.007±0.008
0.08
8
1 the southern −0.007±0.008−0.014±0.008
0.38
0.08
0.03
0.009
0.010−0.002±0.008 at energies 2between
0.006±0.008
0.45−0.002±0.008
0.006±0.008
0.80
0.45
4 and
8EeV, 0.80
while above 8EeV
the
0.013
0.015
16 8EeV, as well
Note.
We show
the results
for theistwothebins previously
reported(The
PierreAuger
Collaboration
2017a), i.e., between 4 and 8EeV and above
where
mean
cosine
of
the
decl.
of
the
 cos
E § obtained
0.781410.3 R.A. close to 0.058−0.008±0.017104±110.059

8
excess 8–16
comes from a region with
100°,0.011
with a0.00816
as dividing the high-energy range into three bins.
15
0.025
events,is forthethemean
sine0.03of the event zenith
4sinR§ the0.6525energiesabove
32EeV.
Table
4reports
theresultsforthecorresponding
deﬁcitintheopposite
direction,inaccordance50
energiesabove
32EeV.
Table
reports
results
16–32
20.2−0.08±0.03
82±20
0.10
0.065
14
0.02
0.018
ApJ, 868:4
E 8 EeV.
The
only
coefﬁcientsare present. They are also consistent with results obtainedin pastℓobsnew
 nonvanishing35 n. 2bins.
angles,andis correlation
the latitudeof the observatory.
corresponding
azimuthanalysisin these
new corresponding
energy
bins.
0.07
28
0.05 andazimuthandanalysisintheseenergywiththeresultsfromtheharmonic
analyses
in
R.A.
32
39.5
−0.08±0.07
115±35
0.11

46
0.08
,0.0326
0.03(2018)S(dQ)
betweenthequantities
reported
in
Table
7
are
analysesin
The
PierreAuger
Collaboration(2015b)and Thexxz
Note that, as is well known, when the coverage of the sky is not
azimuth.
S (d y, Qyz )  0.63 and S (dz, Qzz )  0.91. The nine components
Pierre Auger & Telescope Array Collaborations (2014). Note that

3.2. Reconstruction of the CR Dipole

3.2. Reconstruction of the CR Dipole

Large scale anisotropy above 4 EeV

Prediction of dipole amplitude using
Aab et al.
Auger mixed composition

The Astrophysical Journal, 868:4 (12pp), 2018 November 20

Power law index: 0.79±0.19
Constant amplitude disfavoredwith 3.7σ

015

Aab et al.

Figure 5. Map in Galactic coordinates of the direction of the dipolar component of
the ﬂux for different particle rigidities for CRs coming from Galactic sources andpropagating in the Galactic magnetic ﬁeld model of Jansson & Farrar (2012) (bluepoints) and the bisymmetric model of Pshirkov et al. (2011) (red points). The
points show the results for the following rigidities: 64, 32, 16, 8, 4 and 2 EV (withincreasing distance from the Galactic center). We also show in purple the observed
direction of the dipole for E 8 EeV and the 68% CL region for it. The
background in gray indicates the integrated matter density proﬁle assumed for the
Galactic source distribution (Weber & de Boer 2010).

Figure 6. Comparison of the dipole amplitude as a function of energy withpredictions from models (Harari et al. 2015) with mixed composition and asource density S  104 Mpc3. CRs are propagated in an isotropic turbulent
extragalactic magnetic ﬁeld with rms amplitude of 1nG and a Kolmogorov
spectrum with coherence length equal to 1Mpc (with the results having only
5. Discussion
mild dependence on the magnetic ﬁeld strength adopted). The gray line
indicates the mean value for simulations with uniformly distributed sources,
The most signiﬁcant anisotropy in the distribution of CR
while the blue one shows the mean value for realizations with sourcesin determined
the studies performed
above
4EeV
the largedirection (right panel) of theobserved
3D dipolein different
energybinsisabove
4EeV. Indistributed
the sky as the galaxies in the 2MRS catalog. The bands represent thescaleindipolar
modulation
ofthatthelieﬂux
at energies
above
8EeV.
dispersion for different realizations of the source distribution. The steps
e direction toward the galaxies
the
2MRS
catalog
within100Mpc,
and
the
cross
indicates
Augercollab., ApJ,
(2018)lies in Galactic coordinates atobserved reﬂect the rigidity cutoff of the different mass components.
The maximumof this868:4
modulation

catalog.

(l, b)=(233°,−13°), with an uncertainty of about 15°. This is0

t1/2 azimuthal asymmetry - Results
19

10

E [eV]

8

max

µ
µ
 Xmax Xmax


6

4

2

2

0

0
p
p

0

18

19

10

18

20

10

19

4

2

500-1000 m

X8max ICRCX2015
max ICRC 2015

 Xmax ∆s 1500 m∆arrays 1500 m arraye)array∆s 750 m array∆(secs 750 m
max

6

3

20

Fe

µ
 Xmax 

Fe

Fe

4

N

0

10

Fe 3
2

QGSJetII-04
QGSJetII-04
p

18

19

10

10

30°<  θ  <  45°150°<  |Φ|<  180°

10

18

10

18

10

19

10

19

10

20

p

10

E [eV]
20

10
10E [eV] E [eV]201018.8 eV < E < 1019.2 eV
30°<  θ  <  45°150°<  |Φ|<  180°
proton iron
stat.
stat. + syst.

between 500-1000m
and 1000m-2000m for
QGSJet II-04

19

20

10

10

E [eV]
20 / 24

Data

Fe

Data/MC proton[1910, 2160]

NN
1.72±0.10(stat.)±0.40(syst.)[2760, 3120]

3.14±0.36(stat.)±0.72(syst.)

p
18

10

17

10

0

19

10
17

10

18

10

18

10

EPOS-LHC
EPOS-LHC
p
20

10

19

10

E [eV]20

19

10

4500

2000

R. Takeishi, ISHVECRI
2018, Ph.D. thesis (2017)
Figure 4.25: Correlation between the muon purity and the signal size ratio

Data/MC iron

He

1

0

4500 1500

R (m)

He
0

p

Data/MC ratio• inconsistency

1018.8 eV < E < 1019.2 eV

stat.
QGSJETII-04
mstat. +1000-2000
syst.
Xmax ICRCX2015
2015
58
max
MCICRCQGSJETII-03
proton5 Xmax∆s 1500 m∆array
1500 QGSJETII-03
m array
s MC
iron
(sec
)max
∆s e
750
m array∆s 750 m array
2000
Fe
46
4 µ 1500
 Xmax 

1

17

Lateral distribution Xmax

He

1

0

• Data/MC comparison assumingiron primarycosmic rays
• incompatible
with

2

He

17

10

2
2

0

p

4

2

1

Fe

3

N

Resultsµ Xmax 

18

10

 lnA 

5

4

3

10

QGSJETII-04

 lnA 〈ln A〉〈ln A〉

5

10

1000-2000 m(sece)max

2

10

 Xmax4

h i method - Results E [eV]
E [eV]

E [eV]

8

6

Fe Fe

4

20



(sec
e)max
(sec
e)max

6

10

 XmaxXEPOS-LHC

Data/MC ratio

8

1000-2000m m
500-1000

 lnA QGSJETII-04
EPOS-LHC

 lnA  lnA 1000 m20

10

Average number of particles

1810〈ln A〉

E [eV]

Muon deficits in hadron interaction modelsp〈ln A〉

20

10

and 1000m-2000m for
QGSJet II-04

p

10
1020
E [eV] E [eV]to the MC for 2000 m < R < 4000 m. The black, red, green, blue, yellow an
1.26±0.07(stat.)±0.29(syst.)
points represent |φ| < 30◦ , 30◦ < |φ| < 60◦ , 60◦ < |φ| < 90◦ , 90◦ < |φ120◦ < |φ|16< 150◦ , 150◦ < |φ| < 180◦ , respectively. The open circle, filled
1.74±0.19(stat.)±0.40(syst.)
cross represent θ < 30◦ , 30◦ < θ < 45◦ and 45◦ < θ < 55◦ , respectively. T
thin error bars and shaded thick error bars represent the statistical errors andsum of statistical and systematic errors, respectively.

In higher muon purity, larger deficits in models

Need more muon production in the simulation.
20 / 24

17.5 eV)• incompatible
withXcomposition
(E
>
10
max
P. Sanchez-Lucas, ICRC 2017, Pierre Auger collab., Phys.Rev.D 93,072006 (2016), Pierre Auger collab., Phys.Rev.D 96,122003 (2017)16

Matter distribution in the Universe
50 Mpc:0 - 50 Mpc

Predicted dipole amplitude
0.6

p inj., no cut, PT11
p inj., no cut, JF12
O inj., no cut, PT11
O inj., no cut, JF12
Si inj., w. cut, PT11
Si inj., w. cut, JF12

dipole amplitude

0.5
Matter distribution inGalactic
the Universe
coordinates
50

50 - 100 Mpc

dots are 2MRS galaxies

100 Mpc:

0.4
0.3
0.2

Matter distribution inGalactic
the Universe
coordinates

owers

100 - 250 Mpc

0dots are 2MRS galaxies

250 Mpc:

EM Signal

Muon Signal<Xmax>

Pierre Auger Collaboration / Physics Letters B 762 (2016) 288–295

Auger 2017
Auger + TA 2015
99.9% sensitivity

0.1
100

Correlation measurement292

dots are 2MRS galaxies

20
17

0

20
20

inyakov

Hadron-interaction-model independent measurements

0

10

20

Nuclear Interactions

30
40
Emin [EeV]

50

60

70

May 2018
T. Pierog,flux
ISVHECRI 2018
dRI – UHECR
directional

Pierre Auger Collaboration / Physics Lette∗ for lg( E /eV) = 18.5–19.0. Right: the same distribution for 1000 proton and 1000 iron showers simulated with EPOS-LHC.
∗
Fig. 1. Left: measured X max
vs. S 38

A. di Matteo and P. Tinyakov, https://doi.org/10.1093/mnras/sty277
very similar elongation rate (slope)for allFigure6.models
The magnitude of the dipole
as atofunctionofthe anden-zenith angle. This way weare scaleda referencerateenergyElongation(2MASS
same mass composition evolution
avoid a galaxies)
decorrelation
betweenobservables from combining difergy threshold
Emin forRedshift
the threesurvey
injectionmodels andtwotheGMF
still differences in absolute values
ferent energies and zenith angles in the data set. S (1000) is scaled
models we considered. The points labelled“Auger + TA 2015”
to 38◦ and 10 EeV using the parameterizations from [16]. X max+/- 20g/cm2 is a realistic uncertainty bandand “Auger 2017” show the dipole magnitude
reportedis scaled to 10
EeV using inan Deligny
elongation rate d⟨ X max ⟩/d lg( E /eV) =2
58 g cm−The/decade,an averagewith little variation between(2015) and Taborda (2017) respectively.
dottedlinesvalue
show
primaries and interaction models [9]. Here, these scaled
the 99.9% C.L. detection thresholdsdifferent
using
the current and ∗near∗ . Thus, X ∗∗quantities will be denoted as X max
and S 38
max and S 38
future Auger and TA exposures (seearethe
the text
valuesforof Xdetails).
max and S (1000) one would have observed, had
Galactic coordinates

292

Elongation rate

the shower arrived at 38◦ and 10 EeV. It should be noted that the
specific choice of the reference values is irrelevant, since a transformation to another reference value shifts the data set as a whole,leaving the correlation coeﬃcient invariant.
∗ the∗
As a measure of the correlation between X max
and S 38
∗ , S ∗ ) introduced by Gideon and Holranking coeﬃcient rG ( X max
38
lister [17] is taken. Conclusions are unchanged when using otherdefinitions of correlation coeﬃcients, including the coeﬃcients ofPearson or Spearman, or other ones [18]. As for any ranking coefficient, the rG value is invariant against any modifications leaving
T. Pierog, KIT
J.- 8/30
Bellido et al., ICRC
2017
the ranksof events unchanged (in particular to systematic shifts

Preliminary

fiducial trigger cut (the station with the highest signal should haveat least 5 active neighbor stations). This requirement does not introduce a mass composition bias since in the energy and zenith
ranges considered the surface detector is fully eﬃcient to hadronic
primaries [19,20]. Selecting energies of lg( E /eV) = 18.5–19.0 andzenith angles <65◦ , the final data set contains 1376 events. The
resolution and systematic uncertainties are about 8% and 14% inprimary energy [21], <20 g cm−2 and 10 g cm−2 in X max [9], and<12% and 5% [22] in S (1000), respectively.
The simulations were performed with CORSIKA [23], using
EPOS-LHC, QGSJetII-04 or Sibyll 2.1 as the high-energy hadronicinteraction model, and FLUKA [24] as the low-energy model. All
events passed the full detector simulation and reconstruction [25]with the same cuts as applied to data. For each of the interaction
models the shower library contains at least 10000 showers for proprimaries
and 5000–10000
showers
each for helium, oxygen
ton Fig.
∗ for∗
1. Left: measured X max
vs. S 38lg( E /eV) = 18.5–19.0. Right: the same distribu
and iron nuclei.

Mixed composition at the ankle
are scaled to a reference energy and zenith angle. This way we
4. Results
PierreAuger collab.,
B 762,
288 (2016)
avoid
a decorrelation
betweenPhys.Lett.
the observablesfrom combining
dif-

fidu
at l

17

10

850

90

Preliminary

Auger FD ICRC17 (prel.) ± stat.
Auger SD ICRC17 (prel.) ± stat± sys.

1
TA ICRC 2015

80
70

σ(Xmax) [g/cm2]

800〈Xmax〉 [g/cm2]

Data points fromI. Valino et al.,Proc. of ICRC
2015, D. Ivanov
at al., Proc. ofICRC 2015

Flux × E3/1024 (eV2 m-2 sr-1 s-1)

An overview of energy spectrum, mass composition and anisotropy

750

700
650

Auger ICRC 2015
-11017.5600

18FLUX MAP ABOVE 8 EeV

18.5
19
19.5
log (E (eV))

E > 8 EeV

20

20.5

18

1017

10

10

Equatorial coordinates
Pierre Auger Collab. Science 357, 1266 (2017)18

E > 57 EeV

19

E [eV]

20

10

50
40

30
P. Sanchez-Lucas et al.,20 of ICRC 2017
Proc
Pierre Auger collab.,
10
Phys.Rev.D
96,122003(2017)
0
18
1017
10

E > 100 EeVlines: air shower simulations
using post-LHC hadronic inte

All Sky Survey with TA&PAO

K. Kawata et al., Proc. of ICRC 2015

10

60

Small-scale anisotropy

S. Troitsky et al., Proc. of ICRC 2017

Oversampling with 20°-radius circle

Doublet('T=1.31o)TA

Triplet? orDoublet('T=1.35o)

Galactic center

Increase dipole amplitude above 4 EeV

No correction for
E scale difference
b/w TA and PAO !!

Northern TA : 7 years 109 events (>57EeV)
Southern Auger : 10 years 157 events (>57EeV)Aug

2 doublets above 100 EeV.
Southern hotspot isseenultrahighat Cen A(Pre-trial ~3.6σ)⇒ Needmore
statisticofenergy cosmic
rays(UHECRs)9Æ the probability to have 2 doublets
12

s,
e
2
-

2 ! rms. The results of this analysis applied to both thecontributions
fr
19
HiRes data and to QGSJET-II proton and iron Monte Carlotutions. The coo
calculations are shown in Fig. 4. The HiRes Xmax width
P H Y S I C A L R E V Iand
E WG. LOlsen,
ETTE
PRL
104,
161101(2010)
data are consistent with the predictions of QGSJET-II proGround staff istons. The width of the Xmax distribution of protons within
10
Note added.—unchanging
com
2
HiRes-2 Monocular
the QGSJET01 model tends to be about 5 g=cm broader
M.
Takedaetal.,collaboration
drstrongconstrain
HiRes-1 Monocular2than that of QGSJET-II, while SIBYLL protons are 2–3 g=cm
AGASA
PRL 81, 1163
sition of thehig
resulta trans
narrower than those of QGSJET-II. Both of these shifts are(1998)
lactic cosmic rasmall compared with statistical uncertainties, particularly
Of the hadrat the highest energies.
agreementisIn summary, the HiRes data are consistent with a*CorrespondWithincurrent
R. Abbasi et al.,constant
elongation
rate
of
47:9
"
6:0ðstatÞ
"[1] R. U. Abb1
PRL 100, 101101
sistent
with this
3:2ðsystÞ g=cm2 =decade above 1.6 EeV, and thus with an(2008).
(2008)
Comparison
wit
[2]
J.
Abraham
R. Abbasi etlightelements.
[3]
R.
U. Abbaal., PRL 104,
wouldth
J.W. Cronin / Nuclear Physics B (Proc. Suppl.) 138 (2005) 465–491
485[4] (2010)
K.imply
Greisen161101
17
17.5
18
18.5
19
19.5
20
20.5
21
steadilychangin[5] G. T.
ZatsepFIG. 3. HiRes stereo hXmax i compared with the predictions forlog10(E (eV))
ningenergy[6] theV. Berezins
QGSJET01, QGSJET-II, and SIBYLL protons and iron after full
TheMerida,presentM
detector simulation. The number of events in each energy bin
de Mexico,
spectralresults
is displayed below the data point.
FIG.
3(coloronline).
The
cosmic-ray
energyspectrum
meaE>
40
EeV
[7] W. Heitler,
high-energy
cut
sured by the HiRes detectors operating in monocular mode. The
J. Cronin,
University
in
which
cosmi
spectrum of the HiRes-I and HiRes-II detectors are shown. The
Nucl.Phys.Proc.Suppl. tion of the Xmax width. The uncorrected rms and sample[8] J. Matthew
galactic origin,ofhighest two energy bins for each detector are empty, with the
138:465 (2005)[9] J. H. Boyer
standard deviations are biased estimators of the width [30]
shapedby 457
int
A
482,
II
68% confidence level bounds shown. The spectrum of theand tend to be subject to large fluctuations in distributions
background.
[10] R. U. Abba
AGASA experiment is also shown [7,8].
with broad tails.
ThisJ. Girard
work is[11]
et
In order to focus attention on the center of the Xmax
9321949,
PHY
460, 278(2
distribution
and ofreducesensitivity
tofluctuationsin theFIG.
4.
Results
fitting
HiRes
stereodata
Xdistribution
to[12] R. U. Abba
PHY-0098826,
max
01101-3
tails,
the widthis quantified
asthe width!Superimposed
X of a unbinned
Gaussian
truncated
at
2
!
rms
(blackpoints).
are[13] R. U. Abba
0305516,
PHY✦
likelihood
to a Gaussianof Ia proton
distribution
truncatedat(squares)and iron
expectationsfit basedon QGSJET-II[14] FG03-92E
J. Belz et aNo.
2(triangles)! rms. The
results
this analysis
appliedto points
both theMonte
Carloofcalculations.
Monte
Carlo
are
[15] F. Kakimot
contributions
fr
✦
HiRes
QGSJET-II
protonand ironMonte Carloshown datawith andsmalltooffsets
in energy
to provide
separation.
Sect. A 372
tutions.
The coocalculations are shown in Fig. 4. The HiRes Xmax widthand G. Olsen,
data are consistent with the predictions of QGSJET-II pro✦
161101-4
Ground staff istons. The width of the Xmax distribution of protons within
Note added.—2
the QGSJET01 model tends to be about 5 g=cm broader
3

Flux*E /10

24

2

-2

-1

-1

(eV m s sr )AGASA/HiRes results (2003)

Super-GZK cosmic rays?
Proton dominated?
Small-scale anisotropy?

Due to a significantly improved statistics of ultrahigh-energy cosmic rays
(UHECRs) by the Telescope Array Experiment and Pierre Auger Observatory,
we firmly confirmed a suppression of the energy spectrum at the highest energies
and observed intriguing large-scale anisotropies in arrival directions of UHECRs.
We also encountered a gradually transition to a heavier composition, a deficit of
the number of muons in simulations and a lack of desired small-scale anisotropies.
In this talk, I highlight recent results of the two observatories including on-going
updates and then address scientific goals and requirements for future UHECR
observatories in next decade. I introduce three ideas as a personal decadal survey:① a fine-pixel fluorescence telescope for low-energy extension, ② a layered
Water-Cherenkov detector array for sub-EeV anisotropy and ultrahigh-energy
photon search, and ③ a low-cost fluorescence telescope array suitable formeasuring the properties of UHECRs with an unprecedented aperture.
20

21

Expected results from ongoing upgrades and smoking gun
0.5

Primary (3 Mpc)
3 Mpc
20 Mpc40 Mpc75 Mpc002-2log(Eγ φγ [erg cm ])✦ Expected (optimistic) results in next 5 years-2✦ TA×4: confirmationof the TA hotspot with >5σ-4-6 indication of small scale anisotropies selecting a✦ AugerPrime:

20 eV.
light composition,
protonfractionat
10-8-10 model: < 20 g/cm2 uncertainty on Xmax at 1020 eV✦ Interaction-12log(Detection Probability)

2

3

-0.5
-1
-1.5

w/-2 cosmic radio
background(pessimistic)-2.5
A= 3000 km2

only CMB(optimistic)
A= 3000 km2

-3
6

7

8

9
10
log(Eγ [GeV])

11

12✦ Smoking gun of cosmic ray origin ⇒ γ-rays detection spacialcoincidence with UHECR hotspot

FIG. 2: Energy fluences of UHE photons from a LL GRB-like
iso
50.5
UHECR
burst with
ẼHECR10
erg
for each distance. The✦ γ-ray:limitedsources
in=nearby
universe
primary γ-ray spectrum is also shown (see the caption of Fig.
for the ofsourceparameters).
Thick 1-100
lines show
the non-CRB✦1 bursts
γ-raysabove 10 EeV,
events−13case while thin lines show the CRB case, with BIG = 10G.
2 :expected✦ 3000The
burst kmrates
each distanceare 1/28000,
25 - 40 within
Mpc (TA×4,
Auger)−1
1/94, 1/12, and 1/1.8 yr .
✦ 30000 km2: 40 - 80 Mpc✦ 300 km2

: 5 - 10 Mpc

tons. These boosted photons can create pairs as long as
they are energetic, and the process repeats itself until

10

20

30

40

50 60
D [Mpc]

70

80

90

100

FIG. 3: The comparison of Poisson probabilities to detect
UHE (> 10 EeV) photons and high-energy (> 10 PeV) neutrinos from a LL GRB-like UHECR burst. For UHE photons,
2
2
A = 3000 km without the CRB (solid lines), A = 3000 km
52
with the CRB (dashed lines), A = 3 × 10 km without the
52
CRB (dotted-dashed lines), and A = 3 × 10 km with the2
CRB (dotted lines). For neutrinos, A = 1 km (double-dashedlines), assuming IceCube-like detectors. Thick and thin lines
iso
50.5
iso
51
are for ẼHECR = 10
erg and ẼHECR = 10 erg, respectively.
K. Murase, Phys.Rev.Lett. 103 (2009) 081102

Addressing the intermediate anisotropies

22

E>57 EeVAll SkyTA+Auger
Surveyatwith
TA&PAOAab et al.

SBG flux map, Φ=12.9°

The Astrophysical Journal Letters, 853:L29 (10pp), 2018 February 1

Oversampling with 20°-radius circle

Table 1
Populations Investigated
SBGs

l (n)

b (n)

Distancea (Mpc)

Flux Weight (%)

Attenuated Weight: A/B/C (%)% Contributionb: A/B/C (%)

NGC 253
97.4
−88
2.7
13.6
20.7/18.0/16.6
35.9/32.2/30.2
M82
141.4
40.6
3.6
18.6
24.0/22.3/21.4
0.2/0.1/0.1
NGC 4945
305.3
13.3
4
16
19.2/18.3/17.9
39.0/38.4/38.3
M83
314.6
32
4
6.3
7.6/7.2/7.1
13.1/12.9/12.9
IC 342
138.2
10.6
4
5.5
6.6/6.3/6.1
0.1/0.0/0.0
NGC 6946
95.7
11.7
5.9
3.4
3.2/3.3/3.5
0.1/0.1/0.1
NGC 2903
208.7
44.5
6.6
1.1
0.9/1.0/1.1
0.6/0.7/0.7
NGC 5055
106
74.3
7.8
0.9
0.7/0.8/0.9
0.2/0.2/0.2
NGC 3628
240.9
64.8
8.1
1.3
1.0/1.1/1.2
0.8/0.9/1.1
The Astrophysical Journal Letters, 853:L29 (10pp), 2018 February 1
Aab(>57EeV)
et al.
Northern
TA
:
7
years
109
eventsNGC 3627
242
64.4
8.1
1.1 No correction for 0.8/0.9/1.1 K. Kawata et al., Proc. 0.7/0.8/0.9
of ICRC 2015
4σ
at
39
EeV,
only
9.7%
anisotropic
fraction
Southern Auger : 10 years0.8/0.9/1.1
157 events (>57EeV)
E scale difference 2.1/2.4/2.7
NGC 4631
142.8
84.2
8.7
2.9
Table 1b/w TA and PAO !! 2.3/2.8/3.3
Southern
is seen at
Cen
A(Pre-trial~3.6σ)
M51
104.9
68.6
10.3
3.6
0.3/0.4/0.5
Populations Investigated
PierreAugerhotspot
Collaboration,
ApJL
853:L29(2018)Starburst galaxies detected by Fermi
12
NGC 891
140.4
−17.4
11
1.7
1.1/1.3/1.5
0.2/0.3/0.3
a
b
SBGs
l
)
b
)
Distance(Mpc)
Flux
Weight(%)
AttenuatedWeight:
A/B/C(%)%
Contribution: A/B/C (%)(
n(n
NGC 3556
148.3
56.3
11.4
0.7
0.4/0.6/0.6
0.0/0.0/0.0
NGC 253
660
141.6
−47.4
15
0.9
0.5/0.6/0.8
0.4/0.5/0.6
NGC
97.4
−88
2.7
13.6
20.7/18.0/16.6
35.9/32.2/30.2
NGC 2146
135.7
24.9
16.3
2.6
1.3/1.7/2.0
0.0/0.0/0.0
M82
141.4
40.6
3.6
18.6
24.0/22.3/21.4
0.2/0.1/0.1
NGC 4945
3079
157.8
48.4
17.4
1.0/1.4/1.5
0.1/0.1/0.1
305.3
13.3
4
162.1
19.2/18.3/17.9
39.0/38.4/38.3
NGC 1068
172.1
−51.9
17.9
12.1
5.6/7.9/9.0
6.4/9.4/10.9
M83
314.6
32
4
6.3
7.6/7.2/7.1
13.1/12.9/12.9
NGC
238
−54.6
22.3
1.3
0.5/0.8/0.8
0.9/1.5/1.6
IC3421365
138.2
10.6
4
5.5
6.6/6.3/6.1
0.1/0.0/0.02
Arp 299
141.9
55.4
46
1.6
0.4/0.7/0.6
0.0/0.0/0.0
NGC
6946
95.7
11.7
5.9
3.4
3.2/3.3/3.5
0.1/0.1/0.1

10 EeV γ-ray with 300 km: 5 - 10 Mpc

23

GWTC-1: Gravitational Wave Transient Catalog
1840+8-14 Mpc

credit: Guiseppe Greco, Axel Mellinger

2
γ-ray 10 EeV with 30,000 km : 40 - 80 Mpc

UHE γ-ray burst driven multi-messenger

stimation summary plots V. The contours show 90% and 50% credible regions for the sky locations of all GW events
ction. The probable position of the source is shown in equatorial coordinates (right ascension is measured in hours,
1811.12907
asured in degrees). 50% and 90%arXiv:
credible regionsof posterior probability sky areas for the GW events. Top panel:
O2 GW events [21] (GW170817, GW170104, GW170823, GW170608, GW170809, GW170814) for which alerts were24

Established method for γ-ray detections
Schematic viewÆof ALPACA
83,000 m2
Air
ShowerArray
1m2 plastic
scintillation detector

5,400 m2

15 m300 m

Muon
Detector
Array

Underground water
Cherenkov muon
detector (56 m2 / unit)

佐古＠ICRR共同利用成果発表 2018/12/18

<pe>top/(<pe>top+p

0

0.7

0.1

0.2

0.3

0.4

sin

0.5

0.6

0.7

0.8

r [m]

2

0.6
Surfacedetector array of Layered Observational Water-cherenkov counters

4. Fraction of
0.5photo-electrons (pe) collected in the top layer for the EM (square symbols, a coefﬁcient of the matrix) and muonic (round symbols, b coefﬁcponent for two hadronic models. Left : dependence as a function of zenith angle. Right: dependence as a function of distance to the shower core. The coefﬁcients a a
0.4 the contribution from the photons, e þ e " and
From
left
to
right
ssentially independent of the shower characteristics or detector distance from the shower core, they only depend on the tank geometry.

140
0.1
120

0

100 200
80

400

600

800

r [m]

60

1000319 VEM
1200

1400

Signal [VEM peak]µ± , p, EPOS-LHCµ± , p, QGSJet II.04
A. Letessier-Selvon et al. / Nuclear Instruments and Methods in Physics Research A 767 (2014) 41–49±µ , Fe, EPOS-LHC
120m II.04µ± , Fe,491
QGSJet0.2

Signal [VEM peak]

8

γ

0.3

+
e ,e

100

+
µ ,µ
±
e

e

80

±

43

rec

true

60

40
40
(square symbols, a coefﬁcient of the matrix) and muonic (round symbols, b coefﬁcient)
20
20
topcoefﬁcients a and bngle. Right: dependence as a function of distance to the shower core. The

rom the shower 0core, they only depend on the tank geometry.

ces

Signal [VEM peak]
Signal [VEM peak]uger WCD design.
12050

10040
80
30
iversality 60
2040
y calculate2010
a

-140

-120

-100

-80

-60

-40

-20

0

20

-160

-140

-120

-100

-80

-60

-40

-20

0

20

time [8 ns]

time [8 ns]

Fig. 2. Distribution of the Cherenkov photons production point in a 1.2 m height and 1.8 m radius WCD. From left to right the contribution from the photons, e þ e " and35
muon component of a 30 EeV EAS with 451 zenith angle is shown.

e±
±
e

rec

true

Signal [VEM peak]

-160

0

30
25
20

µ ± rec
Contents lists available at ScienceDirectµ ± true

Nuclear Instruments and Methods in Physics Research A 767 (2014) 41–49

15

Nuclear Instruments and Methods inPhysics Research A

10journal homepage: www.elsevier.com/locate/nimak]and b and to characterize the5bottomof the LSD, 0simulationsof the detector response
0
0
-160
-140
-120
-160have-140been-120 -100-100 -80-80 -60-60 -40-40 -20-20 0 0 2020
20-160 -140 -120 -100 -80
-60
-40
-20
0
20
formed.
Air showers
simulatedwiththeLayered water Cherenkov detector for the study of ultra high energy
time
time[8[8ns]ns] [18]
timens][16], using EPOS-LHC [17] and QGSJetII.04
Fig. 3. Schematic andartistic view
WCD[8design.
cosmicraysof a LSD built from an Auger
interaction models and FLUKA [19] at low energy.
a,n
a
a
a,b
2
AntoineLetessier-Selvon,
PierreBilloir,
MiguelBlanco
,
Ioana
C.
Mariş
, core. Th
5. Simulation35ofNuclear
the signal
collected
in
the
top
(bottom)
part
of
a
10
m
LSD
for
a
11
EeV
showerwith
461
zenith
angle
and
491
m
away
from
the
3.
Performances
construction),aInstruments
and aMethods
in Physics
Research A ±767 (2014) 41–49
25
es have been generated
withuniform
distribuMariangela
Settimo
µ
rec
!
!
!

nstruction of the muon and EM traces compared to the generated
in the right panel. The agreement is striking both in shape (timing information o! ones is shown"r [m]

Surface detector array of Layered Observational Water-cherenkov counters

ructed using the individual signal reconstruction in each of the detectors. Right: individual LSD muonr bars represent the signal resolution.

9.7 ± 0.4)%

15.9 ± 0.3)%

89% p
8% Fe

600

Protons
Irons

Double liner
Design, calibration and simulations

θ = 35

Design options

T1400 [ns]E = 30 EeV
Sfisher = 2.43

400

200
10% p
91% Fe

1

0

0

50

100
S µ (1000) [a.u.]

150

Design options✦ 750 m spacing in triangular arrangement

howers with energy between 8 and 13 EeV for proton (solid line) and iron (dashed line). Right: example
start time and muon signal size. The horizontal axis is the muon size at 1000 m from the shower core in
ed muon signal in each WCD. The vertical axis is the start time at 1400 m (see text).

✦ 10 m2, 800 stations → ~200 km2 ,[12,21]. An
iron mixed
muon signal
al start time
the showernary planar
nt curvature
have used a
ed with theground, a Fisher value better than 2 will be possible. The parameter
T1400 is given as example and it is not meant to be considered as the
optimal variable for the analysis. A detailed study of X max (or age)sensitive parameters and of the mass separation is in progress and is
however out of the scope of this paper.
17.5✦ 16 million USD for detectors✦ 100% efficiency above 10
4. Calibration strategyeV✦ p/Fe separation: 10-1, γ/Hadron separation: 10-3

Current height of the liners:
Ioana C. Mariş for the
doubleliner
proponents
40 cm top liner
80 cm bottom liner
Photomultipliers position
Solution 1:
2PMTs in the top, 1 in the bottom
Implemented in Geant4/Pierre’s code
Solution 2:
3PMTs in the top, 1 in the bottom

26Sub-EeV γ-ray search andanisotropy
Upper limits to the integrated ph
a

b

F
to
Sgr A*Sgr A*

Targeted photon search: galactic centerE2 x flux [TeV cm-2 s-1]

Figure 1: VHE -ray image of the Galactic Centre region. The colour scale indicates counts per 0.02 ⇥0.02 pixel.
Left panel: The black lines outline the regions used to calculate the CR energy density throughout the central molecularzone. A section of 66 is excluded from the annuli (see Methods). White contour lines indicate the density distributionof molecular gas, as traced by its CS line emission30 . The inset shows the simulation of a point-like source. Rightpanel:
Zoomed view of theof
inner⇠ 70 pc PeVatronand the contourresultsof the region
extract the spectrum
the diffuse• Interpretation
H.E.S.S.
forusedtheto galactic
centerofregion[The
H.E.S.S.
Collaboration,
Nature 531 (2016) 476]
emission.
-11
10

E
av

SLOW 10y

E
[Ee

10-12

1

-13

10

2

4

10-14
H.E.S.S. measurement
1σ confidence band of the best-fit spectra
Auger photon GC limit (Γ=2.32 ± 0.11)
H.E.S.S. extrapolation (Γ =2.32 ± 0.11)
H.E.S.S. extrapolation (Γ =2.32 and Ecut = 2.0 EeV)

10-15

10-16
10-1

•
•

1

10

102

3

10

104

Constrain the naive extrapolation to EeV energies
Upper limit on the cutoff energy of 2 EeV

18.07.2017510✦ First GZK γ and γ-ray burst detections.
6✦ Large scale anisotropy using different compositions.

107

10

E [TeV][The Pierre Auger Collaboration, ApJ 837 (2017) L25]

Marcus Niechciol (Pierre Auger Collaboration) | ICRC 2017 (Busan) | CRI183

11 / 12✦ Tuning the hadron interaction model at LHC energies
27

Why sub-EeV?: transition of the equatorial dipole
110

S. Mollerach, E. Roulet / Progress in Particle and Nuclear Physics 98 (2018) 85–118
The Astrophysical Journal, 868:4 (12pp), 2018 November 20

Figure 3. Evolution with energy of the amplitude (left panel) and direction (right panel) of th
map in Galactic coordinates of the right panel the dots represent the direction toward the galax
the direction toward the ﬂux-weighted dipole inferred from that catalog.

Figure 4. Maps in Galactic coordinates of the ratio between the number of observed events
arrival directions, for the four energy bins above 4EeV.

Table 5Three-dimensional Dipole Reconstruction fo
Energy (EeV)

宗像＠ICRR共同利用成果発表 2018/12/18

d⊥

dz

5.0
11.5

0.007
0.0060.003
0.011
0.0600.010−0.024±0.009−0.026±0.015

10.3
20.2
39.5

0.013
0.0580.011
0.025
0.0650.018
0.05
0.080.03−0.008±0.017−0.08±0.03−0.08±0.07

Interval

Median

4–8
8
8–16
16–32
32

Note. We show the results obtained for the two bins previously reported (The Pierre Auger Co
as dividing the high-energy range into three bins.

28

and E 8 EeV. The only nonvanishing correlation coefﬁcients
between the quantities reported in Table 7 are S (dx, Qxz ) S (d y, Qyz )  0.63 and S (dz, Qzz )  0.91. The nine components

are pranaly
PierreTA, preliminary, 〈∆〉=17
g/cm
2
HiRes, 〈∆〉=26 g/cm
HiRes/MIA
CASA-BLANCA
Yakutsk
Tunka
Auger

850
800
750
700
650
600

840

Quite a simple and natural picture!

820

Fe knee

Galactic

pre-knee slope
780

760

p

740

QGSJetII
Sibyll2.1
EPOSv1.99

500

1015

21

knee
800

n
o
t
ro

550

450

〈Xmax〉 [g/cm ]2

2〈Xmax〉 [g/cm2]

Why sub-EeV?: transition of mass composition, galactic/extra-galactic

720

slope x = 3.3

iron
1016

1017

1018

1019

1020

E [eV]

700• same slope!
• natural match
to 18
the knee!
10

3 1015 eV

Kyoto, 11 October 2016

Extra-galactic
19

10

20

10

slope x = 3.3

E [eV]UHECR 2016 / GCR-EGCR transition

E. Parizot (APC, Paris 7)

Figure 8: Measurements of hXmax i with non-imaging Cherenkov detectors (Tunka [118, 126], Yakutsk [127, 128], CASA-BLANCA [123]) and fluorescence
detectors (HiRes/MIA [129], HiRes [130], Auger [131] and TA [132]) compared to air shower
simulations [133] using hadronic interaction modelsE.
[36–38].
HiRes
Parizot
in UHECR 201617
fluorescence
detection
10text). eV
and TA data have been corrected for detector
e↵ects as indicatedby the h iabove
values (seeThe right panel shows a zoom to the ultra-high energy region.

M. Unger,
K. Kampert,
Astropart.Phys.
35 (2012)at660-678The measurement of the longitudinal air shower development
A characteristicfeatureof the laterallight distribution

29

〈Xmax〉 [g/cm ]2

2〈Xmax〉 [g/cm2]

A fine-pixel fluorescence telescopeTA, preliminary, 〈∆〉=17
g/cm
2
HiRes, 〈∆〉=26 g/cm
HiRes/MIA
CASA-BLANCA
Yakutsk
Tunka
Auger

850
800
750
700

840
820
800
780

650
600

n
o
t
ro

760

p

QGSJetII
Sibyll2.1
EPOSv1.99

550
500
450
1015

740

✦ 1720pixel: 3 mm × 3 mm, 8 × 8 pixels per module
700✦ 100 modules ⇒ 25 cm × 25 cm (40 kUSD)

iron
1016

1017

1018

1019

1020

E [eV]

MPPC® (Multi-Pixel Photon Counter) arrays

18

10

19

10

20

10

✦ 80 × 80 pixels covering in 25°
× 25 °
E [eV]✦ 0.3° × 0.3° per pixel FoV with 1 m2 aperture

Figure 8: Measurements of hXmax i with non-imaging Cherenkov detectors (Tunka [118, 126], Yakutsk [127, 128], CASA-BLANCA [123]) and fluorescence
detectors (HiRes/MIA [129], HiRes [130], Auger [131] and TA [132]) compared to air shower simulations [133] using hadronic interaction models [36–38]. HiRes
S13361-3050 series
and TA data have been corrected for detector e↵ects as indicated by the h i values (see text). The✦
right panel shows a zoom to the ultra-high energy region.

Xmax detection using fluorescence technique
15.5
above 10 eV at high altitude of 4000 m (likeMPPC arrays in a chip size package miniaturizedthrough the adoption of TSV structure
ALPACA)A characteristic feature of the lateral light distribution at
The measurementof the longitudinal air shower development30

Future ground array for UHECR detections
Exposure [km yr sr]

HiRes

Auger/TA(hybrid
detector)

2

AGASA

Telescope Array Experiment

106

AGASA/HiRes
5102~ 3000 km

Future groundarraysUpgradephase

2~30,000 km

AugerPrime

K-EUSO

Auger

2~ 100 km

GZK recovery
Mass composition
above 100 EeV

TA×4

104

HiRes
Pierre Auger Observatory

3

10

Anisotropy above
57 EeV with
10×statistics

TA

AGASA

World-onecollaboration

Fly’s Eye
1990

1995

2000

31

2005

2010

2015

2020

2025

2030

2035

2040
Year

Fluorescence detector Array of Single-pixel Telescopes✦ Target : > 1019.5 eV, ultra-high energy cosmic rays (UHECR) and neutral particles✦ Huge target volume ⇒ Fluorescence detector arrayFine pixelated camera

Too expensive to cover a huge area

Segmented mirror telescope
Variable angles of elevation – steps.

Single or few pixels and smaller
15 deg optics

45 degLow-cost and simplified telescope32

ce Detectors

Fluorescence detector Array of Single-pixel Telescopes
Telescope Array:700 km2

Pierre Auger: 3000 km2

ope Array:700 km2

TA (same scale) Auger✦ Each telescope: 4 PMTs, 30°×30° field of view (FoV)

3000 km2

700 km2✦ Reference design: 1 m2 aperture, 15°×15° FoV

per PMTale)✦ Each station: 12 telescopes, 48 PMTs, 30°×360°

3

FoV.

20 km

(not drawn to scale)3✦ Deploy on a triangle grid with 20 km spacing, like

57 EeV“Surface Detector Array”.

Photons at diaphragm✦ With 500 stations, a ground coverage is 150,000 km2.
1

3

✦ ~100 million USD for detectors

1

5 years: 5100 events (E > 57 EeV),
650 events (E > 100 EeV)All Sky Survey with TA&PAO

2

Oversampling with 20°-radius circle

3

Photons at diaphragm

60 stations2
17,000 km

Photons at diaphragm

56 EeVzenith 500
2

Directional studies
on both energy
spectrum and Xmax.
16

Northern TA : 7 years 109 events (>57EeV)
No correction forK.
Kawataetal., Proc.
2015
Southern
Auger of: 10ICRC
years 157
events (>57EeV)
E scale difference

33

FAST fluorescence telescope
Reference: D. Mandat et al., JINST 12, T07001 (2017)100
90

Efficiency [%]

80
70
60
50
40
30

Mirror reflectivity

20

Filter transmission

10

Total efficiency

0
260

280

300

320 340 360 380
Wavelength [nm]

400

420

✦
4 PMTs(20
cm,reflectance
8 dynodes
R5912-03MOD,
base
Figure5. The typical
spectralof the FAST
mirror between 260 nm and 420
nm, along

E7694-01)

spectral transmission of the UV band-pass filter. The resultant total optical efficiency is shown in b

2 Cherenkov telescope of the MAGIC [18] observatory. The filter is construct
filter✦used on thea number of small segments in order to fit the FAST prototype’s octagonal aperture. The ind
segments are fit together using brass “U” and “H” profiles, resulting in an aperture of 1 m2

1 m aperture of the UV band-pass filter(ZWB3), segmented mirror of 1.6 m diameter✦ 3 telescopes has been installed at Utah
6

Telescope support structure

remote operation and automatic shutdown

The✦telescope’s mechanical support structure was built from commercially available alu
profiles. This allows for straightforward assembly/disassembly, and easy packing and transp✦
to their light weight, while also providing an extremely stable and rigid platform for th

425 hours observation by October 2018

Comparison with Simulation - best fit

Reconstruction result of the highest event-2018/05/15

30

Simulation

Data
200
Simulation
150

20
10

100

0
50

− 10
− 20
300

350
400
Time bin [100 ns]

50
40

Data

30

Simulation

200

250

300

Data

30

Simulation

Data
150

Simulation

30

Simulation

10

0

0

− 10

− 10

− 20

− 20
150

200

250

300

350

200

Data
150

− 30
200

400
450
500
Time bin [100 ns]

Simulation

250

300

350
400
Time bin [100 ns]

50
40

Data

30

Simulation

20
100

100

10

10

0

50

0

50

− 10

− 10
0

− 20

4/9

Data

10

20

− 30
200

40

20

Preliminary

200

50

20

− 30

350
400
Time bin [100 ns]

Npe / 100 ns

250

Npe / 100 ns

Npe / 100 ns

− 30
200

0

40

Npe / 100 ns

Data

50

Npe / 100 ns

40

250

Npe / 100 ns

50

Npe / 100 ns

Npe / 100 ns

Work: Justin Albury, Jose Bellido

250

300

350
400
Time bin [100 ns]

200

0
250

300

350
400
Time bin [100 ns]

200

− 20
250

300

350
400
Time bin [100 ns]

− 30
200

ZenithAzimuthCore(X) Core(Y) Xmax
59.8 deg -173.3 deg 7.9 km -9.0 km 842 g/cm2
+0.4
+0.5
+0.2
+0.2
+10
-0.7
-1.5
-0.2
-0.2
-10

250

300

Energy17.3 EeV
+0.3
-0.7

350
400
Time bin [100 ns]

35

Future plan: FAST telescope to be installed in Auger

FD (Los Leones)

MIDAS
Lidar dome
GAP-2018-XXX

Lidar dome

DRAFT

JSPS grant-in-aid for scientific11research, 18H01225

MIDAS

Fig. 4. Solution 2: Pour a new concrete foundation just north of the existing
MIDAS pad.

Staking for hut position
Manpower requirements

The FAST and Auger collaboration members listed above would travel to
Malargüe to install and set up the FAST telescope. Installation is expected totake approximately one week, based on the time required for installation of thefirst three prototypes at TA (the third prototype installation also included
construction of the enclosure).
Following the initial setup campaign, the telescope would be operated remotely
via an ethernet link without any required local manpower. However, in case of apower outage at Los Leones, although our telescope will shut down

536

Summary and future perspectives
Future scientific goals for on-going upgrades:
TA×4: confirmation of the TA hotspot with >5σAugerPrime: indication of small scale anisotropies selecting a light composition, proton fraction at
1020 eVInteraction model: < 20 g/cm2 uncertainty on Xmax at 1020 eV
Smoking gun of cosmic ray origins: γ-rays detection spacial coincidence with UHECR hotspot
As a personal decadal survey,
A fine-pixel fluorescence telescope for low-energy extension
SLOW: A layered Water-Cherenkov detector array for sub-EeV anisotropy and ultrahighenergy photon search
FAST: A low-cost fluorescence telescope array suitable for measuring the properties of UHECRs
with an unprecedented aperture37

A homework from Jim“I hope you can bring the single pixel fluorescence detector to practical application.
While most of my colleagues are pleased with the results of Auger, I am disappointed
we failed to find sources. Instrumentation like yours may make that possible some day”
James Cronin (The 1980 Nobel Prize in Physics)

38