BOSTON UNIVERSITYGRADUATE SCHOOL OF ARTS AND SCIENCES

Dissertation

STUDIES OF NON-STANDARD EFFECTS IN ATMOSPHERIC NEUTRINO
OSCILLATIONS AT SUPER-KAMIOKANDEbyWEI WANG
B.E., Harbin Institute of Technology, 1995
M.S., Harbin Institute of Technology, 1998
M.A., Boston University, 2001

Submitted in partial fulfillment of therequirements for the degree ofDoctor of Philosophy
2007Approved byFirst Reader
Edward T. Kearns, Ph.D.
Associate Professor of Physics

Second Reader
Martin Schmaltz, Ph.D.
Assistant Professor of Physics

Third Reader
Ulrich Heintz, Ph.D.
Associate Professor of Physics

Acknowledgments
I must first express my deep gratitude to the North American graduate school system.
It is only under this system that I was given the chance to pursue a Ph.D. degree in a fieldthat I had little knowledge and background to start with. With its openness, I was able tobuild up my physics knowledge and eventually start advanced research.
Professor Ed Kearns, my academic advisor, is always a supporting mentor. Only withhis patience could I conquer the difficulties and get back onto the right track in the field of
research. His knowledge and understanding of experimental HEP gave me many enlightening moments. Professor Chris Walter who now is at Duke University helped me to start
my research. I learned many things from him and have benefited much from his cheerfulpersonality and high spirit. Thanks to Professor Mark Messier of Indiana University for
his initiation on the studies in this dissertation. His early work provides much guidance
for me in this field.
Working in the SciFi group of the K2K collaboration at the beginning of my Ph.D.
program taught me many things. I would like to thank Professor Sakuda for his kindly
support and guidance.
I would like to express my great appreciation to my other committee members: Professor Kenneth Lane, Professor Sidney Redner, Professor Ulrich Heintz and Professor Martin Schmaltz. I also include Professor James Rohlf, Professor Andrew Cohen, Professor
Sekhar Chivukula and Professor Michael El-Batanouny. I have taken classes and sought
other academic help from all of them. Some of them taught me when I knew little about
advanced physics; some of them entered my physics life later while I started to take on
more advanced topics. They are all great teachers and without their teaching, I could nothave accomplished the research presented in this dissertation.
I want to thank the people in the ATMPD group: Professor Kajita, Professor Suiii

lak, Professor Scholberg, Professor Habig, Professor Moriyama, Professor Okumura, Dr.
Nakayama and Dr. Ishitsuka. I learned many things from them and I also had wonderfultime working together with them.
Thanks to Professor Totsuka, Professor Suzuki, Professor Sobel, and Professor Stone;
they led and are leading a very successful collaboration within which I learned how to
do research and finished my Ph.D. program. Thanks to Professor Nakahata for his effort
toward taking care of the Super-Kamiokande detector.
Thanks to the Kamioka Observatory faculty members and staff for taking care of theexperiment and providing excellent technical support. Thanks to the SciFi group members
for the great experience of working together.
I am very much grateful to Professor Kostelecky of Indiana University, Professor Pakvasa and Professor Learned of Hawaii University for their valuable suggestions and comments during my study.
I would also like to thank the people in the physics department of Carleton University,
especially Professor Armitage and Professor Karlen, for their hospitality, kindness andteaching. I stayed only one year but I learned much.
Thanks to my fellow colleagues at Boston University. They are wonderful and helpful.
The fun moments we shared those years are going to be an important part of my lifelong
treasures.
At last, I would like to thank my parents for their constant encouragement since thevery first day I entered school. Without their support, it would have been impossible forme to come to this point.

iv

STUDIES OF NON-STANDARD EFFECTS IN ATMOSPHERIC NEUTRINO
OSCILLATIONS AT SUPER-KAMIOKANDE(Order No.

)

WEI WANG
Boston University Graduate School of Arts and Sciences, 2007
Major Professor: Edward T. Kearns, Associate Professor of Physics
ABSTRACT
Neutrino oscillation due to mass eigenstate mixing has become the standard theory
accounting for both solar and atmospheric neutrino data. This explanation indicates thatneutrinos have small but non-vanishing masses, which is a sign of new physics beyond theStandard Model. In this dissertation, we will compare the standard explanation with three
types of alternative theories using Super-Kamiokande (SK) atmospheric neutrino data.
The first type of non-standard theory involves sterile neutrinos. By using the neutral
current enhanced data samples of SK and by considering matter effect, we conclude it isunlikely that sterile neutrinos are responsible for SK atmospheric neutrino zenith angle
distributions. Furthermore, we study the allowance of sterile neutrino admixture in atmospheric neutrino mixing and find an admixture of 23% sterile neutrino is allowed at 90%
confidence level based on a 2+2 mass hierarchy model.
The second type of non-standard theory involves neutrino oscillation induced by violations of Lorentz invariance (LIV) and CPT symmetry (CPTV). The neutrino oscillations
induced by the temporal components of the LIV and CPTV terms in the minimal Standard
Model Extension (SME) have different energy and pathlength dependences compared tothe standard oscillation. Our analysis indicates that it is unlikely to explain SK atmospheric
neutrino data with the oscillation effects induced by the temporal components of the minimal SME separately. By treating LIV- and CPTV-induced oscillations as sub-dominant
effects, limits on symmetry-breaking parameters are established.
v

The third category of non-standard theory involves vanishing neutrinos caused by neutrino decoherence and neutrino decay. Our study shows that it is unlikely to explain SK
atmospheric neutrino zenith angle distributions using these two non-oscillatory models.
By treating them as sub-dominant effects, limits on these two types of new physics are set
based on several specific models.
Our study shows that the oscillation between muon neutrinos and tau neutrinos is thebest model explaining SK atmospheric neutrino data among the models we test. In mostcases, limits on new physics established in this study using SK atmospheric neutrino data
are the best currently available.

viContents
Acknowledgementsiii

Abstractv

Table of Contentsvi

List of Tablesx

List of Figuresxi

List of Abbreviations

xvi

1

Introduction

1

1.1

A brief history of the Standard Model neutrinos . . . . . . . . . . . . . .

1

1.2

Neutrino mass . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

1.3

Neutrino mixing and oscillation . . . . . . . . . . . . . . . . . . . . . .

4

1.3.1

Two-Flavor oscillation in vacuum . . . . . . . . . . . . . . . . .

6

1.3.2

Two-flavor oscillation in matter . . . . . . . . . . . . . . . . . .

7

Beyond neutrino oscillation . . . . . . . . . . . . . . . . . . . . . . . . .

9

1.4
2

Atmospheric Neutrinos

11

2.1

11

Primary cosmic rays . . . . . . . . . . . . . . . . . . . . . . . . . . . .
vii

3

2.2

Atmospheric neutrino energy spectrum . . . . . . . . . . . . . . . . . . .

14

2.3

Flavor ratios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

15

2.4

Neutrino pathlength . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

2.5

Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

20

Super-Kamiokande Experiment

22

3.1

Detector setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

23

3.2

Electronics and data acquisition . . . . . . . . . . . . . . . . . . . . . .

25

3.3

Detector calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.3.1

Water properties . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.3.2

The absolute energy calibration . . . . . . . . . . . . . . . . . .

27

Neutrino interaction simulation . . . . . . . . . . . . . . . . . . . . . . .

29

3.4.1

Elastic and quasi-elastic scattering . . . . . . . . . . . . . . . . .

30

3.4.2

Single meson production . . . . . . . . . . . . . . . . . . . . . .

30

3.4.3

Deep inelastic scattering . . . . . . . . . . . . . . . . . . . . . .

31

3.4.4

Nuclear effects . . . . . . . . . . . . . . . . . . . . . . . . . . .

32

3.5

The detector simulation . . . . . . . . . . . . . . . . . . . . . . . . . . .

34

3.6

Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

35

3.4

4

Super-Kamiokande Atmospheric Data

37

4.1

Data reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

37

4.1.1

Reduction of fully-contained events . . . . . . . . . . . . . . . .

38

4.1.2

Reduction of partially-contained events . . . . . . . . . . . . . .

46

4.1.3

Reduction of upward-going muons . . . . . . . . . . . . . . . . .

51

Event reconstruction . . . . . . . . . . . . . . . . . . . . . . . . . . . .

54

4.2.1

Vertex reconstruction . . . . . . . . . . . . . . . . . . . . . . . .

54

4.2.2

Ring counting . . . . . . . . . . . . . . . . . . . . . . . . . . . .

56

4.2

viii

4.3

5

6

Particle identification . . . . . . . . . . . . . . . . . . . . . . . .

57

4.2.4

Precise vertex fitter . . . . . . . . . . . . . . . . . . . . . . . . .

57

4.2.5

Momentum reconstruction . . . . . . . . . . . . . . . . . . . . .

58

4.2.6

Upward-going muon event reconstruction . . . . . . . . . . . . .

59

Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

59

4.3.1

Systematic uncertainties of data reduction . . . . . . . . . . . . .

60

4.3.2

Systematic uncertainties of event reconstruction . . . . . . . . . .

61

Mass-Induced νµ − ντ Oscillation

63

5.1

Data preparation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

64

5.1.1

Data binning . . . . . . . . . . . . . . . . . . . . . . . . . . . .

65

5.2

Oscillation Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

67

5.3

Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . .

69

Constraints on Sterile Neutrinos
6.1

6.2
7

4.2.3

73

Comparison of νµ − ντ and νµ − νs oscillations . . . . . . . . . . . . . .

74

6.1.1

Signatures of νµ − νs oscillation . . . . . . . . . . . . . . . . . .

74

6.1.2

Best-fits of two models . . . . . . . . . . . . . . . . . . . . . . .

76

An admixture model . . . . . . . . . . . . . . . . . . . . . . . . . . . .

80

Violations of Lorentz Invariance and CPT

83

7.1

An overview LIV and CPTV . . . . . . . . . . . . . . . . . . . . . . . .

84

7.2

Allowed limits on LIV parameter . . . . . . . . . . . . . . . . . . . . . .

88

7.3

The violation of CPT symmetry . . . . . . . . . . . . . . . . . . . . . .

90

7.3.1

Testing the LSND mass squared splitting scale . . . . . . . . . .

90

7.3.2

Allowed limits of CPTV . . . . . . . . . . . . . . . . . . . . . .

92

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

93

7.4

ix

8

Neutrino Decay and Neutrino Decoherence

96

8.1

Neutrino decoherence . . . . . . . . . . . . . . . . . . . . . . . . . . . .

97

8.1.1

A pure neutrino decoherence model . . . . . . . . . . . . . . . .

99

8.1.2

Allowed limits of neutrino decoherence . . . . . . . . . . . . . . 102

8.2

8.3
9

Neutrino decay . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106
8.2.1

Comparisons between neutrino decay and neutrino oscillation . . 107

8.2.2

Decay as a sub-dominant effect . . . . . . . . . . . . . . . . . . 109

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112

Conclusion

115

A Outer Detector MC Tuning116

A.1 An overview of OD Monte Carlo simulation . . . . . . . . . . . . . . . . 117
A.2 Tuning MC parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . 121
B PC Reduction Systematic Uncertainty

125

B.1 Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126
B.2 Evaluation of ID cut uncertainties . . . . . . . . . . . . . . . . . . . . . 128
B.3 The Uncertainty of OD cuts . . . . . . . . . . . . . . . . . . . . . . . . . 131
B.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
Bibliography

133xList of Tables
2.1

Systematic uncertainties on the atmospheric neutrino flux . . . . . . . . .

21

3.1

Super-K detector related systematic uncertainties . . . . . . . . . . . . .

35

3.2

Neutrino interaction related systematic uncertainties . . . . . . . . . . . .

36

4.1

Super-K data reduction uncertainties . . . . . . . . . . . . . . . . . . . .

62

4.2

Systematic uncertainties in event reconstruction . . . . . . . . . . . . . .

62

6.1

Portions of NC and CC events in the NC enhanced data samples . . . . .

75

6.2

Chi-square difference breakdown (∆χ2 = χ2νµ −νs − χ2νµ −ντ ) . . . . . . .

79

8.1

Chi-square difference breakdown between the standard oscillation and thedecoherence model (∆χ2 = χ2decoherence − χ2νµ −ντ ) . . . . . . . . . . . . 101

8.2

Chi-square difference breakdown between the standard oscillation and the
decay case II model (∆χ2 = χ2decayII − χ2νµ −ντ ) . . . . . . . . . . . . . . 110

xiList of Figures
2.1

Primary cosmic ray fluxes of protons and helium nuclei [1] . . . . . . . .

12

2.2

Neutron flux fluctuation [2]. . . . . . . . . . . . . . . . . . . . . . . . .

13

2.3

(top) The direction averaged atmospheric neutrino energy spectrum forνµ + ν µ ; (bottom) The ratio of the calculated neutrino flux. . . . . . . . .

2.4

The flux ratio of νµ + ν µ to νe + ν e averaged over all zenith and azimuth
angles versus neutrino energy. . . . . . . . . . . . . . . . . . . . . . . .

2.5

15

16

The flux ratios of νµ to ν µ and νe to ν e versus neutrino energy. The solidline is from [3], the dashed line is from [4] and dotted line is from [5]. . .

17

2.6

Distributions of neutrino pathlength in air . . . . . . . . . . . . . . . . .

18

2.7

The averaged pathlength versus cosine of zenith angle . . . . . . . . . . .

19

2.8

The flux of atmospheric neutrinos versus zenith angle. The solid line isfrom [3], the dashed line is from [4] and dotted line is from [5]. . . . . . .

2.9

20

The flux of upward-going atmospheric neutrinos versus zenith angle. The
solid line is from [3], the dashed line is from [4] and dotted line is from [5]. 21

3.1

3.2

A schematic drawing of the Super-K detector (the number of ID PMTs
during SK-II period was reduced by half) . . . . . . . . . . . . . . . . .

24

The reconstructed π 0 → γγ invariant mass . . . . . . . . . . . . . . . . .

28

xii

3.3

The determination of the absolute energy scale of Super-Kamiokande basedon various techniques . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3.4

28

The mean reconstructed energy of cosmic ray stopping muons divided
by their range (upper) and muon-decay electron (lower) as a function ofelapsed days during SK-I.

3.5. . . . . . . . . . . . . . . . . . . . . . . . .

The scattering angle distribution by neutrino interactions off the H 2 O target from the K2K experiment (data are from Fig.1(b) of Ref. [6]). . . . .

3.6

32

Charged current total cross section divided by E ν for (a) neutrino and (b)anti-neutrino nucleon charged current interactions. . . . . . . . . . . . .

4.1

29

33

The number of hits in the largest outer detector cluster. Solid histogram is
the MC prediction with neutrino oscillation considered. . . . . . . . . . .

39

4.2

Parent neutrino energy spectra . . . . . . . . . . . . . . . . . . . . . . .

61

5.1

Binning strategy for the oscillation analysis . . . . . . . . . . . . . . . .

67

5.2

Chi-square contours on the ∆m2 − sin2 2θ parameter plane for νµ − ντmixing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5.3

70

Chi-square contours on the ∆m2 − sin2 2θ parameter plane for νµ − ντmixing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

70

5.4

1D chi-square projections . . . . . . . . . . . . . . . . . . . . . . . . . .

71

5.5

The distribution of normalized pulls with its Gaussian fit results . . . . .

72

6.1

Zenith angle distributions of neutral current enhanced samples. Left: subGeV; right: multi-GeV . . . . . . . . . . . . . . . . . . . . . . . . . . .

75

6.2

Survival probabilities of muon neutrinos crossing the Earth . . . . . . . .

76

6.3

χ2 contours of νµ − ντ oscillation (with NC enhanced samples) . . . . . .

77

6.4

χ2 contours of the νµ − νs oscillation (with NC enhanced samples) . . . .

78

xiii

6.5

Comparison of the best-fit zenith angle distributions of ν µ − ντ and νµ − νsoscillations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

79

6.6

Distributions of the normalized pull terms of ν µ − ντ and νµ − νs oscillations 80

6.7

A case of 2+2 mass hierarchy and four-flavor mixing . . . . . . . . . . .

81

6.8

The νs admixture allowance . . . . . . . . . . . . . . . . . . . . . . . .

82

7.1

Allowed regions for parameters κ and α . . . . . . . . . . . . . . . . . .

87

7.2

χ2 versus α for a general oscillation formula . . . . . . . . . . . . . . . .

88

7.3

Chi-square contours of LIV. Left: η = 0; Right: η = π . . . . . . . . . .

90

7.4

Contour plot of the ad hoc CPTV model chi-square. Dash line is definedby∆m2 = ∆m̄2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

91

7.5

Chi-square contours considering CPTV as a sub-dominant effect . . . . .

93

8.1

Chi-square contours on the sin2 2θ − γ0 plane . . . . . . . . . . . . . . . 100

8.2

Comparison of pull terms of νµ − ντ and decoherence . . . . . . . . . . . 102

8.3

Comparison of pull terms of νµ − ντ and decoherence . . . . . . . . . . . 103

8.4

Allowed decoherence limits for γ = γ0 case . . . . . . . . . . . . . . . . 104

8.5

Allowed decoherence limits for γ = γ0 /E case . . . . . . . . . . . . . . 104

8.6

Allowed decoherence limits for γ = γ0 (E/GeV )2 case . . . . . . . . . . 105

8.7

Chi-square contours of neutrino decay case I . . . . . . . . . . . . . . . . 108

8.8

Chi-square contours of neutrino decay case II . . . . . . . . . . . . . . . 109

8.9

Comparison of pull terms of the ν µ − ντ oscillation and the neutrino decaycase 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

8.10 Comparison of pull terms of the ν µ − ντ oscillation and the neutrino decaycase 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111
8.11 Chi-square contours of neutrino oscillation + neutrino decay in linear scale 112

xiv

8.12 Chi-square contours of neutrino oscillation + neutrino decay in logarithm
scale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113
8.13 The non-zero local minimum of neutrino decay . . . . . . . . . . . . . . 114
A.1 OD MC performance before this tuning. Upper two: charge distributions
inside/outside of the time window and; lower two: hit time and number ofhits distributions (dots are observation and solid histograms are MC) . . . 118
A.2 SPE pulse height distribution of Hamamatsu R1408 . . . . . . . . . . . . 119
A.3 Measured charge distribution of Hamamatsu R1408 . . . . . . . . . . . . 119
A.4 Distribution of centroids using exponential SPE distribution model (theconversion factor to number of photoelectrons is around 2 in apdetsim)

120

A.5 OD charge distribution inside/outside of time window with the tuned MC
simulation compared with observation . . . . . . . . . . . . . . . . . . . 123
A.6 OD hit distribution with the tuned MC simulation compared with observation124
B.1 Distributions of npe200 of the final PC sample. The circles are the final PC
sample; triangle are MC . . . . . . . . . . . . . . . . . . . . . . . . . . . 128
B.2 cos θfit−first distributions of MC and observed data . . . . . . . . . . . . . 129
B.3 Distribution of mugood with mudistance! 3000 . . . . . . . . . . . . . . 130
B.4 Distribution of dcorner . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130
B.5 Distributions produced by another set of MC parameters . . . . . . . . . 132

xv

List of Abbreviations
ADC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Analog to Digital Converter
CC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Charged Current
CPTV . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . CPT Violation
DAQ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Data Acquisition
FC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Fully Contained
FRP . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Fiber Reinforced Plastic
GPS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Global Positioning System
GRB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Gamma-ray Burst
GUT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Grand Unified Theories
HE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . High Energy
ID . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Inner Detector
IMB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Irvine Michigan Brookhaven
K2K . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . KEK to Kamioka
KEK . . . . . . . . Ken Energy Kenkyusho (High Energy Accelerator Research Organization)
LE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Low Energy
LEP . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Large Electron-Positron Storage Ring
LIV . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Lorentz Invariance Violation
LSND . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Liquid Scintillator Neutrino Detector
MACRO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Monopole and Cosmic-ray Observatoryxvi

MC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Monte CarloMINOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Main Injector Neutrino Oscillation Search
NC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Neutral Current
OD . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Outer Detector
PC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Partially Contained
PMT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Photomultiplier Tube
QTC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Charge to Time Converter
SK . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Super-Kamiokande
SLE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Super low energy
SM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Standard Model
SNO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Sudbury Neutrino Observatory
SPE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Single Photoelectron
TDC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .Time to Digital Converter
TQ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Time Charge
VME . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Versa Module Europe

xviiChapter 1
Introduction
Neutrinos were first postulated by W. Pauli in December 1930 [7] to explain the continuous
energy spectrum of electrons from nuclear β decay. Ever since then, physicists have been
performing different experiments and proposing many kinds of theories to explore and to
explain the behavior of this “ghostly” particle.
In this chapter, we will first briefly go through the history of neutrinos. Then, after
introducing neutrino mass theories and experimental results, we will talk about the phenomenology of neutrino mixing and neutrino oscillation.

1.1 A brief history of the Standard Model neutrinos
Twenty-six years after W. Pauli postulated neutrinos, in 1956, they were detected by
F. Reines and C. Cowan [8]. The neutrinos they detected were anti-neutrinos of electron type from nuclear reactors. In the following year, M. Goldhaber, L. Grodzins andA. W. Sunyar experimentally deduced that neutrinos are always “left-handed” [9]. In 1962,
J. Steinberger, M. Schwartz, L. Lederman and their collaborators found a new kind of neutrino and confirmed that neutrinos actually have different flavors [10]. The new neutrino
12
they found was of muon type since it was always associated with muon creation or decay.
In 1975, a group of scientists led by Martin L. Perl found a new kind of lepton, τ ,at Stanford Linear Accelerator Center (SLAC) [11]. This is the first member of the third
generation of the Standard Model (SM) and the discovery strongly suggested the existenceof the tau type neutrino ντ . Though tau neutrinos were not directly detected until 2000 at
the DONUT experiment in Fermilab [12], people have not doubted their existence.
In late 1980’s, by measuring the invisible decay width of the Z 0 produced in the e+ e−
collider LEP, particle physicists confirmed there are only three generations of light active
neutrinos [13].
In the Standard Model, neutrinos are spin half and left-handed (right-handed for antineutrinos) fermions, and they only participate weak interaction. The results from collider
experiments indicate neutrinos only have three different flavors. The absence of righthanded neutrinos suggests that neutrinos are massless since right-handed neutrinos areneeded by the mass term in the Lagrangian.

1.2 Neutrino mass
Although the Standard Model considers neutrinos as massless particles, the experimental
efforts to measure their masses have never stopped since the postulation of neutrinos. However, all the measurements so far are consistent with zero. The current best limit from direct measurements comes from the tritium β decay experiments at Mainz and Troitsk [14].
Their sensitivity reaches around 2 eV/c2 . The first experimental support for non-zero neutrino masses came from the discovery of neutrino oscillations [15]. We will discuss thephenomenology of neutrino oscillations briefly in the next section.
Due to neutrino mixings, the masses measured in β decay experiments are effective
masses which have contributions from all three mass eigenstates and three mixing angles.

3
The next generation direct measurement experiment Karlsruhe Tritium Neutrino Experiment (KATRIN) expects to enhance the sensitivity by one order of magnitude down to∼ 0.2 eV/c2 [16]. There are also indirect constraints from cosmological large scale structure measurements, which are sensitive to the sum of all the neutrino masses [17, 18]. By
combining the Wilkinson Microwave Anisotropy Probe (WMAP) measurements on theanisotropies of the cosmic microwave background (CMB) and the Sloan Digital Sky Sur!
vey (SDSS) data on galaxy clustering properties, Tegmark et al set a limit ofmν <1.7 eV at 95% confidence level [19]. This limit can be improved further down to 0.42 eV
by including weak lensing, Lyα forest, and supernovae data [20].

Theoretically, we can give neutrinos masses by doing some extensions to the Standard
Model. Basically, there are several ways [21, 22]. The standard approach is that neutrino
masses come from the interaction with Higgs field, just like other fermions in the Standard
Model (SM), the neutrino mass term reads as:
LD = −(ν̄L m†D νR + ν̄R mD νL ).

(1.1)

Where mD = Yν < H > is the Dirac mass generated by interacting with the Higgs field
H and Yν is the coupling constant. This requires the presence of right-handed neutrinos ν R
which must be sterile, i.e. they do not interact with matter. However, it is rather unnaturalthat the coupling of neutrinos with the Higgs field Y ν must be extremely small comparedto other fermions due to the smallness of neutrino masses.
Neutrinos have zero conserved charges (electric and color), which gives the possibility
that neutrinos could be Majorana particles, i.e. neutrinos and antineutrinos are the same
particles. A Majorana mass term can be written as − 12 mψ̄ψ for a Majorana spinor ψ, whereψ is self-conjugate: ψ = ψ c . This possibility gives rise to a very elegant and natural way
to generate the smallness of neutrino masses: “seesaw mechanism” [23, 24, 25, 26, 27].

