The all-particle spectrum of primary cosmic rays in the wide energy range from 10^14 eV to 10^17 eV observed with the Tibet-III air-shower array

0.76 MB·English

by  M. Amenomori

The all-particle spectrum of primary cosmic rays
in the wide energy range from 1014 eV to 1017 eV
observed with the Tibet-III air-shower array
M. Amenomori1, X. J. Bi2, D. Chen3, S. W. Cui4, Danzengluobu5, L. K. Ding2,
8 X. H. Ding5, C. Fan6, C. F. Feng6, Zhaoyang Feng2, Z. Y. Feng7, X. Y. Gao8, Q. X. Geng8,
0 H. W. Guo5, H. H. He2, M. He6, K. Hibino9, N. Hotta10, Haibing Hu5, H. B. Hu2,
0
J. Huang11, Q. Huang7, H. Y. Jia7, F. Kajino12, K. Kasahara13, Y. Katayose3, C. Kato14,
2
K. Kawata11, Labaciren5, G.M. Le15, A. F. Li6, J.Y. Li6, Y.-Q. Lou16, H. Lu2, S. L. Lu2,
n
a X. R. Meng5, K. Mizutani13,17, J. Mu8, K. Munakata14, A. Nagai18, H. Nanjo1,
J M. Nishizawa19, M. Ohnishi11, I. Ohta20, H. Onuma17, T. Ouchi9, S. Ozawa11, J. R. Ren2,

