博士學位 請求論文
指導敎授 崔 榮 一

A Measurement of νµ↔ντ Oscillation Parameters
Using Atmospheric Neutrino Observed in Super
Kamiokande-I and Super Kamiokande-II

成均館大學校 大學院
物理學科
核및粒子物理專攻徐賢寬

博士學位 請求論文
指導敎授 崔 榮 一

A Measurement of νµ↔ντ Oscillation Parameters
Using Atmospheric Neutrino Observed in Super
Kamiokande-I and Super Kamiokande-II이 論文을 理學 博士學位請求論文으로 提出합니다
年 4月 日2006

成均館大學校 大學院
物理學科
核및粒子物理專攻
徐賢寬

이

論文을

博士學位

徐賢寬의

論文으로

理學

認定함

2006年 6月 日

審査委員長 오 필

렬

審査委員

김

수

봉

審査委員

박

차

원

審査委員

최

수

용

審査委員

최

영

일

A Measurement of νµ↔ντ Oscillation Parameters
Using Atmospheric Neutrino Observed in Super
Kamiokande-I and Super Kamiokande-II

Hyunkwan Seo
Department of Physics
Sungkyunkwan University

Abstract
The deficit of atmospheric muon neutrinos observed in Super Kamiokande iswell explained by νµ ↔ ντ oscillation hypothesis. In this study, neutrino oscillation
analysis was done utilizing zenith angle and energy dependent deficit of atmosphericmuon neutrinos. The effect of momentum binning on the constraint of oscillation
parameters was investigated. The best-fit parameters were located at (sin22θ = 1.00,

×∆m2 = 2.5 10-3 eV2) with the minimum χ2 value, χ2min = 767.5 / 737 DOF by the
combined analysis of Super Kamiokande-I and Super Kamiokande-II data in the
optimized momentum binning. The allowed oscillation parameter regions were
sin22θ > 0.93 and 2.1×10-3 < ∆m2 < 3.0×10-3 eV2 at 90% C.L.

Contents
1 Introduction ............................................................................................... 1
1.1 Review of neutrino mass and oscillation .............................................. 1
1.2 Calculation of neutrino oscillation probability ..................................... 2
1.3 Usefulness of atmospheric neutrino for the oscillation study ................. 4
1.4 Detection method of atmospheric neutrino in Super Kamiokande ......... 5
2 Super Kamiokande Detector ...................................................................... 10
2.1 Detector .......................................................................................... 11
2.1.1 Tank ...................................................................................... 11
2.1.2 PMTs ..................................................................................... 11
2.1.3 PMT support structure and others .............................................. 14
2.1.4 Inner detector electronics and data acquisition system .................. 16
2.1.5 Outer Detector Electronics and Data Acquisition ......................... 18
2.2 Detector Calibration ........................................................................ 18
2.2.1 Water Transparency Measurement ............................................. 18
2.2.2 Relative Gain Calibration ......................................................... 34
2.2.3 Relative Timing Calibration ...................................................... 36
2.2.4 Other calibrations.................................................................... 37
2.3 2001 accident and Super Kamiokande-II reconstruction .................... 39
3 Monte Carlo Simulation ............................................................................ 40
3.1 Flux calculation of atmospheric neutrino........................................... 40
3.2 Neutrino Interaction ........................................................................ 47
3.2.1 Elastic and quasi-elastic scattering ............................................ 47

3.2.2. Single meson production .......................................................... 48
3.2.3. Deep inelastic scattering.......................................................... 49
3.2.4. Nuclear effects........................................................................ 49
4 Atmospheric Neutrino Data Reduction and Reconstruction ......................... 52
4.1 Fully Contained Events (FC) ............................................................ 55
4.1.1 FC data reduction.................................................................... 55
4.1.2 FC event reconstruction............................................................ 57
4.2 Partially Contained Events (PC) ....................................................... 61
4.2.1 PC data reduction.................................................................... 61
4.2.2 PC Event Reconstruction .......................................................... 67
4.3 Upward going muon (UPMU) data ................................................... 70
4.3.1 UPMU data reduction .............................................................. 70
4.3.2 Upmu Event Reconstruction ...................................................... 74
4.3.3 Upmu Eye-Scanning................................................................. 75
4.3.4 Background Estimation from Cosmic Ray Muon .......................... 76
5 Atmospheric Neutrino Data ....................................................................... 78
5.1 FC and PC Events ........................................................................... 78
5.2 UPMU events .................................................................................. 81
6 Two Flavor Neutrino Oscillation Analysis................................................... 83
6.1 Overview......................................................................................... 83
6.2 Oscillation Analysis with New Binning .............................................. 84
6.2.1 New binning............................................................................ 84
6.2.2 Oscillation analysis method with new binning .............................. 84
6.2.3 Systematic uncertainty terms used in the oscillation analysis.......... 87
6.2.4 Study of sensitivity to finer binning with SK-I MC......................... 93
6.2.5 Result with SK-I data and SK-II data separately ........................... 93
6.3 Oscillation analysis with SK-I and SK-II combined ...................... 97
6.3.1 Combining strategy .................................................................. 97
6.3.2 Result of the combined analysis ............................................... 97
6.3.3 Oscillation analysis with pseudo data to check the likeliness of theresult obtained with real data ..........................................................103
7 Conclusion ..............................................................................................105

Appendix Tables..............................................................................107

Bibliography..............................................................................................116

List of Tables
4.1 Event rate for each FC reduction step. ...................................................... 58
4.2 The efficiency of each reduction step of FC events. .................................... 58
4.3

The misidentification probabilities for single-ring e-like and µ-like events
estimated with simulated CC quasi-elastic neutrinos. .................................. 63

4.4 Event rate for each PC reduction step. ...................................................... 68
4.5 The detection efficiency of PC events at each reduction step. ....................... 69
4.6

Percentage of events saved or rejected at each fitter. ................................... 74

4.7

Upmu detection efficiency estimated using Monte Carlo. ............................ 75

5.1 The number of observed events in the sub-GeV and multi-GeV samples as wellas the expected number of events in the absence of neutrino oscillations........ 79
5.2

UP/DOWN ratio of fully-contained events. ............................................... 80

5.3 The number of observed events in the upward-going muon data sample andMonte Carlo. The numbers in 100 years SK-I MC were normalized to SK-II MC
live-time of 60 years. ............................................................................. 82
6.1

Summary of the best fit parameters and the allowed range obtained with SK-I
data. .................................................................................................... 96

6.2

Summary of the best fit parameters and the allowed range obtained from the
combined analysis. SK-I result are shown for comparison. .......................... 98

A1. Systematic errors in neutrino flux (1) ............................................108
A2. Systematic errors in neutrino flux (2) ............................................109
A3. Systematic errors in neutrino interaction.......................................110
A4. Systematic errors in event section .................................................111

A5. Systematic errors in event reconstruction.....................................112
A6. Systematic errors in event section .................................................113
A7. Systematic errors in event reconstruction.....................................114
A8. Systematic errors in solar activity..................................................115

List of Figures

1.1

Process of atmospheric neutrino production in the earth’s atmosphere. ............ 5

1.2

A cut away diagram of Super Kamiokande detector. ..................................... 6

1.3

Cherenkov radiation is emitted when a particle travels faster than the speed oflight in the medium. ................................................................................. 7

1.4

A muon neutrino event forming a ring pattern of Cherenkov light. in the Super
Kamiokande detector. .............................................................................. 9

2.1

A cross section view of Super Kamiokande detector. .................................. 12

2.2

Schematic view of a 50 cm PMT used in inner detector of Super Kamiokande.
........................................................................................................... 13

2.3 Quantum efficiency of the photocathode as a function ofwavelength. ............ 13
2.4

Schematic view of PMT support structures. ............................................... 15

2.5 A block diagram of ATM board used for inner detector data acquisition. ....... 16
2.6 The data acquisition system for the inner detector....................................... 17
2.7

OD DAQ block diagram and data flow. .................................................... 19

2.8

Relative photo-sensitivity. ....................................................................... 20

2.9 Effective charge observed as a function of the path length. .......................... 20
2.10 Time variation of the water attenuation length measured by through-goingcosmic ray muon. .................................................................................. 21
2.11 Laser system for water scattering and absorption parameter measurement and its
typical event display. ............................................................................. 22
2.12 The detail view of the arrangement of laser firing setup................................ 23

2.13 The upper plot is the photon arrival time distribution of the top PMTs. Dots isdata and line is Monte Carlo. The lower plot is the ratio, (MC-DATA)/DATA.
........................................................................................................... 24

2.14 Comparison of time variation of light attenuation length measure by laser systemand cosmic ray muon. ............................................................................ 25
2.15 Positions of light injectors at the detector. .................................................. 26
2.16 Diagram of whole light scattering measurement system used in Super
Kamiokande-II. ..................................................................................... 27
2.17 ‘Optical SW (switch)’, which distribute laser light into 8 different fibers
mechanically by step motor. .................................................................... 28
2.18 How to calculate Number of hit (Nhit) for each LI for monitoring water quality.
........................................................................................................... 29

2.19 ‘Nhit/Q’ as a function of light intensity. Horizontal axis is the ratio of lightintensity to normal intensity, which was used in Monte Carlo. ..................... 30
2.20 Upper plot is the spatial distribution of ‘number of hit’ for 3 different intensity
made with Monte Carlo of oldtop LI. Vertical axis is the number of hit PMTs
(Nhit) in the top region of the tank normalized by the intensity and horizontal
axis is the distance from oldtop LI. Bottom plot is the ratio between lower
intensity to normal (filled circle) and larger intensity to normal (open circle). 31
2.21 ‘Nhit/Q’ as a function of light intensity without PMTs within 2 m from LI. .... 32
2.22 Left plot (a) is the distribution of ‘Oldtop LI – newtop LI’ for‘100×DATA(Nhit/Q) / MC(Nhit/Q)’ and Right plot (b) the distribution of‘(Each wavelength) – (Mean of 4 wavelength)’. ......................................... 33
2.23 DATA(Nhit/Q) / MC(Nhit/Q) obtained for each LI for various period of Super
Kamiokande-II to investigate position dependence of the water quality. ........ 34
2.24 The relative gain measurement system. ...................................................... 35
2.25 The timing calibration system. .................................................................. 37
2.26 A typical plot of timing vs. pulse height. This plot is referred as ‘TQ-map’. .... 38
3.1(a) The atmospheric neutrino energy spectrum calculated by several models. (b)The ratio of the calculated flux models. .................................................... 42

3.2 The flux ratio of νµ + anti-νµ / νe + anti-νe as a function of neutrino energy. ... 42
3.3 The flux ratios of νµ to anti-νµ and νe to anti-νe versus neutrino energy. Solid,
dashed and dotted lines show the prediction by Honda, Bartol and Fluka flux. 43
3.4 The flux of atmospheric neutrino versus zenith angle. Solid, dashed and dottedlines show the prediction by Honda, Bartol and Fluka flux. ......................... 44
3.5 The zenith angle distribution of the flux of upward-going atmospheric neutrino
observed as upward-going muon event in Super Kamiokande. Solid, dashed anddotted lines show the prediction by Honda, Bartol and Fluka flux. ................ 46
3.6 The calculated flight length of neutrinos for vertically down-going (0.95 < cosθ
< 1.0) and near horizontal-going (0.05 < cosθ < 0.10) directions................... 46
3.7

Cross section for neutrino (upper) and anti-neutrino (lower). Solid line shows
the calculated total cross section. The dashed, dot and dash-dotted lines show
the calculated quasi-elastic, single-meson and deep-inelastic scatterings,
respectively. Data points are from various experiments. .............................. 51

4.1

A schematic view of the four event categories of the atmospheric neutrinos
observed in Super Kamiokande. .............................................................. 53

4.2 The number of hits in the largest outer detector cluster, which was used toseparate the FC and PC event. The histogram is the MC prediction with neutrino
oscillation (solid line) and no oscillation (dashed line). The upper figure is for
SK-I MC and data. The lower figure is with SK-II, which zoomed in the region
around the separation. In SK-II, the criteria of the separation was defined as 16.
........................................................................................................... 54

4.3 The distribution of the likelihood difference between a single-ring and multiring assumption for sub-GeV (left) and multi-GeV (right) FC events of SK-I.
The points show the data and the histograms show the Monte Carlo.............. 60
4.4 The distribution of the likelihood difference between a single-ring and multi-ringassumption for sub-GeV (left) and multi-GeV (right) FC events of SK-II. The
points show the data and the histograms show the Monte Carlo. ................... 60
4.5

An event display of an observed single-ring e-like event. ............................ 62

4.6

An event display of an observed single-ring μ-like event. ............................ 62

4.7

The distribution of particle identification likelihood for sub-GeV (left) andmulti-GeV (right) FC single-ring events of SK-I. Dots show data and histograms
show Monte Carlo. ................................................................................ 64

4.8

The distribution of particle identification likelihood for sub-GeV (top) andmulti-GeV (bottom) FC single-ring events of SK-II. Dots show data andhistograms show Monte Carlo. ................................................................ 64

4.9 The procedure to isolate upmu events. ...................................................... 71
4.10 The zenith angle distribution of upward muon candidates near the horizon for
two different regions in azimuth made with SK-I data. The white circles and theblack inverted triangle are for thin and thick rock region, respectively. .......... 77
5.1 Event rate as a function of elapsed days of SK-II data. SK-I event rates wereshown as a reference by line. ................................................................... 81
6.1 The schematic view of the binning. Each bin in this figure is divided further by
10 zenith angle bin. The top is the bins used in the old analysis (180 bins) andthe bottom is the new bins (370 bins) used in this analysis. .......................... 86
6.2

Result of sensitivity study to binning effect. 90 % C.L allowed region obtainedwith the virtual data using 180 and 370 bins are shown for various oscillation
parameter region and live-time. ............................................................... 94

6.3 The left plot (a) is the allowed oscillation parameter region corresponding to the
68%, 90% and 99% confidence level (C.L.) obtained with 370 bins. The rightplot (b) is the comparison of the results with 180 bins and 370 bins for 90 % C.L.
........................................................................................................... 95

6.4

90 % C.L allowed region for each sub-sample. The left is the result with 180
bins and the right is 370 bins. .................................................................. 95

6.5 The allowed oscillation parameter region obtained with SK-II data. The left is
the contour obtained with all the samples combined and the right is 90% C.L
allowed region for each sub-sample.......................................................... 96
6.6 The left plot is the allowed oscillation parameters for νµ ↔ ντ with the combined analysis of SK-I and SK-II. Three contours correspond to 68 %,90 % and 99 % C.L. allowed region. The right plot is the comparison withthe result of SK-I only analysis for 90 % C.L.. ........................................... 99
6.7

χ2-χ2min distribution projected to sin22θ and△m axes resulted from SK-I
2

and SK-II combined analysis. ............................................................... 99
6.8

90 % confidence level allowed oscillation parameter regions of each subsample for νµ ↔ ντ oscillations with SK-I and SK-II data combined. ...100

6.9 The zenith angle distribution of fully-contained sub-GeV and multi-GeV sample
obtained by the combined analysis of SK-I and SK-II. The points show data,
box histograms show the non-oscillated Monte Carlo and the lines show the best
fit expectation obtained from the combined analysis. ................................101
6.10 The zenith angle distribution of multi-ring, partially-contained and upward-going
muon

sample obtained by the combined analysis of SK-I and SK-II. The

points show data, box histograms show the non-oscillated Monte Carlo and thelines show the best fit expectation obtained from the combined analysis. .....102
6.11 The width of ∆m2 region at 90% C.L obtained with 20 set of pseudo data and real
data for SK-I (left) and SK-II (right), respectively to check the likeliness of theresult of real data. ................................................................................104

Chapter 1
Introduction
1.1 Review of neutrino mass and oscillation
Neutrinos are still the least known particles among the well-establishedfundamental fermions in the Universe. The most important property of the neutrinos
is whether neutrinos are massive or not. Neutrino interacts very little with other matter.
It can travel through 50 billion miles of water without interacting. On the surface ofthe earth, about a half of a trillion neutrinos from the sun pass through each squarecentimeter of area every second. It is also believed that there are 300 relic neutrinos
per cubic centimeter of a volume in our universe which were produced about one
second after the Bing-bang. Thus, if neutrino has non-zero mass, it can contribute to
the mass and ultimate fate of the universe. The direct mass measurements based on
the decay kinematics of 3H, π and τ have been tried but it appeared that we reached
the end of the road on this effort unless a new generation of the techniques isforthcoming. There is, however, the promising way to probe the neutrino mass, whichis referred to as neutrino oscillation experiments. Oscillation is the changing of aneutrino’s type as it travels through space or matter. This can occur only if theneutrino possesses mass. When a neutrino is created it must be of a specific flavor,that is, it must be an electron neutrino (νe), a muon neutrino (νµ) or a tau neutrino (ντ).
When one speaks of massive neutrinos they are given the labels ν1, ν2 and ν3 whichhave mass m1, m2 and m3, respectively. In general a neutrino of a particular flavor will1

be a linear superposition of three neutrinos, each of a particular mass. Recent
experiments using atmospheric [2-5], solar [41,42], reactor [43] and accelerator
neutrinos [44] have demonstrated that neutrino change flavor as they travel from thesource to the detector, a phenomenon consistent with the hypothesis of neutrino
oscillation. The finding is a major discovery with a far reaching impact in theelementary particle physics, cosmology and astrophysics. The phenomenon ofneutrino oscillation which requires neutrinos to have non-zero mass will alter ourview of the world of elementary particles, and the Standard Model, the currently
prevailing theory of the elementary particles, must be modified. In the Standard
Model the neutrinos have zero mass. The finding will also make the theories of theGrand Unification more viable and attractive, and make the universe a lot heavier
than we currently assume.

1.2 Calculation of neutrino oscillation probability
If neutrinos are massive and mixed, neutrinos are produced and detected in theform of the weak eigenstates whereas when they propagate their motion is dictated by
the mass eigenstates. This leads to the phenomenon of neutrino oscillations [20]. If
neutrinos oscillate, the weak eigenstates and mass eigenstates are related by a unitarymatrix U as

ν W = Uν M ,Where, in the case of the two generations of neutrinos,
 cos θU = 

 − sin θ