4
Using the charge conjugate fields νLc and νRc , ν̄R νL can be rewritten as 12 (ν̄R νL + ν¯c L νRc )and ν̄L νR can be rewritten in the similar way. Then, the Lagrangian term has both Diracmasses and Majorana masses can be written as [28]:


mTD





1
  νL  mL
Lseesaw = − (νLc νR )  + h.c.
2
mD MR
νRc
Where mD is the Dirac mass matrix, mL is the Majorana mass matrix of left-handed neutrinos and MR is the Majorana mass matrix of right-handed neutrinos. Seesaw mechanism
can be best explained in a case with just one family of neutrinos. Further, let us assume
mL = 0 and MR $ mD . We have a 2 × 2 mass matrix

 0 mD 
M =
.
mD MR
This matrix has two eigenvalues, given approximately by M R and −m2D /MR . Though mD
is naturally at the scale of charged lepton masses, a tiny neutrino mass of scale ∼ m 2D /MR
can still be generated since MR $ mD .

1.3 Neutrino mixing and oscillation
Along with the early discoveries of neutrino physics in 1958, Pontecorvo proposed thepossibility that neutrino and antineutrino might oscillate in a way resembling the K 0 − K 0
system [29]. In 1962, after the discovery of the muon neutrino flavor, Maki, Nakagawa
and Sakata proposed that if neutrino mass eigenstates do not coincide with neutrino flavor eigenstates, then mixing will happen, which will cause oscillation phenomenon of thesort we discovered with Super-K [30]. In 1968, Gribov and Pontecorvo in Ref. [31] pro-5
posed neutrino oscillation as a solution to the solar neutrino deficit problem discoveredby R. Davis et al [32, 33]. Super-Kamiokande provided strong evidence of atmospheric
neutrino oscillations in 1998 [15] and long baseline neutrino experiments K2K and MINOS have consistent results with Super-Kamiokande [6, 34]. Solar and reactor neutrino
experiment data are also explained by neutrino oscillations [35, 36].
In the case of three flavors, the mass eigenstates and the flavor eigenstates are connected by an unitary matrix UP M N S , νe  ν  = UP M N S
 µ 


ντ





 ν1  ν ,
 2 


ν3

where the mixing matrix UP M N S is generally expressed in three mixing angles, θ 12 , θ23 ,and θ13 , and a phase term δ,





−iδ

0 sin θ13 e
00
 1
  cos θ13UP M N S = 010 0 cos θ23 sin θ23  
0 − sin θ23 cos θ23
cos θ13− sin θ13 eiδ 0 cos θ12 sin θ12 0 × − sin θ12 cos θ12 0 


0
0
1








Based on the results from Chooz and Palo Verde experiments, θ13 is very small, less than
0.16 at 90% C.L. [37, 38, 39] and the scale of the solar mass-squared splitting is several
orders of magnitudes smaller than the atmospheric one. So, for most experiments, two-6
flavor mixing analysis is sufficient.

1.3.1 Two-Flavor oscillation in vacuum
Let’s study the νµ − ντ mixing case:










 νµ   cos θ sin θ   ν2 
=



ντ
ν3
− sin θ cos θ
In the relativistic limit, which applies to most cases due to the small masses of neutrinos, the dispersion relation could be written as E = p + m 2 /2p. Assuming the two mass
eigenstates have the same momentum, the Hamiltonian of the system is:m22 p + 2p

H =

0

0
m2
p + 2p3




.

(1.2)

What are produced and detected in experiments are weak interaction, i.e. flavor eigenstates. If a neutrino is produced at time t = 0, after time t, the wave packets will be at a!
flavor eigenstate: |ψ(t)& = 3i=2 φi (t)|νi &. And the evolution equation is:
i









d  φ2 (t)  φ2 (t) 

=H
,
dt
φ (t)
φ (t)
3

(1.3)

3

where φi (0) = Uαi , α = µ or τ depending the initial flavor state of the experiment. The
solution to Eq. 1.3 is:









φ2 (0)  φ2 (t) −iHt 
=e.
φ3 (t)
φ3 (0)