7
T. Saito21, T. Y. Saito22, M. Sakata12, T. K. Sako11, M. Shibata3, A. Shiomi9,11, T. Shirai9,
1
H. Sugimoto23, M. Takita11, Y. H. Tan2, N. Tateyama9, S. Torii13, H. Tsuchiya24, S. Udo11,

] B. Wang8, H. Wang2, X. Wang11, Y. Wang2, Y. G. Wang6, H. R. Wu2, L. Xue6,
x
e Y. Yamamoto12, C. T. Yan11, X. C. Yang8, S. Yasue25, Z. H. Ye15, G. C. Yu7, A. F. Yuan5,
-
p T. Yuda9, H. M. Zhang2, J. L. Zhang2, N. J. Zhang6, X. Y. Zhang6, Y. Zhang2, Yi. Zhang2,
e Zhaxisangzhu5, and X. X. Zhou7,
h
[ (The Tibet ASγ Collaboration)

2
v
3
0
8 1DepartmentofPhysics,HirosakiUniversity,Hirosaki036-8561, Japan.
1 2KeyLaboratoryofParticleAstrophysics,Institute ofHighEnergyPhysics,ChineseAcademyofSciences, Beijing100049,
. China.
1
3FacultyofEngineering,YokohamaNationalUniversity,Yokohama240-8501, Japan.
0
8 4DepartmentofPhysics,HebeiNormalUniversity,Shijiazhuang050016, China.
0 5DepartmentofMathematics andPhysics,TibetUniversity,Lhasa850000, China.
: 6DepartmentofPhysics,ShandongUniversity,Jinan250100, China.
v
7Institute ofModernPhysics,SouthWest JiaotongUniversity,Chengdu610031, China.
i
X 8DepartmentofPhysics,YunnanUniversity,Kunming650091, China.
r 9FacultyofEngineering,KanagawaUniversity,Yokohama221-8686, Japan.
a 10FacultyofEducation, UtsunomiyaUniversity,Utsunomiya321-8505, Japan.
11Institute forCosmicRayResearch,UniversityofTokyo,Kashiwa277-8582,Japan.
12DepartmentofPhysics,KonanUniversity,Kobe658-8501, Japan.
13ResearchInstitute forScienceandEngineering,WasedaUniversity,Tokyo169-8555,Japan.
14DepartmentofPhysics,ShinshuUniversity,Matsumoto390-8621,Japan.
15CenterofSpaceScienceandApplicationResearch,ChineseAcademyofSciences, Beijing100080, China.
16PhysicsDepartmentandTsinghuaCenter forAstrophysics,TsinghuaUniversity,Beijing100084, China.
17DepartmentofPhysics,SaitamaUniversity,Saitama338-8570,Japan.
18AdvancedMediaNetworkCenter,UtsunomiyaUniversity,Utsunomiya321-8585, Japan.
19NationalInstitute ofInformatics,Tokyo101-8430, Japan.
20TochigiStudyCenter,UniversityoftheAir,Utsunomiya321-0943,Japan.
21TokyoMetropolitanCollegeofIndustrialTechnology, Tokyo116-8523, Japan.
22Max-Planck-Institutfu¨rPhysik,Mu¨nchenD-80805,Deutschland.
23ShonanInstituteofTechnology, Fujisawa251-8511,Japan.
24RIKEN,Wako351-0198, Japan.
25SchoolofGeneralEducation, ShinshuUniversity,Matsumoto390-8621,Japan.
1
ABSTRACT
We present an updated all-particle energy spectrum of primary cosmic rays in a wide range
from1014eVto1017eVusing5.5×107eventscollectedintheperiodfrom2000Novemberthrough
2004 October by the Tibet-III air-shower array located at 4300 m above sea level (atmospheric
depthof606g/cm2). Thesizespectrumexhibitsasharpkneeatacorrespondingprimaryenergy
around4PeV.Thisworkusesincreasedstatisticsandnewsimulationcalculationsfortheanalysis.
We performed extensive Monte Carlo calculations and discuss the model dependences involved
in the final result assuming interaction models of QGSJET01c and SIBYLL2.1 and primary
composition models of heavy dominant (HD) and proton dominant (PD) ones. Pure proton
and pure iron primary models are also examined as extreme cases. The detector simulation was
also made to improve the accuracy of determining the size of the air showers and the energy of
the primary particle. We confirmed that the all-particle energy spectra obtained under various
plausible model parameters are not significantly different from each other as expected from the
characteristics of the experiment at the high altitude, where the air showers of the primary
energy around the knee reaches near maximum development and their features are dominated
by electromagnetic components leading to the weak dependence on the interaction model or the
primary mass. This is the highest-statistical and the best systematics-controlled measurement
covering the widest energy range around the knee energy region.
Subject headings: cosmic rays—methods: data analysis—stars: supernovae: general
1. Introduction
(Nikolsky & Romachin) during air-shower devel-
opment. Common to all models is the prediction
Although almost 100 years have passed over
of a change of the chemical composition over the
sincethediscoveryofcosmicrays,theirsourceand
knee region. Direct measurementsofprimarycos-
acceleration mechanism are still not fully under-
mic rays on board balloons or satellites are the
stood. The energy spectrum and chemical com-
best ways for the study of the chemical composi-
positionofcosmic rayscanbe keyinformationfor
tion, however, the energy region covered by them
probing their origin, acceleration mechanism and
with sufficient statistics are limited to 1014 eV.
propagation mechanism. The cosmic-ray spec-
Theenergyspectrumandchemicalcompositionof
trum has been observed by many ground-based
primary cosmic rays around the knee, therefore,
experiments to resemble two power laws, having
a form dj/dE ∝ E−γ, with γ = 2.7 below the has to be studied with ground-based air-shower
energy around 4 × 1015 eV, and then steepen- experiments using surface array and/or detectors
for Cerenkov light.
ing to γ = 3.1 abovethis energy (H¨orandel 2003).
The change of the power index at this energy is Many reports have so far been made on the
calledthespectral“knee”. Althoughtheexistence energy spectrum as well as the chemical com-
of the knee has been well established experimen- position of primary cosmic rays. Although the
tally, there are still controversial arguments on global features of the all-particle spectrum agree
its origin. Proposals for its origin range from as- well when we take into account the systematic er-
trophysical scenarios like the change of accelera- rors of about 20% involved in the energy scale
tion mechanisms (Berezhko & Ksenofontov 1999; (H¨orandel 2003), there are still serious disagree-
Stanev et al. 1993;Kobayakawaet al. 2002;Vo¨lk 2004) ments in the chemical composition depending on
at the sources of cosmic rays (supernova rem- experimentalmethods,forexample,thekneecom-
nants, pulsars,etc.), the single source assumption position obtained by the Tibet and KASCADE
(Erlykin & Wolfendale 2005),oreffectsduetothe experiments can be summarized as follows. We
propagation(Ptuskin et al. 1993;Candia et al. 2002) have already reported the energy spectrum of
inside the galaxy(diffusion, drift, escape fromthe protons and helium in the energy range from
Galaxy) to particle-physics models like the inter- 200 TeV to 10,000 TeV (Amenomori et al. 2000;
action with relic neutrinos (Wigmans 2003) dur- Amenomori et al. 2006a)fromair-showercoreob-
ing transport or new processes in the atmosphere servation, suggesting a steep power index of ap-
2
proximately -3.1. This indicates that the power nent with minor contributionof muons, whose in-
index of light component is changed from ap- teractionmodeldependence isknowntoberather
proximately -2.7 as measured by direct observa- largeamongcurrentinteractionmodels leadingto
tions to -3.1 around a few hundred TeV. Hence, alargesystematicerrorintheexperimentscarried
the lightcomponentshouldbecome lessabundant out at the sea level because of the large contribu-
at the knee and the main component responsi- tion of muons. In other words, the air shower ob-
ble to the structure of the knee must be heavier servation at high altitude is sensitive to the most
than helium. Furthermore, the spectral shape of forward region of the hadronic interactions in the
light component seems to keep the power law in- center of momentum system (CMS) where high
stead of the exponential cutoff. On the contrary, energy secondaries are produced, and the electro-
the KASCADE using electron-muon size analysis magneticcomponentasadecayproductofneutral
(Antoni et al. 2005) claims that the knee in the pions dominate the number of shower particles,
all-particle spectrum is due to the steepening of while it is insensitive to the central region of the
thespectraoflightelementswithexponentialtype CMSwherelargenumberofmuonsasdecayprod-
cutoff. uctofchargedpionsareproduced. Thedifferences
The accurate measurement of the all-particle among current interaction models are mainly re-
energy spectrum around the knee is essential to lated to the central region as seen in the prob-
establish the chemical composition at this energy lem of the electron-muon correlation. Hence, the
range. There is no precise measurement of the air-showerexperiment in Tibet can determine the
chemical composition around the knee region yet primary cosmic-rayenergy much less dependently
and it is in fact impossible to discriminate indi- uponthechemicalcompositionandtheinteraction
vidual elements clearly by indirect observations. model than experiments at the sea level.
Therefore, most of the works made so far simply We have already reported the first result on
discussedtheaveragemass<lnA>. Anotherap- the all-particle spectrum around the knee region
proachistheunfoldingoftheall-particlespectrum based on data from 2000 November to 2001 Oc-
usingshowercharacteristicslikeelectron-muonra- tober observed by the Tibet-III air-shower ar-
tio, depth of the shower maximum and so on. In ray (Amenomori et al. 2003a). In this paper, we
thesemethods,thedetailedinformationoftheall- present an updated all-particle energy spectrum
particle spectrum plays an important role in de- using data set collected in the period from 2000
termining the chemical composition. It is also Novemberthrough2004October. Theupdatesare
expected that the specific features of each com- due to (1) increased statistics by approximately
ponent like cutoff energy or source characteristics 2.6 times, (2) use ofnew simulationcodes and (3)
should be reflected in the shape of the all-particle improvementofthelateralstructurefunctionused
spectrum as discussed in the single source model forthesizeestimationofairshowers. Theprevious
(Erlykin & Wolfendale 2005). The important fea-
tures of the all-particle spectrum are the absolute
itsnhhteaernppsoniwteyse,srtiihnnedtpehxoesbisteiizfooenresopfaetnchdteraukfmntee,rew,ththheieckhdniaeffereeraednnedceetphloyef e 111eee+++000789 1 . 0 a<t =Ysaencg(tbhaejtian)g<1.1 101^01^61 e7V eV
z
connected with the acceleration mechanism and n si 1e+06 10^15 eV
theTshoeumrceeroitfocfotshmeicairra-syhso.werexperimentinTibet Electro 100000 10^14 eV
10000
is that the atmospheric depth of the experimen-
tal site (4300 m a.s.l., 606 g/cm2) is close to the 1000 Proton (QGSJET) Proton (SIBYLL)
Iron (QGSJET) Iron (SIBYLL)
maximumdevelopmentoftheairshowerswithen- 100
0 200 400 600 800 1000
ergies around the knee almost independent of the
Atmospheric Depth ( g/cm^2 )
masses of primary cosmic rays as demonstrated
in Fig. 1 for vertically incident cosmic rays. It Fig. 1.— Average transition curves of air-shower
should be also noted that the number of shower size induced by protons and iron nuclei for a ver-
particles is dominated by electromagnetic compo- tical incidence.
3
resultwasobtainedusingalmostthesameanalysis
asusedinTibetI(Amenomori et al. 1996)except
for the parameters which depend on the detector
configuration. Inthepresentpaper,thesimulation
code Cosmos is replaced by Corsika with interac-
tionmodelsQGSJET01candSIBYLL2.1,whichis
now widely used in many analyses by other works
andenablesthecomparisonofthisworkwithoth-
erseasier. Thethirdupdateonthestructurefunc-
tion is made to coverwider energy rangethan be-
fore (see section 4.1.3).
Tibet III Air Shower Array (2003) 36,900 m2
Thus, we obtained the all-particle energy spec-
trumofcosmicraysinawiderangeover3decades
between1014 eVand1017 eVandthe updated re-
sult is compared with previous ones.
2. Tibet experiment
The Tibet air-shower experiment has been op-
◦ ◦
erated at Yangbajing (E90 31’, N30 06’; 4300 m
above sea level) in Tibet, China, since 1990. The
Tibet air-shower array is designed not only for
observation of air showers of nuclear-component
origin but also for that of high energy celestial
gamma rays. Because of such multiple purposes,
the detector is constructed to cover a wide dy-
namic range for particle density covering 0.1 to
FT Detector (512)
5000and a goodangularresolutionfor the arrival FT Detector w/ D−PMT(249)
15 m Density Detector (28)
direction of air showers with energy in excess of a 789 detectors
few TeV being better than 1 degree.
Fig. 3.— Schematic view of the Tibet-III ar-
The Tibet-I surface array was constructed in ray operating at Yangbajing. The Tibet-III array
1990 (Amenomori et al. 1992) using 65 plastic consists of 761 FT detectors and 28 D detectors
scintillation detectors placed on a lattice with 15 around them. In the inner 36,900 m2, FT detec-
m spacing. This array was gradually expanded torsaredeployedat7.5mlatticeintervals,among
761 FT counters, 249 sets of detectors are also
equipped with D-PMT in addition to FT-PMTs.
707 mm
5mm Thick Lead Open-whitesquares: FTdetectorswithFT-PMT;
m
30 m Scintillator Open-black squares: FT detectors with FT-PMT
and D-PMT; Open circles: density detectors only
with D-PMT.
m
m
500
Fast Timing Density PMT
PMT
HV Cable
Signal Cable
Fig. 2.—SchematicviewoftheFTw/D-detector.
4
to the Tibet-II (1994) and Tibet-III (1999) array. the methods of the analysis are developed so that
At present, it consists of 761 fast timing (FT) they canreproducethe inputs ofsimulatedevents
counters and 28 density (D) counters surrounding liketheprimaryenergy,thelocationoftheshower
them. In the inner 36,900 m2, FT counters are axis, the arrival direction and so on. Even the
deployed at 7.5 m lattice intervals. All the FT most basic quantity like the number of particles
counters are equipped with a fast-timing photo- arriving to a detector should be ’defined’ through
multiplier tube (FT-PMT ; Hamamatsu H1161) MC because we do not measure the number of
measuring up to 15 particles. Among the 761 FT particles but the charge of PMT-output which is
counters, 249 sets of detectors (with interval of not simply proportionalto the number of charged
15 m) are also equipped with density photomulti- particles entering into the detector if we take into
plier tube (D-PMT ; Hamamatsu H3178) of wide account the contribution of the electromagnetic
dynamic range measuring up to 5000 particles in processes by photons inside the detectors. An-
addition to FT-PMTs, so that UHE cosmic rays other example of the role of MC is to define the
with energy above the knee can be observed with effective areaof the showerarraywhich should be
a good accuracy. determinedtoavoidthe erroneouscountingofthe
Eachcounterhasaplasticscintillatorplate(BI- eventswhoseshoweraxisesaredroppingoutsideof
CRON BC-408A) of 0.5 m2 in area and 3 cm in the effective area. Therefore detailed MC calcula-
thickness. A lead plate of 0.5 cm thick is put tions are needed on air-shower generation in the
on the top of each counter as shown in Fig. 2 atmosphere and on the detector response. Con-
in order to increase the counter’s sensitivity by sequently, the final result inevitably depends on
converting photons in an electromagnetic shower the interactionmodel and on the primary compo-
into electron-positron pairs (Bloomer et al. 1988; sitionmodelinMC.Thisisthemainsourceofthe
Amenomori et al. 1990). The recording of signals systematicerrorsinvolvedintheair-showerexper-
is made on time and charge information for the iment and we try to show them explicitly in the
FT-PMTs, while only the charge information for present work.
the D-PMTs. The D counterssurroundingthe in- A full Monte Carlo (MC) simulation has been
ner array are also equipped with both FT-PMT carried out on the development of air showers in
and D-PMT, where only the charge information the atmosphere and also on the detector response
of both PMTs are recorded. An event trigger sig- of the Tibet-III array. The simulation code COR-
nalis issuedwhenany fourfoldcoincidence occurs SIKA (version 6.204) including QGSJET01c and
in the FT countersrecordingmore than 0.6 parti- SIBYLL2.1 interaction models (Heck et al. 1998)
cles. Fig. 3 is the schematic view of the Tibet-III is used to generate air-shower events. All shower
array. particlesintheatmospherearetraceddowntothe
Theprimaryenergyofeacheventisdetermined minimumenergyof1MeVwithoutusingthinning
by the shower size N , which is calculated by method.
e
fitting the lateral particle density distribution to Although the chemical composition of the pri-
the modified NKG structure function (see section mary particles around the knee region is not well
4.1.3). The air-shower direction can be estimated established yet, we have to assume it in the sim-
◦
with an inaccuracy smaller than 0.2 at energies ulation. The simplest way to bracket all possi-
above 1014 eV, which is calibrated by observ- bilities is to assume pure proton and pure iron
ingtheMoon’sshadow(Amenomori et al. 2003b). primaries. Since it is almost evident that such
We used the data set obtained during the period assumptions are not realistic and lead to unac-
from 2000 November through 2004 October. The ceptable results showing disagreement with the
effective live time used for the present analysis is direct observations, these results will be men-
805.17 days. tioned as extreme cases. More realistic treatment
of the chemical composition is to extrapolate the
3. Simulation known composition at low energies measured by
direct observations. The uncertainty in extrap-
Monte Carlo simulation (MC) plays an impor-
olating to the high energy range can be treated
tant role in air-shower experiments since most of
by bracketing the reported results on the compo-
5
Table 1: Fractions of the proton(P), helium(He),
CNO(M), NaMgSi(H),SClAr(VH) and iron(Fe)
components in the assumed primary cosmic-
ray spectrum of the HD and PD models
(Amenomori et al. 2000).
HD model 1014-1015 eV 1015-1016 eV 1016-1017 eV
P 22.6% 11.0% 8.1%
HD spectrum
HMe(CNO) 1291..20%% 1212..46%% 187..48%% 2/s)m110043 He All
H(NaMgSi) 9.0% 9.4% 8.1% sr/ (*0.1) Proton
VH(SClAr) 5.6% 6.2% 5.8% 1.5/V102 CNO
Fe 22.2% 39.1% 51.7% Ge10 (*0.01)
PPD model 10143-91.001%5 eV 10153-81.011%6 eV 10163-71.051%7 eV × dJ/dE (10-11 N((**0Va0..0MH00001g1))Si
He 20.4% 19.4% 19.1% 2.5 E10-2
M(CNO) 15.2% 16.1% 16.5% 10-3 TPAIRKBEOENTTOON satellite (*0I.r0o00n01)
H(NaMgSi) 9.4% 9.9% 10.2% 10-4 JRAUCNEJEOB (a)
VH(SClAr) 5.8% 6.2% 6.3% 103 104 105 106 107 108
Fe 9.4% 9.9% 10.2% Energy (GeV)
PD spectrum
2/s)m110043 He All
sr/ (*0.1) Proton
sition study around the knee. In order to exam- 1.5/V102 CNO
ine the composition dependence involved in the Ge10 (*0.01)
all-particle spectrum, we used two kinds of the dE ( 1 N(*a0.M001g)Si
mdoimxeidnacnocmepoofshiteioanvymcoodmeplso.neOnnteairsoubnadsedthoenkntheee × dJ/10-1 (*V0.0H001)
as reported by works (Amenomori et al. 2006a; 2.5 E10-2
TIBET
Apomsietnioonmoisriceatllaeld. 2H0D00;mOogdieol.20A04n)oathnedrthisisbcaosmed- 1100--43 PAJRARKUCEONENJTEOOOBN satellite (*0I.r0o00n01) (b)
on the dominance of light component (p and 103 104 105 106 107 108
Energy (GeV)
He) as reported by works (Antoni et al. 2005;
Raino 2004; Fowler 2001) and called PD model. Fig. 4.— Primary cosmic ray composition for
The energy spectra of individual mass groups in (a) the HD model and (b) the PD model. The
HD model and PD model are shown in Fig. 4(a) all-particle spectrum, which is a sum of each
and Fig. 4(b), respectively. Table 1 shows their component, is normalized to the Tibet data,
fractionalcontentsatgivenenergies. Intotal,four and they are compared with other experiments.
kind of the primary composition, namely, pure PROTON satellite (Grigorov et al. 1971),
proton, pure iron, HD model and PD model are AKENO (Nagano et al. 1984),JACEE
usedinthesimulationwiththe minimumprimary (Asakimori et al. 1998),RUNJOB
energy of 50 TeV. (Apanasenko et al. 2001).
All secondary particles are traced until their
energies become 1 MeV in the atmosphere. The
shower axis was placed on Tibet array at random
withinaradiusof100mfromthecenterofthear-
ray. In order to treat the MC events in the same
way as experimental data analysis, simulated air-
shower events were input to the detector with the
6
samedetectorconfigurationasthe Tibet-III array tribution function which is tuned to reproduce
with use of Epics code (ver. 8.64) (Kasahara) to the above defined number of particles using the
calculatetheenergydepositoftheseshowerparti- Monte Carlo simulation under our detector con-
cles. Experimentally,thenumberofchargedparti- figurations. The number of MC events, typically
clesisdefinedasthePMToutput(charge)divided for QGSJET+HD,is as follows. About 10 million
by that ofthe single particle peak, whichis deter- of air showers are generated with primary energy
mined by a probe calibration using cosmic rays, above 50 TeV. After imposing the selection crite-
typically single muons. For this purpose, a small ria described in sec. 4.2, the surviving number of
scintillatorof25cm×25cm×3.5cmthickwitha events is about 5 million, among which 1.9 mil-
PMT (H1949)is put onthe topofthe eachdetec- lion events belong to the unbiased energy region
tor during the maintenance period. This is called corresponding to E0>100 TeV. Almost the same
a probe detector and is used for making the trig- numberofMC eventsareobtainedforothermod-
ger of the each Tibet-III detector. The response els to compare with each other.
of each detector is calibrated every year through
probe calibration. In the simulation, these events 4. Analysis
triggered by the probe detector was also exam-
4.1. Reconstruction of air-showers
ined by a MC calculation, in which the primary
particles were sampled in the energy range above
An example of the shower profile obtained by
thegeomagneticcutoffenergyatYangbajing(>10
Tibet-III arrayisshowninFig.6(a)andFig.6(b)
GeV), and all secondary particles which pass the
which represent the map of arrival time and par-
probedetectorandtheTibet-IIIdetectorsimulta-
ticle density of shower particles, respectively. Al-
neously were selected for the analysis. Since the
though the Tibet array has quite low energy
valueofPMToutputisproportionaltotheenergy
threshold of a few TeV for the purpose of the
loss of the particles passing through the scintilla-
celestial gamma ray observation, the detection ef-
tor, the peak position of the energy loss distribu-
ficiency for the nuclear component including iron
tion corresponds to the experimental single peak
nuclei is not sufficient at low energy region. Ad-
of the probe calibration.
ditional event selection condition is required for
According to the MC, the peak position of the unbiased detection of all particles and for the
the energy loss in the scintillator is calculated capability of the lateral density fitting. The fol-
as 6.11 MeV (the details are written in the pa- lowing condition was applied on the selection of
per (Amenomori et al. 2007)). We can then cal- the events for the all-particle spectrum analysis.
culate the number of charged particles for each
detector hit as the total energy loss in each scin- N ≥10 with n ≥5 , (1)
D p
tillator divided by 6.11 MeV instead of counting
the number of charged particles arriving to the where ND expresses the number of detectors hit
detector in MC events. We confirmed that the and np the number of particles per a detector.
shape of the energy loss distribution, which is de- This conditionsatisfies the requirements for unbi-
termined by probe-calibration simulation, shows ased analysis in the energy range above 100 TeV
a reasonable agreement with the charge distribu- as described below.
tion of the experimental data as shown in Fig. 5,
4.1.1. Determination of the core position.
wherethe proportionalitybetweenthe energyloss
∆E and the PMT output charge Qi is assumed The core position of each air-shower (Xcore,
as Qi =ki×∆E, where ki is a proportional con- Ycore)isestimatedbyusingthefollowingequation
stant depending on given detector typically being :
around4pC/MeV.Thus,alldetectorresponsesin- Σρiwxi Σρiwyi
(X ,Y )=( , ) , (2)
cluding muons andthe materializationof photons core core Σρ w Σρ w
i i
inside the detector are taken into account. The
whereρ istheparticledensityatthei-thdetector
i
total number of charged particles of each event
and the weight w is an energy dependent param-
(hereafter,we callthis an“ air showersize (Ne)”)
eter varying between 0.8 and 2.0. It is confirmed
wasestimatedusingthemodifiedNKGlateraldis-
7
thatthemeanerrorofthecorepositioncanbees-
Energy Loss (MeV) timated as 5 m by reconstructing MC events (see
2000 5 10 15 20 25 30
Fig. 7). Lower energy event selection condition
Exp (Output Charge)
than eq. (1) leads to worse core resolution which
Events150 MC (Energy Loss) makes lateral density fitting difficult.
mber of 100 4.1.2. Determination of the arrival direction.
Nu50
The arrival direction of the air shower is es-
0
0 20 40 60 80 100 timated using the time signal measured by the
Output Charge (pC)
761 FT (fast-timing) counters. The shape of the
Fig. 5.— Charge distribution in a detector mea- showerfrontisassumedtobe areverse-conictype
sured by probe calibration (see the text). In or- as shown in Fig. 8. Direction cosine of the shower
der to compare the charge distribution with the axis is determined by using a method of least
simulation on the energy loss in a scintillator, the squares in which the difference is minimized be-
MC result is adjusted by multiplying a constant tween the arrival time signals of each detector
tomeetwiththe samepeakpositionasthe exper- andtheexpectedvaluesontheassumedconewith
iment. The fluctuation of the number of photons given direction cosine.
inscintillationlightistakenintoaccountwiththe An experimental check of the angular resolu-
normal distribution in MC. tionby this method is made by the observationof
the Moon’sshadow(Amenomori et al. 2003b)us-
ing large statistics of low energy events (>3 TeV)
andalsoaso-calledeven-oddmethod,whichmany
authorshavebeenusingforestimatingtheangular
resolution (Amenomori et al. 1990). The recon-
(a)
structionof the highenergy MC events assuresus
thatthe meanerrorofair-showerdirectioncanbe
estimated as 0.2◦ at energies above 1014 eV (see
Fig. 9).
4.1.3. Estimation of the shower size.
The lateral density distribution is corrected to
theinclinedplaneperpendiculartotheshoweraxis
and used for the shower size estimation. In this
1400
(b)
1200
mber of events 1 680000000
Nu 400
200
0
0 5 10 15 20 25 30 35 40
R (m)
Fig. 7.— The distribution of core-position error.
The mean error of core position can be estimated
Fig. 6.— An example of the map of (a) arrival
as 5 m. Showerselection criteria: E ≥ 100 TeV,
time of shower particles, (b) particle density, ob- 0
sec(θ) ≤ 1.1 and the core position located at the
tained by Tibet-III array.
inner 135 m × 135 m of the array.
8
1
( a )
0.5
0
a(s)-0.5
-1
-1.5
A cone -2
0.6 0.8 1 1.2 1.4 1.6
S
0
( b )
-1
ARCieorlr aest hilovowece artt iidmoiinrnegctionFitting plane Relative timing )s( b---432
-5
Fig. 8.— Determination of air-shower direction. 0.6 0.8 1 S 1.2 1.4 1.6
Fig. 10.— The numerical values of a and b are
plotted as a function of s, where original defini-
tions of a(s) and b(s) in NKG function are shown
bythedottedlines,andtheopencirclesdenotethe
averagedMC data using QGSJET01c+HDmodel
which are fitted by empirical formulae shown by
red lines. See the text.
work, the determination of the lateral distribu-
tionfunctionofshowerparticlesisveryimportant,
sincethetotalnumberofchargedparticlesineach
3500
event is estimated by fitting this function to the
3000
experimental data. Using the Monte Carlo data
mber of Events 122505000000 oifmubntecantitni,oenwdecuafnondubenredtfithhtetaestdatmwheeelfclootlnloodwtihtiineognlasmtaeorsdatilfihdeeidsetNxrpiKbeuGr--
Nu 1000
tionofshowerparticlesunderaleadplateof5mm
500
thickness:
0
0 0.2 0.4 0.6 0.8 1
Angle Distance (degree)
Fig. 9.— The distribution of opening angle be- f(r,s)= CN(se)(rr′)a(s)(1+ rr′)b(s)/rm′2 (3)
m m
tween true and estimated arrival directions. The
mean error of air-shower direction can be esti-
mated as 0.2◦. Shower selection criteria : E0 ≥ C(s)=2πB(a(s)+2,−b(s)−a(s)−2) (4)
100 TeV, sec(θ) ≤ 1.1 and the core position lo-
cated at the inner 135 m × 135 m of the array. wherer ′ =30m,andthe variablescorresponds
m
to the age parameter, N is the total number of
e
showerparticles and B denotes the beta function.
The original meaning of r in NKG formula is a
m
Moliere unit, being 130 m at Tibet altitude, how-
′
ever,we treatr as a unit scale ofthe lateraldis-
m
tribution suitable to describe the structure of the
9
air showersobserved by Tibet-III array,whose ef-
Table 2: The shower-size resolution are summa-
fective areais 135m× 135m. The functions a(s)
rized for the events induced by primary parti-
and b(s) are determined as follows. In CORSIKA
cles of E ≃ 100 TeV and E ≃ 1000 TeV with
simulation, the shower age parameter s is calcu- 0 0
sec(θ) ≤ 1.1 and the core position located at in-
lated at observation level by fitting the number
ner 135 m × 135 m of the array, based on the
of particles to a function for the one dimensional
QGSJET+HD, QGSJET+PD, SIBYLL+HD and
shower development. It may be possible to as-
SIBYLL+PD model.
sumethatairshowerswiththesameshowerages,
are in almost the same stage of air-shower devel-
E0 = 100 TeV
opment in the atmosphere, i.e. they show almost
Zenith angle QGS. QGS. SIB. SIB.
the same lateral distribution for shower particles
(sec(θ)) +HD +PD +HD +PD
irrespective of their primary energies. The lateral
1.0 - 1.1 9.0% 9.1% 9.2% 9.5%
distributionoftheparticledensityobtainedbythe
E0 = 1000 TeV
simulation with carpet array configuration is nor-
Zenith angle QGS. QGS. SIB. SIB.
malized by the total number of particles which is
(sec(θ)) +HD +PD +HD +PD
derived from the total energy deposit in an in-
finitely wide scintillator. These events are then 1.0 - 1.1 5.0% 5.8% 5.2% 5.9%
classified according to the stage of air-shower de-
velopment using the age parameter and they are
averaged over the classified events. The fitting of
determine the functions a(s) and b(s) and used in
theequation(3)totheaveragedMCdataismade
the analysis described below. Other details are
to obtain the numerical values a and b. Thus,
described in Ref.(Amenomori et al. 2007).
we can obtain the behavior of a and b as a func-
tion of s as shown in Fig. 10(a) and Fig. 10(b), Based on the Monte Carlo simulation, the cor-
whereoriginaldefinitionsofa(s)andb(s)inNKG relation between the true shower size (true size)
function are shown by dotted lines. Although our and the estimated shower size (fit size) is demon-
result shows different dependences of a and b on stratedinFig.12(a)andFig.12(b). Here,thetrue
s from the original NKG function, it is confirmed shower size means the number of particles calcu-
thatthelateraldistributionoftheshowerparticles lated for a carpet array while estimated shower
is better reproduced by our formula (see Fig. 11). size is for the real Tibet-III array using the mod-
This expression is valid in the range of s = 0.6 ∼ ified NKG function mentioned above. As seen in
1.6,secθ < 1.1andr= 5∼ 3000m. Twointerac- these two figures, a good correlation is obtained
tion models of QGSJET01c and SIBYLL2.1 and between the true shower size and the estimated
four primary composition models of pure proton, showersize. Thesystematic deviationoflessthan
pure iron, HD and PD are used independently to 1 % seen around 100 TeV of Fig. 12(a) shows
that we need finer tuning of modified NKG func-
tion at low energies and this error is finally cor-
0.1
rected. The shower size is well reproduced with
0.01 Modified NKG
age = 1.06 standard deviation of 5% around the primary en-
)r( f 0 10.0e.00-000511 OaNgerie g==in 10a..l77 N087K2G x 10^6 Ne = 1.841 x 10^6 eorngythoefQ1G00S0JETTeV+HwDithmo1.d0el<. Tshece(sθh)o≤we1r.s1izbearseesd-
1e-06 olutions estimated are summarized for the events
MC data (Proton)
1e-07 True Ne = 1.833 x 10^6 with different simulation model combinations in
age= 1.03
1e-08 E0 = 2.24 x 10^15 eV Table 2.
Zenith angle =19.3 degree
1e-09
1 10 100 1000
r ( m ) 4.2. Data selection
Fig. 11.— Lateral density distribution of the Thefollowingevent-selectioncriteriaareadopted
charged particle obtained with use of the carpet in the present analysis.
simulation. The shower size Ne is better repro-
1)Morethan10detectorsshoulddetectasignalof
duced by our modified NKG function.
morethanfiveparticlesperdetector,asmentioned
10