sin θ cos θ The equation of motion for the mass eigenstate is

2

•

 E1

iν M = Hν M = 

 0

0
ν
E 2  Mν 

where , ν M =  1 ν 2 

Since neutrinos are expected to be extremely relativistic, if we use the following
approximation,Ei =p 2 + mi2 ≅ p +

mi2
:p≅E
2p

then
•



iν M =  E ⋅1 +Replacing ν M → e

1  m12
2 E  0

0 
ν M
m22 

ν M , the above equation becomesiαt•iν M = ( E + α )ν M +
2

 m1

 0

Where Μ M = 1Μ Mν M2E
0 .

m22 

Taking α = − E ,•iν M =

1
Μ Mν M
2E

Now, using ν M = U (θ )ν W , the equation of motion for weak eigenstates can be+3derived as•iν M =1
U (θ )Μ M U (θ ) +ν W
2E

By solving the above differential equation, we obtain the well–known oscillation
probability for the two generation case.
 1.27∆m 2 (eV 2 ) L( km) Eν (GeV )P (ν α → ν β ) = sin 2 2θ sin 2 

where, θ is the mixing angle, ∆m2 is the difference of the squared mass eigenvalues,
L is the distance from source to detector, and Eν is neutrino energy.

1.3 Usefulness of atmospheric neutrino for the oscillation study
Atmospheric neutrinos are produced from the decays of particles resulting frominteractions of cosmic rays with Earth’s atmosphere. Production of electron and muon
neutrino is dominated by the following processes and its charge conjugate. (Figure1.1)

π + (or K + ) → µ + +ν µ

a e +ν +ν+eµThat gives an expected ratio of the flux of µ to the flux of e of about 2. Of course
this is approximate as the actual ratio depends on energy and to some extent angle andposition. Vertically downward-going neutrinos travel about 15 km while verticallyupward-going neutrinos travel about 13,000 km before interacting in the detector.
Because of good geometrical symmetry of the earth, we can expect up-down
symmetry of neutrino flux. In the atmospheric neutrino analysis [2-5], however, asignificant deficit of atmospheric muon neutrino interactions depending on zenith4angle were observed compared to the expectation and it was interpreted as neutrino
oscillation. Thanks to wide ranges of L and Eν, which atmospheric neutrinos undergo,
we can measure the neutrino oscillation parameters, ∆m2 and sin22θ. This thesis
reports on the atmospheric neutrino analysis observed in Super Kamiokande.

(cosmic rays)

p, He …
π±, K±µ±

νµ(νµ)
νµ(νµ)

Figure 1.1e±
νe(νe)

Process of atmospheric neutrino production in the earth’s atmosphere.

1.4 Detection method
Kamiokandeofatmospheric

neutrino

in

Super

Super-Kamiokande (also called Super-K or SK) is a 50 kiloton water Cherenkov
detector. Figure 1.2 shows a cut-away diagram of the Super-Kamiokande detector(The detail description of Super Kamiokande detector is given in Chapter 2).
Atmospheric neutrinos are detected in two ways in Super-K. The low energy

5neutrinos from 100 MeV to 10 GeV are observed via the following charged current
interactions with nuclei in the water.

ν + N → l + N'
Where, N and N’ are the nucleons, l is the lepton.

Figure 1.2

A cut away diagram of Super Kamiokande detector.

The flavor of neutrino is identified by the flavor of the lepton. This charged lepton
generates Cherenkov light, which is detected by the PMTs installed on the wall of6Super-K. Another way of atmospheric neutrino detection is for high energy muon
neutrinos. High energy νμ can undergo charged current interactions with the rocksurrounding the detector and then produce high energy muons, which enter thedetector. This muon is also observed via Cherenkov light by PMTs in Super-K.
Cherenkov radiation is emitted when charged particles pass through matter with avelocity exceeding the velocity of light in the medium,v > vi = c / nwhere, n = refractive index of the medium, c = velocity of light and vi = thresholdvelocity.

For a particle traveling faster than light, the wave-fronts do overlap, andconstructive interference is possible and this leads to emit visible light known as
Cherenkov radiation (Figure 1.3).

Figure 1.3 Cherenkov radiation is emitted when a particle travels faster than the speed oflight in the medium.

7

The threshold energy corresponding to this threshold velocity is 0.768 MeV forelectron and 158.7 MeV for muon. Cherenkov light is emitted under a constant angle

θ with the particle trajectory, given bycosθ =

1βnWhere β = v/c

When this cone of Cherenkov light reaches the wall of the detector, it forms a ringpattern. The PMTs measure the arrival time of the light and its photon yield. Figure
1.4 shows an event display of neutrino event detected via Cherenkov light in Super
Kamiokande detector.

8

NUM
2
RUN
999999SUBRUN
0EVENT
2DATE **-Dec- 7
TIME
23: 4:10
TOT PE:
7304.2
MAX PE:
42.9
NMHIT : 2238
ANT-PE:
65.4
ANT-MX:
4.4
NMHITA:76

90/00/00:NoYet:NoYet
90/00/00:;R= 0:NoYet
R :Z : PHI : GOOD
0.00: 0.00: 0.00:0.000
CANG : RTOT : AMOM : MS

RunMODE:MonteCarloTRG ID :00000011
T diff.:
0.
us
:
0.
ms
FEVSK :80000003
nOD YK/LW: 4/ 1
Q thr. :0.0
BAD ch.: no mask
SUB EV : 0/ 0
Dec-e: 0( 0/ 0/ 0)
VQ:
0
0
0
0
VT:
0
0
0
0
VQ:
0.
0.
0.
VT:
0.
0.
0.

Comnt;

0.
0.

Figure 1.4 A muon neutrino event forming a ring pattern of Cherenkov light. in the Super
Kamiokande detector.

9Chapter 2
Super Kamiokande Detector
Super-Kamiokande, a 50 kiloton water Cherenkov detector is located in the
Mozumi mine of the Kamioka Mining Company in Gifu prefecture, Japan, under the
peak of Mt. Ikenoyama, providing a rock overburden of 2,700 meter-water-equivalent.
Super-K consists of two concentric, optically separated water Cherenkov detectors.
Super-Kamiokande began to take data in April, 1996. It was shut down formaintenance and upgrade in July, 2001. This 1996-2001 running period is referred toas Super-Kamiokande-I. During filling the water to start Super Kamiokande-II, in
November, 2001, an apparent cascade of implosions triggered by a single
photomultiplier (PMT) implosion destroyed over half of the PMTs installed in thedetector. In 2002, the detector was reconstructed using about 5,000 PMTs (theoriginal number of PMTs used in Super-K is 11,146) encased in acrylic covers to
avoid a similar accident. This partial reconstruction was done quickly in only a year.
Super Kamiokande-II was operated from Dec. 2002 to Oct. 2005 with half theoriginal density of PMTs. After termination of Super Kamiokande-II, the longawaited full reconstruction of the detector began. In June, 2006, the detector’s third
phase, Super Kamiokande-III, will start to take data.

102.1 Detector
2.1.1 Tank
The outer shell of the detector is a cylindrical stainless steel tank, 39 m in diameterand 42 m in height. Figure 2.1 shows a cross section view of Super Kamiokande. The
tank is self-supporting, with concrete backfilled against the rough-hewn stone walls tocounteract water pressure when the tank is filled with water. The capacity of the tankexceeds 50 ktons of water. A cylindrical PMT support structure divides the tank into
two distinct, optically isolated volumes. The structure has inner dimensions, 33.8 m(diameter) by 36.2 m (height), defining the inner detector (ID) which contains 32
ktons of water and was viewed by inward-facing 50 cm PMTs (Hamamatsu R3600).
Approximately 2.5 m distance remains on all outside of the support structure. Thisoutside region defining outer detector (OD), serves as an active veto counter against
incoming particles as well as a passive shield for neutrons and gamma rays from thesurrounding rocks. OD was instrumented with 1,885 outward-facing 20 cm PMTs
(Hamamatsu R1408) in both Super Kamiokande-I and Super Kamiokande-II. The two
detector volumes are isolated from each other by two light-proof sheets on both
surfaces of the PMT support structure. The 55 cm thick support structure comprises
dead space from which light in principle cannot escape. Therefore, any interaction inthat space cannot be detected. The OD PMTs were mounted in water-proof housings
which effectively block light from the dead space.

2.1.2 PMTs
A 50 cm ID PMTs used in this experiment is depicted in Figure 2.2. The bialkali(Sb-K-Cs) photocathode has peak quantum efficiency of about 21% at 360 nm-400
nm as shown in Figure 2.3. The collection efficiency for photoelectrons (pe) at thefirst dynode is over 70%. The transit time spread for a 1 pe signal is 2.2 ns. The

11

average dark noise rate at the 0.25 pe threshold is about 3 kHz in Super Kamiokande-I.
The ID PMTs were operated withgain of 107 at a supply high voltage ranging from 1700 to 2000 V. The neck of each
PMT was coated with a silver reflector to block external light.

Figure 2.1 A cross section view of Super Kamiokande detector.

12Figure 2.2 Schematic view of a 50 cm PMT used in inner detector of Super Kamiokande.

Figure 2.3

Quantum efficiency of the photocathode as a function ofwavelength.

13

2.1.3 PMT support structure and others
Figure 2.4 shows a detail of the PMT support structure. All support structurecomponents are stainless steel. The basic unit for the ID PMTs is a supermodule, aframe which supports a 3x4 array of PMTs. Each supermodule has two OD PMTs
attached on its back side. Opaque black plyethylene telephthalate sheets cover thegaps between the PMTs in the ID surface. These sheets improve the optical separation
between the ID and OD and suppress unwanted low-energy events due to residual
radioactivity occurring behind the PMTs. Light collection efficiency in the OD isenhanced by wavelength shifting (WS) plates attached to each OD PMT. The WS
plates are square acrylic panels, 60 cm on a side and 1.3 cm thick, doped with 50 mg/ℓof bis-MSB. The WS plates function by absorbing UV light, and then re-radiating
photons in the blue-green, better matching the spectral sensitivity of the PMT’s
bialkali photocathode. The light collection of the PMT plus WS unit is improved overthat of the bare PMT by about a factor of 1.5. To further enhance light collection, theOD volume is lined with a reflective layer made from Type 1073B Tyvek
manufactured by DuPont. This inexpensive and very tough paper-like material has
excellent reflectivity in the wavelength range in which PMTs are most sensitive,
especially at short wavelengths. Measured reflectivities are on the order of 90% forwavelengths in excess of 400 nm, falling to 80% at about 340 nm. The presence of
this liner allows multiple reflections of Cherenkov light, which minimizes the effects
of dead PMTs, given their coarse spacing in the OD. Cherenkov light is spread overmany PMTs, reducing pattern resolution but increasing overall detection efficiency.
The average geomagnetic field is about 450 mG and is inclined by about 45 withrespect to the horizon at the detector site. The strength and uniform direction of thegeomagnetic field could systematically bias photoelectron trajectories and timing in
the PMTs. To counteract this, 26 sets of horizontal and vertical Helmholtz coils arearranged around the inner surfaces of the tank. With operation of these coils, theaverage field in the detector is reduced to about 50 mG. The magnetic field at various
PMT locations were measured before the tank was filled with water.

14

Figure 2.4 Schematic view of PMT support structures.

15

2.1.4 Inner detector electronics and data acquisition system
ID PMT signals are processed by custom built TKO modules called AnalogTiming-Modules (ATMs). The TKO (TRISTAN KEK Online) system was originallydeveloped and built by KEK, and is optimized for front-end electronics where a large
number of channels are to be handled. The ATM has the functionality of a combined
ADC (Analog-to-Digital Converter) and TDC (Time-to-Digital Converter), andrecords the integrated charge and arrival time of each PMT signal. Signals from 12
PMTs are fed to an ATM board. Figure 2.5 shows a block diagram of the ATM board.
The PMT signal into the current splitter is divided into four signals. One of them is
fed to the discriminator with the threshold level for each channel set to -1 mV, whichcorresponds to 1/4 pe equivalent. When the PMT signal is above the threshold level, a“hit” signal with 200ns width and 15 mV pulse height is asserted on the ATM front
panel. The HITSUM signal, which is used to generate the global trigger signal, is theanalog sum of all these 15 mV/Channel pulses. At the same time, either of the split
signals, A or B, is held by a QAC (Charge to Analog Converter), and a TAC (Time to

Figure 2.5 A block diagram of ATM board used for inner detector data acquisition.

16

Analog Converter) starts to integrate constant current. If a global trigger is received,
the information in the TAC/QAC is digitized and stored in internal memory buffers.
Since, for each channel, TAC integration is started by the PMT signal itself, thearrival time of the signal can be inferred from the value of the integrated charge. A“self gate” chip generates start/stop signals for the TAC, the gate signal for QAC, andclear signals for both. Output signals from the TAC/QAC are fed to an ADC anddigitized. There are two TACs and QACs for each channel so that events in rapid
succession, such as muon followed by its decay electron, can be processed without
dead-time.
Figure 2.6 shows a schematic view of the inner detector data acquisition (ID DAQ)
system. There are total 946 ATM modules, located in four “quadrant huts”. ATMs
record ADC/TDC data of each PMT above threshold when a global trigger signal isasserted by the VME TRG (TRiGger) module. The global trigger signal and event
number information generated in the TRG module, are distributed to all ATMs via 48
GONG (Go/NoGo) modules.

Figure 2.6 The data acquisition system for the inner detector.

172.1.5 Outer Detector Electronics and Data Acquisition
A block diagram of the OD data acquisition (OD-DAQ) system is shown in Figure
2.7. As with the inner detector, signals from the PMTs are processed and digitized ineach of 4 quadrant electronics huts. From there, the signals are picked off, digitizedand stored. As in the ID, local HITSUM pulses are formed in each quadrant hut andsent to the central hut, where their analog sum creates an OD trigger signal when the
OD HITSUM threshold is exceeded. Figure 2.7 illustrates the data flow in a quadrant
hut as well as the overall data flow. When a trigger from any source occurs, the data
are passed out of the quadrant huts to the central hut for further processing. Once the
PMT signals have been picked off, coaxial ribbon cables feed them to the custom
charge-to-time conversion modules (QTC). The purpose of these modules is tomeasure the hit time and charge of the PMT pulse and convert it to a form which canbe easily read and stored by the time-to-digital converters (TDCs). The output of the
QTC modules is a logic pulse (ECL level) whose leading edge marks the hit arrival
time and whose width represents the integrated charge Q of the PMT pulse.

2.2 Detector Calibration
2.2.1 Water Transparency Measurement
1) Indirect measurement with cosmic rays

The light attenuation length in water can be measured by using through-goingcosmic ray muons. These muons are energetic enough to deposit almost constantionization energy per unit path length (about 2 MeV / cm) independent of particle
energy. This fact makes it possible to use these muons as a “constant” light source.
The advantage of this method is that continuous and abundant samples of muonscome without cost as we take normal data.

18

Figure 2.7 OD DAQ block diagram and data flow.

The disadvantage is that we cannot measure the transparency as a function of
wavelength. However, since what we measure is the Cherenkov spectrum, thisdisadvantage is not necessarily a problem. Under the assumption that light reaching
the PMTs is not scattered, the charge Q observed by a PMT is expressed by

Q = Q0 ⋅f (θ )
l



l



L

⋅ exp −

Where ℓ is the light path length, L the effective attenuation length, Q0 a constantand f(θ) the relative photo-sensitive area, which depends upon the incidence angle θof the light on the PMT, as shown in Figure 2.8. Figure 2.9 is Ql / f (θ ) as afunction of the path length ℓ together with the best fit in the form of the function def-19

Figure 2.8 Relative photo-sensitivity. (a) Measurement result, (b) Definition of the incident
angle.

Figure 2.9 Effective charge observed as a function of the path length.

20

-ined above. The resulting attenuation length is found to be 105.4±0.5 m. As the data
sample for this measurement is accumulated automatically while we take normal data,
it is possible to continuously monitor the attenuation length as a function of time.
Figure 2.10 shows the attenuation length changes with time, which are correlated withwater quality. Such time variations in the attenuation length are corrected in analysis
from the time series of calibration data.

Figure 2.10

Time variation of the water attenuation length measured by through-goingcosmic ray muon.

2) Light Scattering Measurement with laser

To measure the scattering and absorption parameters separately, a combination ofdye and N2 lasers of wavelengths 337, 371, 400, 420 nm are used as light sources. In
Super Kamiokande-I, the laser beam is brought into the top of the detector via opticalfibers. Each laser fires every 6 second during normal data taking. Figure 2.11 shows aschematic view of the laser calibration setup and a typical event of laser light. A
cluster of hits at the bottom of the tank is due to unscattered photons. The remaining
hits, in barrel and top PMTs, are due to photons scattered in the water or reflected by
the bottom PMTs and black sheet. The total charge in the cluster of bottom hits is use-

21

Figure 2.11

Laser system for water scattering and absorption parameter measurement and its

typical event display.

22

-d for normalization factor as a reference of light intensity. Figure 2.12 show thedetail view of the arrangement of laser firing setup. The direction of light fired byeach laser is turned by half mirror and finally go to the fiber bundle where 7 fibers areconnected. One of them is the fiber brought into the top of detector. 4 of them arebrought to the two monitoring PMTs and two monitoring APDs. The other two
remained spare. The intensity of the each laser is adjusted by ND filter. The two
diffusers are put in front of the fiber bundle to diffuse light randomly. Figure 2.13
shows distributions of photon arrival times at the top PMTs. In the upper plot, thepeak near 730 ns is from scattered photons, and the second peak represents reflected
photons at the bottom of the detector. For the Monte Carlo, the total number ofscattered photons is tuned using Rayleigh scattering parameters, and the shape of thearrival time distribution is used to adjust the absorption parameters. The lower plot in

Figure 2.12 The detail view of the arrangement of laser firing setup.

23

Nhit / BottomQ

337nm : TOP
0.01
0.009
0.008

MC
DATA

ABS = 0.46E8
RAY = 1.90E8

0.007
0.006
0.005
0.004
0.003
0.002

(MC-DATA)x100/DATA %

0.001
0
80

700