7
For the atmospheric neutrino case, if a νµ with energy E was created in the upper
atmosphere, after traveling a distance L, the system is in the following state:22|ψ(t)& = cos θei∆m L/2E |ν2 & + sin θe−i∆m L/2E |ν3 &,where ∆m2 = m23 − m22 , assuming m3 > m2 . We have applied approximations E ' p
and L ' ct during the derivation. Then the probability we find it still in the muon neutrino
state, i.e., the survival probability is,Pνµ →νµ = |(νµ |ψ(t)&|2 = 1 − sin2 2θ sin2

∆m2 L
,
4E

(1.4)

1.3.2 Two-flavor oscillation in matterIn the previous section, we derived the neutrino survival probability due to mass eigenstate
mixing in the vacuum. However, when neutrinos propagate in matter, if different interaction states have different forward scattering potentials due to the different interactions withthe medium, then the Hamiltonian is different from Eq. 1.2, thus the survival probability
formula is modified. This phenomenon is called matter effect and it was first pointed out
by L. Wolfenstein in [40, 41].
Among the three SM flavors, the matter effect only exists between the mixing of ν e andνµ or νe and ντ for media containing electrons, which give ν e ’s an extra forward scattering
potential due to the extra charged current interaction. The neutral current interaction is thesame for all three flavors. However, if there are new flavor eigenstates, say sterile neutrinos,
which we will discuss in later chapters, their interactions with matter are different fromthe SM interactions, and there will also be matter effects between them and the three SM
flavors.

8
Without losing generality, we study a two-flavor case here: the mixing between two
interaction eigenstates α and β, the mass eigenvalues are m a and mb (mb > ma and∆m2 ≡ m2b − m2a ), and the mixing angle between them is θ. Let us assume the interaction
eigenstate α gains an extra forward scattering potential ζ in the medium considered, thenthe Hamiltonian becomes (in the previous section, we see that pI 2×2 plays no role in the
final results so we omit it from the Hamiltonian. Also we will apply p ' E and x ' ct
from now on:



H=

m2a
2E

0











0   cos θ − sin θ   ζ 0   cos θ sin θ 
+


.
m2b
sin θ cos θ
0 0− sin θ cos θ
2E

(1.5)

If ζ does not depend on time, by diagonalizing this Hamiltonian, we can solve the evolution
Eq. 1.3.
M2

M2

Suppose we obtain two eigenvalues: 2E2 > 2E1 , by diagonalizing H in Eq. 1.5, we cancalculate the effective mass-squared splitting ∆m2M and the effective mixing angle θM inmatter:

,
2Eζ 2
2

 ∆m2M = M22 − M12 = ∆m2 (cos 2θ − ∆m2 ) + sin 2θ
2q sin 2θM = sin 2θ M∆m
2 −M 2 =2

1

sin 2θ

.

(cos 2θ− 2Eζ2 )2 +sin2 2θ∆mHowever, in the medium with varied density, which means ζ is time dependent, theanalytical results above are not valid any more. In this case, we can divide the neutrino
trajectory into different sections, x1 , x2 , ..., xn , approximating each small section with aconstant densities and calculate the survival probability by evolution method:








φa (0)  φa (x) −iH x −iHx
−iH x  = e n n e n−1 n−1 · · · e 1 1 ,
φb (x)
φb (0)(1.6)

9
which form a flavor state |ψ(x)& =Pα = |(να |ψ(x)&|2 [42].

!bi=a φi (x)|νi &, and the probability finding flavor α is

1.4 Beyond neutrino oscillation
Recently, mass eigenstate mixing induced neutrino oscillation has become the standard
theory accounting for the solar, atmospheric and long baseline neutrino experiments. However, there are other alternative theories which can also induce neutrino oscillations or cancause effects similar to some of the experimental observations.
Mass eigenstate mixing induced neutrino oscillations are basically due to the factthat different mass eigenstates have different energies for the same momentum: E i =pc + p2 /2mi . Thus, in principle, any theories which can modify the dispersion relation
to meet this condition can induce neutrino oscillations, such as Fermi point splitting [43],violations of Lorentz invariance (LIV) and CPT (CPTV) [44, 45].
Besides the non-standard theories which can induce neutrino oscillations, there is also
another category of alternatives which can cause neutrino vanishing instead of oscillation,for example, neutrino decoherence [46, 47, 48, 49, 50] and neutrino decay [51]. Somemodels based on these theories can also explain some experimental data quite well.
Whether these non-standard effects can explain neutrino experiment data is certainly
worth investigating. On the other hand, these alternatives are new physics beyond theStandard Model and testing them is of many theoretical importance. Atmospheric neutrino
data, as we will introduce later, cover a wide range of energies, pathlengths and matter
densities, will prove to be a powerful tool to explore new physics.
In the following chapters, we will introduce atmospheric neutrinos, the Super-Kamiokande
experiment with which we do the observation, present the studies on mass-induced neutrino oscillations, introduce some of the alternative theories and test their plausibility and10allowed limits based on the current data.

Chapter 2
Atmospheric Neutrinos
In the neutrino oscillation section of Chapter 1, we showed that neutrino oscillation probability basically depends on two kinetic variables: the energy and the pathlength. Naturally,
we also need the flavor information. We will briefly introduce the properties of atmospheric neutrino fluxes from these three perspectives. One fact needs to be pointed out isthat the neutrino flux presented here is calculated before considering the neutrino oscillation effect and except for the general principles, the values are for the Super-Kamiokandesite only, since the flux is affected by the geomagnetic field.

2.1 Primary cosmic rays
Atmospheric neutrinos are produced in cosmic ray interactions in the upper atmosphere.
Thus, the atmospheric neutrino flux is a convolution of the primary cosmic ray spectrumat the production point with the yield of neutrinos per primary particle. Roughly speaking,primary cosmic ray consists of 90% protons, 9% alpha particles and a small portion ofother kinds of heavy nuclei [52]. Figure 2.1 shows the measured fluxes of protons andhelium nuclei by different experiments [1].
1112WebberLEAPMASS1
CAPRICE
BESS98
IMAX
AMS
Ryan
JACEE
Ivanenko
Kawamura
Runjob
Fitted curve (p)Fitted curve (He: high)Fitted curve (He: low)

Figure 2.1: Primary cosmic ray fluxes of protons and helium nuclei [1]Cosmic rays are isotropic when they reach the vicinity of the Earth. But, in order toreach the upper atmosphere and interact with the air nucleus, these charged particles haveto pass through the geomagnetic field, which produces local effects — the neutrino flux
is location dependent [1]since the geomagnetic field varies at different locations. The flux
data we present here is for the Super-Kamiokande site only.
Interactions of primary cosmic ray particles with air nuclei produce mesons which aremainly pions and kaons. Subsequently, neutrinos are produced in the decay chains of thesemesons. For example, the decay chain from pions are:
π ± → µ± + νµ /ν̄µµ± → e± + νe /ν̄e + ν̄µ /νµKaon decays follow the similar chains.

.

(2.1)

13

Figure 2.2: Neutron flux fluctuation [2].
Solar activity effects the primary cosmic rays and thus on the atmospheric neutrino
flux. The minimum neutrino flux occur at times of high solar activity.
Figure 2.2 shows the neutron monitor data from Moscow Neutrino Monitor [2]. According to the cosmic ray proton, helium and neutron measurements [53, 54, 2], during the
SK-I period, the cosmic ray flux was near the solar minimum until the summer of 1999,
rapidly decreased during the next year, and was at the minimum value consistent with solar maximum from summer of 2000 until Super-Kamiokande I stopped taking data in July
2001. Therefore, the atmospheric neutrino Monte Carlo is calculated for 3 years of solar
minimum, 1 year of changing activity, and 1 year of solar maximum for Super-K I. During
the SK-II period, which is from Jan 2003 to October 2005, there are about one year solar
minimum in 2003, and two years of changing activity during the 2004 and 2005 period.

14

2.2 Atmospheric neutrino energy spectrum
Roughly speaking, low energy neutrinos around 1 GeV are mainly produced by the primary fluxes around 10 GeV. At this energy range, cosmic ray fluxes are affected by solar
activity and by the geomagnetic field due to their low rigidities. Neutrinos above 100 GeV
are mainly from primary cosmic rays with energies greater than 1000 GeV which are freeof the effects of solar activity and the rigidity cutoff. But details of the higher energy
primary cosmic ray flux are not as well measured.
Following the primary interactions of cosmic ray particles with the air nuclei, the propagation and decay of secondary particles are simulated [1]. The calculated energy spectraof atmospheric neutrinos at Kamioka site are shown in Fig. 2.3 (top). Fig. 2.3(bottom)
shows the comparison of different calculations. The fluxes calculated in Ref. [4] (solidline) and Ref. [5] (dashed line) are normalized by the flux in Ref. [3].
Below 10 GeV, the difference among different calculations is around 10%. This is the
combined effects of the uncertainty of the primary cosmic ray flux measurements [55, 56],
which is about 5% below 100 GeV, and different hadronic interaction models used in
calculations.
Above 100 GeV, the primary cosmic ray spectrum and compositions are less wellknown. Therefore, for neutrino energies much higher than 10 GeV, the uncertainties in
the absolute neutrino flux are much larger. The index of low energy (<100 GeV) proton
spectra is −2.74 ± 0.01 according to Ref. [1]. However, this spectrum does not fit well tothe high energy part. Ref. [3] suggests allowing a different spectral index for the protons
above 100 GeV: the best fit value is −2.71. Also, it is discussed in Ref. [1] that the spectrum index for the He flux can be fit by either −2.64 or −2.74, shown in Fig. 2.1. So, there
could be 0.10 uncertainty in the spectrum index for He. The spectrum indices for heavier
nuclei have uncertainties larger than 0.05 [1].

Flux ra tio

Flu x × Eν2(m -2 s e c -1 s r-1 GeV )

15
10

3

10

2

10
Honda flux
Bartol flux
Fluka flux

1
10

-1

1.4
1.3
1.2
1.1
1
0.9
0.8
0.7
0.6
10

Bartol/Honda
Fluka /Honda
-1

1

10

Eν(GeV)

10

2

10

3

Figure 2.3: (top) The direction averaged atmospheric neutrino energy spectrum for ν µ +ν µ ;(bottom) The ratio of the calculated neutrino flux.
Taking the flux weighted average of these spectrum index uncertainties, we assign
0.03 and 0.05 for the uncertainties in the energy spectrum index in the primary cosmic rayenergy spectrum below and above 100 GeV, respectively.

2.3 Flavor ratios
There are mainly two factors affecting the flavor ratio: neutrino energy and the relative
contributions from kaons and pions. Figure 2.4 shows the ratio of ν µ + ν µ to νe + ν e as afunction of the neutrino energy. Solid, dashed and dotted lines show the prediction by [3],[4] and [5], respectively. The angle dependence has been integrated out. In the low energy
region less than about 5 GeV, most of the neutrinos are produced by the decay chain ofpions in Eq. 2.1. To a very good approximation, the ratio is around 2. The uncertainty

16of this energy range is about 3%, which is estimated by comparing the three calculation
results.
10
98
7

Honda flux
Bartol flux
Fluka flux

Flux ra tio

6
5
4
3

2

10

– /ν +ν
–
νµ+ν
µ e e
-1

1

Eν(GeV)10102

Figure 2.4: The flux ratio of νµ + ν µ to νe + ν e averaged over all zenith and azimuth angles
versus neutrino energy.
However, as neutrino energy increases, the parent muons are more likely to reach theground before they decay and the contribution from K decay increases: at 10 GeV, about
10% of electron type neutrinos and 20% of muon type are from K; it increases to more than
30% at 100 GeV for both flavors. And the uncertainty of K/π production ratio [57, 58] inthis energy region is as big as 20%, which causes a bigger uncertainty on the flavor ratio in
the energy range of 10 to 100 GeV. As seen from Fig. 2.4, the difference of the calculatedνµ + ν µ to νe + ν e ratio is as large as 10% at 100 GeV. Above 5 GeV, we assumed that theuncertainty linearly increases with log Eν from 3% at 5 GeV to 10% at 100 GeV.
Although the Super-K detector does not distinguish the signs of the charges, the observation is still affected by ν/ν̄ ratio since they have different interaction cross sections
and the oppositely charged final state particles produce different responses in the detector.

17
Also, some physics effects are indeed sensitive to the ν/ν̄ ratio, e.g., the matter effect and
CPT violation (CPTV) models which we will discuss in this dissertation.
Figure 2.5 shows the calculated flux ratios νµ to ν µ and νe to ν e . Solid, dashed anddotted lines show the prediction by Honda et al in Ref. [3], Barr et al in Ref. [4] and Battistoni et al in Ref. [5], respectively (same key as Fig. 2.4). Below 10 GeV, for both ratios,
different calculations agree to about 5%. But the disagreement gets larger as neutrino energy increases. The systematic uncertainties on the ν/ν ratio are assumed to be 5% below
10 GeV and linearly increase with log Eν to 10% and 25% at 100 GeV, for the νe to ν e andνµ to ν µ ratios, respectively.
1.6
1.5

Flux ratio

1.4
1.3
1.2
1.1
1
–

0.9
10

–

!e/!e
-1

1

!µ/!µ
10

E! (GeV)

10

-1

1

10

10

2

E! (GeV)Figure 2.5: The flux ratios of νµ to ν µ and νe to ν e versus neutrino energy. The solid line
is from [3], the dashed line is from [4] and dotted line is from [5].

18

2.4 Neutrino pathlength
The pathlength of a neutrino is the function of the production height and the zenith angle.
The typical production height for neutrinos is around 15 km [1]. We use a Monte Carlo
method to simulate a distribution of the neutrino production heights. Figure 2.6 shows
the pathlength distributions in air of different neutrino energies and neutrino flavors forboth vertical and horizontal cases. Including the production height, Fig. 2.4 shows theaveraged pathlength for 1 GeV νµ as a function of zenith angle. It is not hard to imagine
that the uncertainty on pathlength for neutrinos passing a great distance crossing the Earthis relatively small and for horizontal and downward going neutrinos, is relatively larger
since the baselines are shorter.
!µ-vertical

0.12

0.25 - 0.40GeV

0.1

1.0 - 1.6GeV

0.1

0.08

4.0 - 6.3GeV

0.08

0.06

0.06

0.04

0.04

0.02

0.02

0

0

0.12

20

40

0

60

(km)

!µ-horizontal

0.1

0.08

0.08

0.06

0.06

0.04

0.04

0.02

0.02
0

200

400

0

0.12

0.1

0

!e-vertical

0.12

600

800

(km)

0

20

40

60

(km)

!e-horizontal

0

200

400

600

800

(km)

Figure 2.6: Distributions of neutrino pathlength in air
Cosmic rays are isotropic at the vicinity of the Earth. Therefore, geometrically, theneutrino flux should be up-down symmetric. However, due to the geomagnetic effect, lowenergy neutrinos are up-down asymmetric.

19
Figure 2.8 shows the zenith angle dependence of the atmospheric neutrino fluxes forseveral neutrino energies. Solid, dashed and dotted lines show the predictions in Ref. [3],Ref. [4] and Ref. [5], respectively (same key as Fig. 2.4). At low energies, and at the
Kamioka location, the fluxes of downward-going neutrinos are lower than those of upwardgoing neutrinos. This is due to the deflection of primary cosmic rays by the geomagnetic
field, roughly characterized by a minimum rigidity cutoff. For neutrino energies higher
than a few GeV, the calculated fluxes are essentially up-down symmetric, because theprimary particles are more energetic than the rigidity cutoff.
For upward going muons, which are produced by neutrinos interacting with rock beneath the Super-K detector, the typical energies are relatively higher and the decay K
contribution becomes larger. We assume that the K/π production ratio uncertainty is 20%
in the whole energy region [57, 58]. The uncertainties in the zenith angle and energy distributions due to the K/π production uncertainty are included in the systematic errors in the

Figure 2.7: The averaged pathlength versus cosine of zenith angle

Integrated neutrino flux (m-2 sec-1sr-1 )

20
700

200

600

175

300
200
100
0
-1

25

150

500
400

30

125
_

100

νµ+νµ

20

_

νµ+νµ

_

15

75
_

50

νe+νe

25

0.3-0.5 GeV
0

1

0
-1

10

_

νe+νe

5

0.9-1.5 GeV
0

cosΘ

νµ+νµ

1

0
-1

_

νe+νe
3.0-5.0 GeV
01

Figure 2.8: The flux of atmospheric neutrinos versus zenith angle. The solid line is from[3], the dashed line is from [4] and dotted line is from [5].
analysis. Figure 2.9 shows the zenith angle dependence of the atmospheric neutrino fluxes
for higher energy region observed as upward muons in Super-Kamiokande. Again, solid,
dashed and dotted lines show the prediction in Ref. [3], Ref. [4] and Ref. [5], respectively(same key as Fig. 2.4).

2.5 Summary
In summary, the uncertainty on the absolute flux is very large, but we do have pretty
accurate knowledge on the flavor ratio and up-down flux ratio, which will provide themost important information on neutrino oscillation analysis. Table 2.1 has all the systematic uncertainties which show our understanding on atmospheric neutrino flux at SuperKamiokande.

21

0.8

Integrated neutrino flux (m-2sec-1sr-1)

6
5

0.7

_

νµ+νµ

0.006
0.005

_

ν +ν
0.6 µ µ

4

0.004

0.5
0.4

3

0.003

0.3

2

_

νµ+νµ

0.002

0.2
1
0
-1

0.1

9.0-15 GeV
-0.5

0

0
-1

0.001

30-100 GeV
-0.5

cosΘ

0

0
-1

0.3-1.0 TeV-0.5

0Figure 2.9: The flux of upward-going atmospheric neutrinos versus zenith angle. The solidline is from [3], the dashed line is from [4] and dotted line is from [5].

Table 2.1: Systematic uncertainties on the atmospheric neutrino flux
Systematic Uncertainty
σ(%) Index
Absolute Normalizationfree
1(νµ + ν̄µ )/(νe + ν̄e ) (Eν < 5GeV )3.02(νµ + ν̄µ )/(νe + ν̄e ) (Eν > 5GeV )3.03νe /ν̄e (Eν < 10GeV )5.04νe /ν̄e (Eν > 10GeV )5.05νµ /ν̄µ (Eν < 10GeV )5.06νµ /ν̄µ (Eν > 10GeV )5.0
7
Up/down ratio
1.0
8
Horizontal/vertical
1.0
9K/π ratio
20.0
10
Production height
10.0
11Spectral index of primary cosmic ray above 100GeV 3.0
12
Sample-by-sample spectral index (FC Multi-GeV)
5.0
13
Sample-by-sample spectral index (PC and Upµ)5.0
14Solar activity
50.0
48Chapter 3
Super-Kamiokande Experiment
The Super-Kamiokande detector is a 50 kiloton water Cherenkov detector located at the
Kamioka Observatory of the Institute for Cosmic Ray Research, University of Tokyo. Thisfacility is in the Mozumi mine in Gifu prefecture, Japan, under the peak of Mt. Ikenoyama
which provides a 1,000 meters rock overburden (2,700 meters water equivalent). It is amulti-purposed underground experiment.
As the successor of the Kamiokande and IMB experiments [59, 60], the main physics
goals of Super-Kamiokande are to study nucleon decay [61], solar [62] and atmospheric
neutrinos [63]. Its low threshold, low-background environment, large target mass and large
aperture make the detector a promising facility for many other research topics, for example
but not limited to: supernova neutrinos, exotic properties of neutrinos [64] , exotic particle
searches [65], astrophysical [66, 67, 68, 69] and cosmological [70, 71] studies; it served
as the far detector of K2K experiment [72]; it will also served as the far detector of forecoming T2K experiment [73].
The Super-K experiment has three running periods so far. Super-K I (SK-I) started operation from April 1996 and continued until July 2001. In July 2001, the first upgrade wasstarted to refurbish the failed photomultipliers (PMT) and cables and to partially replace22

23
the outer detector (OD) reflector. The upgrade finished in October 2001 but in November
2001, a broken PMT at the bottom of the inner detector (ID) triggered a chain reaction
and destroyed 70% of the PMTs in the detector. More information about the accident canbe found in [74]. The second upgrade started in March 2002 and finished in December
2002, the rebuilt Super-K detector, Super-K II (SK-II), has a fully recovered OD but ID
only has half of the PMTs that SK-I has. SK-II ran from January 2003 to October 2005
and has a livetime of 804 days for atmospheric neutrino observation. The fully recovered
Super-K III (SK-III) started operation in July 2006.
This chapter will be mainly focused on factors relevant to the atmospheric neutrino
sector. More details on all the technical respects of the detector can be found in Ref. [75,62, 63].

3.1 Detector setup
Shown in Fig. 3.1, the Super-K detector consists of two concentric, optically separated
water Cherenkov detectors contained in a stainless steel tank 42 meters high and 39.3
meters in diameter, holding 50 kilotons of pure water. The top of the detector, under thehemisphere, consists of electronics huts and open work area. Super-K water is supplied
from the spring water inside the mine. The water in the tank is continuously circulated
through the water purification system to filter out any microparticle contaminations tokeep the transparency and, more importantly, radioactive materials from air, mainly radon(Rn). The water temperature in SK tank is controlled at ∼ 12 oC.
The ID holds a cylindrical volume of pure water 16.9 m in radius and 36.2 m high.
During SK-I running period, the inner detector (ID) has 11,146 Hamamatsu R3600 50 cm
diameter photomultiplier tubes (PMTs) mounted on the inner wall, comprising a photocathode coverage of about 40%. The number of ID PMTs during SK-II period is reduced

24

Figure 3.1: A schematic drawing of the Super-K detector (the number of ID PMTs during
SK-II period was reduced by half)by half due to the lack of supply. The photocathode coverage thus becomes around 20%.
The 50 cm PMTs were specially designed [76] to have good single photoelectron (p.e.) response, with a root-mean-squared timing resolution of 2.5 nsec. The space between PMTs
is covered by black plastics to prevent reflection.
The ID is surrounded by the outer detector (OD), a cylindrical shell of water 2.62.75 m thick including a dead space 55 cm. The OD is optically isolated from the ID
and mainly functions as a 4π active veto — to identify the incoming and outgoing muons.
The OD is instrumented with 1,885 outward-facing Hamamatsu R1408 20 cm PMTs. A
50cm × 50cm × 1.3cm wavelength shifter [77] is attached to every OD PMT to increase
light collection efficiency. The interior of OD is covered with reflective Tyvek to further

25
increase light collection. We will discuss OD PMTs performance and the OD Monte Carloin Appendix A. Besides functioning as the veto, the OD is also a thick passive radioactivity
shield.
Each PMT is connected to high voltage (HV) supplies and signal processing electronics
via coaxial cables. For the ID, signal and HV cables are separated, while the OD uses asingle cable. The high voltages on the PMTs are set to make sure they have approximately
equal gains within the ID and OD respectively in order to reduce the overall systematic
uncertainty. These cables are all brought up to the open work area on the top of the tank,where they are distributed into four electronic huts. Every hut contains electronic racks
and front-end data acquisition computers serving ID and OD PMTs for one quadrant of
the detector. There are four cable holes on top of the detector to let the four cable bundles
pass into the four electronic huts. Due to the big volumes occupied by the cable bundles in
the OD, the cosmic ray muons passing through those “dead” volumes can not be detectedby the OD so four plastic scintillators are installed on top of them outside the OD to vetothese cosmic ray muons.

3.2 Electronics and data acquisition
Both ID and OD PMT signals are processed by asynchronous, self-triggering circuits thatrecord the time and charge of each PMT hit over a threshold. Each ID PMT signal is digitized with custom Analog Timing Modules (ATMs) [78, 79] which provide 1.2 µsec timing
range at 0.3 nsec resolution and 550 pC charge range at 0.2 pC resolution (∼ 0.1 p.e.). The
ATM has automatically-switched dual channels to provide deadtime-free data acquisition.
The outer PMT signals are processed with custom charge-to-time conversion modules anddigitized with LeCroy 1877 multi-hit TDCs over a −10 µsec to +6 µsec window centeredon the trigger time. More details about the DAQ system of Super-K detector can be found26
in [75, 63].
An event used in the atmospheric neutrino analysis is triggered by the coincidence ofat least 30 PMT hits in a 200 nsec window which is the time for photons crossing thedetector diagonally. The hit threshold for each individual PMT is about 1/4 p.e. This
trigger condition corresponds to the mean number of hit PMTs for a 5.7 MeV electron.
The trigger rate is 10-12 Hz so the total number of events in the raw data is of order ~ 10 6
per day. The trigger rate due to cosmic ray muons is 2.2 Hz. Digitized data are saved at a
total rate of 12 GB per day.

3.3 Detector calibration
3.3.1 Water properties
The optical attenuation length L atten of water is parameterized in the following way,

I = I0 e−l/Latten /l2 ,where I is the light intensity after traveling a distance of l in water; I 0 is the initial intensity. At Super-K, Latten is measured through both direct and indirect ways. The direct
measurement uses a system consisting of a titanium-sapphire laser, a diffuser ball and a
CCD camera. To measure the attenuation length, the diffuser ball which is pumped withlaser is lowered into difference positions in the SK tank and the light intensity on the surface is monitored by the CCD camera. A two inch PMT is used to monitor the stabilityof the laser and also to measure the initial light intensity. Since water quality is expected
to change with time, this kind of measurements is made at regular intervals. The indirect
measurements use through-going cosmic ray muons. Through-going muons lose almost

27
same amount energy per unit distance in SK tank, ∼ 2MeV/cm, which makes them “standard” light sources for calibration purposes. Through-going muons also have advantages
being continuous and plentiful so they can be used to monitor the attenuation length as afunction of time.
Latten is actually the combined effects of absorption and scattering on the intensity oflight: Latten = αabs + αscat . To measure them separately, a combination of dye and N 2
lasers of wavelengths 337, 371, 400 and 420 nm are used as light sources. The laser beamis brought into the SK tank through optical fibers pointing at the bottom of the tank. Eachlaser fires every 6 seconds during normal data taking. From the spatial and the timing
distributions of the photons detected by PMTs, both absorption and scattering coefficientsare studied with the help of Monte Carlo and incorporated into the detector simulation.
Every 5 days, accumulated data are combined and analyzed as a single set.

3.3.2 The absolute energy calibration
Super-K employs several different sources to calibrate the absolute energy scale. For
low energy from ~5 MeV to ~15 MeV, an electron linear accelerator (LINAC) [80] anda deuterium-tritium neutron generator (DTG) [81] are used and a small uncertainty (better
than 1%) in hit counting is achieved [75]; for high energy, which is more relevant to atmospheric neutrino data, the following sources are used: the total number of photo-electrons
as a function of muon track length, where the muon track length is estimated by the reconstructed muon entrance point and the reconstructed vertex point of an electron from themuon-decay; the total number of photo-electrons as a function of Cherenkov angle for low
energy cosmic ray muons; the spectrum of muon-decay electrons; and the invariant massof π 0 s produced by neutrino interactions, see Fig. 3.2. Figure 3.3 summarizes the absoluteenergy scale calibration by these studies.

28

100

Number of events

80

60

40

20

0

0

100

200

300

400

500

600

700

2

Invariant Mass (MeV/c )Figure 3.2: The reconstructed π 0 → γγ invariant mass

(MC-data)/data (%)

10
8
6
4
2
0
-2

decay-e0"
low stop µhigh stop µ

-4
-6
-8
-10
10

2

3

10
10
10
Momentum range (MeV/c)

4

10

5

Figure 3.3: The determination of the absolute energy scale of Super-Kamiokande basedon various techniques
The stability of the energy scale was also monitored continuously using stopping muonsand muon-decay electrons. Figure 3.4 shows the time variation of the mean reconstructedenergy of stopping muons divided by muon range and the mean reconstructed electronenergy from muon-decays. Vertical axes in both figures are normalized to mean values
and each data point corresponds to two month period. The root-mean-square of the energy
scale variation is ±0.9 % over the time of the experiment. From combining the absoluteenergy scale accuracy study (±1.8 %) and the energy scale time variation (±0.9 %), thetotal uncertainty of the energy scale of atmospheric neutrino detection was estimated to be±2.0 %.

µ momentum/track

1.05
1.04
1.03
1.02
1.01
1
0.99
0.98
0.97
0.96
0.95

electron momentum

29

1.05
1.04
1.03
1.02
1.01
1
0.99
0.98
0.97
0.96
0.95

±1%

0

250

500

750

1000 1250 1500 1750 2000

±1%

0

250

500 750 1000 1250 1500 1750 2000
Elapsed days from Apr.-1-1996

Figure 3.4: The mean reconstructed energy of cosmic ray stopping muons divided by theirrange (upper) and muon-decay electron (lower) as a function of elapsed days during SK-I.

3.4 Neutrino interaction simulation
Atmospheric neutrino oscillation analysis relies heavily on the comparison of the experimental data with the theoretical expectation. Thus, we need detailed simulations of bothneutrino interactions and the Super-K detector. For atmospheric neutrino interactions, weneed to simulate the interactions of neutrinos ranging from 10 MeV to 100 TeV with thewater nuclei. In the case of upward muons, we need to simulation the interactions with therock nuclei surrounding the detector. At Super-K, we have two neutrino interaction simulation programs: NEUT [82] and NUANCE [83]. We will only describe NEUT briefly. NEUT
considers the following charged and neutral current neutrino interactions are considered:
1. (quasi-)elastic scattering, νN → lN $ ,2. single meson production, νN → lN $ m,
3. coherent π production, ν 16 O → lπ 16 O,
4. deep inelastic scattering, νN → lN $ hadrons.

30
Here, N and N $ are the nucleons (proton or neutron), l is the lepton, and m is the meson,respectively. For single meson production, K and η production are simulated as well as
the dominant π production processes. If the neutrino interaction occurred in the oxygen
nuclei, generated particles like pions and kaons interact with the nucleus before escaping.

3.4.1 Elastic and quasi-elastic scattering
The formalization of quasi-elastic scattering off a free proton was described by LlewellynSmith [84]. For scattering off nucleons in 16 O, the Fermi motion of the nucleons and Pauli
Exclusion Principle were taken into account. The nucleons are treated as quasi-free particles using the relativistic Fermi gas model of Smith and Moniz [85]. 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.

3.4.2 Single meson production
Rein and Sehgal’s model was used to simulate the resonance productions of single π, Kand η [86, 87]. 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, is restricted to be lessthan 2 GeV/c2 . In addition to the dominant single π production, K and η production isconsidered.

31
To determine the angular distribution of pions in the final state, we also use Rein and
Sehgal’s method for the P33 (1232) resonance. For the other resonances, the directional
distribution of the generated pions is set to be isotropic in the resonance rest frame. The
angular distribution of π + has been measured for νp → µ− pπ + [88] and the results agree
well with the Monte Carlo prediction. We also consider the Pauli blocking effect in the
decay of the baryon resonance by requiring that the momentum of nucleon should be larger
than the Fermi surface momentum. Pion-less delta decay is also considered, where 20%of the events do not have the pion and only the lepton and nucleon are generated [89].
The quasi-elastic and single meson production models have a parameter (axial vector
mass, MA ) that must be determined by experiments. For larger M A values, interactions
with higher Q2 values (and therefore larger scattering angles) are enhanced for these channels. The MA value was tuned using the K2K [6] near detector data. In our atmospheric
neutrino Monte Carlo simulation, M A is set to 1.1 GeV for both the quasi-elastic andsingle-meson production channels, but the uncertainty of the value is estimated to be 10%.
Figure 3.5 shows the K2K 1-kton water Cherenkov detector data on the scattering angle forsingle Cherenkov ring events [6] together with the prediction by the Monte Carlo used inthis analysis. The scattering angle agrees well between the data and Monte Carlo overall.
Coherent single-pion production, the interaction between the neutrino and the entire
oxygen nucleus, is simulated using the formalism developed by Rein and Sehgal [90].

3.4.3 Deep inelastic scattering
In order to calculate the cross-sections of deep inelastic scattering, the GRV94 [91] parton
distribution function is 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 greaterthan or equal to 2 for 1.3 < W < 2.0 GeV/c2 , because single pion production is sepa-

events / 10

o

32

Data
(stat.+syst.)

5000

MC
4000

3000

2000

1000

0

0

50

100

150
[deg.]

Figure 3.5: The scattering angle distribution by neutrino interactions off the H 2 O target
from the K2K experiment (data are from Fig.1(b) of Ref. [6]).
rately simulated as previously described. In order to generate events with multi-hadron
final states, two models are used. For W between 1.3 and 2.0 GeV/c 2 , a custom-madeprogram [92] is used to generate the final state hadrons; only pions are considered in thiscase. For W larger than 2 GeV/c2 , PYTHIA/JETSET [93] is used.
Total charged current cross sections including quasi-elastic scattering, single meson
productions and deep inelastic scattering are shown in Fig. 3.6. The solid line shows thecalculated total cross section. The dashed, dot and dash-dotted lines show the calculatedquasi-elastic, single-meson and deep-inelastic scatterings, respectively. Data points aretaken from different experiments [63].

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

33

1.4
CCQE

1.2

#/E! (10-38 cm-2GeV-1)

CC single "ANL
GGM 77
GGM 79Serpukhov

1

Total #

ANL 82
BNL 86(a)

CCFR 90
CDHSW 87
IHEP-JINR 96
IHEP-ITEP 79
CCFRR 84
BNL 82

0.8
0.6
0.4
0.2
0
-1
10

1

E! (GeV)

10

0.7
CCQE

#/ E! (10-38 cm-2GeV-1)

0.6

(b)

Total #GGM 77
GGM 79Serpukhov

CCFR 90
CDHSW 87
IHEP-JINR 96

IHEP-ITEP 79
CCFRR 84

0.5
0.4
0.3
0.2
0.1
0
-1
10

1

E! (GeV)10

Figure 3.6: Charged current total cross section divided by E ν for (a) neutrino and (b)anti-neutrino nucleon charged current interactions.

34
The interactions of pions are very important because the cross section for pion production
is quite large for neutrino energies above 1 GeV and the interaction cross sections for pions
in nuclear matter is also large.
In our simulation program, we consider the following pion interactions in 16 O: inelastic
scattering, charge exchange and absorption. The initial position of the pion generated
according to the Woods-Saxon nucleon density distribution [94]. The interaction mode isdetermined from the calculated mean free path of each interaction. To calculate the mean
free path, we adopt the model described by Salcedo et al. [95]. The calculated mean freepath depends not only on the momentum of the pion but also on the position of the pion in
the nucleus. If inelastic scattering or charge exchange occurs, the direction and momentumof the pion are determined by using the results of a phase shift analysis obtained from π−N
scattering experiments [96]. When calculating the pion scattering amplitude, the Pauliblocking effect is 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: π 12 C scattering,
π 16 O scattering and pion photo-production (γ+ 12 C → π − + X) [97].

3.5 The detector simulation
A Monte Carlo simulation on how the detector reacts to physical processes is needed tounderstand the systematics of the detector, to do data reduction and to compare our observation with the predictions of different theories or models.
The detector’s responses to the final state particles are simulated with a Monte Carlo
program based on the GEANT3 package [98], in which the propagation of particles, thegeneration and propagation of Cherenkov photons, and the response of the PMTs are considered. The CALOR package [99] was employed to simulate hadronic interactions in35

Table 3.1: Super-K detector related systematic uncertainties
Systematic uncertainty
Value (SK-I) Value (SK-II)
Energy calibration
0.02
0.025
Up/down symmetry of energy calibration
0.006
0.006

Index
37
39

water. This package can reproduce the pion interactions results very well down to low
momentum regions of ∼ 1GeV/c. For still lower momenta (p π ≤ 500MeV/c), a custom
subroutine was made based on experimental data in Ref. [100] on π − 16 O scattering andin Ref. [101] on π − p scattering.
For the propagation of Cherenkov photons in water, Rayleigh scattering, Mie scattering and absorption were considered in our simulation code. The attenuation coefficientsin water were tuned to reproduce the measurements of the laser system described in Section 3.3.1. Light reflection and absorption on detector material, such as the surface of
PMTs and black plastic sheets between the PMTs are simulated based on direct measurements, using probability functions that depend on the photon incident angle.
To give an example how we tune the detector simulation program, an OD simulation
tuning process was presented in Appendix A.

3.6 Summary
Our understanding to the Super-K detector can be summarized into systematic uncertainties. SK-I and SK-II are practically different detectors with similar design so they have thesame systematic uncertainties terms but with different values. Table 3.1 shows the valuesof SK-I and SK-II.
The uncertainties on neutrino reactions are summarized in Table 3.2. They are naturally
the same for SK-I and SK-II.

36

Table 3.2: Neutrino interaction related systematic uncertainties
Systematic uncertainty
Value Index
MA in quasielastic and single-π
0.1015Quasielastic scattering model dependence 0.10
16Quasielastic scattering cross section
0.10
17Single-meson production cross section
0.10
18
DIS model dependence
0.10
19
DIS total cross section
0.05
20
Coherent pion cross section
0.30
21
NC/CC ratio
0.20
22
Nuclear effect in 16 O
0.30
23
Pion energy spectrum
0.10
24
CC ντ cross section
0.30
25
Hadron simulation
0.1030

Chapter 4
Super-Kamiokande Atmospheric Data
Raw Super-K data consist mainly of downward-going cosmic ray muons and low energy
radioactivity events. The the ID and OD double structure can remove cosmic ray muonseasily with high efficiency. The atmospheric neutrino energy range, E vis > 30MeV, essentially makes the analysis free of any low energy background. Visible energy is definedas the energy of an electromagnetic 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.
After background rejection, Super-K observes ∼ 20 atmospheric neutrino events per day.
In this chapter, we will briefly review data reduction and event reconstruction processes.

4.1 Data reduction
Based on the topology of the SK detector, neutrino events are divided into three categories: fully-contained (FC) events, partially-contained (PC) events and upward-going
muons (Upµ). Both FC and PC events are neutrino events which interact inside the ID. FC
events are the ones whose final state particles are completely contained inside the ID while
PC means there are some final state particles entering the OD, typically the muon from a37

38
CC νµ interaction.
FC and PC events share the same set of good run selection criteria so they have identical
livetime. SK-I has a livetime of 1489.2 days and SK-II 803.919 days. Upµ event’s parent
neutrinos interact within 5 km underneath the SK detector in the Earth and the final stateparticle, nearly 100% muons, enter the SK detector. Upµ events mainly rely on fitting thelong muon tracks in ID so they are less susceptible to detector effects thus have looser
data quality cuts so longer livetimes are obtained. During the SK-I period, the Upµ sample
accumulated a livetime of 1645.91 days and for SK-II, 827.744 days. For both SK-I and
SK-II running periods, same data reduction schemes and algorithms are used but the cut
values are modified accordingly since SK-I and SK-II are, in principle, two different butvery “similar” detectors. Three different data reduction paths are used to separate FC, PC
and Upµ events.
Starting from raw data, to separate FC and PC events, a fast spatial clustering algorithm
is applied to the OD hits: if the number of hits in the largest OD cluster is less than 10,
the event is defined as FC; otherwise, it is defined as PC. Figure 4.1 shows the comparison
between MC and observation. Since the OD is essentially the same configuration, this cut
value is same for SK-I and SK-II.

4.1.1 Reduction of fully-contained events
The major background events in FC sample come from cosmic ray muons, low energy
events from radioisotopes and electronic noise. The FC event reduction process has five
steps.

39

10

10

4

10 4

FC/PC separation

number of events

10

3

2

10 3

10

2

10
0

1

2

3

4

5

6

7

8

9

10

1
0

20

40

60

80

100

120

number of outer detector hits

Figure 4.1: The number of hits in the largest outer detector cluster. Solid histogram is theMC prediction with neutrino oscillation considered.
FC first reduction
The first step of FC reduction uses simple and efficient criteria to reduce the amount ofdata to a manageable size. It has 3 cuts:• pe300 > 200 p.e.: pe300 is the total charge collected in ID within a 300 nsec time
window. By requiring more than 200 p.e.’s, low energy background events are rejected. This cut is reduced to pe300 > 100 p.e. for SK-II due to the fact that the ID
PMTs are reduced by 50%.
• nODhit800 < 50: nODhit800 is the total hits in OD within a 800 nsec time window.
This cut removes big amount of cosmic muons.
• tdif f > 100 µsec: tdif f is the time interval between the current event and the previous40one. By requiring the time interval, decay electrons entering ID from cosmic muons
are removed.
Those cuts reduce event rate from ∼ 106 events/day down to ~3000 events/day.
FC second reduction
The second step of FC is to further remove cosmic muon events and low energy backgrounds. It has 2 cuts:• nODhit800 ≤ 25 if petot < 100, 000 p.e.: petot is the total number of p.e.’s in ID.
This cut is similar to the second cut of step 1 but tighter to remove the relatively low
energy cosmic muons.
• pemax /pe300 < 0.5: pemax is the maximum number of p.e.’s recorded by any single
ID PMT. By requiring this ratio, low energy background events and flasher eventswith one single big hit are removed. A PMT sometimes flashes due to a discharge
around the dynode structure. The noise event cause by such a problematic PMT iscall a flasher event. Under this circumstance, one tube tends to record most of thephotons so this cut can remove some flasher events. More elaborate cuts will beapplied to remove the remaining flasher events in later steps.
After this reduction step, the event rate goes down to ~200 events/day.
FC third reduction
After background events are significantly reduced in step 1 and 2, more complex criteria(thus more computation time is needed) can be applied to further reject cosmic ray muons,
flasher events and low energy events which escaped previous simple cuts. At this stage,

41
the goal is to remove several specific kinds of background events and event reconstruction
tools are used to help the reduction.
• Remove through-going muon events which do not produce enough hits in OD. To
identify these through-going muons, following criteria are required:– pemax > 250 p.e.
– mugood> 0.75: mugood is the goodness being a through-going muon and isdefined as follow:

mugood =!

.
(ti −Ti )2
exp
−
/σi2
i
2(1.5σi )2!
2
i 1/σi,(4.1)

where ti and σi are the hit time of the i-th PMT and its uncertainty, T i is the
expected hit time based on the time when the muon enters ID and its reconstructed track.
– nODhit,in ≥ 10 or nODhit,out ≥ 10: nODhit,in(nODhit,out ) is the number of hit
PMTs in OD within 8m around the entrance(exit) point in a 800 nsec time
window. The entrance point is defined as the position of the earliest hit ID
PMT with neighbor hits; the exit point is defined as the center of the saturated
ID PMTs.
• Remove stopping muon events which escape cuts in steps 1 and 2. The following
criteria are used to identify stopping muons:– If mugood> 0.5, then nODhit,in ≥ 5 or nODhit,in ≥ 10: mugood of stoppingmuons is also defined as Eq. 4.1. The entrance point of a stopping muon isreconstructed in the same way for the through-going muon described above.

42• Remove cosmic ray muons entering from cable holes. As we described before inSec. 3.1, cosmic ray muons passing through the cable bundles in OD can not be
detected by OD and they will be mistakened as neutrino events. Topologically, thestopping ones will be considered as FC events. These events can be identified by thefollowing cuts:– One veto counter (described in Sec. 3.1) hit.
– lveto < 4m: lveto is the distance from the reconstructed vertex to the hit vetocounter.
• Remove accidental events. Occasionally, a low energy event and the following cosmic ray muon happen in the same trigger gate. These events can not be rejected by
the cuts in steps 1 and 2 because the absence of OD activities during the gate and theenough number of ID p.e.’s due to the muon. By combining the two following cuts,these accidental events can be removed:– nODhit,of f 500 > 20: nODhit,of f 500 is the number of hit OD PMTs in the fixed
time window between +400nsec and +900nsec after the trigger.
– peof f > 5000p.e.: peof f 500 is the number p.e.’s in ID in a fixed time window
from +400 nsec to +900 nsec after the trigger• Remove flasher events. Most flasher events have wider time distribution comparedto neutrino events. Flasher events can be identified by the following criteria:– nmin100 ≥ 14 or nmin100 ≥ 10 if nIDhit < 800: nmin100 is the minimum numberof ID hits in a 100 nsec sliding time window; n IDhit is the total number of ID
hits.

43
• Remove low energy events. The main source of low energy events at this stage areelectronic noise and decay of radioisotopes. They are distributive sources so these
events can be identified by the following cut:– nIDhit50 < 50: nIDhit50 is the number of ID hits within a 50 nsec sliding time
window after subtracting the each hit time by the time of flight assuming all thephotons come from the same vertex. The vertex is determined as the position
at which the time residual distribution is peaked.
After those cuts, the number of events is around 45/day.
FC fourth reduction
The forth reduction is dedicated to the remaining flasher events. As described previously,
flasher events are due to the discharges of PMT dynodes. Such “hot” PMT intend to discharge repeatedly before they are identified and shut off. Thus, the flasher events producedby such a tube have similar spatial hit patterns and can be identified by calculating the
spatial correlations between events. The number of ID PMTs is of order 10 4 so the computation power needed is very high. In order to reduce the amount of computation, ID
PMTs are divided into 1450 patches and each patch represents a spatial point in the calculation of the correlation coefficients. The correlation coefficient between events a and b isdefined as,Nr ab =

1 / (Qai − (Qa &)(Qbi − (Qb &)
N i=1
σaσb

,

where N is the number of spatial points (number of patches in our case), Q ai (Qbi ) is the
charge at each patch and σ a (σ b ) is its standard deviation of event a(b), and (Qa &((Qb &) isthe average charge of event a(b).

44
If the distance between the PMTs with the highest number of p.e.’s of two events issmall, which means the highest charge is produced around the same area, then they are
more likely to be a pair of flasher events. Considering the fact of the radius of ID tube isaround 50cm, an extra value of 0.15 is given to the correlation coefficient to account thisextra information if the distance is less than 75cm. The threshold of tagging similarity is
set atab= 0.168 log10
rth

peatot + pebtot+ 0.1302,where peatot (pebtot ) is the total number of p.e.’s in ID of event a(b). Above this threshold,
events are identified as similar events.
The calculation is made by pairing with the 10,000 events before and after the target
one, and the number of matching up two similar events are counted. If the number ofmatches is beyond certain values which depend on the maximum correlation coefficientamong the 20,000 pairs, the target event is identified as flasher.
After this reduction step, the event rate is 18 events/day.
FC fifth reduction
The last reduction step is to remove the remaining background events coming from different sources. This step is similar to the third reduction step but with more delicate reconstruction algorithms.
• Stopping muons. The cut to identify stopping muons is similar to the one in step3, nIDhit,in ≥ 5, but the entrance point by extrapolating the fitted track of the event
backward instead of using the position of the earliest hit PMT.

45
• Invisible muons. An invisible muon event is cause by the decay electron of a cosmic
ray muon which is below Cherenkov threshold in ID due to the energy loss along
the way entering ID. Such events have signals in OD before the triggers and can be
identified by the following criteria:– petot < 1000p.e.: petot is the total number of p.e.’s in IDmax– nmaxODhit,early > 4: nODhit,early is the maximum number hit PMTs in a 200 nsec

sliding time window from -8,900 nsec to -100 nsec
– nmaxODhit,early + nODhit,500 > 9 if lclust < 500cm or nODhit,early > 9 otherwise:nODhit,500 is the number of hit PMTs in OD in a fixed time window from 400 nsec to 100 nsec; lclust is the distance between two OD clusters used during
the calculation of nODhit,early and nODhit,500 respectively.
• Accidental events. More detailed calculations are employed to remove the remaining
accidental events:– pe500 < 300p.e.: pe500 is the number of p.e.’s in a fixed time window between-100 nsec and 400 nsec
max
– nmax
OD,late > 20: nOD,late is the maximum number of hit OD PMTs in a 200 nsec

sliding time window from +400 nsec to 1,600 nsec
• Long tail flasher events. These events are identified by the following cut:min
– nminIDhit,100 > 5 if the goodness of point fit < 0.4: n IDhit,100 is the minimum

number of hit ID PMTs in a 100 nsec sliding window from +300 nsec to+800 nsec; point fit will be explained later in the event reconstruction section
After the fifth step, event rate is 16 events/day.

46

4.1.2 Reduction of partially-contained events
PC events require OD activities, which makes the reduction rather hard because it gives
cosmic ray muons chances to get mixed in. PC reduction, like FC reduction, also has five
steps which are arranged in the same strategic way — simple and efficient reduction cuts
first then more delicate ones.
PC first reduction
Due to the high event rate, the strategy of the first PC reduction is to use fast and effective
cuts to reject through-going cosmic ray muons and low energy backgrounds. PC eventsshould meet the following criteria:• petot ≥ 1000 p.e. : petot is the total number of p.e.’s observed in ID. This cut is toreject the low energy backgrounds. A signal strength of 1000 p.e.’s corresponds to amuon with momentum of 310 MeV/c, while a PC muon event should have traveled
at least 2m thus its momentum is ∼ 500 MeV/c, well beyond this cut.
• Twid < 260 nsec: Twid is the width of the hit timing distribution among OD PMTs.
This cut is targeted at through-going cosmic ray muons — through-going muons
enter OD twice thus have a wide Twid distribution.
• nODcluster < 2: nODcluster is the number of OD hit PMT clusters. A hit cluster isformed by grouping PMTs with more than 8 p.e.’s within 8m. This cut is also to
reject through-going muons since they, in principle, leave two hit clusters in OD.
Obviously, the first step has no attempt to reject stopping cosmic ray muons. After the firststep, the event rates goes down to 14,000 events/day.

47
PC second reduction
The second PC reduction is designed to further remove the remaining through-going cosmic ray muons. A different clustering algorithm is applied in this step to help identifying
through-going muons. The OD (ID) walls are divided into 11 × 11 (21 × 21) patches andthe total charge in each patch is counted. Clusters are formed by calculating the charge
gradient. Following cuts are required for PC events:• nODcluster2 < 2: nODcluster2 is the number of clusters in OD using this new clustering
algorithm and require a hit PMT must have at least 6 p.e.’s.
min• nminODhit < 7: nODhit is minimum number of PMTs in the hit clusters of the top,bottom, or barrel region. This cut is to remove corner clipping through-going muon
events which leave hits in both top or bottom and barrel regions.
• pe200cm > 1000 p.e.: pe200cm is the number of p.e.’s in ID collected within 2maround the highest charged OD PMT in the OD cluster. This cut is to• n2nd−cluster : number of hits in the 2nd (ordered by number of hits in the cluster) OD
cluster• nmin : number of hits in the first cluster or the rest of OD hits depending on which
one is smaller• nouter : number of hits in the most charged OD cluster
PC third reduction
This step is to remove flasher events and cosmic ray stopping muons. To remove flasher
events which have wider time distribution compared to neutrino events, the same cut usedin FC third step is used:48• nmin100 ≥ 14 or nmin100 ≥ 10 if nIDhit < 800: nmin100 is the minimum number of
ID hits in a 100 nsec sliding time window; n IDhit is the total number of ID hits.
To remove stopping muons, a simple fitter call “point fit” is used to determine the directionof the ring. Then through backward extrapolating, the “entrance point” of the track in OD
is identified. For a true neutrino event, there should not be hit PMTs around this point
except noise hits; while for a real cosmic ray muon, it is the true entrance point and we
expect some OD hits in this region. So the stopping muons are identified by the following
cut:• nehit8m > 10: number of hits in OD within 8 meters of the back extrapolated entrance point on OD wall.
After this reduction step, the event rate is around 100 events/day.
PC fourth reduction
PC 4th reduction is targeted at the remaining cosmic ray muons. In this step, combined
information of point fit and through-going muon fit is used. The vertex position and track
direction are from point fit and mugood and mudistance are from through-going muon
fitter. PC events must pass the following cuts:• cosθfit−first > −0.8: θfit−first is the angle between the fitted direction and the line
connecting the vertex and the first hit tube in the ID cluster. These two vectors
should be back-to-back for cosmic ray muons.
• mudistance> 30 m if mugood> 0.85: mugood is the goodness value being a throughgoing muon which is defined in Eq. 4.1. Real cosmic ray muons have better goodness.

49• dcorner > 150 cm : dcorner is the distance from the fit vertex to the nearest ID corner.
This cut is to remove corner clipping through-going muons.
The event rate after the forth reduction drops to ~20 events/day.
PC fifth reduction
This last reduction step is to reject all the remaining background events. The cuts aredesigned specifically for each kind of distinctive background source.
• Remove low energy events. petot ≥ 3000: petot is the total number of p.e.’s in ID.
• Remove through-going muons. There are several different cuts targeting at the different signatures of through-going muons and we will explain them one by one:
– d1st−2nd > 20 cm and pe2nd > 10: d1st−2nd is the distance between the highest
charged cluster and the 2nd highest charged cluster; pe2nd is the number ofp.e.’s in the second highest charged cluster. These clustering information isfrom the algorithm in PC second.
– ncluster5 > 1: ncluster5 is the number of clusters in OD using a modified versionof the clustering algorithm in the 2nd PC reduction. Instead of using 11 × 11
patch dividing scheme for OD, a 6 × 6 one is used and the threshold of hit
PMTs changed from 2.0 p.e.’s to 0.5 p.e.
– nODhit,top > 7 and peOD,top > 10, nODhit,bottom > 7 and peOD,bottom > 10, and< Tdif f < 1.5 × 40m: nODhit,top (nODhit,bottom ) is the number of hit
0.75 × 40mccon the top (bottom) of OD; pe OD,top (peOD,bottom ) is the number of p.e.’s on the
top (bottom) of OD; Tdif f is the time difference between the average hit timesof the top and the bottom hits. This cut is targeted at the through-going muons
passing through the detector along the ID wall.

50< Tdif f < 1.5 × 40m: this– nODhit,in ≥ 5 and nODhit,out ≥ 5 and 0.75 × 40mcc
cut uses the information gained from a precise fitter, MS-fit, to remove theremaining through-going muons. n ODhit,in (nODhit,out ) is the number of hits
in OD within 8 m of entrance (exit) point using the precise fitter; T dif f is theaverage hit time difference between the entrance and the exit point hit PMTs.
• Remove stopping muons. There are three criteria in series to identify stoppingmuons with different features. A precise fitter, MS-fit, is used to provide the information removing stopping muons for the first two criteria. The third signature
utilizes stopping muon fitter.
– nODhit,in ≥ 10: this is the number of hit PMTs in OD within 8 m of theentrance point.
– θT DC−f it > 90◦ or θM S−f it > 90◦ : θT DC−f it (θM S−f it ) is the angle between
the ring direction using TDC- (MS-) fitter and the direction of OD cluster.
≥ 0.6 and nODhit,in > 6: mugood is the goodness being– mugood> 0 and pepecone
tot
a stopping muon; pecone is the number of p.e.’s inside the 42 ◦ cone and petot is
the total number of p.e.’s in ID.
• Remove cable hole events. Unlike FC events, the veto counters on top the cablebundles can not identify all the cosmic muons since there are PC events exiting thedetector and hit the veto counters. Two signatures are used to identify cosmic raymuons for PC reduction:– one veto counter hit– d0ring · d0vertex−veto > −0.8: d0ring is the ring direction reconstructed by TDCfitter and d0vertex−veto is the unit vector from the hit veto counter to the reconstructed vertex.

51

4.1.3 Reduction of upward-going muons
The goal of the upward-going muon (UPµ) reduction is to save the muons produced
through neutrino interactions in Earth. The main background are downward-going andhorizontally-going cosmic ray muons. The reduction of UPµs consists of two steps, first,
a automated reduction program is applied to reduce the background as much as possible; then, physicists will scan these events one by one using an event displayer to pickout the neutrino events. However careful, it is impossible to eliminate the cosmic raymuons around the horizontal direction. Thus, the cosmic ray muon contamination near thehorizontal direction is estimated by extrapolating the cosmic ray muon zenith angle distributions above the horizon. Events with only high-E triggers are considered as candidatesof UPµs.
Energy cut
We only use the events with pathlength greater than 7 meters in order to reduce the π ±background events produced by downward-going cosmic ray muons around the detector [74]. Thus, we can safely cut away the low energy events. On the other hand, it isdifficult to reconstruct ultra high energy events which saturate the Super-K ID. So, thefollowing cut on the ID charge is used,• 8000 ≤ petot < 1, 750, 000 : petot is the total number of p.e.’s in the ID
Background event rejection
Several different muon fitters specialized in different kind of UPµ events are used to reconstruct the pathlengths and the directions of the muons. The zenith angle cuts have twostages:52
1. In this stage, muon fitters stopmu1st and muboy are used and the following events
are rejected:• Corner clipper events: either of the following cut is met– muboy classifies the event as corner-clipper AND goodness value is greaterthan 0.35– muboy classifies the event as stopping or through-going AND the muon
pathlength is less than 4 meters AND goodness value is greater than 0.35
AND downward• Double muon events: both muboy and stopmu1st classify the event as multiple
muon• Downward-going cosmic ray muons: muboy classifies the event as downwardgoing and goodness-of-fit is greater than 0.45
The following events are classified as good UPµs and saved for the final eye scan:• Muboy output direction is upward and f cone > 0.8
According to the output information of muboy, the remaining events are separated
into through-going muons, which will be passed onto thrumu1st in the 2nd stage,
and stopping muons, which will be passed onto fitter stopmu2nd in the 2nd stage.
2. The second stage uses fitters stopmu2nd, thrumu1st, fstmu, thrumu2nd and NNfit tofurther reject the downward-going events. The logic in this stage is very simple —all the events from stage 1 will pass through all the fitters, if identified as upwardgoing events with goodness-of-fit above certain thresholds, they will be saved as
UPµ candidates for the final eye scan; if identified as downward-going events withgoodness-of-fit above certain thresholds, they will be rejected as cosmic ray muons;53
if identified neither as downward-going nor as upward-going but as horizontal by atleast one of the fitters or upward by thrumu1st, they will be saved as UPµ candidates
for the final eye scan. The rest are rejected.
Eye scanThe final eye scan over the UPµ candidate events are performed by experienced physicists
using event displayers. The scanning is performed by two physicists independently first.
Then they come together to make final decision.
Event classification
The following quantities are calculated for each event:• nOD,entry : number of OD hits between 800 ns and 1,300 ns and 8 m around the entry
point of OD• nOD,exit: number of OD hits between 800 ns and 1,300 ns and 8 m around the
projected exit point of OD
An event is considered as a stopping muon event if n OD,entry ≥ 10 and nOD,exit < 10;otherwise, a through-going muon event.
Through-going muon events are further divided into non-showering events and showering events by comparing the observed p.e.s with the expected number of p.e.s assuming
pure ionization energy loss. Detailed procedure is recorded in Chapter 7 of Dr. Desai’s
Ph.D. dissertation [74].

54

4.2 Event reconstruction
After the data reduction, event reconstruction algorithms are applied on all the neutrino
events. This procedure is same for both Monte Carlo and data. As we mentioned in theprevious section, some of the reconstruction algorithms are used during the data reduction
also. For the Super-K water Cherenkov detector, we have the following reconstruction
algorithms: vertex fitting, ring counting, particle identification, precise vertex fitting, andmomentum determination.

4.2.1 Vertex reconstruction
Vertex fitting process has three parts, point fitter, ring edge search and TDC fitter.
Point fitter
Vertex position is first roughly estimated using the timing information of ID hits assuming
all the photons are emitted at the same time from a point source at position (x, y, z). A
variable called “time residual” is defined in this step as follow,
ti = t0i −0(x − xi )2 + (y − yi )2 + (z − zi )2 /vwater,where t0i is the hit time of the ith PMT, whose position is (xi , yi, zi ) , in ID; vwater is speedof light in water. Finding the point (x, y, z) which maximizes the following goodness gives
us the rough position of the vertex,

N

“

”

2
i −t̄
1 / − t1.5σ
/2
e
Gp =
N i=1

,

55
where, N is the number of hit PMTs, t̄ is the average of ti ’s, and σ is the PMT timing
resolution (2.5 ns).
The particle direction is also roughly estimated by adding all vectors which are weighted
by the corresponding PMT charges connecting the vertex and all the hit PMTs:

d0 =

N
/

qi d0i,i=1where, d0i is the unit vector connecting the vertex and the i th PMT and qi is the number of
p.e.s deposited in the ith PMT. Certainly, d0i is not normalized here yet.
Ring edge search
In this step, the direction and the outer edge of the dominant ring is determined. First, anestimator is defined as follow,

Q(θedge ) =1 θedge0dθpe(θ)sin θedge2dpe(θ)dθ32θ

−( edge
σ

e

−θC 2
)

θ

,

where pe(θ) is the angular distribution of the observed charge as the function of opening
angle θ with respect to d0 from the point fitter, the azimuthal direction is integrated out; θ Cand σθ are the Cherenkov radiation angle and its resolution.
θedge is found by requiring:
1. θedge > θpeak , where θpeak is the position where pe(θ) takes its maximum;
2. d2 pe(θ)/dθ2 |θedge = 0, where θedge is the first point away from θpeak which satisfies
this equation.

56
TDC fitter
After the ring edge is found, TDC fitter is called to determine the vertex more precisely.
TDC fitter is based on the same principles used by point fitter: finding the point which reproduces the observed photon distributions the best. Using the information from previous
two step, by considering the following factors:
1. Photons are emitted along the track of charged particles;
2. Photons observed outside the ring edge are treated as scattering light,
TDC fitter is able to provide more precise estimation for the vertex position.

4.2.2 Ring counting
After the vertex, direction and the Cherenkov opening angle of the most energetic ring isconstructed using the vertex fitters, further attempts are tried to see whether the event is asingle- or multi-ring one. The algorithm consists of two steps:
1. Search for ring candidates using Hough transformation,
2. Test whether the candidates are real Cherenkov rings.
Searching ring candidates
Hough transformation is used to search for ring candidates. In order to perform Houghtransformation on our cylindrical surface, the ID PMTs are remapped onto a sphere whose
center is the vertex found by the vertex fitter and the surface of the sphere is divided into
72 × 36 bins in the θ − φ plane. To find the extra rings, expected Cherenkov charges
are subtracted from the PMTs and the remaining charges are redistributed using Houghtransformation with a Cherenkov opening angle of 42 ◦ . The charge peaks after the transformation are the directions of ring candidates.

57
Testing ring candidates
The candidates found by Hough transformation are tested using likelihood method. The
basic idea is to compare the likelihood function values of the case of n rings and the casen + 1 rings, staring from n = 1 till the n-ring case is preferred over the assumption ofn + 1 rings. The details are provided in Dr. Ishitsuka’s dissertation [102] or Dr. Clark’s
dissertation [103].

4.2.3 Particle identification
Super-K detector can identify two types of Cherenkov rings: e-like and µ-like, according
to the ring pattern and the opening angle. Tracks of electrons, particularly low-energy
electrons, are not straight due to multiple scatterings thus electrons produce rather diffused
Cherenkov rings. However, muons, as they are heavier, have more straight tracks thus
intend to produce much sharper Cherenkov rings. This difference is the basis of the particle
identification (PID).
The PID likelihood functions being e-like, L e , or µ-like, Lµ , for a Cherenkov ringare constructed using the expected charges (assuming the charged particle producing thisring as electron or muon) and the observed charges of the PMTs inside the opening angleof 1.5θC . By comparing the values of those two likelihoods, the charged particle type isdecided.

4.2.4 Precise vertex fitter
TDC fitter is not sensitive to the shift along the particle track due to the fact it uses thetime residuals. To overcome this short-coming, MS fitter is made to improve the vertex
resolution in the longitudinal direction.
In MS fitter, the likelihood function used in PID is employed with a fixed Cherenkov

58
angle θC which is calculated under the assumption of the decided particle type. The algorithm iterates in the following way:
1. Change the vertex position in the transverse direction with respect to the particle
track direction to maximize the goodness defined in TDC fitter.
2. Move the vertex position along the track direction to find the position where the likelihood function, Le or Lµ depending on the particle type, defined in PID algorithm
is maximized.
3. Adjust the direction to maximize L e or Lµ .
4. repeat till the shift on the vertex is less than 5 cm and the change on the direction isless than 0.5◦ .

4.2.5 Momentum reconstruction
The momentum of the charged particle producing a certain Cherenkov ring is reconstructedwith the information of the vertex, the direction and the PID. For Cherenkov radiation,
there is a very good correlation between the radiated photons and the particle momentum.
In order to get the corrected number of photons radiated by the particle, the detected photoelectrons are corrected for by the effective photo-sensitive area of the PMT, the pathlengthand the light attenuation and scattering in water. The total number of corrected photoelectrons, Rtot , is defined as the total number of photoelectrons in a 70 ◦ opening angle
and between the time window of -50 ns and +250 ns around the time residual peak. The
relation between Rtot and the particle momentum is obtained by Monte Carlo simulation.

59

4.2.6 Upward-going muon event reconstruction
The reconstruction of UPµ events are different from the FC and PC events. There arethree fitters used: OD fitter, TDC fitter, and Upmu fitter. According to the track length,
total p.e.s in ID, the number OD clusters and the direction differences of these three fitters,
the reconstruction process choose different fitters as the results. Details are recorded inRef. [74].
Due to fact that the primary neutrino interactions of Upµ events are in Earth, the correlation between the muon momentum and parent neutrino energy is relatively weaker.
In addition, the momenta of upward through-going muons can not be reconstructed since
they penetrate the Super-K detector. So the most valuable information from the UPµ reconstruction is the direction, which has better correlation with the parent neutrino directions
than FC and PC events.

4.3 Summary
Super-K events are categorized into following categories according to the locations ofneutrino interactions and the topology of the events:• Contained events: primary neutrino interactions are in ID– Fully contained events (FC): final state particles are completely contained within
ID∗ Single-ring e-like events∗ Single-ring µ-like events∗ Multi-ring e-like events: the most energetic ring is e-like∗ Multi-ring µ-like events: the most energetic ring is µ-like60
Each of these FC event samples are further divided into sub-GeV and multiGeV samples according to their visible energy,∗ Sub-GeV: Evis ≤ 1.33 GeV∗ Multi-GeV: Evis > 1.33 GeV– Partially-contained events (PC): final state particles enter OD∗ PC stopping: the muon penetrate ID but stop in OD∗ PC through-going: the muon penetrate both ID and OD• Upward-going muons (UPµ): primary neutrino interactions are in Earth and themuons enter Super-K detector– Stopping muons: muons that stop in ID– Non-showering muons: muons penetrate the detector and energy loss is mainly
through ionization– Showering muons: muons penetrate the detector and losing energy mainly
through radiation
Figure 4.2 shows the parent neutrino energy spectra for different data samples based on
the Monte Carlo simulation. We will talk how we get those neutral current enhanced
samples in Chapter 6. As we see from the spectra, Super-K atmospheric neutrino covers awide energy range from ∼ 100 MeV to ∼ 1 TeV. This is a big advantage for testing new
physics at high energies.

4.3.1 Systematic uncertainties of data reduction
The data reduction algorithms for SK-I and SK-I are same except the cuts related to numberof ID PMTs. So, the same kinds of systematic uncertainties are introduced in. However,61

Figure 4.2: Parent neutrino energy spectra
the values are different between SK-I and SK-II. Since SK-II is not as good as SK-I due
to less ID PMTs, systematic uncertainties which depend on the ID PMT coverage becomelarger.
Table 4.1 shows the values for SK-I and SK-II. Appendix B gives an example on how
to evaluate the systematic uncertainty in PC reduction.

4.3.2 Systematic uncertainties of event reconstruction
The reconstruction algorithms are same for SK-I and SK-II except some parameters. So,
the same set of systematic uncertainties are evaluated for SK-I and SK-II which are shownin Table 4.2.

62

Table 4.1: Super-K data reduction uncertainties
Systematic Uncertainty
Value (SK-I) Value (SK-II)
FC reduction
0.002
0.002
PC reduction
0.026
0.040
UPµ reduction
0.01
0.01
FC/PC separation
0.009
0.005
Non-ν background (e-like)
0.01
0.01
Non-ν background (µ-like)
0.01
0.01
Stopping/through-going µ separation
0.004
0.003
Non-showering µ background
0.015
0.030
Stopping µ background
0.17
0.24
Non-ν background (Multi-GeV single-ring e-like)
0.14
0.326
Non-ν background (Multi-GeV multi-ring e-like)
0.20
0.554
Non-showering/showering separation
0.01180.017
Showering µ background
0.13
0.24
PC stopping/through-going separation
0.12
0.166

Table 4.2: Systematic uncertainties in event reconstruction
Systematic uncertainty
Value (SK-I) Value (SK-II) Index
Ring separation
0.10
0.10
34
PID of single-ring events
0.01
0.01
35
PID of multi-ring events
0.10
0.10
36
Energy cut of stopping µ
0.011
0.011
38

Index
26
27
28
29
31
32
33
40
41
42
43
45
46
47

Chapter 5
Mass-Induced νµ − ντ Oscillation
In the framework of the Standard Model, there are three types of neutrinos, νe , νµ and ντ .
That there are no more active neutrino flavors beyond these three is strongly supported by
the experimental evidence from LEP experiments [13]. However, the LSND experimentmeasured a mass squared splitting which is different from both the solar and the atmospheric scales [104]. If the LSND result is confirmed, it could suggest there is a fourth
flavor. Due to the constraints from LEP experiments, the fourth flavor must be sterile
neutrino. We will disscuss issues related to sterile neutrinos in the next chapter.
As we described in the first chapter, previous analysis [105] and results from experiments Chooz [38] and Palo Verde [39] indicate that the mixing between muon neutrinos
and electron neutrinos are extremely small. Thus, a two-flavor analysis of ν µ − ντ oscillation for atmospheric neutrino is quite sufficient at the current experimental precision to
extract mixing parameters between the second the third generations.
In this chapter, we will carry out the oscillation analysis assuming that muon neutrinos
mix with tau neutrinos using both SK-I and SK-II atmospheric neutrino data. In order toprovide more precise information on the mixing parameters, a finer binning of the reconstructed event energy is applied in the data analysis.
63

64

5.1 Data preparation
Mixing between muon neutrinos and tau neutrinos implies that the missing atmosphericmuon neutrinos are turning into the tau type. For a water Cherenkov detector, it is not
so easy to identify tau events [106] because tau lepton decays almost instantaneously andthe final state particles consist of electrons, muons, or one or more pions which produce
multi-ring signals without an easily identified leading charged lepton. In addition, due tothe energy threshold of tau lepton production which is ∼ 3.5 GeV, only a small portion of
tau neutrinos interact in SK detector produce tau leptons. Thus, there is no much statistics
available either. For the whole SK-I period, we only expect 78 tau events [106] assuming∆m2 = 2.4 × 10−3 eV2 . Super-K collaboration studied the tau appearance study using
both likelihood method and neural network approach, and the results are consistant withνµ − ντ oscillation assumption [106].
Due to the difficulties of identifying tau events and their low statistics, SK is basically a atmospheric neutrino disappearance experiment for the νµ − ντ oscillation. These
small amount of tau events (78 for SK-I and 43 for SK-II normalizing the livetimes) are
categorized into e−like events and they have very little effects on the oscillation analysis
compared to the total number of events. The survival probability of muon neutrinos is asfollow:Pνµ →νµ = 1 − sin2 2θ sin2

∆m2 L
∆m2 (eV 2 )L(km)= 1 − sin2 2θ sin2 1.27.
4E
E(GeV )

(5.1)

As we will show later, the zenith angle oscillation analysis is based on the comparison
between the observed number of events and the expected ones. In order to know the expected number of neutrino events, we need to know the neutrino flux. As we mentioned

65
in Chapter 2, the systematic uncertainty on the absolute flux intensity of atmospheric neutrinos is as large as ~20%. However, we do have better understanding on the flavor ratio.
Thus, though electron neutrinos do not participate oscillations in the two generation framework for the atmospheric neutrino sector, e−like events can help to cancel the systematic
uncertainty of the absolute neutrino flux. In order to increase the statistics of e−like events,
we also extract charge current enhanced samples from FC multi-ring e-like events. The algorithm is based on a likelihood analysis using PID likelihood, momentum fraction of themost energetic ring, number of muon decay electrons, and distance between muon decay
electron and primary neutrino interaction position [105].
The following data samples are used in our oscillation analysis:• Fully-contained single ring e-like events• Fully-contained multi-ring CC e-like events• Fully-contained single ring µ-like events• Fully-contained multi-ring µ-like events• Partially contained events• Upward-going muons
Neutral current events are not included since they do not provide any information for aνµ − ντ oscillation analysis. But we will see that neutral current enhanced samples arehelpful in the sterile neutrino analysis which we will discuss in next chapter.

5.1.1 Data binning
Examining the survival probability formula Eq. ??, the constraint on the mixing angle termsin2 2θ is mainly decided by available statistics. The mass-squared difference is coupled

66
with L/E, which means the better resolution on L/E, the better constraint we gain on∆m2 . Zenith angle oscillation analysis has the advantage that it can utilize all the neutrino
events to provide the best constraints on the mixing angle. In this analysis, pathlengths
and energies of parent neutrinos are not reconstructed because the reconstruction can not
be applied to some of the data samples, like Upµ, or the reconstructed L and E are notcorrelated to the true values, like low energy multi-ring FC events.
Pathlength of a neutrino is determined by the production height and the zenith angle.
The experimental variable which is the most closely correlated to the pathlength of theparent neutrino is the reconstructed zenith angle (θrec ) of the event. The reconstructed
direction of an event is decided based on the event type and reconstruction information:
For FC single ring events, the reconstructed ring direction is the event direction; For FC
multi-ring events, the total reconstructed momentum defines the event direction; For PC
and Upµ events, directions of the leading muons define the event directions. We bin all
neutrino events into 10 zenith angle bins based on cos θ rec value.
For contained neutrino events, i.e. FC and PC events, the experimental variables which
are the most closely correlated to the parent neutrino energy are the visible energy or the
total reconstructed momentum depending on the event types. Binning is carried out basedon the logarithm of the reconstructed momentum or the visible energy. For PC events, inorder to further improve the energy resolution, we separate the PC stopping events andthe PC through-going events by comparing the expected OD hits and the observed OD
hits [102]. The PC stopping and PC through-going events separation was made possible
after the OD Monte Carlo was improved by tuning some of the simulation parameters. The
tuning process is recorded in Appendix A.
For Upµ events, because the neutrino interactions happen in the Earth, there is no
good correlation between the parent neutrino energy and the energy deposited in SuperK detector. However, the three categories of Upµ samples do statistically correspond to

67

Figure 5.1: Binning strategy for the oscillation analysis
neutrinos with different typical energies. No further binning is done for Upµ events interms of energies.
SK-I and SK-II data are binned in the same fashion. Figure 5.1 shows the energy
binning for different SK samples. In total, we have 760 bins combining SK-I and SK-II
data samples.

5.2 Oscillation Analysis
A least chi-square method called “pull method” [107] is used to extract the mixing parameters out of the data. The chi-square is defined as follows:
2

χ

n
/

m
/
Niobs
1j
exp
obs
obs
=
2(Ni − Ni + Ni ln exp ) +( sys )2 .
Ni
σj
i=1
j=1

(5.2)

Here, Niobs is the number of observed events in the ith bin, σjsys is the j th systematic
uncertainty, 1j is the pull value on the j th systematic uncertainty, Niexp is the expected

68
number of events in the ith bin, and n and m are number of data bins and number ofsystematic uncertainties respectively. In total, for ν µ −ντ oscillation analysis, there are 760
data bins as shown in the previous section and 70 systematic terms recorded in previous
chapters.
Among the systematic uncertainties, as we shown in previous chapters, there are 14
neutrino flux uncertainties recorded in Table 2.1; there are 12 neutrino interaction related
ones which are the same for SK-I and SK-II recorded in Table 3.2; there are 12 detector
related systematic uncertainties which are different between SK-I and SK-II recorded inTable 3.1, Table 4.1 and Table 4.2; and there is one solar activity term recored in Table 2.1
which has the same uncertaity value but different pulls for SK-I and SK-II. The effects ofthe systematic uncertainties on the expected number of events in each bin are taken into
account in the following way:
Niexp = (1 +m/

fji 1j )Niexp0

j=1Here, Niexp0 is the predicted of number of events in the ith bin without considering thesystematic uncertainties and fji is the response coefficient of the number of the eventsin the ith bin with respect to the j th systematic uncertainty. Niexp0 is calculated in thefollowing way:

Niexp0 = Ninosc Psurvival ,(5.3)

where Ninosc is the predicted number of events without considering oscillation and P survivalis the neutrino survival probability. And f ij is evaluated using MC method.
To find the parameters minimizing the χ 2 , the solutions to the pull terms can be ob-

69tained by taking ∂χ2 /∂1k = 0:!nn/
Niobs fki1ki=1!m
− sys =
fki Niexp0.
iσkj=1 (1j fj + 1)
i=1
Keeping only the linear terms of 1/

(5.4)

!m

i
j=1 (1j fj + 1), since we expect 1i ’s are small, thisnon-linear equation is transfered into a linear equation and can be solved by inverting am × m matrix. However, such a procedure does not apply to the mixing parameters. So,
we prepare a grid in the mixing parameter space and solve the linear equation set at eachgrid point for the pulls 1i . The grid point with the minimum chi-square values gives thebest-fit values for the mixing parameters.

5.3 Results and Discussion
Based on the previous experimental results of Super-K, the parameter space we choose forthe νµ − ντ oscillation analysis is sin 2 2θ ∈ [0.7, 1.2] and log ∆m2 ∈ [−3, −2] with theunit of ∆m2 as eV2 . The parameter space is divided into 201 × 201 points uniformly. The
best-fit values are sin2 2θ = 1 and ∆m2 = 2.5 × 10−3eV2 .
Figure 5.2 shows the best-fit zenith angle distributions. We see that the ν µ − ντ oscillation accounts well for the missing muon neutrinos. Figure 5.3 shows the chi-squarecontours corresponding to 68%, 90% and 99% confidence levels respectively and Fig. 5.4
shows the χ2 behaviors around the best-fit values of sin2 2θ and ∆m2 .
At 68% confidence level, the allowed parameter space is sin2 2θ > 0.95 and 2.2 ×
10−3 eV 2 < ∆m2 < 2.7 × 10−3 eV 2 . The chi-square value at the best-fit point is 848 with
the number of degrees-of-freedom (dof ) 755. Due to the fine binning we adopted, manyof the bins have very small number of entries. Thus, the p − value corresponding to this70

Figure 5.2: Chi-square contours on the ∆m2 − sin2 2θ parameter plane for νµ − ντ mixing

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-30.7

0.8

0.9
1
sin22$

1.1

1.2

Figure 5.3: Chi-square contours on the ∆m2 − sin2 2θ parameter plane for νµ − ντ mixing

71

868
866
864
862

&2

860
858
856
854
852
850
848
0.7

0.8

0.9

1

sin22$

1.1

1.2

-3

10-2

10

%m2(eV2)

Figure 5.4: 1D chi-square projections
chi-square value should be calculated via Monte Carlo instead of assuming a chi-square
distribution with dof = 755. A custom toy MC program is designed by assuming eachbin follows a Poisson distribution with the Poisson parameter λ equals to the expected
number of entries. Our simulation shows that the p − value is 18% for the least χ 2 value
we obtained.
Figure 5.5 shows the pull term distributions. It follows a Gaussian distribution, which is
expected [107]. A Gaussian fit gives the mean µ = 0.04 and the standard deviationσ = 0.4.
In Ref. [107], a Gaussian distribution with a standard deviation of 1 is expected in principle. However, due to the fact that we are generally conservative on some systematic
uncertainties, the Gaussian distribution for the pulls we obtained here has a smaller standard deviation.

72

20

pulls

Number of Pulls

18
16&2 / ndf 14.96 / 1714

Constant

12

Mean0.03935

10

Sigma

0.4251

5.518

8
6
4
2
0
-3

-2

-1

0

1

2

3

Normalized Pulls
Figure 5.5: The distribution of normalized pulls with its Gaussian fit results

Chapter 6
Constraints on Sterile Neutrinos
Mass induced neutrino oscillations have been widely accepted in recent years as experimental evidence accumulates [15, 35, 108, 36, 6]. The analysis we did in the last chapter
is based on this theory, and indeed, the model fits the observation very well. However,whether there are sterile neutrinos (νs ) involved in neutrino oscillations is of both experimental and theoretical interests. There are observations which are in favor of the existenceof sterile neutrinos, for example: if there are sterile neutrinos, the LSND anomaly could be
explained assuming the existence of sterile neutrinos [104, 109]; supernova rapid neutron
capture process (r-process) problem can be solved [110, 111]. Theoretically, a simple way
to give neutrino mass is to have Dirac right-handed neutrinos which are SU(2) singlet thus
must be sterile [21, 22].
Previous studies using 1100 live days of SK-I data rejected the ν µ − νs oscillation
at 99% confidence level [112]. In this chapter, using the pull method described in theprevious chapter and the complete SK-I and SK-II data sets, we design a systematic way
to do model comparison to test which model explains SK data better. We will also study
how much admixture of sterile neutrinos is allowed in the atmospheric neutrino oscillation.

73

74

6.1 Comparison of νµ − ντ and νµ − νs oscillations
6.1.1 Signatures of νµ − νs oscillation
For the two competing models, ν µ −ντ and νµ −νs , though the 2nd neutrino is not identified,there are some signatures for us to tell two oscillations apart. Firstly, sterile neutrinos do
not interact with matter at all so they do not even make any neutral current signals in SK
detector; secondly, because the potentials experienced by ν µ and νs in matter are different,
matter effect modifies the survival probabilities of muon neutrinos, as we described in
Chapter 1.
To utilize the first signature, we need to identify the neutral current events. For a water
Cherenkov detector, by identifying the π 0 particles which are produced in neutral current
single-pion neutrino reactions, neutral current events can be enriched quite well [113].
However, in order to keep as many events as we can, we do not adopt the tight cuts in
Ref. [113]. As we will see later, even without the tight cuts in Ref. [113], the neutral
current enhanced events based on the less restricted cuts defined below can differentiate
two models quite well already. For sub-GeV events, the cuts we use to enhance neutral
current events are:• Multi-ring events: π 0 particles produced in neutral current reactions produce two
e-like rings.
• The most energetic ring is e-like: this cut is chosen to get rid of CC ν µ events.
• 400 MeV < Evis < 1330 MeV: this cut is chosen to guarantee certain directional
information of the parent neutrinos.
As for the multi-GeV samples, we simply use the left-over events after enhancing thecharged current enhanced events using the likelihood method developed in the previous

75
225

200

200

180
160
140

150

Number of Events

Number of Events

175

125
100
75

100
80
60

50

40

25
0

120

20
-1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

cos$

0.8

1

0

-1

-0.8 -0.6 -0.4 -0.2

0

0.2

0.4

0.6

0.8

1

cos$

Figure 6.1: Zenith angle distributions of neutral current enhanced samples. Left: sub-GeV;
right: multi-GeV
Table 6.1: Portions of NC and CC events in the NC enhanced data samples
NC events CC νµ CC νe
Sub-GeV sample
37%
22%
41%
Multi-GeV sample
24%
35%
41%
chapter.
After applying these cuts, we have the zenith distributions of neutral current enhanced
samples shown in Fig. 6.1. The hatched areas are neutral current events based on theMonte Carlo simulation. The percentages of neutral current events and charge current ν µcontaminations are shown in Table 6.1.
The second signature of νµ − νs oscillation is the matter effect. To have a quantitative idea about the matter effect experienced by νµ − νs oscillations for atmospheric
neutrinos, we calculate the survival probability of muon neutrinos crossing the Earth using the method described in Chapter 1 and the Earth density profile provided by PREM
model [114]. As we can see from Fig. 6.2, for atmospheric neutrinos, the oscillation probability is indeed suppressed by the matter effect around and above 10 GeV.

76

Figure 6.2: Survival probabilities of muon neutrinos crossing the Earth

6.1.2 Best-fits of two models
We now have both charged current and neutral current enhanced data samples to do theoscillation analysis. Considering that the uncertainty on total number of multi-ring e-likeevents is relatively big, we give neutral current enhanced an independent flux normalization, just like we did for the charged current enhanced e-like events. As for the binningof these neutral current enhanced samples, the sub-GeV neutral current events are treatedas one energy bin. The multi-GeV sample is divided into 4 energy bins based on log E vis
values since the parent neutrinos have wider energy spread. Both samples are divided into
10 zenith angle bins based on the reconstructed event direction.
Using the same least chi-square method, we perform a νµ − ντ oscillation using bothCC and NC data samples. The same best-fit values for the mixing parameters are obtained.
However, the goodness-of-fit is worse compared to the results based on CC samples only.
The minimum chi-square is 971.2 with a number of degree-of-freedom 853. Using thesame toy Monte Carlo procedure, we found the p−value is 7.3%. Figure 6.3 shows the chisquare confidence contours. As for the νµ − νs oscillation, we simply replace the P survival
for muon neutrinos in Eq. 5.3 with a routine which has the matter effect implemented.
Maximal mixing is still preferred by ν µ − νs oscillation as shown in Fig. 6.4. However,the ∆m2 value is driven up a bit, 3.5 × 10−3 eV2 . This is because, as we argued in theprevious sub-section, the survival probabilities of ν µ − νs oscillations are suppressed, thus,

77

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-30.7 0.75 0.8 0.85 0.9 0.95
sin22$

1

Figure 6.3: χ2 contours of νµ − ντ oscillation (with NC enhanced samples)in order to fit the Super-K data, the best-fit value of ∆m2 is increased to compensate forthis suppression, which is also why the constraint on mixing angle is tighter than the ν µ −ντcase. The chi-square value at the best-fit point is 1023.6 with the same number of degreeof-freedom, which is much worse than νµ − ντ model. The difference 52.4 corresponds toa 7.2σ exclusion level for the pure νµ − νs model.
Figure 6.5 shows the comparison of best-fit zenith angle distributions. For many binsof the data sample, the νµ − νs model reproduces SK observation as well as the νµ − νsmodel. However, for sub-GeV neutral current enhanced events, PC through-going events,
upµ stopping events and non-showering µ events, the ν µ − ντ model fits the data better
than the νµ − νs model. Table 6.2 shows the detailed chi-square differences for different
categories of events. The ∆χ2 from data bins is 38.3. And Fig. 6.6 shows the comparisonof the distributions of the pull terms. The pull term distribution of the ν µ − νs oscillation is

78

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-30.7 0.75 0.8 0.85 0.9 0.95
sin22$

1

Figure 6.4: χ2 contours of the νµ − νs oscillation (with NC enhanced samples)wider than the νµ − ντ one, which is another sign that νµ − ντ oscillation is a better modelaccounting for the atmospheric neutrino data. Pull terms contribute 14.1 to the chi-squaredifference out of the total difference 52.4.
The νµ − νs model can not reproduce sub-GeV neutral current enhanced events as wellas the νµ − ντ model. The νµ − ντ model expects the data fairly flat which is what we
observe. However, for the νµ − νs model, the loss of neutral current events makes the
expected distribution tilted toward upward going bins where muon neutrinos oscillate into
sterile ones.
As we showed in Fig. 6.2, around the atmospheric ∆m 2 scale, the strongest matter
effect happen around 10 GeV. The typical energy of PC through-going events, Upµ stopping events and non-showering µ events is around this scale, which is why these events
can differentiate the two models better than other types of neutrino events. We also see

79

Figure 6.5: Comparison of the best-fit zenith angle distributions of ν µ − ντ and νµ − νsoscillations
Table 6.2: Chi-square difference breakdown (∆χ2 = χ2νµ −νs − χ2νµ −ντ )Data Samples
Bins (SK-I+SK-II) χ2νµ −νs χ2νµ −ντ ∆χ2Single ring sub-GeV e-like
50+50
104.8 104.0
0.8Single ring multi-GeV e-like
50+50
108.6 110.7 -2.1
Multi-ring multi-GeV CC e-like
50+50
86.6
85.8
0.8
Single ring sub-GeV µ-like
50+50
104.9 106.1 -1.2
Single ring multi-GeV µ-like
30+30
64.8
66.8
-2.0
Multi-ring µ-like
40+40
79.3
75.5
3.8
NC-enhanced sub-GeV events
10+10
19.5
14.5
5.0
NC-enhanced multi-GeV events
40+40
105.7 104.5
1.2
PC stopping µ40+40128.7 125.8
2.9
PC through-going µ40+40
114.4 102.1 12.3
Upward stopping µ
10+10
21.1
14.1
7.0
Upward non-showering µ
10+10
28.1
16.9
11.2
Upward showering µ
10+10
24.5
25.0
-1.5
TOTAL
430+430
991.1 952.8 38.3

80

20
18
Standard

Number of Pulls

16

Sterile

14
12
10
8
6
4
2
0
-3

-2

-1

0

1

2

3

Normalized Pulls
Figure 6.6: Distributions of the normalized pull terms of ν µ − ντ and νµ − νs oscillations
that though neutral current enhanced events help to tell two models apart, it is mainly thematter effect that distinguish two models due to the advantages such as long baselines,
high energies and high matter densities that atmospheric neutrinos possess compared tolong baseline experiments.

6.2 An admixture model
Compared to νµ − ντ oscillation, the fact that the νµ − νs oscillation is excluded at 7.2σ
level does not mean that there is no involvement of sterile neutrinos in Super-K atmospheric neutrino oscillations. It is still possible that a certain portion of muon neutrinos are
oscillating into sterile neutrinos besides the dominating ν µ − ντ oscillation. In this section,

81

m4
m3

νs
ντ
νµ
νe

m2
m1

Figure 6.7: A case of 2+2 mass hierarchy and four-flavor mixing
we will study a case of the admixture of sterile neutrinos.
In the case of 2+2 mass hierarchy shown as the diagram in Fig. 6.7, according to
G. L. Fogli et. al. [115], based on the current experimental results, the four flavor mixing
can be transformed into two 2-flavor mixings by constructing two superposition states ofντ and νs in the following way:










 ν+   cos ξ sin ξ   ντ 
=
.


ν−
νs
− sin ξ cos ξ
For the atmospheric sector, the oscillation now is between ν µ and ν+ which is a superposition state: cos ξ|ντ &+sin ξ|νs &. Thus, the portion of sterile neutrinos is sin 2 ξ. Accordingly,
the matter effect strength is weaken by a factor of sin2 ξ. Now, we have 3 parameters in ouroscillation: the mixing angle between |ν m3 & and |νm4 &, ∆m2 = m24 − m23 and the portionof sterile neutrinos sin2 ξ. Since the maximal mixing is a very strong constraint, we aregoing to assume maximal mixing in our analysis, i.e. sin 2 2θ = 1.
Using the same least chi-square framework, the best-fit is sin 2 ξ = 0 which means
Super-K data prefer zero amount of sterile neutrino involvement, see Fig. 6.8. However, at
90% C.L., we do allow 23% of sterile neutrino admixture. Compared to the limit of 67%
at 90% C.L. set by Fogli et al in Ref. [115] using the published Super-K data, our limit
using the complete SK-I and SK-II data sets is much more stringent.

82

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-3

0

0.2

0.4

0.6

0.8

sin2'

Figure 6.8: The νs admixture allowance

1

Chapter 7
Violations of Lorentz Invariance and
CPT
The mixing of mass eigenstates is not the only physical origin of neutrino oscillation as
pointed out in by Coleman and Glashow [44] and Kostelecky et. al. [45], Lorentz Invariance Violation (LIV) and CPT symmetry violation (CPTV) can also cause neutrinos tooscillate even they are massless.
In addition to the possibility of explaining current observations, it is also of theoreticalimportance to test the limits of those two fundamental symmetries [44, 116]. Lorentz
Invariance and CPT symmetry have been long considered as fundamental laws of nature [117, 118, 119, 120] and serve as the corner stones of the quantum field theory. It
is not completely clear that those symmetries will remain perfectly unbroken at all energies. Some quantum gravity theories predict that LIV and CPTV might happen at Planckscale, MP ∼ 1019 GeV. Thus testing the limits of their validity is a way to explore newphysics beyond the Standard Model [44, 121, 116, 122, 123]. Neutrino oscillation, as anatural interferometer, is a promising method to test the low energy effects of LIV and
CPTV.
83

84

7.1 An overview LIV and CPTV
LIV and CPTV are expected in some Planck scale theories through spontaneous symmetry breaking mechanism [44, 124, 125] or topological defects [126]. A. V. Kostelecky
et. al. set up a general framework of LIV and CPTV called the Standard Model Extension (SME) [124] to study the possible observational effects at low energies. By treating
the SM as a low energy effective theory of underlining fundamental physics, all the possible Lorentz symmetry breaking terms are added into the SM Lagrangian using the SM
fields and gravitational fields.
Lorentz invariance is intimately related to CPT symmetry as stated in the CPT theorem [117, 118, 119, 120]: for a local field theory, Lorentz invariance leads to CPT symmetry. While breaking CPT symmetry naturally causes Lorentz invariance violation, which
is proven in [127]. To generate neutrino oscillation effects, it turns out that the minimal
renormalizable version of SME [44, 123] is sufficient. The “minimal” SME (mSME) only
has a CPT-violating term and a CPT-conserving but LIV term in addition to the originalStandard Model terms:↔
1
LmSM E = −(aL )µAB L̄A γ µ LB + i(cL )µνAB L̄A γ µ D ν LB ,
2

(7.1)

where the first term is CPT-odd and the second term is CPT-even. LA/B , A and B runover the three generations e, µ, and τ , are the left-handed SU(2) doublet fields, D is thecovariant derivative, and the CPTV coefficients aL and the LIV-only coefficients cL arise
from spontaneous symmetry breaking of more fundamental theories.
From the perspective of experimental observation, the LIV phenomena can be categorized into rotationally symmetric and rotationally asymmetric cases. We will study therotational symmetric cases by keeping only the temporal components of those two inter-85
action terms. Also, in this study, we always confine ourselves within two generations, ν µand ντ .
Under these considerations, the effective Hamiltonian in the flavor basis is:(Hef f )AB = |0p|δAB +

5
1 4
(aL )0AB p0 − (cL )00AB p0 p0 .
|0p|

Two 2 × 2 matrices (aL )0AB and (cL )00AB , where A and B are generation indices, are
not necessarily diagonal in the flavor basis. First, let us consider the matrix (c L )00AB .
Diagonalizing it by rotating the neutrino fields by θ c gives two eigenvalues, ca and cb , andtwo eigenstates, |νa & and |νb &. In this new basis, if neutrinos are massless, the Hamiltonian
then becomes: ca 0 (Hef f )ab = pc + 
 p.
0 cb
The energies of the two eigenstates are

Ei = pc + pci = p(c + ci ), i = a, b.

This can be interpreted as that two eigenstates have different speeds of light “c”, i.e. neutrinos travel at slightly different speeds in vacua. Those two eigenstates are defined as
Maximal Attainable Velocities (MAV) [128, 44] eigenstates, which can be viewed as theenergy eigenstates at infinite momentum for massive particles. Flavor eigenstates are superpositions of MAV eigenstates and vice versa:

86



 |νµ & = cos θc |νa & + sin θc |νb & |ντ & = cos θc |νa & − sin θc |νb &,where θc is the mixing angle between two different MAV eigenstates, |ν a & and |νb &. We
define cT T = ca − cb ; without loss of generality, we can always assume c T T > 0. The
survival probability of muon neutrinos after traveling a distance L is:Pνµ →νµ = 1 − sin2 2θc sin2 (cT T LE/2) .

(7.2)

A more general formula of this form is discussed in Ref. [129], which considers the energy
dependence of the LIV effect. Ignoring neutrino mass, the survival probabilities of muon
neutrinos are as follow:

Pνµ →νµ = 1 − sin2 2θc sin2 (κLE −α ) .

(7.3)

Where real number α is from the energy dependence effect of LIV, and the Planck scale,
MP , is already absorbed into κ which has dimension of energy to the power of α − 1. This
formula will also address the case of CPTV-induced oscillation when α = 0, as we will
discussion later.
With the analysis framework established in the previous chapter, we can easily test
how well this model fits our observation using the Super-K atmospheric neutrino data.
The parameter space here is spanned by sin2 2θc , κ and α. Using the charged current data
samples, the best-fit parameters are: sin2 2θc = 1, κ = 4.8 × 10−3 eV0.16 , and α = 1.16.

87)(a.u.)

10-2

10-3
-1

-0.5

0

0.5
(

1

1.5

2

Figure 7.1: Allowed regions for parameters κ and α
Projecting the chi-square onto κ − α plane, Figure 7.1 shows the allowed regions of theparameters at confidence levels of 68%, 90% and 99%. To better examine the effect ofenergy dependence, the chi-square vs. α plot is shown in Fig. 7.2.
For the oscillation of the format L × E, the chi-square value is 1060. Compared tostandard L/E oscillation, L × E type oscillation is strongly disfavored. At 68% C.L.,
0.95 < α< 1.30. Clearly, Super-K observation favors the standard L/E oscillation.
Now let us examine the effect of matrix (aL )0AB . The effective Hamiltonian in thebasis of diagonalizing (aL )0AB is: aa 0 (Hef f )ab = pc ± .
0 abThus, the CPTV-breaking term will produce neutrino oscillation of the form ∆aL, where

88

1050

&2

1000
950
900
850
-1 -0.5

0

0.5

1

1.5

2

(
Figure 7.2: χ2 versus α for a general oscillation formula∆a = aa − ab and aa and ab are the two eigenvalues of (aL )0AB . This oscillation isenergy independent. Obviously, this is not the case we observe at Super-K. As shownin Fig. 5.2, the difference between the null oscillation prediction and the observation has
a clear energy dependence. Also shown in Fig. 7.2, the α = 0 case, which is the ∆aL
oscillation, is strongly disfavored.
We then can ask the question: how much LIV and CPTV could be allowed in the
Super-K atmospheric observation? We will address this issue in the following sections.

7.2 Allowed limits on LIV parameter
From the analysis in the last section, we see that it’s unlikely for LIV alone to explain theatmospheric neutrino observation. However, if we consider the LIV-induced oscillation as

89a sub-dominant effect in addition to the mass-induce oscillation, the survival probability
becomes:Pνµ →νµ = 1 − sin2 2Θ sin2 (ΩL/4),
where, Ω sin 2Θ = |∆m2 sin 2θ/E + 2cT T eiη E sin 2θc | Ω cos 2Θ = ∆m2 cos 2θ/E + 2cT T E cos 2θc.

In this expression, ∆m2 and θ are the mass eigenstate mixing parameters as usual and η
is the phase difference between two mixing matrices of the mass eigenstates and the MAV
eigenstates.
Considering the two cases, η = 0 and η = π, assuming the best-fit mass eigenstate mixing parameter values based on Super-K atmospheric neutrino data, ∆m2 = 2.5 × 10−3 eV2and sin2 2θ = 1, we can test the limits allowed by Super-K atmospheric neutrino observation on the Lorentz invariance violating parameter c T T . Using the chi-square analysis
framework set up in Chapter 5, we can find the best-fit values in the parameter spacespanned by sin 2θc and cT T .
In the η = 0 case, the best-fit values are sin2 θc = −0.12 and cT T = 5 × 10−25 .
Chi-square contours corresponding confidence levels 68%, 90% and 99% are shown in
Fig. 7.3. We note that we obtain non-zero values for LIV parameters, sin 2θv = −0.12 andcT T = 0.5 × 10−24 . However, the zero value is within the first sigma level. At 90% C.L.,the allowed LIV limit is: c T T < 1.2 × 10−24 . In the η = π case, chi-square contours
are shown in the right panel of Fig. 7.3. The best fit point is at sin 2θ v = −0.02 and

90

-23

-23

cTT(10 )1

cTT(10 )

1

10-1

10-2
-1

10-1

10-2
-0.5

0
sin2$v

0.5

1

-1

-0.5

0
sin2$v

0.5

1

Figure 7.3: Chi-square contours of LIV. Left: η = 0; Right: η = π
cT T = 0.6 × 10−24 . At 90% C.L., the allowed LIV limit is: c T T < 1.3 × 10−24 .

7.3 The violation of CPT symmetry
7.3.1 Testing the LSND mass squared splitting scale
One ad hoc way to reconcile the anomalous LSND result [104] with other experiments is
to allow neutrinos and anti-neutrinos have different mass squared splittings [130], which
will violate CPT symmetry. As the conversed of the CPT theorem [117, 118, 119, 120],CPTV naturally leads to LIV and Ref. [127] provides a detailed proof.
Super-K atmospheric neutrinos, as we described in Chapter 2, have both neutrino andanti-neutrino components. It is not possible for the Super-K detector to identify the incident neutrino is neutrino or anti-neutrino event by event. However, due to the facts thatneutrinos and anti-neutrinos have different fluxes, different interactions, and different final state particles, Super-K is able to see a statistical effect if neutrinos and anti-neutrinos
have different mass squared splittings. Early work done by M. Messier on this topic canbe found in Ref. [131].
Assuming maximal mixing for both neutrinos and anti-neutrinos, MC simulated zenith
angle distributions are generated for different mass squared splittings for neutrinos and91

2% m (eV2)10-2

10-3 -31022%m (eV )10-2

Figure 7.4: Contour plot of the ad hoc CPTV model chi-square. Dash line is definedby∆m2 = ∆m̄2 .
anti-neutrinos and compared with the observation using the least chi-square method described in the last chapter. Figure 7.4 shows the contour plot of the least chi-square on theplane of ∆m2 − ∆m̄2 .
The best-fit values of the minimum chi-square point is ∆m 2 = 3.7 × 10−3 eV 2 and∆m̄2 = 1.5 × 10−3 eV 2 . Shown in Fig. 7.4, the 90% confidence level contour covers
part of the equal line of neutrino and anti-neutrino mass squared splittings. It is clear that
Super-K atmospheric observation has no indication of the LSND scale for anti-neutrinomass squared splitting.

92

7.3.2 Allowed limits of CPTV
As we see in Sec. 7.1, it is unlikely that Super-K atmospheric observation can be explainedby the massless neutrino oscillation induced by the temporal component of the CPTV
term in the minimal Standard Model. The scale of the LSND mass squared difference is
excluded by Super-K observation. Seeing these facts, even if CPT had been violated, itwould have be a sub-dominant effect.
Considering the eigenstates defined by the CPT-odd interaction (a L )0AB L̄A γ 0 LB sameas the mass eigenstates [132], the dispersion relation including this CPTV interaction is asfollows:

Ei = p +m2i± ai ,
2p(7.4)

here ai , i = 1, 2, is the eigenvalue of the matrix (aL )0AB and “±” signs corresponds toneutrinos and anti-neutrinos. Without losing generality, we can define ∆a ≡ a 1 − a2 ,assuming a1 > a2 . The survival probability of muon neutrinos (anti-neutrinos) becomes:2pνµ /ν̄µ →νµ /ν̄µ = 1 − sin 2θ sin

22∆m2 L ∆aL
±
4E
2

3

.

(7.5)

Where “+” and “-” correspond to neutrino and anti-neutrino oscillations respectively. ThisCPTV model allows neutrinos and anti-neutrinos to have different mass squared split2), which is the ad hoc model we studied in Sec. 7.3.1
tings (by absorbing ±∆a into ∆m
4E
by requiring that the deviations of the mass squared splitting of of neutrinos and antineutrinos from the best-fit value are same. Equation 7.4 enables us to study the allowed
limit of CPTV by treating it as a sub-dominant effect in addition to the dominant mass-93%m2(eV2)

10-2

10-3 00.5

1% a(10 GeV)
-23

1.5

2

Figure 7.5: Chi-square contours considering CPTV as a sub-dominant effectinduced oscillations. We can assume the mixing is still maximal for this model and theparameter space is spanned by ∆m2 and ∆a. Using Eq. 7.5 to produce MC zenith angle
distributions and to compare with Super-K atmospheric data, the chi-square contour plot
on the ∆m2 − ∆a plane is shown in Fig. 7.5.
The best-fit value of ∆a is 0, which means the standard ν µ − ντ is recovered. At
90% C.L., the limit on the temporal component of the CPTV violating interaction is ∆a <1.05 × 10−23 GeV.

7.4 Discussion
The Standard Model Extension has all the possible terms breaking Lorentz invariance and
CPT symmetry, which provides an ideal model to test the low energy effects of LIV and
CPTV. Different types of experiments are sensitive to different types of effects. For example, due to the fact that atmospheric neutrinos come from all directions, they are notso advantageous on testing directional LIV effects. Baseline experiments with neutrino

94
beams with well-defined directions are suitable for testing such effects [133]. However,
the wide energy coverage of atmospheric neutrinos gives Super-K advantages for testing
effects sensitive to high energies.
As comparisons, the LSND experiments limit on a L -type SME coefficients is of order10−19 GeV and the limit on cL -type is of order 10−17 [133]. By using the complete data
sets of SK-I and SK-II and treating LIV as sub-dominant effect of mass-induced neutrino
oscillations, we are able to provide a better limit on c T T than the one set in Ref. [134],
which uses published Super-K oscillation analysis results and sets a limit c T T < 6 × 10−24
at 90% C.L. The limit on cT T set by Gonzalez-Garcia and Maltoni in Ref. [135] using bothSuper-K and K2K data is 0.81 × 10−24 , which is slightly better than our limit. However,
our limit on ∆a is better than the one set by Gonzalez-Garcia and Maltoni, which is ∆a <3.2 × 10−23 GeV at 90% C.L. Besides the different data sets, the main difference between
the analysis in Ref. [135] and ours is that our analysis has more systematic uncertainty
terms and their effects on the data bins are simulated using Monte Carlo method, while thestudy by Gonzalez-Garcia and Maltoni is based on a smaller set of systematic uncertainties
and the “coupling” factors, i.e. the response factors we defined in Chapter 5, are calculatedwith a more coarse binning. The Super-K limit on c T T is of the same order of the limit set
by the MACRO experiment, which is 6 × 10−24 at 90% C.L.
Compared to other types of experiments, limits on LIV and CPTV from neutrino experiments are quite comparable and sometimes even better. Based on the cosmic ray spectrum, limits on the deviation of MAV from the speed of light, i.e. |1 − c|, can be set.
<∼ 10−24 − 10−15 depending on the particles species and details of the analysis [44, 136];
based on nuclear magnetic resonance frequencies, authors of Ref. [44] interpret the resultsfrom Ref. [137] as a limit of 3 × 10−22 on this deviation. As for the limit on CPTV, theK 0 − K̄ 0 system gives a limit of m K 0 − mK¯0 < 0.44 × 10−18 GeV [138]. According toBarger et al [132], data from g − 2 experiment can be analyzed to probe a0L to the orders

95of ∼ 10−22 − 10−25 GeV.

Chapter 8
Neutrino Decay and Neutrino
Decoherence
Due to the limited energy and angular resolutions of the Super-K experiment, models predicting “vanishing neutrinos” caused either by neutrino decoherence or by neutrino decay
can produce effects which are similar to neutrino oscillations. Testing whether such models could explain what we observe at Super-K is obviously importrant. On the other hand,both neutrino decoherence and neutrino decay require new physics beyond the Standard
Model. With the wide energy and pathlength coverage, Super-K atmospheric neutrino data
samples are in an advantagous position to test the effects of these new physics.
The neutrino decoherence phenomenon is due to the dissipative interaction of a neutrino system with an environment [139, 50, 48, 49, 46]. While the neutrino decay models
we consider here are based on the assumption that the neutrino mass eigenstates could beunstable [140, 141, 51]. Although the physics behind the two theories are quite different,at the experimental level for Super-K, the expected observations and the techniques to test
their validity and limits are very similar — they both cause the disappearance of neutrinos
instead of oscillation.
96

97

8.1 Neutrino decoherence
Neutrino decoherence may arise when we consider a neutrino system couples to an environment. It could be due to the small violations of quantum mechanics caused by quantum
gravity effect [139, 50, 46, 142], a very large neutrino background [48, 49], or flavor dependent interactions into extra dimensions [48]. If the coherence lengths of muon neutrinos
are less than the diameter of the Earth due to the presence of those possible dissipativeinteractions, similar effects to neutrino oscillation would be observed at Super-K — theamount of observed muon type neutrinos is less than predicted. In this section, we test
whether neutrino decoherence could be accounted for the Super-K atmospheric neutrino
zenith angle distributions.
The coupling to the environment produces mixed quantum states so a Schödinger equation description for the system is not possible any more, and we are forced to use the Liouville equation instead. Let ρ(t) be the density matrix which is Hermitian and tr[ρ] = 1 ,then the Liouville equation is:dρ= −i[Hm , ρ].
dt

(8.1)

In the case of a two neutrino flavor system, i.e., a two-level system, the Hamiltonian of1diag{m22 , m23 }. Solving the evolution equation 8.1 and computing
the system is Hm = 2E

tr(|νµ &(νµ |ρ(t)) gives the survival probability of the muon neutrinos.
The decoherence effect can be put into the evolution equation in the following fashion [46]:dρ= −i[Hm , ρ] − D[ρ].
dt

(8.2)

98
The dissipative interaction term D[ρ] leads the transitions of the system from pure to mixed
states. Operator D has the dimension of an energy, and its inverse defines the typicalcoherence length of the quantum system.
Although the physics behind the term D[ρ] could have different origins, there are certain phenomenological requirements it has to meet in order to make the model physical.
The first general one is the so-called complete positivity — the mapping from ρ(0) to ρ(t)must be linear, Markovian, and trace-preserving. With this requirement, D[ρ] must be of
Lindblad form [46]:D[ρ] =

n
/
4
a=1

5
{ρ, Da† Da } − 2Da ρDa† .

(8.3)

Where Da ’s are Lindblad operators arising from averaging away environment dynamic.
It is also natural to require that the von Neumann entropy S = −tr(ρ log ρ) increases
with time monotonically and the conservation of the statistical average of the energy ofthe system, i.e., dtd tr(Hm ρ) = 0, which requires Da ’s are Hermitian and [Hm , Da ] = 0
respectively. For a two-level system, Da can be expanded in terms of Pauli matrices,
Da = d0a · 0σ , where da ’s are real and the dimension is the square root of an energy. With
those three requirements, after some algebra, [46, 47], the survival probability of muon
neutrinos with energy E after traveling a distance L can written as:

Pνµ →νµ = 1 −Where γ =

1 2
∆m2 L
sin 2θ(1 − e−γL cos
).
2
2E

(8.4)

!n

2
a=1 da has the dimension of energy and it’s the coherence length of thesystem under the influence of the environment dynamics. When the coherence length
1/γ goes to infinity, the survival probability reduces to the pure mass-induced expression

99
Eq. 1.4.
However, the exact functional form of the parameter γ is not known but considering
the fact that Super-K atmospheric neutrino energy spectrum covers several magnitudes, it’scertain that our data would be very sensitive to its dependence on the neutrino energy [46].
So, we will try three cases, n = −1, 0 and 2 for its power-law dependence on E, γ =γ0 (E/GeV )n suggested in [46].

8.1.1 A pure neutrino decoherence model
As pointed out in Ref. [46], in the case of γ ∝ 1/E, even without the contribution fromneutrino oscillation, neutrino decoherence can explain atmospheric neutrino data very
well. We will first test this pure neutrino decoherence model using SK-I and SK-II atmospheric neutrino data.
The survival probability in this pure decoherence case is:Pνµ →νµ = 1 −1 2
sin 2θ(1 − e−γ0 L/E ).
2

(8.5)

This situation may happen when two masses are very close to each other or degenerate.
Beside the different energy and pathlength dependence, it is obvious that neutrino decoherence will have less neutral current events also due to the fact that neutrinos are just simply“disappearing”. Thus, we will include the NC-enhanced data samples in our analysis also.
Using the same chi-square analysis formalism, replacing the survival probability formula with Eq. 8.5, we are able to find the best values for parameters γ 0 and sin2 2θ for thispure decoherence model to reproduce the SK atmospheric neutrino zenith distributions.
Chi-square contours of 68%, 90% and 99% on the sin 2 2θ − γ0 are shown in Fig. 8.1. The
best fit value for this model is: γ 0 = 1.3 × 10−21 GeV , sin2 2θ = 1, and the minimum
chi-square at the best-fit point is 987 with the number of degree-of-freedom 853.

100

1

sin22$

0.98
0.96
0.94
68% C.L.
90% C.L.

0.92

99% C.L.

0.9
0

0.05

0.1

0.15

0.2

0.25

0.3

* (10-21GeV2)
Figure 8.1: Chi-square contours on the sin2 2θ − γ0 plane
Compared to the minimum chi-square of the ν µ − ντ oscillation model, which is 971
with the same number of degree-of-freedom, this pure neutrino decoherence model is excluded at confidence level of 4 σ.
Let us check which data samples are contributing to the chi-square difference. Table 8.1 shows the breakdown of the chi-square difference, χ2decoherence − χ2νµ −ντ . The
difference of the two models comes from the decoherence model being non-oscillatory.
Precisely speaking, we can not simply isolate the contributions from individual data samples since they are being analyzed together and thus correlated with one another. But still,
the samples that have the strongest oscillation behavior should prefer the ν µ − ντ model,
while it would be harder for other samples to tell the difference. And, this is what we
see in the table. Multi-ring µ-like events, PC events, upward stopping muons and upward
non-showering muons provide the strongest preference to oscillation. Because of the badenergy and angular resolutions of NC-enhanced multi-GeV events, the standard oscillation
does not account for their zenith angle distribution as well as the decoherence. Overall, the101Table 8.1: Chi-square difference breakdown between the standard oscillation and the decoherence model (∆χ2 = χ2decoherence − χ2νµ −ντ )Data Samples
Bins (SK-I+SK-II) χ2decoherence χ2νµ −ντ ∆χ2Single ring sub-GeV e-like
50+50
104.50104.06 0.44
Single ring multi-GeV e-like
50+50
110.15110.75 -0.61
Multi-ring multi-GeV CC e-like
50+50
85.25
85.76 -0.51
Single ring sub-GeV µ-like
50+50
103.66106.11 -2.44
Single ring multi-GeV µ-like
30+30
67.51
66.79 0.72
Multi-ring µ-like
40+40
81.34
75.49 5.85
NC-enhanced sub-GeV events
10+10
15.22
14.47 0.75
NC-enhanced multi-GeV events
40+40
101.44104.51 -3.07
PC stopping µ40+40127.15125.77 1.38
PC through-going µ40+40
104.09102.14 1.94
Upward stopping µ
10+10
17.79
14.08 3.71
Upward non-showering µ
10+10
18.37
16.86 1.51
Upward showering µ
10+10
25.90
25.97 -0.07
TOTAL
430+430
962.37952.77 9.6
standard oscillation model explains the SK data better than the decoherence model. The
total chi-square difference from the data bins is 9.6.
The chi-square difference from pull terms is 7.4. Figure 8.2 shows the term by termcomparison of the pull values between the standard oscillation and the decoherence case.
For most of the systematic uncertainty terms, the decoherence model needs to pull themmore to fit the data. To see the overall effect, Fig. 8.3 shows their statistical comparison.
The decoherence pull terms follow a wider distribution than the standard ν µ −ντ oscillation
model. This is another sign that the decoherence model is disfavored by our data. From
Fig. 8.2, we see the largest pull for the decoherence model is the 14th one, which is the
PC+Upµ sample-by-sample normalization factor. The decoherence model needs to “pull”
the PC+Upµ normalization factor down more than the standard oscillation model to fit thedata, which is because the exponential dependence of the survival probability on the L/Evalue.

102
3! µ-! +

2

Decoherence

Pull

1
0
-1
-2
-3
0

10

20

30

40

50

60

70

Uncertainty Index

Figure 8.2: Comparison of pull terms of ν µ − ντ and decoherence

8.1.2 Allowed limits of neutrino decoherence
It is unlikely for neutrino decoherence alone to explain the SK atmospheric neutrino data
based on the model comparison study in Sec. 8.1.1. However, we can not eliminate thepossibility that neutrino decoherence indeed happens in atmospheric neutrino oscillations
as a sub-dominant effect. Then, survival probabilities of neutrinos will depend on bothneutrino oscillation and neutrino decoherence, shown in Eq. 8.4. By fitting Eq. 8.4 against
SK atmospheric neutrino data, the contribution of neutrino decoherence, assuming its existence, can be separated from neutrino oscillation. We will test three different kinds ofenergy dependence, n = 0, 1, 2, of the decoherence factor γ = γ0 (E/GeV )n as suggested
in Ref. [46].
The most trivial case, n = 0, assumes that decoherence is independent of neutrinoenergy. The best-fit recovers the standard neutrino oscillation, i.e. γ = γ 0 = 0 GeV2 ,sin2 2θ = 1 and ∆m2 = 2.4 × 10−3 eV 2 . Projecting the chi-square onto ∆m2 − γ0plane, confidence contours of 68%, 90% and 99% are shown in Fig. 8.4. At 90% C.L.,
γ0 can be as large as 0.76 × 10−23 GeV. Compared to the limit set by Lisi et al is
γ0 < 3.5 × 10−23 GeV [46] based on the published SK-I zenith angle distributions, our103

Number of Pulls

20
18! µ -! + oscillation

16

Decoherence

14
12
10
8
6
4
2
0
-3

-2

-1

0

1

2

3

Normalized Pulls
Figure 8.3: Comparison of pull terms of ν µ − ντ and decoherence
limit is approximately 5 times better due to the fine detail allowed by our analysis and thecomplete data sets of SK-I and SK-II.
The inverse energy dependence of γ, i.e. γ = γ0 /(E/GeV ) is because we do not
expect to see different decoherence effects in different frames and L/E is Lorentz invariant [46].
In this case, the best-fit gives a non-zero decoherence: γ0 = 8 × 10−23 GeV . As shownin Fig. 8.5, this non-zero value of γ0 is consistent with zero. The best-fit values for the
mixing parameters are sin2 2θ = 1 and ∆m2 = 2.4 × 10−3 eV 2 which are consistent withthe νµ − ντ best-fit values. At 90% C.L., γ0 < 0.61 × 10−21 GeV. The limit set by
Lisi et al in Ref. [46] is γ0 < 2 × 10−21 GeV. Again, our analysis based on the finely
binned experimental data gives a much better limit.
The third assumption for energy dependence of γ is inspired by the Planck scale
MP [46]. We generally believe that the low energy effects of Planck scale physics aresuppressed by MP . Considering that γ should have the dimension of an energy, the Planck104%m2(eV2)

10-2

10-3

0

0.5

1

-23

1.5

2

* (10 GeV)
0

Figure 8.4: Allowed decoherence limits for γ = γ 0 case%m2(eV2)

10-2

10-3

0

0.5

1

-21

1.5

2

* (10 GeV)
0

Figure 8.5: Allowed decoherence limits for γ = γ 0 /E case105

%m2(eV2)

10-2

-3

10

0

0.05

0.1
-27

0.15

0.2

* (10 GeV)
0

Figure 8.6: Allowed decoherence limits for γ = γ 0 (E/GeV )2 casescale physics should be suppressed in the form of γ = γ 0 E 2 /MP . By absorbing MP into
γ0 , we have γ = γ0 (E/GeV )2 and we still keep [γ0 ] = 1 for convenience.
In this case, the best-fit also gives a non-zero decoherence effect: γ0 = 4 × 10−29 GeV .
As shown in Fig. 8.5, this non-zero value of γ 0 is consistent with zero. The best-fit values
for the mixing parameters are sin2 2θ = 1 and ∆m2 = 2.4 × 10−3 eV 2 which are consistent
with the νµ − ντ best-fit values. At 90% C.L., γ0 can be as large as 4.8 × 10−28 GeV, whichis twice better than the limit of 9 × 10 −28 GeV set by Lisi et al in Ref. [46].
In conclusion, we do not observe any significant non-zero decoherence effects. However, Super-K data can provide very stringent limits which can restrict possible new physics
predicting neutrino decoherence effects. One thing we need to point out is that the decoherence models here are purely phenomenological results without the details of fundamental
theories involved. They are very general, but it also means that these decoherence param-

106
eters can only be interpreted given specific theories.

8.2 Neutrino decay
According to the mass-induced neutrino oscillation theory, flavor eigenstates are superpositions of mass eigenstates. However, there is no particular reason to consider the mass
eigenstates are stable [51, 141] and the decay of the mass eigenstates can be realized
through different theories described in Ref. [140, 141] and references therein.
In the atmospheric neutrino sector, the muon neutrino is a superposition of two mass
eigenstates: |νµ & = cos θ|ν2 & + sin θ|ν3 &, where |ν2 & and |ν3 & are mass eigenstates. Without
losing generality, we can always assume m2 > m3 . If the |ν2 & state is not stable and its
lifetime is τ , the survival probability of ν µ is:Pνµ →νµ = sin4 θ + cos4 θe−αL/E + 2 sin2 θ cos2 θe−αL/2E cos

∆m2 L
,
2E

(8.6)

where α = m2 /τ is the experimental parameter we can test since we do not know thevalue of either m2 or τ . This survival probability formula can have two simplified cases
based on two particular models:
1. If we assume ν2 → ν¯3 + J and J is a massless scalar, then it can be shown that
∆m2 must be greater than 0.73 eV 2 in order to satisfy the constraints from K →2

Laverages to zero andµ+neutral particles decay [140]. Then, the term cos ∆m2E

Eq. 8.6 is simplified into:Pνµ →νµ = sin4 θ + cos4 θe−αL/E .

(8.7)

2. If there exists a sterile neutrino and it does not mix with other neutrinos (or mix very

107
little),










 |νµ &   cos θ sin θ 0   |ν2 & 

 


 |ν &  =  − sin θ cos θ 0   |ν &  ,
τ
3

 



 


|νs &
|ν4 &
0
0
1
and ν2 → ν¯4 + J, then the ∆m2 in Eq. 8.6 is not constrained any more and can bevery small. The survival probability becomes:
Pνµ →νµ = (sin2 θ + cos2 θe−αL/2E )2 .

(8.8)

We will first test how the two simplified cases explain the Super-K zenith angle distributions.

8.2.1 Comparisons between neutrino decay and neutrino oscillation
Using the same least chi-square framework described in Chapter 5, replacing the survival
probability Psurvival with Eq. 8.7 and Eq. 8.8, we can find the best decay parameters reproducing Super-K data. The decay parameters we use are sin2 θ and α, here [α] = 2.
For pure decay case 1, Pνµ →νµ = sin4 θ + cos4 θe−αL/E , the best-fit values are sin2 θ =0.04 and α = 1.5 × 10−21 GeV2 . The minimum chi-square is 1256.5 with dof = 853.
Figure 8.7 shows the chi-square contours on sin 2 θ − α plane. Compared to νµ − ντ
oscillation, the minimum chi-square is larger by 285.3, which corresponds to an exclusionlevel 16.9 σ. This type of decay is very unlikely to happen.
For pure decay case 2, Pνµ →νµ = (sin2 θ + cos2 θe−αL/2E )2 , the best-fit values are
sin2 θ = 0.68 and α = 2.2×10−22GeV2 . The minimum chi-square is 983 with dof = 853.
Figure 8.7 shows the chi-square contours on α − sin 2 θ plane. Compared to νµ − ντ

108

0.2
68% C.L.
90% C.L.

sin2$

0.15

99% C.L.

0.1

0.05

0
0

0.1

0.2

0.3

0.4

0.5

((10-20GeV2)
Figure 8.7: Chi-square contours of neutrino decay case I
oscillation, the minimum chi-square is larger by 12, which corresponds to an exclusionlevel 3.5 σ. This type of decay fits Super-K zenith angle distributions much better thancase I, however, compared to the standard oscillation, it is still unlikely to happen.
Table 8.2 summarizes the contributions to the chi-square differences from different categories of events. They more or less contribute to the chi-square difference in the same
way as the decoherence case. This is because both model produce similar observation effect although the physics behind them is very different. However, the case 2 of neutrino
decay fits SK data much better than the neutrino decoherence case. The chi-square difference from data bins is 7.2. The rest of the chi-square difference comes from the pull terms,
which is 4.8. Shown in Fig. 8.9, the decay model (case 2) has a wider pull term distribution. Figure 8.10 shows the term-by-term comparisons between the pulls of the standard
neutrino oscillation and the neutrino decay case 2. Again, the largest pull in the decay case

109

0.8
68% C.L.
90% C.L.

sin2$

0.75

99% C.L.

0.7

0.65

0.6
0

0.1

0.2

0.3

0.4

0.5

0.6

((10-21GeV2)
Figure 8.8: Chi-square contours of neutrino decay case II
is the 14th one, which is the sample-by-sample normalization for PC+Upµ events.

8.2.2 Decay as a sub-dominant effect
As we see from the comparisons between neutrino decay and neutrino oscillation, it isunlikely that neutrino decay can account for the Super-K zenith angle distributions. However, this does not mean that neutrino decay is not happening at all — it may coexist withneutrino oscillation. Treating neutrino decay as a sub-dominant effect, we can study theallowed limit on neutrino decay which can provide quantitative guidance to models thatpredict neutrino decay.
Equation 8.6 has both neutrino oscillation and neutrino decay effects included. Fitting
this formula against the Super-K zenith angle distributions, we can find the best parameters for both neutrino oscillation and neutrino decay parameters. In the 3-dimensional

110

Table 8.2: Chi-square difference breakdown between the standard oscillation and the decaycase II model (∆χ2 = χ2decayII − χ2νµ −ντ )Data Samples
Bins (SK-I+SK-II) χ2decayII χ2νµ −ντ ∆χ2Single ring sub-GeV e-like
50+50
104.11 104.06 0.05
Single ring multi-GeV e-like
50+50
110.27 110.75 -0.49
Multi-ring multi-GeV CC e-like
50+50
85.28
85.76 -0.47
Single ring sub-GeV µ-like
50+50
105.28 106.11 -0.83
Single ring multi-GeV µ-like
30+30
65.40
66.79 -1.39
Multi-ring µ-like
40+40
79.06
75.49 3.57
NC-enhanced sub-GeV events
10+10
15.44
14.47 0.97
NC-enhanced multi-GeV events
40+40
102.00 104.51 -2.51
PC stopping µ40+40
126.75 125.77 0.98
PC through-going µ40+40
104.97 102.14 2.82
Upward stopping µ
10+10
16.92
14.08 2.84
Upward non-showering µ
10+10
18.30
16.86 1.44
Upward showering µ
10+10
26.18
25.97 0.20
TOTAL
430+430
959.96 952.77 7.19

20
18

! µ -! + oscillation

Decay II

Number of Pulls

16
14
12
10
8
6
4
2
0
-3

-2

-1

0

1

2

3

Normalized Pulls
Figure 8.9: Comparison of pull terms of the ν µ −ντ oscillation and the neutrino decay case
2

111
3
! µ-! +

2

Decay (II)

Pull

1
0
-1
-2
-3
0

10

20

30

40

50

60

70

Uncertainty Index

Figure 8.10: Comparison of pull terms of the ν µ − ντ oscillation and the neutrino decaycase 2parameter space spanned by α (0 < α < 3 × 10−23 GeV2 ), sin2 θ and ∆m2 , we find thebest-fit values are: α = 0, sin2 θ = 0.45 and ∆m2 = 2.4 × 10−3 eV2 . Super-K obviously
prefers zero amount of neutrino decay. Figure 8.11 shows the 2-dimension projection of
the 3-dimension chi-square function.
However, the 99% C.L. contour is not completely covered by the parameter space we
explored. Since we already know that α = 0 is the global minimum in the exploredparameter space using linear scale, we can now switch to logarithmic scale to explore
larger parameter space. After trying several rounds, we find the proper parameter space tocover the complete 99% C.L. contour. Figure 8.12 shows the complete chi-square contours
in a much larger parameter space.
As we can see from the plot, there is actually a local minimum corresponding non-zeroneutrino decay. The local minimum is shown more clearly in the one dimension projectionof the chi-square function in Fig. 8.13. The parameter values at the local minimum are:α = 4.5 × 10−22 GeV2 , sin2 θ = 0.7 and ∆m2 = 2.5 × 10−3 eV2 . We see that α = 0 is
indeed the global minimum and the ∆m 2 at the local minimum is of the Super-K atmospheric mass squared splitting scale. The existence of the non-zero local minimum is due to

112

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-3
0

0.2

0.4
0.6
-22
((10 GeV2)

0.8

1

Figure 8.11: Chi-square contours of neutrino oscillation + neutrino decay in linear scale
the fact that Super-K has limited resolutions on neutrino energy and direction and limited
statistics around the first oscillation maximum. For models whose oscillation probabilities
have L/E dependence, the Super-K L/E analysis has more power to tell the details ofthese L/E dependence, thus differentiate various models apart more clearly [143]. In theframework of L/E analysis, this local minimum only appears at 99% C.L.
Despite the existence of the local minimum, we can set a limit on α at 99% C.L.:α < 2.1 × 10−21 GeV2 .

8.3 Discussion
Although neutrino decoherence and neutrino decay have completely different physical origins, they produce similar effects — neutrinos vanishing instead of oscillating. In the case

113

10-2
68% C.L.
90% C.L.

%m2(eV2)

99% C.L.

10-3 -2
10

10-1

1
10
((10-22GeV2)

102

3

10

Figure 8.12: Chi-square contours of neutrino oscillation + neutrino decay in logarithm
scaleof neutrino decoherence, neutrinos indeed disappear due to the dissipative term D[ρ] in
Eq. 8.2. However, in the case of neutrino decay, neutrinos are not “disappearing”. They
decay into either particles which are detectable in our detector or another mass eigenstate
which is dominantly a sterile state and does not mix. The power to distinguish neutrinodecoherence and neutrino decay from neutrino oscillation mainly comes from the high energy and the long baseline events, which is understandable since both exotic effects havedependences such as e−γL and e−αL/E .

114

990

68% C.L.

988

90% C.L.

986

99% C.L.

&2

984
982
980
978
976
974
972
10-2

10-1

1

10

((10-22GeV2)

102

103

Figure 8.13: The non-zero local minimum of neutrino decay

Chapter 9
Conclusion
Using Super-K atmospheric neutrino data, we compare the standard ν µ − ντ oscillation
model with three types of non-standard ones: the ν µ − νs oscillaton which is due to themass eigenstate mixing, neutrino oscillations induced by LIV and CPTV, and neutrino
vanishing caused by neutrino decoherence and neutrino decay. We find that the ν µ − ντ
oscillation induced by mass eigenstate mixing is the best model accounting for the Super-K
atmospheric zenith angle distributions.
Furthermore, by treating these non-standard phenomena as sub-dominant effects inaddition to the standard νµ − ντ oscillation, we are able to constrain the limits of new
physics thus provide quantitive guidance to future studies in these fields. Owing to thewide coverage of neutrino energy, neutrino pathlength and the matter density of Super-K
atmopheric neutrino events, most of the limits obtained in this dissertation are the best
currently available. Our study shows that neutrinon oscillation, as a natural interferometer,
is a powerful tool to probe new physics beyond the Standard Model.

115

Appendix A
Outer Detector MC Tuning
Neutrino oscillation analysis of Super-K heavily depends on the comparsion between observation and Monte Carlo (MC) simulation. It is essential for us to have a reliable MC
which simulates the detector response well. As we mentioned in Chapter 3, the OD
functions mainly as a veto detector for cosmic ray muons so we can identify the events
happening inside or outside the detector and incoming cosmic ray muons. The accurate
performance of the OD simulation has not been an issue until we had the need to further
separate PC events into samples with different energies, driven mainly by the requirementsof L/E analysis [144, 102]. Statistically speaking, parent neutrinos of PC events which
penetrate the OD (through-going) intend to have higher energy than ones that stop in theOD (stopping). This separation means better energy resolution for PC events.
In order to have a comparable MC simulation for the stopping and through-going
events, the OD MC performance needed to be improved. It is under this situation thatthe work described in this chapter was carried out. The separation of PC through-goingand stopping events also helps the zenith angle analysis give tighter constraints on masssquared difference comparing with the binning adopted in the previous atmospheric neutrino analysis of Super-K [63].
116

117

A.1 An overview of OD Monte Carlo simulation
As we described in Chapter 3, Super-K MC simulation consists of two parts: first neutrino
interactions are generated by NEUT then the final state particles are passed onto the detector
simulation apdetsim which simulates the detector responses. The signals are recorded
as time and charge information for every single PMT in the detector. Ideally, the simulated
signals should reproduce the observed signal. But before the tuning, shown in Fig. A.1 ,it is clear that the OD charge distribution of actual data spreads much wider than MC andfollows a different shape. The total number of hits of MC in OD are significantly lower
than observation. Data and MC used here are cosmic ray muons stop in ID. The criteria
selecting those events will be presented in the next section.
One of the reasons that MC doesn’t agree with observation we find is due to the Gaussian single photoelectron (SPE) distribution originally assumed for the OD PMTs. In
principle, the SPE distribution of a PMT with good signal-to-noise (S/N) level follows a
Poisson distribution and the Poisson distribution can be approximated by a Gaussian if themean value is large enough. The PMT used in the SK-I OD is Hamamatsu 8 inch R1408.
Figure A.2 shows the SPE pulse height distribution provided by Hamamatsu company.
Its peak-to-valley (P/V) ratio is not very good (note: the temperature of Super-K water isaround 12◦ C). Another factor worth pointing out is that the OD tubes are also very muchaged — they are recycled from IMB experiment. The actual SPE spectrum is even worse,
see Fig. A.3. It shows the charge distribution of OD PMT #1101 from dark noise, which
is a good approximation of the SPE pulse distribution. Obviously, a Gaussian distribution
is not a good model for our OD tubes any more.
Realizing this problematic modeling of the SPE distribution used in the MC simulation,
we decide to replace it with an exponential distribution:

118

Comparison between Data and MC
10 4ID
Entries
Mean
RMS

10 4

400
21029
0.52390.1237

ID
Entries
Mean
RMS

410
38798
1.092
0.8823

10 3
10 3
10 2
10 2

10

1

0

1

2
3
4
OD Charge Dist(t<0)
ID
Entries
Mean
RMS

10 4

5

420
80724
1009.
214.3

0

1
2
3
4
OD Charge Dist(900<t<1300)

50

ID
Entries
Mean
RMS

5

430
582
139.3
30.03

40
10 3

30
20

10

2

10

0

500

1000
1500
Time Dist

2000

0

0

100

200
300
400
Number of Hits

500

Figure A.1: OD MC performance before this tuning. Upper two: charge distributions
inside/outside of the time window and; lower two: hit time and number of hits distributions(dots are observation and solid histograms are MC)p(x) = e−x/g /g, x > 0 .

Where x is the pulse height and g is the centroid of the distribution. The value of thecentroid is proportional to the gain of the PMT. Figure A.4 shows the distribution of OD
PMT centroids obtained by fitting all their charge distributions outside the signal window.
It follows a nice Gaussian distribution with mean (g 0 ) of 0.5 and variance (σg ) of 0.135,
which means the gains of the OD PMTs are well tuned to be close to one another. The g 0

119Figure A.2: SPE pulse height distribution of Hamamatsu R1408

Figure A.3: Measured charge distribution of Hamamatsu R1408

120

Figure A.4: Distribution of centroids using exponential SPE distribution model (the conversion factor to number of photoelectrons is around 2 in apdetsim)
value for OD tubes is fixed at 0.5 but its variance σg is a parameter in our simulation.
An exponential distribution are generated by:

x = −g ln(1 − r).
Where r is a random number between 0 and 1. Due to plateau near 0 in Fig., wewill make a small modification: the SPE spectrum is designed as a joint distribution of an
uniform distribution and an exponential. The joint point position is related to the thresholdand it will be an parameter in the MC simulation and we call it plateau.

121

A.2 Tuning MC parameters
In order to tune the MC, we need some standard events which we can rely on. Cosmic raymuons stopping in ID are “standard rulers” because we can measure their energies very
well using their pathlengths. Following criteria are used to select stopping muons events:
1. Using muboy [74] muon fitter to identify stopping muons and further require tightcuts on goodness-of-fits:• G.O.F. > 0.8: G.O.F. is the goodness-of-fit value which describes how wellthe expected muon signal reproduce the observed one. Its detailed definition
can be found in Ref [74]• fcone > 0.8: fcone is the fraction of charges inside the reconstructed muon
Cherenkov cone in ID. Its definition can be found in Ref. [74]2. Pathlength in ID greater than 7m
3. Qtot (total photoelectron in ID) > 8000
4. dwall > 2m: stopping point from the nearest wall greater than 2m.
After reconstructing the stopping muon events, we simulate identical events with the sameentry points, directions, and stopping distances.
Now, we have a set of “standard rulers”. The next step is to choose the parameterswhich need to be tuned. There are many OD MC parameters and most of them are correlated. We decide to tune the following parameters:
1. plateau: the joint point of the flat distribution and the exponential distribution for
OD PMTs SPE spectrum, confined on [0, 0.5]

122
2. σg : it describes the width of the OD PMT gain distribution scaled to an unit Gaussian
distribution. We start with 0.27 which is the current value 1
3. σpedestal : the variance on the OD discriminator pedestal value. We assume
pedestal follows a Gaussian distribution. Setting the central value at zero for
MC is convenient but the width of the Gaussian depends on how well the discriminators are set to the same values.
4. threshold: the threshold value of OD discriminators, initially set at 0.3 p.e.
5. σsmear : noise smearing factor which is to simulate the noise between neighboring
channels. We assume it follow a Gaussian distribution of (1, σ smear ).
After some trials of comparing MC and observation, we find the following set of parameters with which MC simulation agrees with actual data to the best:
1. plateau=0.3
2. σg = 0.27
3. σpedestal = 0.12
4. threshold=0.55
5. σsmear = 0.5
With this set of parameters, together with the ones we keep the original values, the new
apdetsim simulation is shown in Fig. A.5 and Fig. A.6. The MC reproduces actual observation much better. The agreement on OD hits could be tuned better but that will make thecharge distribution agreement bad.
For the purpose of convenience, we use a unit Gaussian distribution with a scale factor 0.5 to model thecentroid distribution of OD PMTs.
1

123

Comparison between Data and MC(sim26)
ID
Entries
Mean
RMS
UDFLW
OVFLW

1000

1200

400
24042
0.83210.7155
0.
79.00

1000

800

410
38882
1.174
1.033
0.
4884.

800

600

600

400

400

200

200

0

ID
Entries
Mean
RMS
UDFLW
OVFLW

0

1

2
3
4
OD Charge Dist(t<0)
ID
Entries
Mean
RMS
UDFLW
OVFLW

10 3

5

400
24042
0.83210.7155
0.
79.00

0

0

10 3

1
2
3
4
OD Charge Dist(900<t<1300)
ID
Entries
Mean
RMS
UDFLW
OVFLW

5

410
38882
1.174
1.033
0.
4884.

10 2
10 2
10

1

0

1

2
3
4
OD Charge Dist(t<0)

5

0

1
2
3
4
OD Charge Dist(900<t<1300)

5Figure A.5: OD charge distribution inside/outside of time window with the tuned MC
simulation compared with observation

124

Comparison between Data and MC(sim26)
450
582
149.1
25.31
0.
0.

ID
Entries
Mean
RMS
UDFLW
OVFLW

50
40
30
20
10
0

50

75

100

125

ID
Entries
Mean
RMS
UDFLW
OVFLW

40
35

150
175
200
Number of Hits(all)
430
284
137.9
21.24
0.
0.

30

30

225

250

275

ID
Entries
Mean
RMS
UDFLW
OVFLW

25

300

440
298
159.8
24.26
0.
0.

20

25
15

20
15

10

10

5

5
0

50

100
150
200
250
Number of Hits(top entry)

300

0

50

100
150
200
250
Number of Hits(side entry)

300

Figure A.6: OD hit distribution with the tuned MC simulation compared with observation

Appendix B
PC Reduction Systematic Uncertainty
Systematic uncertainties of experiments are important in two respects: first, they are thequantitative evaluations of the quality of the experiment; second, they play essential rolesin the accuracy of the measured quantities. As an example, I will describe and summerize
how we evaluate the systematic uncertainty on the PC reduction efficiency (1PC ).
PC reduction is a set of cuts in a computer program used to separate the partiallycontained neutrino events from the rest. Conceptually, PC reduction efficiency is defined
as follows:

1P C =

Nkept.
NtotalWhere Nkept is the number of PC events after PC reduction and Ntotal is the total numberof PC events in a Monte Carlo sample that should survive the reduction, namely, neutrino
interactions with particles that exit into the OD.
In this chapter, I will describe part of the work I did for Super-K I in detail: how
we evaluate the systematic uncertainty on reduction steps 1 (PC1) through 4 (PC4), i.e.
125

126
δ+PC1−4 . The total systematic uncertainty (δ +PC ) on the PC efficiency is the quadratic sumof the uncertainties on all the reduction steps, PC1-5:

δ+PC =,
δ+2PC1−4 + δ+2PC5 .

(B.1)

B.1 Overview
PC reduction consists of 5 steps for Super-K I. The time and space information of hits in
ID and OD decide whether an event is PC or not. We record the cuts for the first 4 steps,
PC1-4, briefly here:
PC1:• twidth : the distribution the time between two consequential hits
• nOD−cluster : number of OD hit clusters
PC2:• n2nd−cluster : number of hits in the 2nd (ordered by number of hits in the cluster) OD
cluster• nmin : number of hits in the first cluster or the rest of OD hits depending on which
one is smaller• nouter : number of hits in the most charged OD cluster• npe200 :number of photo-electrons within 2 meters in the corresponding ID area ofthe most charged OD cluster

127
PC3:
• nehit8m : number of hits in OD within 8 meters of the back extrapolated entrance
point on OD wallPC4:
• cosθfit−first : θfit−first is the angle between the fitted direction and the line connecting
the vertex and the first hit tube in the ID cluster• mugood: the goodness value being a through-going muon
• mudistance: the distance it should travel in ID assuming as through-going muon
• dcorner :the distance from the fit vertex to the nearest ID corner
These cuts can be grouped into two categories: ID cuts and OD cuts. For ID cuts, weevaluate the uncertainties one by one. The method for each cut is different and details
will be described later. Generally speaking, the common part is to check the cut value
distribution to see how much uncertainty could have been introduced due to non-perfectagreement between MC and data. Data used are the PC events after the complete PC
reduction. There are totally 882 PC events in the sample we use here. For OD cuts, wechoose to shift the MC parameters to certain extremes and then exam the effect on the PC
reduction. The change on the number of kept PC events is considered as the uncertainty
due to OD-related cuts. Then, we have:δ+PC1−4

6
7 n
7/
=8
δ2

2
cut,i + δOD

.

(B.2)

i=1

Where δcut,i is the uncertainty on each ID cut and δOD is the uncertainty on all the OD cuts.

128

Figure B.1: Distributions of npe200 of the final PC sample. The circles are the final PC
sample; triangle are MC

B.2 Evaluation of ID cut uncertainties
There are 5 ID cuts: npe200 , cos θfit−first , mugood, mudistance, and dcorner .
Cut npe200 is coupled with the OD cut nouter in the following way: if npe200 ! 1000 ornouter " 6, then the event is considered as PC. Figure B.1 shows the distributions of final
PC samples of both data and MC.
Observing the distribution, we assume that below 2000 hits, n pe200 follows an exponential distribution. So fitting both MC and observed data below 2000 using an exponentialfunction and checking the difference of those two integrals on [0, 1000], we conclude thatthe systematic uncertainty of cut npe200 ! 1000 is 1%.
Figure B.2 shows the distributions of cut cos θ fit−first ! −0.8. Assuming flat distributions below −0.5 for both MC and observed data, the difference between them in therejected region, i.e. [−1, −0.8], is 5, thus the systematic uncertainty is 5/(5+882) is 0.6%.

129

Figure B.2: cos θfit−first distributions of MC and observed data
Cuts mugood and mudistance are coupled together: if mugood< 0.85 or mudistance<
3000, events are kept as PC events. To evaluate the uncertainty, we examine the changeof number of events in the region defined by mugood! 0.85 and mudistance! 3000 when
we adjust the cut values to compensate for the difference between MC and observed data.
Figure B.3 shows the distribution of mugood with cut mudistance! 3000 of the final PC
data samples.
The average values of the cut mugood are 0.63 and 0.65 for observed data and MC
respectively. The difference is 3%. If we adjust the cut according to this difference, thecut value will shift from 0.85 to 0.82. The shift on the cut value will change the number of
events in the final observed PC sample by 7. We thus assign the systematic uncertainty ofthis combined cut of mugood and mudistance 7/882 = 0.8%.
Now, let us examine the cut dcorner . Figure B.4 shows the distributions of the final
PC samples of both MC and observation. The cut is d corner > 150cm. We assume the130

Number of Events

Figure B.3: Distribution of mugood with mudistance! 3000

50
45
40
35
30
25
20
15
10
5
0

0

250

500

750

1000

1250

1500

1750

2000

2250

2500

Dist. to the Nearst Corner

Figure B.4: Distribution of dcorner

131
distribution is linear over the range [0, 500] for both MC and observation. Fitting the distributions from 150cm to 500 cm using and extrapolating to rejected region, the differenceof MC and observation over the range [0, 150] is 5. Thus, the systematic uncertainty is
5/(882 + 5) = 0.6%.

B.3 The Uncertainty of OD cuts
We adopt a different approach for the evaluation of the OD cut uncertainty. The argument
is as follow: cut values are decided mainly based on MC simulations, so the uncertainties
are introduced in due to the imperfectness of the MC. If we explore the imperfectness of
the MC and the change produced by such imperfectness, we can evaluate the uncertainty
due to the MC. This uncertainty is the systematic uncertainty of the OD cuts.
Recall the chapter on OD MC simulation, see Fig. , we were aware of the fact that the
charge and the hit distributions are not simulated well simultaneously. We decided to pay
more attention to the reproduction of the charge distribution chose the parameters corresponding to the better agreement on charge distributions. However, choosing to reproduce
the hit distribution better is an equally good set of MC parameters. This uncertainty canbe considered as the OD cut uncertainty.
The parameters we choose to tune are:
• Collection efficiency of OD PMTs: 0.22(0.005)• Threshold of OD discriminator: 0.55(0.01)• Reflectivity of OD Tyvek: 0.80(0.01)
The values in parentheses are the steps each time we choose a new set of parameters.

132

Figure B.5: Distributions produced by another set of MC parameters
After trying several trials, we achieved the distributions shown in Fig. B.5. This set of
MC parameters reproduces the OD hit distribution better comparing with the official set ofparameters which reproduces the OD charge distribution better.
Using this set of parameters, we regenerate 10 years of MC events and let them gothrough PC reduction steps 1 through 4. The PC efficiency based on this new set of parameters is 85.9%. Comparing with the efficiency using the official parameter in the OD
tuning chapter, 87.3%, we assign 1.4% as the systematic uncertainty produced by all the
OD cuts from reduction steps from 1 to 4.

B.4 Summary
Systematic uncertainties are important knowledge of experiment systems. Though its concept is simple, there is no standard ways evaluating them. Using different approaches, weevaluated the systematic uncertainties of different PC reduction cuts for PC1 through 4.
According to Eq. B.2, the uncertainty introduced into the PC reduction efficiency 1 PC1−4
during steps PC1 through 4 is:δ+PC1−4

6
7 n
7/
=8
δ2

2
cut,i + δOD = 1.5%

i=1

.

133
The cuts in PC5 are all ID-related and following similar methods described in Sec. B.2,
the uncertainty is found to be 2.4%.
In summary, according to Eq. B.1, the systematic uncertainty of PC reduction efficiency is,

δ+PC =,

δ+2PC1−4 + δ+2PC5 = 3.2% .

Bibliography[1] T. K. Gaisser and M. Honda. Flux of atmospheric neutrinos. Annual Review of
Nuclear and Particle Science, 52:153–199, 2002.
[2] Moscow Neutron Monitor. http://helios.izmiran.troitsk.ru/cosray/main.htm.
[3] M. Honda, T. Kajita, K. Kasahara, and S. Midorikawa. A new calculation of the atmospheric neutrino flux in a 3-dimensional scheme. Physical Review D, 70:043008,
2004.
[4] G. D. Barr, T. K. Gaisser, P. Lipari, S. Robbins, and T. Stanev. A three-dimensional
calculation of atmospheric neutrinos. Physical Review D, 70: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.).
[5] G. Battistoni et al. High energy extension of the FLUKA atmospheric neutrino flux.
2003. (The discontinuity seen in Fig. 2.4 has been addressed by later calculations.
See: http://www.mi.infn.it/ battist/neutrino.html.).
[6] M. H. Ahn et al. Indications of neutrino oscillation in a 250-km long- baseline
experiment. Physical Review Letters, 90:041801, 2003.
[7] K. Winter, editor. Neutrino Physics. Cambridge University Press, 1991.

134

135
[8] F. Reines and C. Cowan. The neutrino. Nature, 178:446, 1956.
[9] M. Goldharber, L. Grodzins, and A. W. Sunyar. Helicity of neutrinos. PhysicalReview, 109:1015, 1958.
[10] G. Danby et al. Observation of high-energy neutrino reactions and the existence of
two kinds of neutrinos. Physical Review Letters, page 36, 1962.
[11] M. L. Perl et al. Evidence for anomalous lepton production in e+-e- annihilation.
Physical Review Letters, page 1489, 1975.
[12] K. Kodama et al. Observation of tau neutrino interactions. Physics Letters B, page
218, 2001.
[13] The LEP Collaborations ALEPH DELPHI L3 OPAL the LEP Electroweak Working Group the SLD Heavy Flavour and Electroweak Groups. A combination of preliminary electroweak measurements and constraints on the standard model. CERNEP/99-15, 1999.
[14] C. Weinheimer. The neutrino mass direct measurements, 2003. The 10th International Workshop on Neutrino Telescopes, Venice, Italy.
[15] Y. Fukuda et al. Evidence for oscillation of atmospheric neutrinos. Physical Review
Letters, 81:1562–1567, 1998.
[16] KATRIN collaboration. KATRIN: A next generation tritium beta decay experimentwith sub-eV sensitivity for the electron neutrino mass, 2001.
[17] S. Wang et al. Weighing neutrinos with galaxy cluster surveys. Physical Review
Letters, 95:011302, 2005.
[18] A. D. Dolgov. Neutrinos in cosmology. Physics Reports, 370:333, 2002.

136[19] M. Tegmark et al. Cosmological parameters from sdss and wmap. Physical Review
D, 69:103501, 2004.
[20] U. Seljak et al. Cosmological parameter analysis including sdss ly-alpha forest andgalaxy bias: constraints on the primordial spectrum of fluctuations, neutrino mass,and dark energy. Physical Review D, 71:103515, 2005.
[21] H. Murayama. The origin of neutrino mass. Physics World, May, 2002.
[22] R. N. Mohapatra and A. Y. Smirnov. Neutrino mass and new physics. Annual
Review of Nuclear and Particle Science, 56:569, 2006.
[23] S. L. Glashow. Quarks and Leptons, Cargése 1979, page 707, 1980.
[24] Y. Yanagida. Proceeding of the Workshop on Unified Theories, Tsukuba, 1979.
[25] T. Yanagida. Horizontal symmetry and masses of neutrinos. Progress of Theoretical
Physics, 64:1103, 1980.
[26] R. N. Mohapatra and G. Senjanović. Neutrino mass and spontaneous parity nonconservation. Physical Review Letters, 44(14):912–915, Apr 1980.
[27] M. Gell-Mann, P. Ramond, and R. Slansky. Supergravity edited by D. Freedman etal. North Holland, Amsterdam, 1979.
[28] R. D. Peccei. Neutrino physics. American Institute of Physics Conference Proceeding, 490:80, 1999.
[29] B. Pontecorvo. Soviet JETP, 7:172, 1958.
[30] Z. Maki, M. Nakagawa, and S. Sakata. Remarks on the unified model of elementary
particles. Progress of Theoretical Physics, page 870, 1962.

137[31] V. N. Gribov and B. Pontecorvo. Neutrino astronomy and lepton charge. PhysicalLetters B, 28:493, 1969.
[32] J. C. Evans, R. Davis, and J. N. Bahcall. Brookhaven solar neutrino detector andcollapsing stars. Nature, 251:486–488, 1974.
[33] R. Davis. Nobel lecture: A half-century with solar neutrinos. Reviews of Modern
Physics, 75:985–994, 2003.
[34] D. G. Michael et al. Observation of muon neutrino disappearance with the minos
detectors and the numi neutrino beam. Physical Review Letters, 97:191801, 2006.
[35] Q. R. Ahmad et al. Direct evidence for neutrino flavor transformation from neutralcurrent interactions in the Sudbury neutrino observatory. Physical Review Letters,
89:011301, 2002.
[36] K. Eguchi et al. First results from KamLAND: Evidence for reactor anti- neutrino
disappearance. Physical Review Letters, 90:021802, 2003.
[37] M. Apollonio et al. Limits on neutrino oscillations from the CHOOZ experiment.
Physics Letters B, 466:415–430, 1999.
[38] M. Apollonio et al. Search for neutrino oscillations on a long base-line at the
CHOOZ nuclear power station. European Physics Journal C, 27:331–374, 2003.
[39] A. Piepke et al. Final results from the Palo Verde neutrino oscillation experiment.
Progress in Partical and Nuclear Physics, 48:113–121, 2002.
[40] L. Wolfenstein. Neutrino oscillations in matter. Physical Review D, 17(9):2369–
2374, May 1978.

138[41] L. Wolfenstein. Neutrino oscillations and stellar collapse. Physical Review D,
20:2634–2635, 1979.
[42] V. D. Barger, K. Whisnant, S. Pakvasa, and R. J. N. Phillips. Matter effects on
three-neutrino oscillations. Physical Review D, 22:2718, 1980.
[43] F. R. Klinkhamer. Neutrino oscillations from the splitting of Fermi points. JETP
Letters, 79:451–456, 2004.
[44] S. Coleman and S. L. Glashow. Cosmic ray and neutrino tests of special relativity.
Physics Letters B, 405:249, 1997.
[45] V. A. Kostelecky. Neutrino oscillations and Lorentz violation. 2004. Invited talk at
the Third International Symposium on Quantum Theory and Symmetries.
[46] E. Lisi, A. Marrone, and D. Montanino. Probing possible decoherence effects in
atmospheric neurino oscillations. Physical Review Letters, 85:1166, 2000.
[47] T. Ohlsson. Equivalence between neutrino oscillations and neutrino decoherence.
Physics Letters B, 502:159–166, 2001.
[48] Y. Grossman and M. P. Worah. Atmospheric νµ deficit from decoherence. 1998.
[49] C. P. Sun and D. L. Zhou. Quantum decoherence effect and neutrino oscillation.
1998.
[50] Y. Liu, L.Z. Hu, and M.L. Ge. Effect of violation of quantum mechanics on neutrino
oscillation. Physical Review D, 56(10):6648–6652, Nov 1997.
[51] V. D. Barger et al. Neutrino decay and atmospheric neutrinos. Physics Letters B,
462:109–114, 1999.

139[52] T. K. Gaisser. Cosmic rays and particle physics. Cambridge University Press, 1990.
[53] Y. Shikaze et al. Solar modulation effect on the cosmic-ray proton spectra measured
by BESS. 2003. Proceeding of the 28th International Cosmic Ray Conferences(ICRC 2003), Tsukuba, Japan, 31 Jul - 7 Aug 2003, Vol.7, p.4027.
[54] Climax and Haleakala neutron monitor. http://ulysses.sr.unh.edu/NeutronMonitor/.
[55] J. Alcaraz et al. Cosmic protons. Physics Letters, B490:27–35, 2000.
[56] T. Sanuki et al. Precise measurement of cosmic-ray proton and helium spectra with
the BESS spectrometer. Astrophysical Journal, 545:1135, 2000.
[57] G. Battistoni, A. Ferrari, T. Montaruli, and P. R. Sala. The FLUKA atmospheric
neutrino flux calculation. Astroparticle Physics, 19:269–290, 2003.
[58] V. Agrawal, T. K. Gaisser, P. Lipari, and T. Stanev. Atmospheric neutrino flux above
1 GeV. Physical Review D, 53:1314–1323, 1996.
[59] M. Fukugita and A. Suzuki. Physics and Astrophysics of Neutrinos. SpringerVerlag, 1994.
[60] M. Fukugita and A. Suzuki. Physics and Astrophysics of Neutrinos. SpringerVerlag, 1994.
[61] K. Kobayashi et al. Search for nucleon decay via modes favored by supersymmetric
grand unification models in super- kamiokande-i. Physical Review D, 72:052007,
2005.
[62] J. Hosaka et al. Solar neutrino measurements in super-kamiokande-i. PhysicalReview D, 73:112001, 2006.

140[63] Y. Ashie et al. A measurement of atmospheric neutrino oscillation parameters bysuper-kamiokande i. Physical Review, D71:112005, 2005.
[64] D. W. Liu et al. Limit on the neutrino magnetic moment using 1496 days of SuperKamiokande-I solar neutrino data. Physical Review Letters, 93:021802, 2004.
[65] Y. Takenaga et al. Search for neutral Q-balls in Super-Kamiokande II. 2006. to be
published on Physical Letters B.
[66] T. Futagami et al. Observation of the east-west anisotropy of the atmospheric neutrino flux. Physical Review Letters, 82:5194–5197, 1999.
[67] S. Fukuda et al.

Search for neutrinos from gamma-ray bursts using Super-

Kamiokande. Astrophysical Journal, 578:317–324, 2002.
[68] K. Abe et al. High energy neutrino astronomy using upward-going muons in SuperKamiokande-I. The Astrophysical Journal, 652:198, 2006.
[69] M. E. C. Swanson et al. Search for diffuse astrophysical neutrino flux using ultrahigh energy upward-going muons in Super-Kamiokande I. The Astrophysical Journal, 652:206, 2006.
[70] M. Malek et al. Search for supernova relic neutrinos at Super-Kamiokande. PhysicalReview Letters, 90:061101, 2003.
[71] S. Desai et al. Search for dark matter WIMPs using upward through-going muonsin Super-Kamiokande. Physical Review D, 70:083523, 2004.
[72] M. H. Ahn et al. Measurement of neutrino oscillation by the K2K experiment.
Physical Review D, 74:072003, 2006.
[73] Y. Itow et al. The JHF-Kamioka neutrino project. 2001. hep-ex/0106019.

141[74] Shantanu Desai. High energy neutrino astrophysics with Super-Kamiokande, 2004.
Boston University doctoral thesis.
[75] Y. Fukuda et al. The super-kamiokande detector. Nuclear Instruments and Methods
A, 501:418–462, 2003.
[76] A. Suzuki et al. Improvement of 20-inch diameter photomultiplier tubes. Nuclear
Instruments and Method A, 329:299–313, 1993.
[77] R. Claus et al. A wave shifter light collector for a water cherenkov detector. Nuclear
Instruments and Method A, 261:540–542, 1987.
[78] T. Tanimori et al. Design and performance of semicustom analog ic including twotacs and two current integrators for ’super- kamiokande’. IEEE Trans. Nucl. Sci.,
36:497–501, 1989.
[79] Hirokazu Ikeda et al. Front end hybrid circuit for super-kamiokande. Nuclear
Instruments and Methods A, 320:310–316, 1992.
[80] M. Nakahata et al. Calibration of Super-Kamiokande using an electron linac. Nuclear Instruments and Methods A, 421:113–129, 1999.
[81] E. Blaufuss et al.

1

6n as a calibration source for Super-Kamiokande. Nuclear

Instruments and Methods A, 458:638–649, 2001.
[82] Y. Hayato. NEUT. Nuclear Physics Proceeding Supplements, 112:171–176, 2002.
[83] D. Casper. The nuance neutrino physics simulation, and the future. Nuclear Physics
Proceeding Supplements, 112:161–170, 2002.
[84] C. H. Llewellyn Smith. Neutrino reactions at accelerator energies. Physics Reports,
3:261, 1972.

142[85] R. A. Smith and E. J. Moniz. Neutrino reactions on nuclear targets. Nuclear Physics,
B43:605, 1972.
[86] D. Rein and L. M. Sehgal. Neutrino excitation of baryon resonances and single pion
production. Annual Physics, 133:79, 1981.
[87] D. Rein. Angular distribution in neutrino induced single pion production processes.
Zeitschrift fur Physik, C35:43–64, 1987.
[88] T. Kitagaki et al. Charged current exclusive pion production in neutrino deuterium
interactions. Physical Review D, 34:2554–2565, 1986.
[89] S. K. Singh, M. J. Vicente-Vacas, and E. Oset. Nuclear effects in neutrinoproductionof delta at intermediate energies. Physics Letters B, 416:23–28, 1998.
[90] D. Rein and L. M. Sehgal. Coherent π 0 production in neutrino reactions. Nuclear
Physics, B223:29, 1983.
[91] M. Gluck, E. Reya, and A. Vogt. Dynamical parton distributions of the proton andsmall x physics. Zeitschrift fur Physik C, 67:433–448, 1995.
[92] Masayuki Nakahata et al. Atmospheric neutrino background and pion nuclear effect for kamioka nucleon decay experiment. Journal of Physical Society of Japan,
55:3786, 1986.
[93] Torbjorn Sjostrand. High-energy physics event generation with PYTHIA 5.7 andJETSET 7.4. Computer Physics Communications, 82:74–90, 1994.
[94] R.D. Woods and D.S. Saxon. Physical Review, 95:577–578, 1954.
[95] L. L. Salcedo, E. Oset, M. J. Vicente-Vacas, and C. Garcia-Recio. Computer simulation of inclusive pion nuclear reactions. Nuclear Physics A, 484:557, 1988.

143[96] Glenn Rowe, Martin Salomon, and Rubin H. Landau. An energy dependent phaseshift analysis of pion - nucleon scattering below 400-mev. Physical Review C,
18:584–589, 1978.
[97] D. Ashery et al. True absorption and scattering of pions on nuclei. Physical Review
C, 23:2173–2185, 1981.
[98] R. Brun and F. Carminati. GEANT detector description and simulation tool. CERN
Programming Library Long Writeup W5013, 1993.
[99] C. Zeitnitz and T. A. Gabriel. The GEANT - CALOR interface and benchmark calculations of ZEUS test calorimeters. Nuclear Instruments and Methods A, 349:106–
111, 1994.
[100] E. Bracci.

Compilation of cross sections : II k − and k + induced reactions.

CERN/HERA, 72-1, 1972.
[101] A. S. Carroll et al. Pion - nucleus total cross-sections in the (3,3) resonance region.
Physical Review C, 14:635–638, 1976.
[102] Masaki Ishitsuka.

L/E analysis of the atmospheric neutrino data from super-

kamiokande, 2004. University of Tokyo doctoral thesis.
[103] Scott Clark. Searches for proton decay with Super-Kamiokande detector, 2006.
Boston University doctoral thesis.
[104] A. Aguilar et al. Evidence for neutrino oscillations from the observation of ν̄ e appearance in a ν̄µ beam. Physical Review D., 64:112007, 2001.
[105] J. Hosaka et al. Three flavor neutrino oscillation analysis of atmospheric neutrinos
in super-kamiokande. Physical Review D, 74:032002, 2006.

144[106] K. Abe et al. A measurement of atmospheric neutrino flux consistent with tau
neutrino appearance. Physical Review Letters, 97:171801, 2006.
[107] G. L. Fogli, E. Lisi, A. Marrone, D. Montanino, and A. Palazzo. Getting themost from the statistical analysis of solar neutrino oscillations. Physical Review
D, 66:053010, 2002.
[108] S. Fukuda et al. Determination of solar neutrino oscillation parameters using 1496
days of Super-Kamiokande-I data. Physics Letters B, 539:179–187, 2002.
[109] R. N. Mohapatra, S. Nasri, and H.-B. Yu. Seesaw right handed neutrino as the sterile
neutrino for LSND. Physical Review D, 72:033007, 2005.
[110] G. C. McLaughlin et al. An active-sterile neutrino transformation solution for rprocess nucleosynthesis. Physical Review C, 59:2873, 1999.
[111] D. O. Caldwell, G. M. Fuller, and Y. Z. Qian. Sterile neutrinos and supernova
nucleosynthesis. Physical Review D, 61:123005, 2000.
[112] S. Fukuda et al. Tau neutrinos favored over sterile neutrinos in atmospheric muon
neutrino oscillations. Physical Review Letters, 85:3999–4003, 2000.
[113] S. Nakayama et al. Measurement of single π 0 production in neutral current neutrino
interactions with water by a 1.3-gev wide band muon neutrino beam. Physics Letters
B, 619:255–262, 2005.
[114] A.M. Dziewonski and D.L. Anderson. Preliminary reference earth model. Physics
Earth and Planet. Inter., 25:297–356, 1981.
[115] G. L. Fogli, E. Lisi, and A. Marrone. Four-neutrino oscillation solutions of theatmospheric neutrino anomaly. Physical Review D, 63(5):053008, Feb 2001.

145[116] Irina Mocioiu, Maxim Pospelov, and Radu Roiban. Breaking CPT by mixed noncommutativity. Physical Review D, 65:107702, 2002.
[117] J. D. Bjorken and S. D. Drell. Relativistic Quantum Fields. New York, McGrawHill, 1965.
[118] R. P. Feynman and S. Weinberg. Elementary particles and the laws of physics.
Cambridge University Press, 1987.
[119] M. E. Peskin and D. V. Schroeder. An Introduction to Quantum Field Theory. Westview Press, 1995.
[120] A. Zee. Quantum Field Theory in a Nutshell. Princeton University Press, 2003.
[121] S. Coleman and S. L. Glashow. High-energy tests of lorentz invariance. PhysicalReview D, 59:116008, 1999.
[122] V. A. Kostelecky and M. Mewes. Lorentz and CPT violation in neutrinos. PhysicalReview D, 69:016005, 2004.
[123] V. A. Kostelecky and M. Mewes. Lorentz and CPT violation in the neutrino sector.
Physical Review D, 70:031902, 2004.
[124] D. Colladay and V. A. Kostelecky. Lorentz-violating extension of the standard
model. Physical Review D, page 116002, 1998.
[125] Y. Grossman, C. Kilic, J. Thaler, and D. G. E. Walker. Neutrino constraints on
spontaneous Lorentz violation. Physical Review D, 72:125001, 2005.
[126] F. R. Klinkhamer and G. E. Volovik. Emergent CPT violation from the splitting of
Fermi points. International Journal of Modern Physics A, 20:2795, 2005.

146[127] O. W. Greenberg. CPT violation implies violation of Lorentz invariance. PhysicalReview Letters, 89:231602, 2002.
[128] S. L. Glashow, A. Halprin, P. I. Krastev, C. N. Leung, and J. Pantaleone. Commentson neutrino tests of special relativity. Physical Review D, 56:2433, 1997.
[129] D. Hooper, D. Morgan, and E. Winstanley. Lorentz and CPT invariance violation in
high-energy nuetrinos. Physical Review D, 72:065009, 2005.
[130] G. Barenboim, L. Borissov, and J. Lykken. Neutrinos that violate CPT, and theexperiments that love them. Physics Letters B, 534:106–113, 2002.
[131] M. D. Messier. Testing CPT conservation using atmospheric neutrinos. The 3rd
Meeting on CPT and Lorentz Symmetry (CPT 04), Bloomington, Indiana, 2004.
[132] V. D. Barger, S. Pakvasa, T. J. Weiler, and K. Whisnant. CPT-odd resonances inneutrino oscillations. Physical Review Letters, 85:5055, 2000.
[133] L. B. Auerbach et al. Tests of Lorentz violation in muon antineutrino to electron
antineutrino oscillations. Physical Review D, 72:076004, 2005.
[134] G. L. Fogli, E. Lisi, A. Marrone, and G. Scioscia. Testing violations of special
and general relativity through the energy dependence of ν µ ↔ ν τ oscillationsin the Super-Kamiokande atmospheric neutrino experiment. Physical Review D,
60(5):053006, Jul 1999.
[135] M. C. Gonzalez-Garcia and Michele Maltoni. Atmospheric neutrino oscillations
and new physics. Physical Review D, 70:033010, 2004.
[136] O. Gagnon and G. D. Moore. Limits on Lorentz violation from the highest energy
cosmic rays. Physical Review D, 70:065002, 2004.

147[137] S. K. Lamoreaux et al. New limits on spatial anisotropy from optically-pumped
201

Hg and 199 Hg. Physical Review Letters, 57(25):3125–3128, Dec 1986.

[138] W-M Yao et al. Review of Particle Physics, volume 33. Institute of Physics, 2006.
[139] S. W. Hawking. The unpredictability of quantum gravity. Communications in Mathematical Physics, 87:395 – 415, 1982.
[140] V. D. Barger, W. Y. Keung, and S. Pakvasa. Majoron emission by neutrinos. Physical Review D, 25(3):907–910, Feb 1982.
[141] V. D. Barger, J. G. Learned, S. Pakvasa, and T. J. Weiler. Neutrino decay as an explanation of atmospheric neutrino observations. Physical Review Letters, 82(13):2640–
2643, Mar 1999.
[142] G. L. Fogli, E. Lisi, A. Marrone, and D. Montanino. Status of atmospheric ν µ → ντoscillations and decoherence after the first K2K spectral data. Physical Review D,
67:093006, 2003.
[143] Super-K Collaboration. Studies of exotic models using L/E analysis. Internal reports and research in progress.
[144] Y. Ashie et al. Evidence for an oscillatory signature in atmospheric neutrino oscillation. Physical Review Letters, 93:101801, 2004.