750

800

850

900

950

700

750

800

850

900

950

1000

60
40
20
0
-20
-40
-60
-80
1000
time (nsec)

Figure 2.13 The upper plot is the photon arrival time distribution of the top PMTs. Dots isdata and line is Monte Carlo. The lower plot is the ratio, (MC-DATA)/DATA.

Figure 2.13 shows the time dependence of the ratio between Monte Carlo and data,
which agree within ~2%.
Using the measured absorption and scattering parameters, the light attenuation
length of water is calculated using the following equation.

24

Latten =1
α abs + α scatt

Figure 2.14 compares the time variation of the light attenuation length measured by
the 420 nm laser and cosmic ray muon during Super Kamiokande-I period. The
result from cosmic ray muon is with average Cherenkov spectrum so it does notagree with the one from the laser system. However, we can see the tendency of time
variation is consistent between the two measurements.

Figure 2.14 Comparison of time variation of light attenuation length measure by laser systemand cosmic ray muon.

253) Upgraded Light Scattering Measurement system in Super Kamiokande-II

Super Kamiokande is a large water tank. So it was questioned if the water quality
depends on positions in the tank. Actually a clue for the position dependence isobserved in the detail analysis with the laser system in Super Kamiokande. However
we were not sure whether it is due to position dependence of the water quality or toother effects. To investigate the position dependence of the water quality, light
injectors (LI) similar to the one used in Super Kamiokande-I were installed at thedifferent 7 positions of the detector and started to inject laser with the beginning ofSuper Kamiokande-II. Figure 2.15 shows the position of light injector (LI) installed in
the tank. The LI used in Super Kamiokande-I was named as Oldtop LI, the new
injector installed at the top as Newtop LI and other LIs were named as B1, B2, B3, B4,
B5 and Bottom LI according to their positions. Figure 2.16 shows the diagram of thewhole light scattering measurement system used in Super Kamiokande-II. To input
laser beam to different 7 fibers connected to different 7 positions in the tank, one laser
beam source is moved by step motor in ‘optical SW’ depicted in Figure 2.16. The step
motor is controlled by software through serial communication (RS-232). Figure 2.17
shows the inside of the ‘optical SW’.

Figure 2.15

Positions of light injectors at the detector.

26

Figure 2.16

Diagram of whole light scattering measurement system used in Super

Kamiokande-II.

27

Figure 2.17

‘Optical SW (switch)’, which distribute laser light into 8 different fibers

mechanically by step motor.

To monitor the water quality, the value, ‘Nhit/Q’ was used, where ‘Nhit’ means thenumber of hit PMT within a region defined around the LI and ‘Q’ is the total charge
deposited at Super Kamiokande used as normalization factor. In SK-I, the charge inthe bottom cluster was used as the normalization factor but in SK-II, total charge wasused because total charge is less affected by possible change of beam shape or beam
direction and their imperfect reproduction in Monte Carlo (MC). The region to
calculate ‘Nhit’ was defined radius 10 m around LI in case of barrel LI and whole top(bottom) region is used for top (bottom) LI. Figure 2.18 explain how to calculate‘Nhit’ visually. Another systematic effect to the value,‘Nhit/Q’ is from light intensity.
We normalized the number of hit by total charge to compensate light intensity depen-28

Number of hit

Number of hit (Nhit)（Nhit)

Radius = 10 m

Top & Bottom LI
Figure 2.18

Barrel LI

How to calculate Number of hit (Nhit) for each LI for monitoring water quality.

-dence of ‘Nhit’. However, even ‘Nhit/Q’ has light intensity dependence slightly.
Figure 2.19 shows ‘Nhit/Q’ as a function of light intensity estimated by Monte Carlowith oldtop LI for 337 nm laser. The tendency is that low intensity gives high‘Nhit/Q’ and high intensity low ‘Nhit/Q’. The reason why ‘Nhit/Q’ depends on lightintensity was investigated and it was found that the multiple hits near LI cause thistrouble. Figure 2.20 shows spatial distribution of ‘number of hit’ for three different
intensity (that is, normal intensity, three times larger intensity than normal intensity
and three times lower intensity) made with Monte Carlo of oldtop LI. Vertical axis is
the number of hit PMTs (Nhit) in the top region of the tank normalized by theintensity and horizontal axis is the distance from oldtop LI. We can see thediscrepancy near LI, which is due to multiple hits. To reduce this multiple hits effect,
we decided not to use PMTs within 2 m from LI in the analysis. Figure 2.21 shows‘Nhit/Q’ as a function of light intensity similar to Figure 2.19 but without PMTs
within 2 m from LI. We can see the light intensity dependence was reduced byremoving PMTs near LI.

29

Nhit/Q

Oldtop LI (337nm, MC)
0.07
0.069

+/- 3%

0.068
0.067
0.066
0.065
0.064
0.063

3.0

2.6

2.2

1.8

1.4

1.0

0.71

0.56

0.46

0.39

0.33

0.062
Intensity/Normal intensity

Figure 2.19 ‘Nhit/Q’ as a function of light intensity. Horizontal axis is the ratio of lightintensity to normal intensity, which was used in Monte Carlo.

30

Normal / 3
NormalNormal × 3

Multiple hits causes thisdiscrepancy

This causes light intensity
dependence

Figure 2.20 Upper plot is the spatial distribution of ‘number of hit’ for 3 different intensity
made with Monte Carlo of oldtop LI. Vertical axis is the number of hit PMTs (Nhit) in the top
region of the tank normalized by the intensity and horizontal axis is the distance from oldtop LI.
Bottom plot is the ratio between lower intensity to normal (filled circle) and larger intensity tonormal (open circle).

31

Nhit/Q

Oldtop LI (337nm, MC)

0.058+/- 3%

0.057

0.056

0.055

0.054

0.053

0.052

3.0

2.6

2.2

1.8

1.4

1.0

0.71

0.56

0.46

0.39

0.33

PMTs within 2m from LI were removed

Intensity/Normal intensity

Figure 2.21 ‘Nhit/Q’ as a function of light intensity without PMTs within 2 m from LI.

In the upgraded system of SK-II, each laser fires every 1 sec during normal datataking. This is more frequent than 6 sec of SK-I because we should inject light in total
8 positions in SK-II. To analyze the laser data, every 10 days data were combined and‘DATA(Nhit/Q) / MC(Nhit/Q)’ was compared for each 7 LI to investigate Z (Z is the
direction along the height of the detector) dependence of water quality. To estimate
systematic error between LIs, ‘100×DATA(Nhit/Q) / MC(Nhit/Q)’ was compared
between oldtop LI and newtop LI. Systematic error between different wavelength oflaser was also estimated by comparing‘100×DATA(Nhit/Q) / MC(Nhit/Q)’

between four wavelength lasers for one LI. Figure 2.22 (a) shows the distribution of‘Oldtop LI – newtop LI’ for ‘100×DATA(Nhit/Q) / MC(Nhit/Q)’ and (b) shows the
distribution of ‘(Each wavelength) – (Mean of 4 wavelength)’. The 1 σ denoted in the32All LIs
Entries137

of eventsof events

All wavelength25

Entries

70

1096

60

20
50

15

σ= 2.6 %σ= 1.7 %

40

30

10
20

5
10

0

-10

-8

-6

-4

-2

0

2

4

6

8

10

0

-10

-8

-6

-4

-2

0

Old top – New top

4

6

8

10

(Each wavelength) - (Mean of 4 wavelengths)(a)Figure 2.22

2Mean - each wave(%)

Oldtop - Newtop(%)(b)

Left plot (a) is the distribution of ‘Oldtop LI – newtop LI’ for×‘100 DATA(Nhit/Q) / MC(Nhit/Q)’ and Right plot (b) the distribution of ‘(Each wavelength)– (Mean of 4 wavelength)’.

figure is the one estimated from the distribution in the plot, which means it is fromcombination of a couple of variables. If we calculate 1 σ for one variable, they are
1.9 % and 2.1 % for ‘Oldtop LI – newtop LI’ and ‘(Each wavelength) – (Mean of 4
wavelength)’, respectively.
Figure 2.23 shows DATA(Nhit/Q) / MC(Nhit/Q) obtained for each LI for various
period of Super Kamiokande-II. According to the result in Figure 2.23, Z-dependenceof water quality is not clearly seen in Super Kamiokande. However, we are trying to
reduce systematic errors of the laser system to see Z-dependence, which may be
hidden under the systematic effect.

33

Bottom
BottomB3

B5

B5

Bottom

85

B4

85

B3

85

B2

90

B1

90

Newtop90

Bottom

95

B5

95

B4

95

B3

100

B2

1.Feb.2004

105

100

B1

Z-dependence_(run040201)

110

100

Figure 2.23

B4

105

B2

B5

BottomB4

B3

1.Dec.2003

B5

110

85

115

B4

105

Z-dependence_(run031201)

B3

1.Oct.2003

115

Data(N/Q) / MC(N/Q)

110

85
Data(N/Q) / MC(N/Q)

B1

Z-dependence_(run031001)

B2

90

B1

90

Newtop90

Bottom

95

B5

95

B4

95

B3

100

B2

100

Newtop100

1151.Aug.2003

105

B2

105

110

B1

1.Jun.2003

Z-dependence_(run030801)

Newtop

110

115

B1

105

Z-dependence_(run030601)

Newtop1.Apr.2003

115

Data(N/Q) / MC(N/Q)

Data(N/Q) / MC(N/Q)

110

85
Data(N/Q) / MC(N/Q)

Z-dependence_(run030401)

NewtopData(N/Q) / MC(N/Q)

115

DATA(Nhit/Q) / MC(Nhit/Q) obtained for each LI for various period of Super

Kamiokande-II to investigate position dependence of the water quality.

2.2.2 Relative Gain Calibration
The high voltage of each PMT is set to provide approximately equal gain for all
PMTs in the Super-Kamiokande detector. Each PMT came with a nominal highvoltage value, which was specified by Hamamatsu Phonics just after its production,using three calibration light sources : a DC light source, a Xe-lamp pulsed light source,
and a light source at the single pe level. However, because of possible long-term drift
in the factory calibration system, which may have caused a systematic difference
between tubes manufactured at different times, all PMTs had their gain re-calibrated.
Normalized Relative gain was measured by the system shown in Figure 2.24. Light
generated by a Xe lamp is passed through an ultraviolet (UV) filter and injected into a
scintillator ball via an optical fiber. The scintillator ball is an acrylic ball with BBOT

34

Figure 2.24

The relative gain measurement system.

(2, 5-bis(5'-tert-butyl-2-benzoxazolyl)thiophene) wavelength shifter and MgO powder
diffuser. BBOT absorbs UV light and emits light with a peak at 440 nm, which is in
the sensitive region of the PMTs used for Super-Kamiokande for Cherenkov lightdetection. In this measurement, each PMT detects a few tens of pe. The intensity ofthe primary UV light is monitored by two photodiodes and one PMT. The output of
the monitor PMT is also used for triggering. The relative gain Gi of the i-th PMT is
obtained by:

Gi =

Qi
l 
⋅ l i ⋅ exp i 
Q0 f (θ ) L

35

where Qi is the observed charge by the PMT, Q0 a constant, L the effective lightattenuation length, ℓi the distance from the light source to the PMT and f(θ) the
relative photo-sensitivity as a function of the incident angle of light a on the PMT
defined in Figure 2.8. The high voltage value for each PMT is set so that the“corrected Q” of each PMT is approximately the same as for all the others. Here,“corrected Q” is the pulse height corrected for light attenuation, acceptance of the
PMT, and uniformity of the scintillator ball. It is further normalized by the Xe
monitor pulse height. This measurement is done for various positions of thescintillator ball and settings of the high voltage value. After the initial gain adjustment,
long term stability of PMT gain was monitored using the same system. The absolutevalue of PMT gain cannot be measured while the experiment is running.
The relative gain spread σgain is defined as the standard deviation obtained byfitting the relative gain distribution of all PMTs to a Gaussian. It was 7.0% at the
beginning of Super Kamiokande-I and 1.8% at the beginning of Super Kamiokande-II.

2.2.3 Relative Timing Calibration
The relative timing of PMT hits is important for event reconstruction. It depends on
the length of the signal cable between the PMT and the ATM, and also depends on
observed charge because of the discriminator slewing effect. The timing difference ineach individual PMT has to be measured precisely to get better timing resolution.
Figure 2.25 shows the system for measuring the relative timing of hit PMTs. The N2
laser emits intense light of wavelength 337 nm within a very short time (less than 3ns). The wavelength is shifted to 384 nm, which is near the lower edge of sensitivityof a PMT to the Cherenkov light, by a dye laser module. The light intensity ischanged by an optical filter, and the PMT timing is measured at various pulse heights.
After passing through the optical filter, the light goes to the diffuser ball in thedetector through an optical fiber. The schematic view of the diffuser ball is also
shown in Figure 2.25. The diffuser tip located at the center of the ball is made from36

Figure 2.25 The timing calibration system.

TiO2 suspended in optical cement. The light emitted from the tip is further diffused by
LUDOX® manufactured by Grace Davison, silica gel made of 20 nm glass fragments.
The combination of diffuser tip and LUDOX® make modestly diffused light without
introducing significant timing spread. A typical scatter plot of the timing and pulse
height is shown in Figure 2.26, which is called a ‘TQ-map’. Each PMT has its own
TQ-map, which is applied in data reduction. The timing resolution of PMTs as afunction of pulse height is estimated from the TQ maps. Typical timing resolution at
the single pe level is better than 3 ns.

2.2.4 Other calibrations
Beside the calibrations which were introduced so far, there are other calibrations

37

using stopping muon, neutral pion, LINAC, etc. The detail of the other calibrations
can be found in [36].

Figure 2.26 A typical plot of timing vs. pulse height. This plot is referred as ‘TQ-map’.

38

2.3 2001 accident and Super Kamiokande-II reconstruction
Super Kamiokande underwent upgrade work to replace PMTs that began in midJuly, 2001 and ended in mid-September. It had begun to be filled with water since 18,
September, 2001. On the morning of 12, November, when the purified water filled thetank to the level 31.7 m above the bottom (about 3/4 of the tank), about 60 % of the
50cm photomultipliers (PMTs) used for the Super Kamiokande detector weredestroyed in a few seconds. Because data had been being taken continuously since thepurified water started to fill the tank, we could analyze the data just after the accident
to investigate what happen at the instant of accident. From the analysis of the datataken during a few msec period after the accident, it is most likely that the accident
was caused by an implosion of one of the bottom PMTs. A lot of investigation about
the possible cause of the implosion was done and the most possible cause wasconcluded like following.

1. Stress accumulated during the upgrade work caused the implosion of a PMT.
2. During the installation or transporting process a slight crack was made at the
ASSY of a PMT caused the implosion.

And the test with the real PMTs and simulation proved that shock waves generated
by an implosion of a single PMT trigger a chain reaction which destroys other PMTs.
To prevent this kind of accident, it was necessary to improve PMT assembling andinstallation procedures and to encase PMTs to prevent onset of shock waves and chain
reaction.
In 2002, the detector was reconstructed with reduced PMT coverage (47% of Super
Kamiokande-I). All the PMTs were encased in acrylic covers. Super Kamiokande II
were operated for about three years and terminated at October, 2005 to begin a fullreconstruction.

39

Chapter 3
Monte Carlo Simulation
3.1 Flux calculation of atmospheric neutrino
To carry out detailed studies of neutrino oscillations using atmospheric neutrinos, it
is important to know the expected flux without neutrino oscillations. The difficulties
and the uncertainties in the calculation of atmospheric neutrino fluxes differ between
high and low energies. For low energy neutrinos around 1 GeV, the primary fluxes of
cosmic ray components are relatively well known. The low energy cosmic ray fluxesof less than about 10 GeV are modulated by solar activity, with the minimum flux
occurring at times of high solar activity. At these energies, the primary cosmic rays
are also affected by the geomagnetic field through a rigidity (momentum/charge)
cutoff. For high energy neutrinos, above 100 GeV, primary cosmic rays with energies
greater than 1000 GeV are relevant. At these energies, solar activity and the rigidity
cutoff do not affect the cosmic rays, but details of the higher energy primary cosmic
ray flux are not measured well. There are several flux calculations [6–14]. Unlike
older calculations [15, 16], in which the secondary particles were assumed to travel inthe direction of the primary cosmic ray (1-dimensional calculations), the current
calculations employ three dimensional Monte Carlo methods.
In this section, we outline the methods of the calculation and results from threeatmospheric neutrino flux calculations [10, 13, 14], which cover the energy range

40

relevant to the present analysis, are compared. They are used in the analysis of theSuper-Kamiokande data. The measured flux of primary cosmic rays is used, including
solar modulation and geomagnetic field effects. The interaction of cosmic ray
particles with the air nucleus, the propagation and decay of secondary particles are
simulated. A neutrino flux was calculated specifically for the Kamioka site. The
calculated energy spectra of atmospheric neutrinos at Kamioka are shown in Figure
3.1 (a). Also shown in Figure 3.1 (b) is the comparison of the calculated fluxes as afunction of neutrino energy. The agreement among the calculations is about 10%
below 10 GeV. This can be understood because the accuracy in recent primary cosmic
ray flux measurements [17, 18] below 100 GeV is about 5% and the data fromindependent experiments agree within the quoted uncertainties. However the primarycosmic ray flux is much less accurate above 100 GeV. Therefore, for neutrinoenergies much higher than 10 GeV, the uncertainties in the absolute neutrino flux
could be substantially larger than the disagreement level among the calculations. 0.05
was assigned for the uncertainty in the energy spectrum index above 100 GeV in the
primary cosmic ray energy spectrum. This causes approximate 0.05 uncertainties in
the neutrino spectrum index above 5 GeV. Figure 3.2 shows the calculated flux ratioof νµ + anti-νµ to νe + anti-νe as a function of the neutrino energy, integrated oversolid angle. This ratio is essentially independent of the primary cosmic ray spectrum.
Especially, in the energy region of less than about 5 GeV neutrino energies, most of
the neutrinos are produced by the decay chain of pions and the uncertainty of thisratio is about 3%, which is estimated by comparing the three calculation results. In thehigher energy region (10 GeV for νµ + anti-νµ, and 100 GeV for νe + anti-νe), thecontribution of K decay in the neutrino production is more important. There, the ratio
depends more on the K production cross sections and the uncertainty of the ratio is
expected to be larger. A 20% uncertainty in the K/π production ratio causes at least a
few % uncertainties in the νµ + anti-νµ to νe + anti-νe ratio in the energy range of 10
to100 GeV. However, as seen from Figure 3.2, the difference in the calculated νµ +anti-νµ to νe + anti-νe ratio is substantially larger than that expected from theK/π uncertainty in the high energy region. The reason for this difference has not been

41

Figure 3.1 (a) The atmospheric neutrino energy spectrum calculated by several models. (b)The ratio of the calculated flux models.

Figure 3.2 The flux ratio of νµ + anti-νµ / νe + anti-νe as a function of neutrino energy.

42understood yet. As a consequence, an energy dependent uncertainty for the νµ + antiνµ to νe + anti-νe ratio above 5 GeV neutrino energy is assumed based on thecomparison of the three flux calculations. It is taken to be 3 % at 5 GeV and 10% at
100 GeV. Figure 3.3 shows the calculated flux ratios of νµ to anti-νµ and νe to anti-νe.
The calculations agree to about 5% for both of these ratios below 10 GeV. However,the disagreement gets larger above 10 GeV as a function of neutrino energies. The
systematic errors in the neutrino over anti-neutrino ratio are assumed to be 5 % below
10 GeV and 10% at 100 GeV for the νe to anti-νe ratio, and 5 % below 10 GeV and
25% at 100 GeV for the νµ to anti-νµ. Figure 3.4 shows the zenith angle dependenceof the atmospheric neutrino fluxes of downward-going neutrinos are lower than those

Figure 3.3 The flux ratios of νµ to anti-νµ and νe to anti-νe versus neutrino energy. Solid,
dashed and dotted lines show the prediction by Honda, Bartol and Fluka flux.

43Figure 3.4

The flux of atmospheric neutrino versus zenith angle. Solid, dashed and dottedlines show the prediction by Honda, Bartol and Fluka flux.

of upward-going neutrinos. This is due to the deflection of primary cosmic rays by thegeomagnetic field, roughly characterized by a minimum rigidity cutoff. For neutrinoenergies higher than a few GeV, the calculated fluxes are essentially up-down
symmetric, because the primary particles are more energetic than the rigidity cutoff.
The enhancement of the flux near horizon for low energy neutrinos is a feature
characteristic of full three dimensional flux calculations [6–14]. However, in SuperKamiokande, the horizontal enhancement cannot be seen in the lepton zenith angle
distribution, due to the relatively poor angular correlation between neutrinos andleptons below 1 GeV. The uncertainties in the up-down and vertical-horizontal ratiosof the number of events are estimated by comparing the predicted ratios of various

44flux models. These uncertainties generally depend on the energy and the neutrino
flavor. The uncertainty in the up-down event ratio is about 1 to 2% in the sub-GeV
energy region and is about 1% in the multi-GeV energy region. The main source of
the uncertainty in the vertical-horizontal ratio around a GeV is the size of thehorizontal enhancement of the flux due to the three dimensional effect; theuncertainty is estimated to be less than a few percent. In the higher energy region,
where upward through-going muons are relevant, the largest source of the uncertainty
in the vertical-horizontal ratio is the K production cross section. We assume that theπ /K production ratio uncertainty is 20% in the whole energy region. The uncertainties
in the zenith angle and energy distributions due to the π /K production uncertainty areincluded in the systematic errors in the analysis. This error is most important forhigher energy neutrinos. For example, the vertical-horizontal uncertainty for upward
through-going muons due to the π/K production uncertainty is estimated to be 3%[19]. Figure 3.5 shows the zenith angle dependence of the atmospheric neutrino fluxes
for higher energy region observed as upward muons in Super-Kamiokande. The flight
length of neutrinos is an important ingredient in the analysis of neutrino oscillation.
For neutrinos passing a great distance through the Earth, the flight length can be
accurately estimated. However, for horizontal and downward going neutrinos, theheight of production in the upper atmosphere must be accurately distributed by theMonte Carlo method. Figure 3.6 shows the distribution of calculated flight length ofneutrinos for vertically down-going and horizontal directions. In summary of theatmospheric neutrino flux, we remark that, the νµ + anti-νµ to νe + anti-νe flux ratio ispredicted to an accuracy of about 3% in the energy region relevant to the data analysis
discussed in this thesis. The zenith angle dependence of the flux is well understood,and especially, above a few GeV neutrino energies, the flux is predicted to be updown symmetric.

45

Figure 3.5

The zenith angle distribution of the flux of upward-going atmospheric neutrino

observed as upward-going muon event in Super Kamiokande. Solid, dashed and dotted lines
show the prediction by Honda, Bartol and Fluka flux.

Figure 3.6

The calculated flight length of neutrinos for vertically down-going (0.95 < cosθ <1.0) and near horizontal-going (0.05 < cosθ < 0.10) directions.

46

3.2 Neutrino Interaction
The neutrino oscillation analysis relies heavily on detailed comparison of theexperimental data with the theoretical expectation. An important element of this is tosimulate the interaction of neutrinos from 10 MeV to 100 TeV with the nuclei ofwater, or in the case of upward going muon, the nuclei of the rock surrounding thedetector, assumed to be SiO2. In the simulation program, the following interactions
were considered.
- quasi-elastic scattering, νN → lN '- single meson production, νN → lN ' m- coherent π production, ν O → lπ O
1616

- deep inelastic scattering, νN → lN ' hadrons
Where, N and N’ are the nucleons, l is the lepton, and m is the meson, respectively.

3.2.1 Elastic and quasi-elastic scattering
The formalization of quasi-elastic scattering off a free proton, which was used in
the simulation programs, was described by Llewellyn-Smith[21]. For scattering off
nucleons in 16O, the Fermi motion of the nucleons and Pauli Exclusion Principle weretaken into account. The nucleons are treated as quasi-free particles using therelativistic Fermi gas model of Smith and Moniz[22]. The momentum distribution of
the nucleons were assumed to be flat up to the fixed Fermi surface momentum of 225
MeV/c. This Fermi momentum distribution was also used for other nuclear
interactions. The nuclear potential was set to 27 MeV/c.

47

3.2.2. Single meson production
Rein and Sehgal's model was used to simulate the resonance productions of single

π, K and η [23, 24]. In this method, the interaction is separated into two parts:

ν + N → l + N *,N* → m + N 'where m is a meson, N and N’ are nucleons, and N* is a baryon resonance. The
hadronic invariant mass, W, the mass of the intermediate baryon resonance, isrestricted to be less than 2 GeV/c2. In addition to the dominant single π production,
K and η production is considered.
To determine the angular distribution of pions in the final state, Rein and Sehgal's
method for the P33 (1232) resonance were used. For the other resonances, thedirectional distribution of the generated pions is set to be isotropic in the resonancerest frame. The angular distribution of π+ has been measured for νp → µ pπ−+mode[25] and the results agree well with the Monte Carlo prediction. The Pauliblocking effect in the decay of the baryon resonance was considered by requiring thatthe momentum of nucleon should be larger than the Fermi surface momentum. Pionless delta decay is also considered, where 20 % of the events do not have the pion andonly the lepton and nucleon are generated [26]. The quasi-elastic and single meson
production models have a parameter (axial vector mass, MA) that must be determinedby experiments. For larger MA values, interactions with higher Q2 values (andtherefore larger scattering angles) are enhanced for these channels. The MA value wastuned using the K2K [27] near detector data. In our atmospheric neutrino interaction
Monte Carlo simulation, MA was set to 1.1 GeV for both the quasi-elastic and single
meson production channels, but the uncertainty of the value is estimated to be 10%.
Coherent single pion production, the interaction between the neutrino and the entire
oxygen nucleus, was simulated using the formalism developed by Rein and Sehgal
[28].

48

3.2.3. Deep inelastic scattering
In order to calculate the cross sections of deep inelastic scattering, the GRV94[29] parton distribution function was used. In the calculation, the hadronic invariant
mass, W, is required to be greater than 1.3 GeV/c2. However, the multiplicity of pions
is restricted to be greater than or equal to 2 for 1.3 < W < 2.0 GeV/c2, because single
pion production is separately simulated as previously described. In order to generate
events with multi-hadron final states, two models were used. For W between 1.3 and
2.0 GeV/c2, a custom made program [30] was used to generate the final state hadrons;
only pions are considered in this case. For W larger than 2 GeV/c2, PYTHIA/JETSET[31] was used.
Total charged current cross sections including quasi-elastic scattering, single meson
productions and deep inelastic scattering are shown in Figure 3.7.

3.2.4. Nuclear effects
The interactions of mesons within the 16O nucleus are also important for theatmospheric neutrino analysis. Basically, all of the interactions were treated by using
a cascade model.

The interactions of pions are very important because the crosssection for pion production is quite large for neutrino energies above 1 GeV and theinteraction cross sections for pions in nuclear matter is also large. In the simulation ofpion interactions in 16O, inelastic scattering, charge exchange and absorption wereconsidered. The procedure to simulate these interactions is as follows. The initial
position of the pion was generated according to the Woods-Saxon nucleon density
distribution [32]. The interaction mode is determined from the calculated mean freepath of each interaction. To calculate the mean free path, the model described bySalcedo et al. [33] was adopted. The calculated mean free path depends not only on
the momentum of pion but also on the position of the pion in the nucleus. If inelastic
scattering or charge exchange occurs, the direction and momentum of the pion are49

determined by using the results of a phase shift analysis obtained from π-N scatteringexperiments[34]. When calculating the pion scattering amplitude, the Pauli blocking
effect is also taken into account by requiring the nucleon momentum after interaction
to be larger than the Fermi surface momentum at the interaction point. The pion
interaction simulation was checked using data for the following three interactions:π

C scattering, π16O scattering and pion photo-production ( γ + C → π + X )
12

12

[37].

50

−

Figure 3.7

Cross section for neutrino (upper) and anti-neutrino (lower). Solid line shows thecalculated total cross section. The dashed, dot and dash-dotted lines show the calculated quasielastic, single-meson and deep-inelastic scatterings, respectively. Data points are from various
experiments.

51

Chapter 4
Atmospheric Neutrino Data Reduction
and Reconstruction
Atmospheric neutrinos observed in Super Kamiokande are classified into four
categories. The neutrinos, which were observed via charged current interactions with
nuclei in the 22.5 kiloton water fiducial volume, are classified into fully-contained(FC) events if all of the energy is deposited inside the inner detector of Super-K and
into partially-contained (PC) events if high energy muons produced by neutrino
interactions exit the inner detector. The neutrino energies of partially-contained events
are typically 10 times higher than that of fully-contained. Neutrinos can also undergocharged current interaction with the rock surrounding the detector and the outgoinghigh energy muons can intersect the detector. These neutrino induced muons in therock cannot be distinguished from the cosmic ray muons traveling downward.
However, the muons traveling in upward direction must be neutrino induced because
cosmic ray muon cannot traverse the depth of earth in upward direction. These
neutrino induced upward-going muons (UPMU) are classified into upward stopping
muon if the muons come to rest in the inner detector, and into upward through-going
muon if the muons pass through the entire inner detector volume. The neutrinoenergies of upward stopping muon events are about 10 GeV, which is roughly sameas that of partially-contained events. The upward through-going muon events are the

52

µ, e

µ

µ

ν

µ

ν

Fully Contained

Partially Contained

UpwardUpward through-

(FC)

(PC)

stopping µ

going µFigure 4.1

A schematic view of the four event categories of the atmospheric neutrinos

observed in Super Kamiokande.

most energetic. Their parent neutrino energy is about 100 GeV on average. A
schematic view of the four event categories are shown in Figure 4.1.
Three different reduction schemes were used to isolate fully contained events,
partially contained events, and upward going muons, respectively. The fully contained
and partially contained data sets shared a common set of good run selection criteria tochoose good quality of data and therefore have identical live-time. The upward going
muon data set relies mostly on fitting long muon track length. So it was less affected
to detector conditions and had looser data quality cuts and therefore has different livetime from fully contained or partially contained data. To separate fully contained andpartially contained events, a fast spatial clustering algorithm was applied to the outerdetector hits. If the number of hits in the largest OD cluster was less than 10 (16 in
SK-II), the event was classified as fully contained a (FC), otherwise, it was classified
as partially contained (PC). Figure 4.2 shows the distribution of the number of theouter detector hits in the largest OD cluster for data (dot) and MC (histogram). A
clear separation of FC and PC events is seen at 10 hits in SK-I. The criterion of theseparation was tuned to 16 in SK-II. Systematic uncertainty of the FC and PC separat-

53

FC/PC separation in SK-I

Number of outer detector hits

FC/PC separation in SK-II

SK-I

SK-II

Number of outer detector hits

Figure 4.2 The number of hits in the largest outer detector cluster, which was used toseparate the FC and PC event. The histogram is the MC prediction with neutrino oscillation(solid line) and no oscillation (dashed line). The upper figure is for SK-I MC and data. The
lower figure is with SK-II, which zoomed in the region around the separation. In SK-II, thecriteria of the separation was defined as 16.

54

-ion was estimated by scaling the number of outerdetector hits to match the

distribution among data and MC.

4.1 Fully Contained Events (FC)

4.1.1 FC data reduction
The Super Kamiokande event sample consists mainly of downward going cosmic
ray muons and low energy radioactivity from parents such as radon. Cosmic raymuons are easily discriminated by requiring little or no ring activity in the outerdetector. For atmospheric neutrino analysis, we consider the events with visibleenergy above 30 MeV, where visible energy (Evis) is defined as the energy of anelectromagnetic shower that gives a certain amount of Cherenkov light (for example,
a muon of momentum 300 MeV/c yields a visible energy of about 110 MeV). To
make the final FC data sample, the following five steps of data reduction criteria wereused.

(1) First reduction and (2) Second reduction
Simple and efficient criteria were applied in the first and second reduction steps
 PE300 > 200 p.e (100 p.e for SK-II)

Where, PE300 is the maximum total photo electrons collected in the inner detector
within 300 nsec sliding time window.
 PEmax/PE300 < 0.5

Where, PEmax is the maximum p.e. in any single PMT in inner detector.
 NHITA800 < 25

Where, NHITA800 is the number of hit PMT in outer detector (OD) within 800 nsec

55

sliding time window.
 TDIFF > 100 μsec

Where TDIFF is the time interval from the preceding event. This is to reject
electrons from stopping muon decays.

(3) Third reductionMore complex criteria were applied in the third reductions step with the help of
event reconstruction tools, for further rejection of cosmic ray muons and low energy
events.
 NHITAin < 10 and NHITAout < 10 (In SK-II, NHITAin < 16 and NHITAout < 16)

Where, NHITAin (NHITAout) is the number of hit OD PMTs within 8 m from theentrance (exit) point of muon track fitted using the inner detector light pattern.
 NHIT50 < 50 (25 for SK-II)

Where, NHIT50 is the number of hit ID PMTs in 50 nsec sliding time window.

(4) Fourth reduction
Additional selection criteria were used to eliminate spurious events, such as“flashing” PMTs that emit light from internal corona discharges. Flasher events were
removed by two different methods.
 Typical flasher events have broader PMT timing distributions than the neutrino

events. Events with broader timing distributions were eliminated.
 Since flasher events have a tendency to be repeated with similar spatial hit

distribution, the pattern information of observed charge was used on to eliminate
these events. A correlation parameter based on the charge pattern was calculatedwith other data events and a “matched” tag was assigned for highly correlated events.

56

A cut was applied based on maximum correlation value and the number of“matched” with other events.

(5) Fifth reduction
Two further event types were eliminated in the fifth reduction step.
 NHITAoff < 10

Where NHITAoff is the number of OD hits in 200 nsec sliding time window
preceding the trigger time (-8900 ~ -100 nsec). This cut eliminates decay electrons
from invisible cosmic ray muons that are below Cherenkov threshold in the innerdetector.
 Cosmic ray muons are removed using a more precise fitter and the same criteria as

the first cut of (3).

Finally, the vertex was required to be within a fiducial volume, 2 meters from thewall of the inner detector, and the visible energy was required to be greater than 30
MeV. Table 4.1 and 4.2 shows the event rate for each reduction step and theefficiency.

4.1.2 FC event reconstruction
The fully contained events, which passed through the data reduction, underwent a
series of reconstruction.

(1) Vertex FittingFirst, the vertex position was determined using PMT hit times. The point which
best fit the distribution of PMT time was obtained by adjusting it for the time of flightof the Cherenkov light and was defined as the initial vertex position. The vertex wasreconstructed again after particle identification was established, to correct for particle

57

Table 4.1 Event rate for each FC reduction step

events/day~ 106

Raw data
1

2

3

 Low energy cut (total p.e.)~ 3000 Decay-e cut (time interval from previous trigger) Cosmic ray µ cut (OD hits in 800 ns time window)~ 200

 Low E & Flasher cut (Max p.e/Total p.e.) Cosmic ray µ cut (OD hits around reconstructed entry or exit)
 Low E cut (ID hits in 50 ns time window)

4  Flasher cut (broad time distribution and pattern algorithm)
5

 Invisible cosmic ray µ cut (OD hits in earlier time) Cosmic ray µ cut (OD hits around precisely reconstructed entry or exit)

Fiducial volume and visible energy cut

~ 18
~ 16
~8

Table 4.2 The efficiency of each reduction step of FC events.

Reduction step

~ 45

SK-II (%)

SK-I (%)

1st reduction

99.92

99.95

2nd reduction

99.89

99.94

3rd reduction

99.71

99.85

4th reduction

99.39

99.17

5th reduction

99.32

99.15

99.17

97.59

Fiducial volume andvisible energy cut

58track length. This final vertex fitter worked well for both SK-I and SK-II. However, it
was found that the reconstruction accuracy of the initial vertex fit was deteriorated in
SK-II. So in order to keep comparable performance with SK-I, new algorithm to fit
Cherenkov edge was adopted and finally the performance of the initial vertex fit of
SK-I and SK-II had become comparable. The final vertex resolution was estimated to
be about 30 cm for single ring fully contained events in both SK-I and SK-II.

(2) Ring Counting
After an initial ring direction and vertex were found, other possible rings weresearched using a Hough transform [38] based technique. The technique is iterative. A
second ring was searched for by choosing possible ring directions based on the Houghmap, and a likelihood technique determined if a second ring from the list of thepossible rings is more consistent with the data than just one ring. If a second ring is
judged to be more consistent, then this procedure was repeated as many as necessary
until finally no further rings were necessary to fit the data. Figure 4.3 and 4.4 show
the likelihood difference between the 2 ring assumption and a 1 ring assumption for
SK-I and SK-II, respectively. In these figures, the oscillation at the best fit parameters(sin 22θ = 1.0, Δm2 = 2.1x10-3 eV2) was taken into account for the Monte Carloprediction. Hereafter, if there is no mention, Monte Carlo prediction is withoscillation effect. A cut was made at likelihood difference of 0 to separate single andmulti ring events. The efficiency for identifying charged current (CC) quasi-elastic νe(νμ) events as single-ring was 93.2 (95.8) % for SK-I and 92.6 (95.8) % for SK-II.
(3) Particle Identification
To determine the identity of the final state particles, a particle identification
algorithm was applied which exploited systematic differences between Cherenkov
rings produced from different particles. Cherenkov rings from electron or gamma-ray
exhibit diffuse light distribution because of electromagnetic cascades and multiple
scattering and are called e-like events. On the other hand, Cherenkov rings from muonor charged pion give sharper ring edge and are called µ-like events. Figures 4.5 and

59

450

Super-Kamiokande

Super-Kamiokande

1201489.2 days

400
Single-ring
6918

350

1489.2 days
Multi-ring
2361

Single-ring
1397

100

300

Multi-ring
1504

80

250
60
200
150

40

100
20
50
0

-50

-40

-30

-20

-10

0

10

20

30

40

50

0

-50

-40

Ring counting likelihood, Sub-GeV

-30

-20

-10

0

10

20

30

40

50

Ring counting likelihood, Multi-GeV

Figure 4.3 The distribution of the likelihood difference between a single-ring and multi-ringassumption for sub-GeV (left) and multi-GeV (right) FC events of SK-I. The points show thedata and the histograms show the Monte Carlo.
200

Super-Kamiokande II Preliminary

225

Super-Kamiokande II Preliminary

627.1 days

627.1 days

175

200

Single-ring
2979

Multi-ring
1011

Single-ring
558

150

Multi-ring
642

175

150

125

125

100

100
75

75
50

50
25

25

0

0

-50

-40

-30

-20

-10

0

10

20

30

40

50

-50

-40

-30

-20

-10

0

10

20

30

40

Ring counting likelihood, Multi-GeV

Ring counting likelihood, Sub-GeV

Figure 4.4 The distribution of the likelihood difference between a single-ring and multi-ringassumption for sub-GeV (left) and multi-GeV (right) FC events of SK-II. The points show thedata and the histograms show the Monte Carlo.

60

50

4.6 show observed single-ring e-ike and μ-like events, respectively. The opening
angle of the Cherenkov cone, which depends on β(≡v/c), was also used to separate
e-like and µ-like events at low momenta. The validity of the method was confirmed
by a beam test experiment at KEK [39]. The misidentification probabilities for singlering e-like and µ-like events were estimated using simulated charged CC (Charged
Current) quasi-elastic neutrinos and summarized in the Table 4.3. The distribution of
the likelihood variable used to discriminate single-ring e-like and µ-like events areshown for both the data and MC for the sub-GeV and multi-GeV samples in Figure
4.7 and 4.8 for SK-I and SK-II, respectively. In both of these cases, there is a clear
separation of the likelihood variable.

(4) Precise Vertex Fitting & Momentum Determination
Next, the Cherenkov rings were re-fit taking into account the expected light pattern
given by the particle identification, and in the case of single-ring events, a specialized
event fitter was applied. After the rings were re-fit, the total photo-electrons wereapportioned to all of the rings. The momentum of a particle was determined from thetotal number of p.e.s within a 70° half-angle of the cone relative to the track direction,with corrections for light attenuation and PMT angular acceptance.
As a final procedure, the final energy and the angle information of the ring werechecked to remove the rings which are not real most likely.

4.2 Partially Contained Events (PC)
4.2.1 PC data reduction
The data reduction for PC events is different from the reduction for FC events in
respect of being required the presence of OD activity. Because these extra OD hits,
which result from the exiting particle, a simple criterion based on the number of hits
in the OD tubes could not be applied to reject cosmic ray background. The reduction

61

NUM
22
RUN
21916
EVENT
286316DATE 2003-Mar-12
TIME
23:50:46
TOT PE: 37351.1
MAX PE:
91.7
NMHIT :
4717
ANT-PE:
89.6
ANT-MX:
20.6
NMHITA:
43

RunMODE:NORMALTRG ID :00000011
T diff.:0.342E+05us: 34.2ms
FSCC:
A0027F90
TDC0: 8896.2
Q thr. :0.0
BAD ch.: masked
SUB EV : 0/ 2

Figure 4.5

An event display of an observed single-ring e-like event.
NUM
18
RUN
21804
EVENT
5052534
DATE 2003-Feb-20
TIME
12:46:22
TOT PE:
7087.7
MAX PE:
73.5
NMHIT :
1632
ANT-PE:
35.8
ANT-MX:
5.5
NMHITA:
40

RunMODE:NORMALTRG ID :00000011
T diff.:0.301E+05us: 30.1ms
FSCC:
E8027F90
TDC0: 8910.0
Q thr. :0.0
BAD ch.: masked
SUB EV : 0/ 3

Figure 4.6 An event display of an observed single-ring μ-like event.

62

Table 4.3

The misidentification probabilities for single-ring e-like and µ-like events

estimated with simulated CC quasi-elastic neutrinos

SK-II

SK-I

e μ

0.9 %

0.8 %

e

1.5 %

0.9 %

e μ

0.4 %

0.4 %

e

1.1 %

1.4 %

Sub-GeVμMulti-GeV

μused to identify partially contained events is described in this section. The description
is based on SK-I data reduction. Though SK-II has basically similar reduction step,some critical values were tuned to keep similar event rate and efficiency with SK-I.

(1) First reduction
 PEtot > 1000 p.e.

Where, PEtot is the total photo electrons observed in the ID. This criterion correspondto muons (electrons) with momentum less than 310 (110) MeV/c. By definition, an
exiting particle in the PC sample must have reached the OD from the inner fiducial
volume, and so must have had a minimum track length of about 2.5 m, whichcorrespond to muons greater than 700 MeV/c in momentum.

(2) Second reduction
 TWIDA < 260 nsec

63

300

200

Super-Kamiokande

Super-Kamiokande
1489.2 days
250

1489.2 days

175µ-like3227

e-like3353µ-like651

e-like
746

150

200
125
100

150

75

100
50

50
25

0

-10

-8

-6

-4

-2

0

2

4

6

8

0

10

-30

-20

-10

0

10

20

30

PID likelihood, Multi-GeV, 1-ring event

PID likelihood, Sub-GeV, 1-ring event

Figure 4.7 The distribution of particle identification likelihood for sub-GeV (left) and multiGeV (right) FC single-ring events of SK-I. Dots show data and histograms show Monte Carlo.

Super-Kamiokande II PreliminarySuper-Kamiokande II Preliminary

225

627.1 days

627.1 days

100

200

µ-like1370

e-like1447µ-like278

e-like
280

175
80

150

60

125

100
40

75

50
20

25

0

0

-10

-8

-6

-4

-2

0

2

4

6

8

10

-30

-20

-10

0

10

20

PID likelihood, Multi-GeV, 1-ring event

PID likelihood, Sub-GeV, 1-ring event

Figure 4.8 The distribution of particle identification likelihood for sub-GeV (top) and multiGeV (bottom) FC single-ring events of SK-II. Dots show data and histograms show Monte
Carlo.

64

30

Where, TWIDA is the width of the time distribution of hits in the OD.
 NCLSTA < 2

Where NCLSTA is the number of hit clusters in the OD.
These cuts eliminated through going muons, which left two well separated clusters in
the OD.
 NHITACmin < 7

Where, NHITACmin is the minimum number of hit PMTs among top (or bottom) andside regions in the OD hit cluster. Muons which clipped the edges of the detector
were eliminated by this cut.
 PE200 > 1000 p.e.s

Where, PE200 is the number of photoelectrons within 200 cm from OD cluster.
Cosmic ray muons which entered and stopped in the inner volume of the detector
were eliminated by this cut. This cut did not remove PC neutrino events because PC
events produce large numbers of photoelectrons (typically 3500 p.e.'s in SK-I) in theregion where the particle exited.

(3) Third reduction
In third step, a simple vertex fit and weighted direction estimation were used.
 NHITAin≤ 10 (16 for SK-II)

Where, NHITAin is the number of OD hits within 8 m of the back-projected entrance
point.
 Flasher events were removed by using their broader timing distribution feature.

(4) Fourth reduction65 dpfitdPMT > -0.8

Where, dpfit is the reconstructed direction by point fit and dPMT is the direction fromthe reconstructed vertex to the earliest hit ID PMT. This cut eliminate backgroundmuons which left few or no entrance hits in the OD.
 DCORN > 1.5 m

Where, DCORN is the distance from the reconstructed vertex to the corner of ID.
This cut removes the remaining corner clipping muons.
 Fitted muon track length > 30 m (if goodness of through-going muon fit > 0.85)(5) Fifth reduction
In the last reduction step, various remaining background events were eliminated byseveral selection criteria.
 NHITAclust ≥ 10 (16 for SK-II)

Where, NHITAclust is the OD hits in the most highly charged cluster in OD. Fullycontained events were eliminated by this cut.
 PEtot > 3000 p.e

This cut is to get rid of low energy background events. 3000 total p.e.s in the innerdetector correspond to 350 MeV of visible energy, well below that of any exitingmuon.
 DISTclust > 20 m

Where, DISTclust is the distance between the highest charge OD hit cluster and thesecond highest one. Clusters in the OD were searched for again with the sameclustering algorithm used in the 2nd reduction step but with different clusteringparameters. Some obvious through-going muons were removed by this cut.

66 The most remaining background events so far are due to the imprecision of the fast

fitters we used. A precise fitting algorithm was then applied to obtain more accurate
information on ring direction and vertex position. With much more accurate
information of the event, most remaining through-going and stopping events wereable to be eliminated based on their distinct geometry and OD signatures.
 Some through-going muons have a very special geometry. They passed through thetank vertically along the wall of ID. These events were eliminated by counting thenumber of p.e.s. and hits within the region defined by an 8 m sphere around the top
and bottom fringes and checking the time interval between the average timings ofthose top and bottom hits.
 Cable hole muon cut

There are four holes on the top of OD through which cables run and which make OD
inefficient to the penetration of cosmic ray muon. To compensate the inefficiency forthe holes, scintillation counters, which are called veto counters, were installed overthe holes on the top of the detector. Events with a veto counter hit and consistency
with cosmic ray muon were eliminated.

After the final reduction step, events were scanned by physicists to check the data
quality. However, no event is rejected based on the scanning. Finally, the vertex was
required to be within a fiducial volume, 2 meters from the wall of the inner detector.
The final event sample is an almost 100 % pure sample. The backgroundcontamination has been estimated to be about 0.2 %. Table 4.4 and 4.5 shows theevent rate for each PC reduction step and the detection efficiency for both SK-I and
SK-II.

4.2.2 PC Event Reconstruction

67

The partially-contained events were reconstructed using inner detector PMT
information by similar vertex, direction fit and ring-counting algorithms to the one offully contained events. For some PC events, however, the direction fit was slightly
modified and some criteria were tuned for SK-II PC event. The description of thecriteria in this section is based on SK-I event reconstruction.

Table 4.4 Event rate for each PC reduction step

events/day
Raw data~ 106

1  Low energy cut (total p.e.)~ 20000

 Cosmic ray µ cut (OD time distribution and # of OD hit cluster)2  Edge clipper cut (OD cluster topology)~ 3500 Stop µ cut (low ID p.e. near OD cluster)3

 Cosmic ray µ cut (OD hits around back-projected entry)
 Flasher cut (broad time distribution)~ 250

 Low OD activity cosmic ray µ cut (compare earliest ID hit & back-4

projected entry)
 Edge clipper cut (fitted vertex position near ID corner)~ 36

 Through-going µ cut (long fitted track length) FC cut (OD hits in highly charged cluster) Low E cut (total p.e.) Through-going µ cut (# of OD hit cluster with different clustering

5parameter from 2nd reduction)~1

 Stop or through-going µ cut (precise fitter info)
 Cosmic ray µ cut passing through vertically along ID wall
 Cable hole µ cut (veto counter)

Fiducial volume and total p.e. cut

~ 0.6

68

Table 4.5 The detection efficiency of PC events at each reduction step.

SK-II (%)

SK-I (%)

98.7

99.0

2 reduction

94.2

94.2

3rd reduction

93.1

93.1

4th reduction

85.8

87.9

5th reduction (Junk cut)

83.2

85.8

Reduction step
1st reductionnd

In case that the number of tubes in the largest OD cluster was greater or equal to 20
and there is no clear ID exit point or ID is judged to be saturated, outer detector
spatial information was used rather than the results of the precise fitting algorithm.
‘No clear ID exit point’ and ‘Saturated ID’ were defined like following.
 No clear ID exit point

If the ID PMT nearest to the projected ID exit point of the fitted track with more than
200 p.e.'s was more than 200 cm away, then it was regarded there is no clear ID exit
point.
 Saturated ID

If there were more than 800 ID PMT hits with more than 200 p.e.'s in each of them, it
was judged to be saturated.

The estimated vertex position resolution for PC events was 64 cm. The angularresolution for the penetrating particle in a PC event was estimated to be 2.8°. Finally,
the fiducial volume cut was applied. The event rate in the fiducial volume was 0.61
events/day and 0.52 events/day in SK-I and SK-II, respectively without taking into
account Monte Carlo efficiency and background.

694.3 Upward going muon (UPMU) data
The upward-going muons observed in Super-K are classified into two categories,that is, upward stopping muon events having only an entrance signal in the OD and
upward through-going muon events having both entrance and exit signals in the OD
as explained at the beginning of this chapter. The criterion used to separate FC andPC events was used to determine an event entry or exit signal. 10 (16 for SK-II) OD
hits within 8 m of the reconstructed entry or exit point were defined as a muon entry
or exit signal. We required the geometrical trajectory of through-going muons to begreater than 7 m in the inner detector, and we imposed a momentum threshold of 1.6
GeV to stopping muon, which is equivalent to 7 m path length. The purpose of thedata reduction is to isolate the upward muon events and the horizontal muon events,
which are needed for background estimation, to provide a classification of stopping or
through-going muon type, and to reject the background from cosmic ray muons andnoise such as flashing PMTs. Figure 4.9 shows the overall procedure to isolate upmu
events with rough event rate for each step. First, the very clear background was
rejected in the upmu reduction step and then the output of the reduction was passed tothe precise fitter. With the fitting information the precise fitter gives, downward
direction cut, 7 m track length cut for through-going muon and 1.6 GeV momentum
cut for stopping muon were applied. Finally, the events passed through the precise
fitter cuts were scanned by physicists and the remaining backgrounds were rejected.

4.3.1 UPMU data reduction
A charge cut of 8,000 ≤ Q (total ID charge) < 1,750,000 (3,000 ≤ Q < 800,000
for SK-II) p.e.s was applied. At very high ID charge corresponding to ≈ 1, 750,000(800,000 for SK-II) p.e.s the ID electronics becomes saturated causing the muon
fitters fail. To isolate the rate of about one neutrino induced upward going muon per
day from the remaining background of about a million cosmic ray muons, we used70

Raw Data~ 106 events / day

Upmu Reduction~ 50 events / day

Precise Fitter~ 3 events / day

Eye Scan~ 1.4 events / day

Final Sample

Figure 4.9 The procedure to isolate upmu events.

seven different muon fitters. Some of these fitters were specialized to fit stoppingmuons, others were specialized for through-going muon events, and some of themwere specialized to fit muon events with Bremsstrahlung. The 7 fitters used are listed
below in the same order that an event pass through.

1.

Stopmu1st‘Stopmu1st’ is the first fitter used in the upmu reduction. It classifies topology ofmuon events into multiple muons, stopping muon and through-going muon using
only OD information and passes the events to the next fitter without rejecting any
event.

2.

Muboy‘Muboy’ is the only fitter which provides fitting results (entrance position anddirection, goodness of fit, track length, number of track, etc) as well as classifying

71

muon topology (through-going muon, stopping muon, multiple muon and cornerclipper). Most of the cosmic ray muon backgrounds were removed by this fitter.

3.

Stopmu2nd‘Stopmu2nd’ was specially designed to fit stopping muons and run on only those
event s which were classified as stopping muons by ‘Muboy’ or ‘Stopmu1st’.

4.

Thrumu1st‘Thrumu1st’ is a fast muon fitter specialized for through-going muon events.

5.

Fstmu

‘Fstmu’ is another fitter specialized for through-going muons.

6.

Thrumu2nd‘Thrumu2nd’ is a fitter specialized to fit multiple muons and Bremsstrahlung muon
events.

7.

NNfit‘NNfit’ is a fitter specialized for fitting Bremsstrahlung muon events.

The basic thing what the fitters do is to decide whether an event should be saved,
rejected or passed to the next fitter.
 The events to be saved:

if an event is classified as through-going or stopping muonby either ‘Stopmu1st’ or ‘Muboy’ in upward direction above fitter’s goodness
threshold, it is saved. If an event, which ‘Muboy’ passed to the next fitter,‘Stopmu2nd’, is the one classified as stopping muon by either ‘Stopmu1st’ or‘Muboy’, ‘Stopmu2nd’ saves it if it is classified as horizontal/upward or if it gives
no-fit but ‘Muboy’ classifies it as horizontal. If an event is classified through-gong72

or stopping by either ‘Stopmu1st’ or ‘Muboy’ and its direction is horizontal above
goodness threshold, the event is tagged as ‘horizontal’ and saved for backgroundcontamination study.
 The events to be rejected:

if an event is classified as downward direction above

fitter’s goodness threshold by any fitter, then the event is rejected. Beside, cornerclipper muon cut and double muon cut are applied, which are described below.

Conditions for corner clipper cut :-‘Muboy’ classifies the event as corner clipper above goodness
threshold.

-‘Muboy’ classifies the event as through-going or stopping and track
length is less than 4 m and direction is downward above goodness
threshold.

Conditions for double muon cut:-

Both ‘Muboy’ and ‘Stopmu1st’ classify an event as multiple muon.

 The events to be passed to the next fitter:

if an event is fitted with low goodness

value or the fitter cannot be run on the event, then it is passed to the next fitter.

Finally, if all fitter cannot fit it, then the event is rejected. However, if at least one
fitter classify it as horizontal or upward, then it is saved. All events after the upwardmuon reduction, which contain many background events, were then passed to theprecise fitter, which gives more accurate fitting results (direction, track length,
momentum etc) about muon. The precise fitter is described in the next section. Table
4.6 show the percentage of saved events and rejected events for each step in upmu
reduction estimated with SK-I Monte Carlo.

73

Table 4.6 Percentage of events saved or rejected at each fitter.

Fitter% of events saved% of events rejected

Stopmu1st

0%

3.4 %

Muboy

9.0 %

93.5 %

Stopmu2nd

20.8 %

0.1 %Thrumu1st

29.4 %

0.9 %

Fstmu

2.8 %

1.0 %

Thrumu2nd

2.9 %

0.9 %

NNfit

6.9 %

0.1 %

No fit by all fitter or not

28.2 %

0.1 %

4.3.2 Upmu Event Reconstruction
All events from the output of the upward-going muon reduction were passed to theprecise fitter. The accurate muon vertex and direction are determined by the precise
fitter and they are used for the neutrino oscillation analysis. Though the fitters in the
upmu reduction give fit direction and entrance point, they are just designed to rejectclear background and to reduce number of events to pass to precise fitter because
precise fitter take long time to fit events.
The basic algorithm is identical to that used for the vertex and direction
determination for single ring fully-contained and partially-contained events. The fitterassumes that the particle is a muon and the vertex position of the event is at the innerdetector surface. However, when the muon produces an energetic electro-magnetic
shower, the assumption of single non-showering muon does not give an accurate
direction. For these events, the information of OD hit is used to determine the particle
direction. The angular resolution of the fitter was estimated by comparing true muon
direction and fitted direction using Monte Carlo and it was about 1.0° for boththrough-going and stopping muons.

74Based on the information the precise fitter gives, we selected the events, whichwere fitted in upward direction, longer track length than 7 m in case of through-going
muon and higher momentum than 1.6 GeV in case of stopping muon. Table 4.7 shows
the upward-going muon detection efficiency after upmu reduction and precise fittercut estimated using Monte Carlo.

Table 4.7 Upmu detection efficiency estimated using Monte Carlo.

SK-II (%)

SK-I (%)

Total

100

100

Reduction

97.7

99.3

96.4

99.2

Reduction step

Precise fitter(with direction and track length cut)

4.3.3 Upmu Eye-Scanning
The events, which passed through the cuts from the precise fitter, still contain somejunk events not rejected by misfitting. To eliminate this kind of events, two
independent physicists scanned the events by event display program. Finally, anotherupmu expert compared the judgment from the two physicists and make final decision
about rejecting events. The eye-scanning is used just to reject the remaining
backgrounds and never change the fitting results like direction, vertex or thestop/through-going judgment decided by the fitter.

75

4.3.4 Background Estimation from Cosmic Ray Muon
The main background for upward going muon is cosmic ray muons. Most of these
events were rejected in the upmu reduction. After the reduction, the precise fitter gave
more accurate fitting information and rejected downward going muon. The remaining
backgrounds after the selection from the precise fitter were eliminated by eyescanning. Nevertheless, there may still exist some background in near horizontal
direction because of the finite fitter resolution and multiple scattering of muons in the
nearby rock. Figure 4.10 shows the zenith angle distribution of upward muoncandidates near the horizon (-0.1 < cosθ < 0.1 for stopping muon and -0.1 < cosθ <0.08 for through-going muon) for two different region (that is, think rock region andthick rock region) in azimuth. The region, which shows negligible cosmic ray
contamination in near horizontal direction, is the thick rock region. However, we cansee the other region, that is, thin rock region have non-negligible contamination. The
shape of zenith angle distribution of thin rock region above horizon was extrapolatedby the following fitting function

f [cos(θ )] = p1 + e[ p 2+ p 3 cos(θ )]The cosmic ray contamination was estimated form the second term in the right of theabove function.

76

10

Number of Events

Normalized Number of Events / φ=360°

103
Region (1) φ<60°,φ>310°
Region (2) 60°<φ<310°

10

2

4
Region (1) φ<60°,φ>240°
Region (2) 60°<φ<240°

10

3

10

2

10

10

Upwardthrough-going µ

1

upward stopping µ

1
-0.1

-0.05

0

0.05

-1

0.1

cosΘ

10

-0.1

-0.05

0

0.05

cosΘ

Figure 4.10 The zenith angle distribution of upward muon candidates near the horizon for
two different regions in azimuth made with SK-I data. The white circles and the black inverted
triangle are for thin and thick rock region, respectively.

77

Chapter 5
Atmospheric Neutrino Data

5.1 FC and PC Events
We have accumulated 1489.2 days of FC and PC data from May 17, 1996 to July16, 2001. These are referred as SK-I data and compared with statistically larger
samples of Monte Carlo events equivalent to 100 years of the SK-I detector exposure.
After 2001 accident and reconstruction, we have accumulated 627.1 days of FC and
PC data from Jan 16, 2003 to Mar 2, 2005. These are referred as SK-II data andcompared with Monte Carlo equivalent to 60 years of the SK-II detector exposure.
The Monte Carlo samples were processed by the same event selection and event
reconstruction steps as the real events.

Fully-contained events were divided into two

sub-samples according to the reconstructed visible energy. We refer to the event
sample below 1.33 GeV as sub-GeV, and above 1.33 GeV as multi-GeV. Fullycontained events were further divided into the events with single reconstructed
Cherenkov ring and events with more than one ring. All single-sing events wereclassified as either e-like or μ-like based on the PID result. Lower energy cuts wereapplied only to the single-ring sample, Pe > 100MeV/c for a e-like and Pμ >200MeV/c forμ-like.

In addition, multi-ring events were used to study theatmospheric neutrino flux. Table 5.1 summarizes the number of observed events in th-78

Table 5.1

The number of observed events in the sub-GeV and multi-GeV samples as well as

the expected number of events in the absence of neutrino oscillations.

SK-I

SK-II

Data

MC

Data

MC

Single-ringSub-e-like3353

2879.81447

1256.2GeV

μ-like3227

4212.81370

1782.8Multi-ring

2361

2791.71011

1199.1

Single-ringMulti-

e-like

746

680.5

280

296.4

GeV

μ-like

651

899.9

278

390.1

Multi-ring

1504

1891.6

642

767.0

911

1129.6

325

461.2

PC

-e sub-GeV, multi-GeV and PC samples as well as the expected number of events inthe absence of neutrino oscillations.
A double ratio defined as R≡(μ/e)DATA/(μ/e)MC was measured like following without
oscillation effect in Monte Carlo prediction for sub-GeV events.
SK-I : Rsub −GeV = 0.658 ± 0.016 ± 0.035+0.026

SK-II : Rsub −GeV = 0.667 −0.025 ± 0.038

A substantial fraction of muons in the multi-GeV energy range exit from the innerdetector and are detected as PC events. The partially-contained event sample isestimated to be 97 % pure CC νμ interactions, even without requiring any particle
identification or ring-number cuts. Therefore, we add FC single-ring and PC event
totals when calculating R in the multi-GeV range. The double ratio R in the multi-

79

GeV energy range was measured like following without oscillation in Monte Carlo.

+0.032

SK-I : Rmulti−GeV + PC = 0.702 −0.030 ± 0.101+0.057

SK-II : Rmulti−GeV + PC = 0.750 −0.052 ± 0.065

The double ratios of SK-I and SK-II agree well.
The up-down ratio, where U is the number of upward-going events (-1 < cosθ < 0.2) and D is the number of downward-going events (0.2 < cosθ < 1), is summarizedin Table 5.2, where θ is zenith angle. While the ratio for e-like events is the consistent
with 1, the μ-like up-down ratio for the multi-GeV data differs from 1 by more than
12 standard deviations.

Table 5.2 UP/DOWN ratio of fully-contained events.

UP/DOWN ratio

e-likeSubGeVμ-likee-likeMultiGeV

μ-like

SK-I

SK-II

SK-II / SK-I

1.11+0.04-

0.97+0.06-

0.88+/-0.06 (2.0σ)

0.04

0.06

0.77+0.03-

0.84+0.05-

0.03

0.05

0.96+0.09-

1.07+0.16-

0.08

0.14

0.48+0.05-

0.40+0.05-

0.04

0.04

80

1.09+/-0.08 (1.1σ)

1.11+/-0.18 (0.6σ)

0.87+/-0.15 (0.9σ)

5.2 UPMU events
The SK-I upward-going muon data used in this analysis were taken from May 1996
to July 2001. The detector live-time is 1645.9 days. The SK-II upward-going muon
data were taken from Jan 5, 2003 to Mar 2, 2005. The detector live-time was 609 days.
Figure 5.1 shows the event rates as a function of the elapsed days for SK-II upward
going muons. The event rates for these samples are stable. Table 5.3 summarizes thenumber of observed events in the upward-going muon data sample and Monte Carlo.
Upmu event rate are consistent between SK-I and SK-II for both data and MC.
Fewer upward stopping muons were observed than predicted, while the observed
number of upward through-going muons was consistent with the theoreticalprediction within the errors.

Figure 5.1

Event rate as a function of elapsed days of SK-II data. SK-I event rates wereshown as a reference by line.

81

Table 5.3 The number of observed events in the upward-going muon data sample and Monte
Carlo. The numbers in 100 years SK-I MC were normalized to SK-II MC live-time of 60 years.

Fitted event
class

Data
SK-I

MC
SK-II

Stop mu

417.7

170

(events / day)

(0.28)(0.28)

Thru mu

1841.6662(events / day)

(1.13)(1.09)

82

SK-I

SK-II

9504

9285

22374

22022

Chapter 6
Two

FlavorNeutrino

Oscillation

Analysis

6.1 Overview
In the previous atmospheric neutrino analysis [2-5], a significant deficit ofatmospheric muon neutrino interactions was observed compared to the expectationand it was interpreted as neutrino oscillation. Oscillation between electron neutrinos
and muon neutrinos cannot explain the data because the number of observed upwardgoing electron neutrinos is consistent with null-oscillation expectation. The other
hand, atmospheric νµ oscillation into ντ can explain the νµ disappearance well, as the
majority of the flux is below the 3.5 GeV neutrino energy, which is the threshold forcharged current τ production. So previously the atmospheric neutrino data wereanalyzed based on hypothesis of pure 2 flavor neutrino oscillation between muon
neutrino and tau neutrino and the best fit parameters of νµ ↔ ντ oscillation were
established using the whole SK-I atmospheric neutrino data [35]. Also in this study,
the pure 2 flavor oscillation was assumed but to measure the oscillation parameters
more precisely, the new strategy of binning, which is one of the oscillation analysis
procedures, was applied and the combined analysis of Super-Kamiokande-I and83

Super-Kamiokande-II data was accomplished.

6.2 Oscillation Analysis with New Binning
6.2.1 New binning
In the oscillation analysis, the data and Monte Carlo are classified by reconstructed
event type, momentum, and zenith angle. The fully-contained (FC) sample is divided
into sub-GeV and multi-GeV according to the visible energy (E vis < 1.33 GeV is for
sub-GeV and Evis > 1.33 GeV is for Multi-GeV), into single-ring and multi-ring by
the number of reconstructed Cherenkov rings, and into e-like and μ-like by particle
identification (PID) of the most energetic ring. Upward going muon events are
divided into upward stopping and upward through-going according to the energy
deposited in OD. Each category is divided further by momentum and zenith angle. In
previous oscillation analysis [35], total 180 bins were used: 150 for the FC sample, 10
for PC sample, 20 for upward going muon sample.
In new binning, the PC sample is divided further into OD stopping and OD
through-going according to energy deposited in OD, multi-ring e-like sample wereadded newly, and momentum of multi-GeV samples were divided more finely toincrease sensitivity to oscillation parameters. The number of total bins has become
370 in this new binning strategy: 270 for FC sample, 80 for PC sample and 20 for
upward-going muon sample. Figure 6.1 shows the schematic view of the old bins (180bins) and new bins (370 bins).

6.2.2 Oscillation analysis method with new binning
The number of events of data in each bin is compared with expectation and theagreement between data and expectation is estimated by the following χ2 statistic

84

using Poisson probability distribution.

χ =2∑ 2(N − N )+ 2N
370 



i =1 


exp
i

obs
i

obs
i

N iexp = N i0 ⋅ P(ν µ → ν τ ) ⋅ (1 +

∑

obs
44 
N

εj
ln  iexp  +   N i  j =1  σ j 2∑ f ⋅ε )
46

j =1

i
j

j

Where,
Niobs is the number of observed events in the ith bin,
Ni exp is the expected number of events taking into account the oscillation effect,

σi is the combined statistical uncertainties of the data and Monte Carlo simulation,
Ni 0 is the expected number of events predicted without neutrino oscillation,

εj is the systematic variation,fi j is the fractional change in the predicted event rate in the ith bin due to a variation of
the parameter εj.
In the estimation of Ni exp, systematic variations due to uncertainties in the neutrino
flux model, neutrino cross section model, and detector response were taken into
account. These uncertainties are represented by 46 parameters, εj, which are describedin the next section.
The second sum in the χ2 definition collects the contributions from variables whichparameterize the systematic uncertainties in the expected neutrino rates but among the46 parameters, only 44 parameters contribute because the two parameters (‘theabsolute normalization’ and ‘sample normalization of multi-ring e-like events’) are
regarded as free. During the χ2 calculation, the systematic parameters, εj are varied tominimize χ2 for each choice of oscillation parameters, sin22θ and ∆m2.
A global scan was made on a (sin22θ, log∆m2) grid searching the point minimizingχ2. Assuming a linear dependence of Ni exp on the each systematic parameter, εj which
minimize χ2 is obtained by calculating ∂χ2 / ∂εj = 0 for each of the parameters εj. As aresult, the minimization of χ2 is equivalent to solving the following equation for k=1,

85

46

∑

∑ ((1 + ( fε ) ...)N − N )f

370

 370
1

N iobs f ji f ki + 2 δ jk ε k =


σj
i =1

 i =1

2

obs
i

exp
i

i
j

where σj is the estimated uncertainty in the parameter εj .

Figure 6.1 The schematic view of the binning. Each bin in this figure is divided further by 10
zenith angle bin. The top is the bins used in the old analysis (180 bins) and the bottom is the
new bins (370 bins) used in this analysis.

86Since this equation has non-linear terms for our χ2 definition, an iteration method wasused to obtain the approximate solution. The following data set were used in thisanalysis as mention at the chapter 5.

SK-I

FC, PC

SK-II

Data

MC

Data

MC

1489.2 days

100

627.1 days

60 years

608.6 days

60 years

years
UPMU

1645.9 days

100
years

6.2.3 Systematic uncertainty terms used in the oscillation analysis
In this section, the 46 systematic parameter terms are described briefly. The
estimated 1 σ for each parameter are shown in Appendices with the best fit values
obtained from the combined analysis described in section 6.3.

Systematic uncertainties in atmospheric neutrino flux

 Absolute normalization ----- (1)Absolute normalization uncertainty comes from the uncertainty in the primarycosmic ray flux and hadronic interaction models. It was roughly estimated to be
20% below 100 GeV and 30% above 100 GeV. However, this parameter wastreated as a free parameter in the analysis because of its large uncertainty.

 Flavor ratioνµ / νe ratio (Eν < 5 GeV) ----- (2)νµ / νe ratio (Eν > 5 GeV) ----- (3)The flavor ratio is predicted from the independent flux calculations. At the low
energy, the uncertainty mainly comes from the estimation of pion spectrum in the87primary hadronic interactions, which affects decay probabilities of secondary
muons in atmosphere. At high energy, the uncertainty mainly comes from the K/πproduction ratio.

 Anti-neutrino / neutrino ratio

anti-νe / νe (Eν < 10 GeV) ----- (4)anti-νe / νe (Eν > 10 GeV) ----- (5)anti-νµ / νµ (Eν < 10 GeV) ----- (6)anti-νµ / νµ (Eν > 10 GeV) ----- (7)Systematic uncertainty of anti-neutrino to neutrino ratio comes from π+/π- ratio in
hadronic interaction in the flux calculation.

 up / down ratio ----- (8)The neutrino flux below a few GeV has up/down asymmetry because of therigidity cutoff by the geomagnetic field. The systematic uncertainty in theup/down ratio of atmospheric neutrino flux is derived from the treatment of thegeomagnetic field in the calculation.

 horizontal / vertical ratio ----- (9)The source of the uncertainties in horizontal/vertical ratio of atmospheric
neutrino flux depends on the neutrino energy. For relatively lower energy
neutrinos below 10 GeV, the uncertainties come form the difference in 3dimensional calculation method for each calculation. For higher energy neutrinos,
the uncertainties are derived from the generated K/π ratio in hadronic interactionsin the atmosphere. However, these uncertainties were treated as the correlated.

 K/π ratio ----- (10) Flight length of neutrino ----- (11)The uncertainty of flight length of neutrino is not related to flux but it affect on
oscillation probabilities of neutrinos. The systematic uncertainty for upward-

88

going neutrinos is negligible but downward-going and horizontal-going neutrinos
are affected largely by the uncertainty of the production height of the neutrinos.
The systematic uncertainty in the calculation of the production height andneutrino flight length is measured by changing the density structure of
atmosphere.

 Energy spectrum of primary cosmic ray ----- (12)The spectrum of primary cosmic ray is well fit to the form Eγ. The spectral index forthe primary proton flux used is γ=2.74. We assigned 0.05 for the uncertainty in theenergy spectrum index above 100 GeV in the primary cosmic ray energy spectrum.

 Sample-by-sample normalization(FC Multi-GeV) ----- (13)(PC + Upward-stopping muon) ----- (14)Actually, the uncertainty of flux estimated from independent flux model depends
on neutrino energy. Therefore, the uncertainties of absolute normalization andenergy spectrum were not enough. So the uncertainties of the flux for Multi-GeV
and PC+upward-stopping muon samples were introduced as sample-by-sample
normalization.

 Normalization of multi-GeV multi-ring e-like events ----- (15) Uncertainty from solar activity ----- (16)Primary cosmic ray flux is affected by solar activity changing by about 11 years
period. The effect of solar activity on the atmospheric neutrino flux calculation
was taken into account. +/- 1 year of uncertainty was introduced to themodulation of solar activity.

Systematic uncertainties in neutrino interaction

 MA in quasi-elastic scattering and single meson production ----- (17)89The value of MA is set to be 1.11 based on experiment. The uncertainty of thisparameter makes effect on angular correlation between incident neutrino and outgoing lepton.

 Model dependence

quasi-elastic scattering ----- (18)deep inelastic scattering ----- (19)

 Cross sectionquasi-elastic scattering ----- (20)single meson production ----- (21)deep inelastic scattering ----- (22)coherent π cross section ----- (23) Neutral current / charged current ratio ----- (24) Nuclear effect in 16O ----- (25)
In neutrino interaction with nucleons in 16O, the Fermi motion of the nucleons
and Pauli Exclusion Principle are taken into account. The nucleons are treated as
quasi-free particles using the relativistic Fermi gas model of Smith and Moniz.
The momentum distribution of the nucleons is assumed to be flat up to the fixedFermi surface momentum of 225 MeV/c and the nucleon potential is set to 27
MeV/c. These configurations and the calculation of kinematics in neutrino
interaction are thought to have uncertainty.

 Nuclear effect (π spectrum) ----- (26)The uncertainty is related to the interaction of π and nucleons in 16O.

 Charged current ντ interaction ----- (27)
The number of charged current ντ interactions in the FC single-ring e-likesamples was considered. The systematic uncertainty in the cross section of
charged current ντ interaction is estimated from τ analysis in Super Kamiokande.

90

Systematic uncertainties in event section

 Reduction efficiency

FC ----- (28)
PC ----- (29)
Upward-going muon ----- (30)

Systematic uncertainty of reduction efficiency was estimated by comparing the
distribution of each cut variable of the data with that of the Monte Carlo.

 FC/PC relative normalization ----- (31)This uncertainty is mainly from FC/PC separation made by OD hits. Fiducial
volume, non-neutrino background and hadron simulation also affect on thisuncertainty.

 Separation of upward-stopping and through-going muon ----- (32)This separation is done by the number of hits in OD PMTs within 8 m from thereconstructed exit point (NHITEX). The uncertainty was estimated by comparing
NHITEX distribution of the data with that of Monte Carlo.

 Separation PC stop/through events ----- (33) Hadron simulation (contamination of NC in single-ring µ-like events) ----- (34) Non-neutrino background

flasher for e-like ----- (35)cosmic ray muon in µ-like events ----- (36) Background contamination

multi-GeV single-ring e-like events ----- (37)multi-GeV multi-ring e-like events ----- (38)

 Background subtraction

upward through-going muon ----- (39)
upward stopping-going muon ----- (40)91The background subtraction of upward-going muon events in the most horizontal
zenith angle bin was described in 4.3.4.

Systematic uncertainties in event reconstruction

 Single-ring / multi-ring separation ----- (41)This uncertainty was estimated by comparing the evaluation function distributionof data with that of Monte Carlo.

 Particle identification

single-ring events ----- (42)multi-ring events ----- (43)The error for the particle identification was estimated by comparing the
likelihood distribution of the data with that of Monte Carlo.

 Energy calibration ----- (44)The absolute energy scale was estimated with several sources over a wide energy
range. The uncertainty of the absolute energy reconstruction was less than 2 %.

 Energy cut for upward-stopping muon ----- (45)The systematic uncertainty of an energy cut (P > 1.6 GeV), which is applied to
upward-stopping muon samples to isolate well reconstructed events, wasestimated by varying the energy scale by 2 %.

 Up/down symmetry of energy calibration ----- (46)The difference of the energy scale for upward-going and downward-going event
was measured using decay electrons from cosmic ray stopping muons.

926.2.4 Study of sensitivity to finer binning with SK-I MC
The effect of the new binning on the 2 flavor analysis was tested using Monte Carlobefore applying it to the real data. To do it, the number of events in each bin was
calculated using 100 years SK-I MC and was scaled to have the live-time of 5 years,
20 years and 80 years respectively. The numbers of real data events in each bin were
replaced with that of the scaled MC. To investigate sensitivity for various oscillation
parameters regions, four different set of parameters, that is, (sin22θ, Δm2) = (0.85,
5.6×10-3), (0.85, 1.0×10-3), (0.85, 2.5×10-3), (1.0, 2.5×10-3) were applied when
making the virtual data, which was made from the scaled MC. With these virtual data,
the 2 flavor oscillation analysis was done in the old binning (180 bins) and the new
binning (370 bins).
Figure 6.2 shows the 90% C.L. allowed regions obtained with the virtual data in
180 bins and 370 bins. It is clear that the analysis in the new binning constrain Δm2
more strongly for all oscillation parameters regions and live-time. With this clear
improvement in MC, we decided to apply this new binning to real data analysis.

6.2.5 Result with SK-I data and SK-II data separately
SK-I real data were analyzed in the new binning with the method described in 6.2.2.
The minimum χ2 value, χ2min = 376.5 / 367 DOF, was obtained at (sin22θ = 1.00, ∆m2= 2.5×10-3 eV2). Including the unphysical region (sin22θ > 1) in the scan, theminimum χ2 value, χ2min = 376.1 / 367 DOF is located at (sin22θ = 1.03, ∆m2 =2.5×10-3 eV2). The minimum χ2 in the unphysical region is lower than that in thephysical region by 0.44. The allowed oscillation parameter region corresponding to
the 68%, 90% and 99% confidence level are shown in Figure 6.3 (a) and the contourof 90% confidence level was compared with that of the old analysis using 180 bins inFigure 6.3 (b). Assuming no oscillation, that is, (sin22θ = 0, ∆m2 = 0), χ2 value is795.7 / 369 DOF. The 90 % C.L allowed region of the oscillation parameters is also

93

-2

-2

10

2

2

2

2

∆m (eV )∆m (eV )

10

-3

-3

10

0.7

0.8

0.9
1
2
sin 2θ

1.1

1.2(a) 180 bins

10

0.7

0.8

0.9
1
2
sin 2θ

1.1

1.2(b) 370 bins

Figure 6.2 Result of sensitivity study to binning effect. 90 % C.L allowed region obtainedwith the virtual data using 180 and 370 bins are shown for various oscillation parameter region
and live-time.

estimated by using FC single-ring sub-GeV below 400 MeV/c, FC single-ring subGeV above 400 MeV/c, FC single-ring multi-GeV, PC, FC multi-ring and UPMU
events separately. The results are shown in Figure 6.4, which compare the old andnew binning analysis. The finely binned samples, ‘Multi-GeV’,’PC’ and ‘Multi-ring’
make contribution to constrain ∆m2 region tightly. FC single-ring and UPMU events
give no effect because their binning was defined same as the old. Table 6.1
summarized the best fit parameters and the allowed range of the oscillation
parameters at 90% confidence level in old and new binning analysis. The same
analysis was also done with SK-II data and 60 year SK-II MC. Figure 6.5 show thecontours plot with all the samples combined and the 90% C.L. contour plot for eachsub-sample.

94

SK-I

-2

10

2

2

2

2

SK-I

-2

∆m (eV )∆m (eV )10

90 % C.L.

99% C.L.
90% C.L.
68% C.L.

370 bins
180 bins

-3

10

-3

0.7

0.75

0.8

0.85

0.9

0.95

10

1

0.7

0.75

0.8

2

0.9

0.95

1

sin 2θ

(a) 370 bins
Figure 6.3

0.852

sin 2θ

(b) 180 vs. 370 bins

The left plot (a) is the allowed oscillation parameter region corresponding to the

68%, 90% and 99% confidence level (C.L.) obtained with 370 bins. The right plot (b) is thecomparison of the results with 180 bins and 370 bins for 90 % C.L.

SK-I
2∆m (eV )222∆m (eV )

10

-2

-2

10

10

Combined

-3

-3

10

10

Sub-GeV low
Sub-GeV highMulti-GeV
PC
Multi-ring
Upward-going muon

Combined
Sub-GeV low
Sub-GeV highMulti-GeV
PC
Multi-ring
Upward-going muon

180 bins

370 bins

-4

-4

10

SK-I

-1

-1

10

0.4

0.5

0.6

0.7

0.8

0.9

1

1.1

1.2

1.3

2

1.4

sin 2θ

Figure 6.4

10

0.4 0.5 0.6 0.7 0.8 0.91

1.1 1.2 1.3 1.4
2
sin 2θ

90 % C.L allowed region for each sub-sample. The left is the result with 180 bins

and the right is 370 bins.

95

Table 6.1 Summary of the best fit parameters and the allowed range obtained with SK-I data.

χ2 / dof

sin22θ

∆m2 (eV2)

90% C.L

Physical

174.9 / 177

1.00

2.11×10-3

sin22θ> 0.92

Unphysical174.7 / 177

1.02

2.11×10-3

1.5×10-3 < ∆m2 < 3.3×10-3

Physical

376.5 / 367

1.00

2.48×10-3

sin22θ> 0.93

Unphysical376.1 / 367

1.03

2.48×10-3

2.0×10-3 < ∆m2 < 3.0×10-3

180 bins

370 bins

SK-II

-2

SK-II
-1

10

2

2

∆m (eV )10-2∆m2 (eV2)

10

-3

10

Sub-GeV low
Sub-GeV highMulti-GeV
PC
Multi-ring
Upward-going muon

99% C.L.
90% C.L.
68% C.L.
-3

10

0.7

Combined

-4

0.75

0.8

0.85

0.9

0.95

1

10

2

sin 2θ

Figure 6.5

0.4 0.5 0.6 0.7 0.8 0.91

1.1 1.2 1.3 1.4
2
sin 2θ

The allowed oscillation parameter region obtained with SK-II data. The left is thecontour obtained with all the samples combined and the right is 90% C.L allowed region foreach sub-sample.

966.3 Oscillation analysis with SK-I and SK-II combined
6.3.1 Combining strategy
SK-II data, which started to be taken from January, 2003, have different systematiceffect in detector response. Therefore we do not mix SK-I and SK-II data directly. To
combine SK-I and SK-II data in oscillation analysis, the number of bins wereincreased double for SK-II data. That is, total 740 bins were used in this analysis: 370
bins for SK-I and 370 bins for SK-II. The definition of 370 bins for SK-II is same as
the one of SK-I. Though the systematic errors related to detector response aredifferent between SK-I and SK-II, the systematic errors related to neutrino flux andneutrino interaction are identical for SK-I and SK-II. So the systematic errors ofneutrino flux and interaction were treated in common in the combined analysis and incase of systematic errors of detector response and the uncertainty of solar activity
were regarded as independent between SK-I and SK-II. In total, 64 systematic error
terms were used in the combined analysis: 14 for neutrino flux, 12 for neutrino
interaction, 19 for SK-I detector response, 19 for SK-II detector response and 2 forsolar activity of SK-I and SK-II. The same oscillation analysis as section 6.2.2 wasrepeated with increased bins and systematic error terms.

6.3.2. Result of the combined analysis
In the combined analysis of SK-I and SK-II, the minimum χ2 value, χ2min = 767.5 /737 DOF, is located at (sin22θ = 1.00, ∆m2 = 2.5×10-3 eV2). Including the unphysicalregion (sin22θ > 1) in the scan, the minimum χ2 value is obtained at (sin22θ = 1.01,∆m2 = 2.5×10-3 eV2). Contours corresponding to the 68%, 90% and 99% confidence
level are shown in Figure 6.6 with comparison of SK-I result for 90% confidence
level. Figure 6.7 shows the χ2 - χ2min distributions projected to sin22θ and ∆m2 axes, inwhich the minimum χ2 - χ2min values for each sin22θ and ∆m2 are plotted. Table 6.297summarized the best fit parameters and the allowed range of the oscillation
parameters at 90% confidence level in SK-I and the combined analysis. Assuming nooscillation, (sin22θ = 0, ∆m2 = 0), a χ2 value is 1306.8 / 739. The allowed neutrino
oscillation parameters obtained by using each six sub-samples are shown in Figure 6.8.
The allowed region contours from combined data show a little bit stronger constraint
than SK-I. The best fit systematic parameters, which were fitted during χ2 calculation,
are shown in Table A1~10 in appendix with their 1σ uncertainty.
The zenith angle distributions of the FC, PC and upward-going muon samples areshown in Figure 6.9 and 6.10. The data (dot) are compared with the Monte Carlo
expectation assuming no oscillations (box) and the best-fit expectation for νµ ↔
ντ oscillations (line). The νµ ↔ ντ oscillation hypothesis provides a consistent
explanation of the full data set.

Table 6.2 Summary of the best fit parameters and the allowed range obtained from the
combined analysis. SK-I result are shown for comparison.

χ2 / dof

sin22θ

∆m2 (eV2)

90% C.L

Physical

376.5 / 367

1.00

2.48×10-3

sin22θ> 0.93

Unphysical376.1 / 367

1.03

2.48×10-3

2.0×10-3 < ∆m2 < 3.0×10-3

Physical

767.5 / 737

1.00

2.48×10-3

sin22θ> 0.93

Unphysical

767.5 / 737

1.01

2.48×10-3

2.1×10-3 < ∆m2 < 3.0×10-3

SK-I

SK-I+SK-II

98

SK-I + SK-II-2

90 % C.L.

2
2

99% C.L.
90% C.L.
68% C.L.

SK-I + SK-II
SK-I

-3

10

SK-I + SK-II-2∆m2 (eV2)10∆m (eV )

10

-3

0.7

0.75

0.8

0.85

0.9

0.95

10

1

0.7

0.75

0.8

2

Figure 6.6

0.85

0.9

0.95

1

2

sin 2θ

sin 2θ

The left plot is the allowed oscillation parameters for νµ↔ ν with the combinedτ

analysis of SK-I and SK-II. Three contours correspond to 68 %, 90 % and 99 % C.L. allowed
region. The right plot is the comparison with the result of SK-I only analysis for 90 % C.L..

SK-I + SK-II

SK-I + SK-II

20

20
99% C.L.
90% C.L.
68% C.L.

17.5
15

15

χ2 - χ2min

12.5

2

χ - χ min

12.5

2

10
7.5

10
7.5

5

5

2.5

2.5

0
10

99% C.L.
90% C.L.
68% C.L.

17.5-3∆m (eV )
2

2

10

-2

0
0.8 0.85 0.9 0.95

11.05 1.1 1.15 1.2

sin22θ

Figure 6.7 χ2-χ2min distribution projected to sin22θ and
SK-II combined analysis.

99△m axes resulted from SK-I and2

SK-I + SK-II
-1

2

2

∆m (eV )

10

-2

10

-3

10

Combined
Sub-GeV low
Sub-GeV highMulti-GeV
PC
Multi-ring
Upward-going muon

-410

0.4 0.5 0.6 0.7 0.8 0.9

Figure 6.8
for νµ1

1.1 1.2 1.3 1.4
2
sin 2θ

90 % confidence level allowed oscillation parameter regions of each sub-sample↔ ν oscillations with SK-I and SK-II data combined.
τ

100

SK-I + SK-II
400
300

300

200

200

100

100

0
-1

Number of Events

400

Sub-GeV e-likeP < 400 MeV/c

400

-0.5

0

0.5

1

0
-1

Sub-GeV µ-likeP < 400 MeV/c

-0.5

0

0.5

1

Sub-GeV µ-likeP > 400 MeV/c

Sub-GeV e-likeP > 400 MeV/c

400

300
200

200
100
0
-1
200

-0.5

0

0.5

1

Multi-GeV e-like

0
-1
200

150

150

100

100

50

50

0
-1

-0.5

0

0.5

1

cosΘ

0
-1

-0.5

0

0.5

1

Multi-GeV µ-like

-0.5

0

0.5

1

cosΘ

Figure 6.9 The zenith angle distribution of fully-contained sub-GeV and multi-GeV sample
obtained by the combined analysis of SK-I and SK-II. The points show data, box histograms
show the non-oscillated Monte Carlo and the lines show the best fit expectation obtained fromthe combined analysis.

101

SK-I + SK-II
multi-ring µ-like

multi-ring e-like

150

200

100
100

50
0
-1

-0.5

0

0.5

1

PC stop

0
-1
300

-0.5

0

0.5

1

PC through

40
200
20

0
-1
150

100

-0.5

0

0.5

0
-1

1

Upward stopping µ

-0.5

0

0.5

1

Upward through-going µ

400
100
200

50
0
-1 -0.8 -0.6 -0.4 -0.2 0

0-1 -0.8 -0.6 -0.4 -0.2 0

cosΘ

cosΘ

Figure 6.10 The zenith angle distribution of multi-ring, partially-contained and upward-going
muon

sample obtained by the combined analysis of SK-I and SK-II. The points show data,

box histograms show the non-oscillated Monte Carlo and the lines show the best fit expectation
obtained from the combined analysis.

102

6.3.3. Oscillation analysis with pseudo data to check the likeliness of
the result obtained with real data
When we compared the allowed oscillation parameter region obtained from the
combined analysis with that of just SK-I (refer to Figure 6.6), there was almost no
improvement with the combined analysis. It was inferred that that was because theresult with just SK-II data was on rather bad side unlikely so it couldn’t contribute forthe improvement by the combined analysis.
To check this inference, 20 set of pseudo data, which have the same live-time withreal data and have the best-fit oscillation parameters (sin22θ = 1.0, ∆m2 = 2.5×10-3
eV2), were made using Monte Carlo for SK-I and SK-II, respectively. The oscillation
analysis was done with this 20 set of pseudo data and the allowed regions of ∆m2 at
90% C.L were compared with that of real data.
The Figure 6.11 shows the width of ∆m2 region at 90% C.L obtained with 20 set ofpseudo data and real data for SK-I and SK-II, respectively. The first bins in the graphof Figure 6.11 indicate the result of real data. We can see that the result of SK-I real
data is on the lucky side and SK-II is rather unlucky. Thus we conclude that becauseof the unlucky situation of SK-II, the combined analysis of SK-I and SK-II madealmost no improvement on the result of just SK-I.

103

SK-I
2
-03

2.25

Width of ∆m at 90% C.L. (× 10

2
1.75

8
7

2

2

Width of ∆m at 90% C.L. (× 10

-03

2

eV )

eV )

SK-II

1.5

6
5

1.25
4
1
3

0.75

Pseudo data

0.5

2
1

0.25
0

Pseudo data

00 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

Real SKSK-I data

Real SKSK-II data

Figure 6.11

The width of ∆m2 region at 90% C.L obtained with 20 set of pseudo data and real

data for SK-I (left) and SK-II (right), respectively to check the likeliness of the result of real
data.

104

Chapter 7
Conclusion

Atmospheric neutrino events observed in Super-Kamiokande have the energy range
approximately from 100 MeV to 10 TeV, and the neutrino flight-length from about 10
km to 13,000 km. These wide ranges of energy and flight length and high statistics
enable us to study neutrino oscillations phenomena. In this study, the 2 flavor neutrino
oscillation (νµ↔ ν ) analysis was done using 1489.2 days exposure of SK-I fullyτ

contained and partially-contained events, 1645.9 days of SK-I upward-going muon
events, 627.1 days of SK-II fully-contained and partially contained events, and 608.6
days of SK-II upward-going muon events.
SK-II detector has half density of ID (inner detector) PMTs compared to SK-I and
the PMTs are covered with acrylic covers to prevent chain reaction of PMT implosion.
As SK-II has different detector configuration and responses from SK-I, some cut
variables in data reduction were tuned and some reconstruction algorithm weremodified in SK-II. The systematic uncertainties related to SK-II data reduction andreconstruction were estimated again and SK-II yielded consistent detection efficiency
and oscillation analysis result with SK-I. To investigate any possible position
dependence of SK water quality, the additional light injectors (LIs) were installed in 7
different positions (one on the top of the tank, five on the barrel and one on thebottom) and started to inject light from 8 different positions of the tank since thebeginning of SK-II while only one injector on the top was used in SK-I. However, the105

position dependence of the water quality was not seen in the result obtained from theupgraded LI system. Therefore the position dependence effect of the water quality
was not applied in the neutrino data analysis.
The oscillation parameters (sin22θ, ∆m2) between νµ and ντ were constrained
utilizing the zenith angle and energy dependent deficit observed in atmospheric muon
neutrino data. Especially it was shown that the determination of ∆m2 allowed region is
improved much by dividing momentum bins of fully-contained multi-GeV and PC
samples further. Oscillation analysis was done with SK-I data, SK-II data andcombined data of SK-I and SK-II in the new finer binning, which was proved moresensitive to constraining the oscillation parameter, ∆m2. The results with the above
three kinds of data set were consistent each other. The best fit oscillation parameters
obtained with the combined data set are sin22θ = 1.00 and ∆m2 = 2.5×10-3 eV2 with
χ2/DOF = 767.5/737. The oscillation parameters are constrained as sin22θ> 0.93 and2.1×10-3 < ∆m2 < 3.0×10-3 at 90% C.L. with the combined data set. Zenith angle
distribution of the observed atmospheric neutrinos in Super Kamiokande is explained
well at the best fit expectation obtained from pure νµ↔ ν oscillation analysis.
τThe oscillation parameters can also be constrained by L/E dependent oscillation
analysis of FC and PC samples. This independent study yielded consistent result [40].
Though 2 flavor νµ↔ ν oscillation is dominant, other sub-dominant effect such asτθ13 and solar terms (θ12 and ∆m122) might be observable in the atmospheric neutrinos
in Super Kamiokande. We are making effort to detect this small effect. To achieve
that, the neutrino flux and interaction should be understood more accurately.

106Appendix Tables107

Table A1. Systematic errors in neutrino flux (1)(common to SK-I and SK-II)

Absolute normalization(νµ+anti-νµ)/(νe+anti-Eν < 5 GeV

σ (%)

Best-fit

Free

7.2

3.0

-2.7

ν e)
anti−νe/νeanti−νµ/νµ

Up/down

Eν > 5 GeV
Eν < 10 GeV
Eν > 10 GeV
Eν < 10 GeV
Eν > 10 GeV

a

5.0
5.0

b

5.0
5.0

c

2.7
2.0
-0.9
-0.8
2.3

e-like0.5

0.1µ-like0.8

0.2

e-like2.1

0.6µ-like1.8

0.5

e-like1.5

0.4µ-like

0.8

0.2

PC

0.4

0.1

Sub-GeV multi-ring µ

0.8

0.2

Multi-GeV multi-ring µ

0.7

0.2

< 400 MeV

> 400 MeV

Multi-GeVa

3.0

Energy dependent error linearly increasing with logEν from 3.0% at 5 GeV to 15% at 100

GeV
Energy dependent error linearly increasing with logEν from 5.0% at 10 GeV to 10% atb

100 GeVc

Energy dependent error linearly increasing with logEν from 5.0% at 10 GeV to 10% at

100 GeV108

Table A2
A2. Systematic errors in neutrino flux (2)(common to SK-I and SK-II)

Horizontal/verticalσ (%)

Best-fite-like0.3

0.1µ-like0.3

0.1

e-like1.2

0.5µ-like1.2

0.5

e-like2.8

1.1µ-like

1.9

0.8

PC

1.4

0.6

Sub-GeV multi-ring µ

1.5

0.6

Multi-GeV multi-ring µ

1.3

0.5

< 400 MeV

> 400 MeV

Multi-GeV

K/π ratio

20Lν (production height)

10

Energy spectrum
Sample-by-sample

a

-2.7

b

0.05

1.1

c

0.045

d

FC Multi-GeV

5.0

PC+upward stopping µ

5.0

20% uncertainty in K/p ratio in hadronic interaction in the atmosphere.

a

10% uncertainty in the atmospheric density structure.

b

0.05 uncertainty in the spectral index of the primary cosmic ray above 100 GeV.

c

109

-7.3
-5.8

Table A3. Systematic errors in neutrino interaction(common to SK-I and SK-II)σ (%)
MA in quasi-elastic and single-π

10

a
b

Best-fit
3.4

Quasi-elastic scattering (model dependence)

1.0

Quasi-elastic scattering (cross-section)

10

9.4

Single-meson production (cross-section)

10

-3.1

c

-0.8

Multi-π production (model dependence)

1.0

Multi-π production (cross-section)

5.0

-1.6

Coherent π production

30

1.0

NC/CC ratio

20

-7.8

Nuclear effect in

16

O

30

Energy spectrum of pions

1.0

Charged current ντ interaction

30

Hadron simulation

a

d
e

f

1.0

h

1.46

-6.9
0.36
1.7
-0.25

10% uncertainty in the axial vector mass, MA value.

Difference from the model, S.K.Singh and E.Oset, is set to 1.0.

bcDifference from the model, A.Bodek and U.K.Yang, is set to 1.0

d30% uncertainty in the mean free path of hadron in the 16O nucleus.

eDifference in the predicted pion energy spectrum by two interaction models (Neut andNuance) was taken.
30% uncertainty assumed for the CC ντ interaction cross-section.

fDifference from Fluka modelh110

Table A4
A4. Systematic errors in event section(For only SK-I)σ (%)

Best-fit

FC reduction

0.2

0.0

PC reduction

2.6

1.0

upward stopping µ

1.3

0.0

upward through-going µ

0.5

0.0

0.9

-0.4

e-like0.4

0.0µ-like0.1

0.0

e-like0.2

0.0µ-like

0.1

0.0

0.2

0.0

Upward stopping/through-going µ separation

0.4

0.0

Relative normalization of PC stop/through

9.6

4.7

BG subtraction of upward through-going µ

3.0

-0.2

BG subtraction of upward stopping µ

17

4.6

e

14

-0.6

e

20

-2.8

Detection efficiency

FC/PC separation
Non-ν BG

Sub-GeV

Multi-GeV

PC

Non-νe contamination in multi-GeV single-ring
Non-νe contamination in multi-GeV multi-ring111

Table A5. Systematic errors in event reconstruction(For only SK-I)σ (%)

Best-fite-like6.3

3.2µ-like2.4

1.2

e-like3.4

1.7µ-like1.3

0.7

e-like15.9

8.1µ-like6.2

3.2

Sub-GeV multi-ring µ

-3.7

-1.9

Multi-GeV multi-ring µ

-7.2

-3.7

e-like-0.6

0.1µ-like0.6-0.1

e-like-0.4

0.1µ-like0.4-1.0

Sub-GeV multi-ring µ

3.4

0.3

Multi-GeV multi-ring µ

4.7

0.4

Energy calibration for FC event

2.0

-0.2

Energy cut of upward stopping µ

1.1

-0.2

Up/down symmetry of energy calibration

0.6

0.0

Normalization of multi-GeV single-ring e

10

-0.05

Ring separation< 400 MeV

> 400 MeV

Multi-GeV

Particle identification

Sub-GeV

Multi-GeV

112

Table A6
A6. Systematic errors in event section(For only SK-II)σ (%)

Best-fit

FC reduction

0.2

0.0

PC reduction

2.6

-1.9

upward stopping µ

1.9

0.1

upward through-going µ

0.9

0.0

0.9

0.2

e-like0.4

0.0µ-like0.1

0.0

e-like0.2

0.0µ-like

0.1

0.0

0.2

0.0

Upward stopping/through-going µ separation

0.3

0.0

Relative normalization of PC stop/through

9.6

13.8

BG subtraction of upward through-going µ

3.0

1.1

BG subtraction of upward stopping µ

17

-11.4

Non-νe contamination in multi-GeV single-ring e

14

0.4

Non-νe contamination in multi-GeV multi-ring e

20

-1.8

Detection efficiency

FC/PC separation
Non-ν BG

Sub-GeV

Multi-GeV

PC

113

Table A7. Systematic errors in event reconstruction(For only SK-II)σ (%)

Best-fite-like7.3

0.4µ-like2.5

0.1

e-like3.0

0.1µ-like1.2

0.1

e-like7.4

0.4µ-like2.6

0.1

Sub-GeV multi-ring µ

-3.3

-0.2

Multi-GeV multi-ring µ

-3.5

-0.2

e-like-0.9

0.0µ-like0.9

0.0

e-like-0.3

0.0µ-like0.3

0.0

Sub-GeV multi-ring µ

16.8

10.2

Multi-GeV multi-ring µ

6.7

4.1

Energy calibration for FC event

2.1

-2.0

Energy cut of upward stopping µ

1.1

0.0

Up/down symmetry of energy calibration

0.5

-0.1

e

10

0.08

Ring separation< 400 MeV

> 400 MeV

Multi-GeV

Particle identification

Sub-GeV

Multi-GeV

Normalization of multi-GeV single-ring114

Table A8. Systematic errors in solar activity

σ (%)
Solar activity (SK-I)

0.2

Solar activity (SK-II)

0.5+/- 1 year of solar activity uncertainty was taken into account.

a,b

115

a

b

Best-fit
0.02
0.02

Bibliography[1] G. Barr et al., Phys. Rev. D39, 3532 (1989); V. Agrawal, et al., Phys. Rev. D53,
1313 (1996); T. K. Gaisser and T.Stanev, Proc. 24th Int. Cosmic Ray Conf.
(Rome) Vol. 1 694 (1995)[2] Y. Fukuda et al. (Super-Kamiokande), Phys. Lett. B436, 33 (1998), hepex/9805006.
[3] Y. Fukuda et al. (Super-Kamiokande), Phys. Lett. B436, 33 (1998), hepex/9805006.
[4] Y. Fukuda et al. (Super-Kamiokande), Phys. Rev. Lett. 82, 2644 (1999), hepex/9812014.
[5] Y. Fukuda et al. (Super-Kamiokande), Phys. Lett. B467, 185 (1999), hepex/9908049.
[6] G. Battistoni et al., Astropart. Phys. 12, 315 (2000), hepph/9907408.
[7] M. Honda, T. Kajita, K. Kasahara, and S. Midorikawa, Phys. Rev. D64, 053011(2001), hep-ph/0103328.
[8] Y. Tserkovnyak, R. Komar, C. Nally, and C. Waltham, Astropart. Phys. 18, 449(2003).
[9] Y. Liu, L. Derome, and M. Buenerd, Phys. Rev. D67, 073022 (2003), astroph/0211632.
[10] G. Battistoni, A. Ferrari, T. Montaruli, and P. R. Sala (2003), hep-ph/0305208.

116

[11] J. Wentz et al., Phys. Rev. D67, 073020 (2003), hepph/0301199.
[12] J. Favier, R. Kossakowski, and J. P. Vialle, Phys. Rev. D68, 093006 (2003),
astro-ph/0305460.
[13] M. Honda, T. Kajita, K. Kasahara, and S. Midorikawa, Phys. Rev. D70, 043008(2004), astro-ph/0404457.
[14] G. D. Barr, T. K. Gaisser, P. Lipari, S. Robbins, and T. Stanev, Phys. Rev. D70,023006 (2004), (Also, private communication with the authors. A slightly
different flux from the one discussed in the paper is used in the present analysis.),astroph/0403630.
[15] M. Honda, T. Kajita, K. Kasahara, and S. Midorikawa, Phys. Rev. D52, 4985
(1995), hep-ph/9503439.
[16] V. Agrawal, T. K. Gaisser, P. Lipari, and T. Stanev, Phys. Rev. D53, 1314
(1996), hep-ph/9509423.
[17] J. Alcaraz et al. (AMS), Phys. Lett. B490, 27 (2000).
[18] T. Sanuki et al., Astrophys. J. 545, 1135 (2000), astroph/0002481.
[19] P. Lipari, Nucl. Phys. B (Proc. Suppl.) 81, 159 (2001).
[20] B. Pontecorvo, JETP, 6, 429 (1958); Z. Maki, M. Nakagawa and S. Sakata, Prog.
Theor. Phys. 28, 870 (1962).
[21] C. H. Llewellyn Smith, Phys. Rept. 3, 261 (1972).
[22] R. A. Smith and E. J. Moniz, Nucl. Phys. B43, 605 (1972).
[23] D. Rein and L. M. Sehgal, Ann. Phys. 133, 79 (1981).
[24] D. Rein, Z. Phys. C35, 43 (1987).
[25] T. Kitagaki et al., Phys. Rev. D34, 2554 (1986).
[26] S. K. Singh, M. J. Vicente-Vacas, and E. Oset, Phys. Lett. B416, 23 (1998).
[27] M. H. Ahn et al. (K2K), Phys. Rev. Lett. 90, 041801 (2003), hep-ex/0212007.

117

[28] D. Rein and L. M. Sehgal, Nucl. Phys. B223, 29 (1983).
[29] M. Gluck, E. Reya, and A. Vogt, Z. Phys. C67, 433 (1995).
[30] M. Nakahata et al. (KAMIOKANDE), J. Phys. Soc. Jap. 55, 0 3786 (1986).
[31] T. Sjostrand, Comput. Phys. Commun. 82, 74 (1994).
[32] R. Woods and D. Saxon, Phys, Rev. 95, 577 (1954).
[33] L. L. Salcedo, E. Oset, M. J. Vicente-vacas, and C. Garcia-Recio, Nucl. Phys.
A484, 557 (1988).
[34] G. Rowe, M. Salomon, and R. H. Landau, Phys. Rev. C18, 584 (1978).
[35] Y.Asie et al. (Super-Kamiokande), Phys. Rev. D 71, 112005 (2005), hepex/0501064.
[36] S.Fukuda et al. The Super-Kamiokande Detector, Nucl. Instrum. Meth.
A501(2003) 418-462.
[37] D. Ashery et al., Phys. Rev. C23, 2173 (1981).
[38] E.Davis, Machine Vision : Theory, Algorithms, Practicalities (Academic Press,
San Diego, 1997).
[39] S.Kasuga et al., Phys. Lett. B374, 238 (1996)[40] Y.Ashie, et al., Phys.Rev.Lett. 93 (2004).
[41] M.B.Smy et al. (Super-Kamiokande), Phys. Rev. D69, 011104 (2004), hepex/0309011.
[42] S.N.Ahmed et al. (SNO), Phys.Rev.Lett. 92, 181301 (2004), nucl-ex/0309004.
[43] K.Eguchi et al. (KamLAND), Phys.Rev.Lett. 90, 021802 (2003), hep-ex/0212021.
[44] E.Aliu et al. (K2K), Phys.Rev.Lett. 94, 081802 (2005), hep-ex/0411038.

118

A Measurement of νµ↔ντ Oscillation Parameters
Using Atmospheric Neutrino Observed in Super
Kamiokande-I and Super Kamiokande-II서현관성균관대학교 물리학과Abstract수퍼카미오칸데에서ντ 진동을

가정하면

잘

검출된

대기

설명된다.

이

zenith angle과

에너지에

진동

분석하였다.

현상을끼치는 영향을
과

수퍼

따라서

감소되는

진동변환상수를

데이터를

중성미자의

연구에서는

연구하였고, 최적의 운동량카미오칸데-II진동변환상수는

뮤온

대기

관측사실을

부족은

뮤온

νµ

↔

중성미자가

이용하여

중성미자

운동량

binning이결정하는데,

binning에서 수퍼 카미오칸데-I

합하여

분석한

결과,

최적의

sin22θ = 1.00 과 ∆m2 = 2.5×10-3 eV2 에 위치하고, 그 때의

χ2/DOF는 767.5/737 이다. 90% C.L.에서 진동변환상수의 범위는 sin22θ > 0.93
과

2.1×10-3 < ∆m2 < 3.0×10-3 eV2 이다.