Is it possible to detect the η' → e+e- decay?: A simulation of the η' decay from e+e- collisions

Is it possible to detect the decay?

A simulation of the decay from collisions

Daniel Hamnevik

Supervisor: Andrzej Kupsc

Uppsala University, Department of Physics and Astronomy, Nuclear Physics

Subject reader: Tord Johnsson

Uppsala University, Department of Physics and Astronomy, Nuclear Physics

Examensarbete 15 hp

Juni 2014

Abstract

Decays of light, pseudoscalar mesons into lepton pairs are one possible way of detecting physics beyond the Standard Model. This thesis will focus on the η

→ e

e

decay, because only few experiments have studied this decay.

By the use of a basic simulation, the feasibility of conducting a search for the η

→ e

e

decay at DAΦNE, a electron-positron collider in Frascati, Italy, was investigated. The simulation implements an idea from previous experiments, where the η

decay was investigated through the e

e

→ η

→ X process.

The results of the simulation show that experiments at DAΦNE could produce an observable signal of the process when the produced η

decays into one of the most probable decay modes. However, in order to draw a definite conclusion, a more detailed study is needed.

Sammanfattning

S¨ onderfallet av l¨ atta, pseudoskal¨ ara mesoner till ett par av leptoner

¨

ar ett m¨ ojligt s¨ att att detektera fysik bortom Standardmodellen. Detta examensarbete fokuserar p˚ a η

→ e

e

s¨ onderfallet, p˚ a grund av att f˚ a experiment har studerat detta s¨ onderfall.

Genom att anv¨ anda en enkel simulering, unders¨ oktes m¨ ojligheten att studera η

→ e

e

s¨ onderfallet vid DAΦNE, en elektron-positron accelerator i Frascati, Italien. Simuleringen f¨ oljde samma genomf¨ orande som tidigare experiment, d¨ ar η

s¨ onderfallet unders¨ oktes genom e

e

→ η

→ X processen.

Resultaten av simuleringen visar p˚ a att experiment vid DAΦNE

kan producera en detekterbar signal av processen d˚ a de producer-

ade η

s¨ onderfaller till de mest troliga s¨ onderfallen. P˚ a grund av

begr¨ ansningarna i den utf¨ orda simulationen, beh¨ over dock en mer

detaljerad studie g¨ oras f¨ or att dra definitiva slutsatser.

Contents

1 Introduction 1

2 Theory of the η

decay 2

3 The Simulation and it’s Limits 5

4 Pictures from the simulation 7

5 Results of the Simulation 19

6 Discussion and Conclusion 21

References 22

Appendix A : Code used in Simulation 24

1 Introduction

The Standard Model (SM) of particle physics is today the most tested and verified explanation of how matter in our universe is made and how this matter interacts, however it fails to explain everything. For example, the hypothesized dark matter, which has not been detected beyond its interaction by gravity, and the asymmetry between the amount of matter and antimatter cannot be explained within the SM. Since observation has been made of physics unexplained by the SM, it stands to reason to look for physics beyond it.

The SM describes the fundamental matter as different particles, together with three forces (Gravity is not included in the SM) with corresponding force carriers [1]. However, the fact that the SM is insufficient for explaining everything means that one should look for new undiscovered particles. One way is to look for experimental data which differ from the theoretical values calculated from the SM. Large differences might point to the influence of exotic particle interactions, not included in the SM. The problem is that no such discrepancies have, with a high degree of certainty, been detected.

There are indications for physics beyond the SM. One experiment is the π

→ e

e

decay by the KTeV collaboration in 2007 [2] which showed a possible existence of a process beyond the SM. The experiment gave a result 3.3 standard deviations higher than the calculations using known particles from the SM. This result points to the decays of a light, pseudoscalar meson to a lepton pair as a process where physics beyond the SM is at work.

Due to this result, further studies and experiments are needed on light, pseudoscaler mesons. Especially on the decays where no actual event have been reported, e.g. the η

→ e

e

decay. Previous decay experiments on π

and η mesons have been carried out, among others, at WASA-at-COSY [3], Celsius/WASA [4] and at HADES [5]. However, decay experiments are not effective in the search for the lepton pair decay of the η

, due to the low branching ratio of the decay. Instead formation experiments are used where the η

mesons are created by collisions of e

e

leptons and then decays [6, 7].

This thesis will focus on the η

→ e

e

decay. This is motivated by

the fact that the only measurements of the decay comes from a formation

experiment at VEPP-2M [6] and at VEPP-2000 [7], where no η

was detected

and the experimental limits are much higher than the predicted values from

the SM. I will present the theory I used of the decays and try to argue about

the needed experimental equipments to make feasible experiments.

2 Theory of the η ⁰ decay

Decays of neutral pseudoscalar mesons into a lepton-antilepton pair, P →l

l

, represent a potentially important process to look for effects of physics beyond the SM [8, 9]. The lowest order mechanism for this decay within the SM is a second order, electromagnetic, one-loop process involving two virtual photons P →γ

γ

(see Figure 1).

Figure 1: Dominating SM mechanism for pseudoscalar meson decay into lepton- antilepton pai.

Due to the loop appearing in the process the decay is sensitive to values of the meson transition form factors

(represented by the circle in the figure) for any q

, where q is the momentum of a photon in the loop. For lighter mesons, the imaginary part of the decay amplitude can be uniquely related to the decay width of the P →γγ decays leading to a lower Unitary Bound [10]; the limit where the two photons are real.

For heavier mesons η

and η

there are more intermediate states, like ργ,

Figure 2: Feynman diagram of the formation of the η

meson and the subsequent decay into the most probable decay products.

The resonance cross section for the peak of the e

e

→ η

process with η

decaying into some final state out is given by [11]:

σ(E) = 2J + 1 (2S

+ 1)(2S

+ 1)

4π k

"

Γ

/4 (E − E

)

+ Γ

/4

#

B

B

(1) where E is the center-of-mass (c.m.) energy, J is the spin of the resonance, 2S

+ 1 and 2S

+ 1 are the number of polarization states of the two incident particles. The c.m. momentum in the initial state is k, E

is the c.m. energy at the resonance, and Γ is the full width at half maximum height of the resonance [11]. B

is the branching ratio of η

→ e

e

and B

is the branching fraction for η

decaying into some final state. Eq. (1) can be reduced to:

σ

= 4π m

B

B

(2)

if S

= S

= S(e) = 1/2, J (η

) = 0, E ≈ E

= m

and k = m

/2.

Here

= 13.7 (GeV)

= 5.3 mb, by the conversion factor (¯ hc)

= 0.389 (GeV)

mb [12].

From Eq. (2) we have the number of particles created as:

N = Lσ

= L 4π m

B

B

(3)

where  is the acceptance for the final state and L is luminosity.  will be

assumed to be  = 1 at the start of the simulation and will be calculated

from the simulated results.

By using Eq. (3), it is simple to roughly calculate the expected number of η

mesons. By adding different B

for a specific decay branch we get the expected number of particles for that decay.

I limited myself to the decays shown in the first column of Table 1 to reduce the work, but it still accounts for approximately 60% of the total number of η

mesons created. Furthermore, those decay products are relative easy to detect, making them interesting to investigate.

Branching ratio B

B(η

→ π

π

η) 42.9(±0.7)%

B(η

→ π

π

η)B(η → γγ) 16.9(±0.7)%

B(η

→ π

π

η)B(η → π

π

π

) 9.8(±0.8)%

B(η

→ π

π

η)B(η → π

π

π

) 14.0(±0.7)%

B(η

→ π

π

η) 22.2(±0.8)%

B(η

→ π

π

η)B(η → γγ) 8.7(±0.8)%

B(η

→ π

π

η)B(η → π

π

π

) 5.1(±0.8)%

B(η

→ π

π

η)B(η → π

π

π

) 7.3(±0.8)%

Table 1: Branching ratios used as B

, according to P.D.G [12].

Two different values are used for B

in this thesis, see Table 2. The experimental limit, or Upper Limit, is from the latest experiment done [7]

and the theoretical value is from SM calculations [13].

Branching ratio B

3 The Simulation and it’s Limits

This work focus on simulating the formation and decay of η

. The theory of the decay, from the SM, is well understood but we do not have sufficient experimental data to test the predictions. My simulation (See Appendix A) was written using ROOT

, a ”Data Analysis Framework” developed at CERN, designed to handle large amount of data from high energy physics and provide pre-written scripts for common processes.

Using established theories (Section 2) and both experimental and theoret- ical values (Table 2), I created a simulation of the creation of η

mesons from e

e

collisions and the subsequent decays of the η

mesons into the most likely decay products (Table 1). I focused on those final decay products because they are around 60% of the total possible decays of η

and because they consist of π

, π

, π

and γ, which are easy to detect.

Figure 3: KLOE detector [16]. The detector is incapable of detecting particles at θ < 15

and θ > 165

; the blue regions in the picture (The green regions are of no importence for this work).

2Can be found at: http://root.cern.ch/

I formed my simulation based on the capabilities of DAΦNE, a electron- positron collider at the INFN Laboratory in Frascati, Italy, operating since 1999. The main goal of the DAΦNE project is a ”study of neutral and charged kaons, intensively produced at the energy corresponding to the maximum of φ(1020) resonance” [14].

The reasons for choosing DAΦNE is because of the high integrated lumi- nosity which can be obtained at the collider and because of the energy used in the production. The DAΦNE energy is only 62M eV higher then the energy at E = √

s = m

= 958M eV (Section 2) that I chose in my simulation.

Using the luminosity from DAΦNE and assuming runtime of one year, I estimated an integrated luminosity of 1f b

, based on a daily integrated lu- minosity of 8pb

[15, 16]. This is a rough estimate of DAΦNE’s performance during a year, but it will give a picture of the available integrated luminosity for making this experiment.

According to theory expectations, only a tiny fraction of events would lead to the η

state. Most of the electron-positron collisions will lead to formation of other particles and even the same final decay products searched for without creating η

, i.e. e

e

→ π

π

η, which results into a large background one needs to account for. Very good discrimination of the background is needed, but this is outside the scope of this work.

Instead, I focused on phase space to generate the possible number of end decay products, such that one can quantitative see if one can detect any η

event. Cuts of the resulting data was done according to the dimensions and limitation of the main detector at DAΦNE, the KLOE detector (Fig. 3).

Due the beam tubes leading the electron-positron beams into the detector,

any particle which exit with an scattering angle θ < 15

and θ > 165

from

the center will not be detected. However, in actual experiment the angle

limit was θ < 22

and θ > 158

[16] and is what I followed. Furthermore, the

4 Pictures from the simulation

This section consists of figures that show results from the simulation. In every figure was B

= 3.0 ∗ 10

used, and every picture is from one run of the simulation. It can be seen from the figure that phase space coverage is quite good after the cuts. This means that the cuts can well be accounted for in the efficiency corrections.

It should be noted all particles among the end products must be detected for the event to be recorded. This means the difference between the number of generated particles and the number of recorded particles can be quite large, which could be a possible problem in executing an experiment on the decay.

Figure 4: Picture of the difference between generated and recorded [π

] respective

[γ] for the η

→ π

π

η → [π

] π

[γ] γ events. Full lines are the total number of

π

respective γ generated and dashed lines are those particles recorded after the

angle and energy cuts. The two higher, blue peaks are from π

and the two lower,

red ones are from γ.

Figure 5: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η → [π

] π

[γ] γ events. The upper image is of [π

] and the lower

one is of [γ]. Red lines correspond to the π

and γ particles which can be detected,

while the blue, shaded lines are of those particles which won’t be detected due to

the detector geometry (Sec. 3).

Figure 6: Picture of the difference between generated and recorded [π

] respective



π

for the η

→ π

π

η → [π

] π

π

π

π

events. Full lines are the total num-

ber of π

respective π

generated and dashed lines are those particles recorded after

the angle and energy cuts. The two higher, blue peaks are from π

and the two

lower, red ones are from π

.

Figure 7: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η → [π

] π

π

π

π

events. The upper image is of [π

] and the lower

one is of

π

. Red lines correspond to the π

and π

particles which can be de-

tected, while the blue, shaded lines are of those particles which won’t be detected

due to the detector geometry (Sec. 3).

Figure 8: Picture of the difference between generated and recorded [π

] respective



π

for the η

→ π

π

η → [π

] π

π

π

π

events. Full lines are the total num-

ber of π

respective π

generated and dashed lines are those particles recorded after

the angle and energy cuts. The two higher, blue peaks are from π

and the two

lower, red ones are from π

.

Figure 9: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η → [π

] π

π

π

π

events. The upper image is of [π

] and the lower

one is of

π

. Red lines correspond to the π

and π

particles which can be de-

tected, while the blue, shaded lines are of those particles which won’t be detected

due to the detector geometry (Sec. 3).

Figure 10: Picture of the difference between generated and recorded

π

respective

[γ] for the η

→ π

π

η →

π

π

[γ] γ events. Full lines are the total number of π

respective γ generated and dashed lines are those particles recorded after the angle

and energy cuts. The two higher, blue peaks are from π

and the two lower, red

ones are from γ.

Figure 11: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η →

π

π

[γ] γ events. The upper image is of

π

and the lower one

is of [γ]. Red lines correspond to the π

and γ particles which can be detected,

while the blue, shaded lines are of those particles which won’t be detected due to

the detector geometry (Sec. 3).

Figure 12: Picture of the difference between generated and recorded

π

respective

[π

] for the η

→ π

π

η →

π

π

[π

] π

π

events. Full lines are the total number

of π

respective π

generated and dashed lines are those particles recorded after

the angle and energy cuts. The two higher, blue peaks are from π

and the two

lower, red ones are from π

.

Figure 13: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η →

π

π

[π

] π

π

events. The upper image is of

π

and the lower

one is of [π

]. Red lines corresponds to the π

and π

particles which can be

detected, while the blue lines are of those particles which won’t be detected due to

the detector geometry (Sec. 3).

Figure 14: Picture of the difference between generated and recorded

π

respective



π

for the η

→ π

π

η →

π

π

π

π

π

events.

π

is one of the π

created

in the decay of η

and

π

is from the decay of the η. The full lines are the

total number of π

respective π

generated and the dashed lines are those particles

recorded after the angle and energy cuts. The two higher, blue peaks are from π

and the two lower, red ones are from π

.

Figure 15: Pictures depicting the angle cuts used in the simulation for the

η

→ π

π

η →

π

π

π

π

π

events. The upper image is of

π

and the lower

one is of

π

. The red lines correspond to the π

and π

particles which can be

detected, while the blue, shaded lines is of those particles which won’t be detected

due to the detector geometry (Sec. 3).

5 Results of the Simulation

By using the experimental limit and the theoretical value (Table 2) and Equation (3), I calculated for the expected number of events; the expected number of η

mesons created (Table 3 and 4). The first parts of Table 3 and 4, the η

→ π

π

η and the η

→ π

π

η part, are the sum of expected events, respectively, and the the rest are fraction of the corresponding sum, according to Table 1.

Branch Number of event (N

)

η

→ π

π

η (

) 68200(±450)

η

→ π

π

η → π

π

γγ 26900(±190) η

→ π

π

η → π

π

π

π

π

15600(±120) η

→ π

π

η → π

π

π

π

π

22300(±160)

η

→ π

π

η (

) 35300(±280)

η

→ π

π

η → π

π

γγ 13900(±110) η

→ π

π

η → π

π

π

π

π

8090(±65) η

→ π

π

η → π

π

π

π

π

11500(±92)

Table 3: Expected number of decay events according to Equation (3) with B

= 3.0 ∗ 10

and B

according to Table 1.

Branch Number of event N

η

→ π

π

η (

) 432(±3.0)

η

→ π

π

η → π

π

γγ 170(±1.2) η

→ π

π

η → π

π

π

π

π

98(±0.78) η

→ π

π

η → π

π

π

π

π

144(±1.0)

η

→ π

π

η (

) 223(±1.8)

η

→ π

π

η → π

π

γγ 88(±0.70) η

→ π

π

η → π

π

π

π

π

51(±0.41) η

→ π

π

η → π

π

π

π

π

73(±0.58)

Table 4: Expected number of decay events according to Equation (3) with

B

= 1.9 ∗ 10

and B

according to Table 1.

The simulation generates a number of particles determined by the branch- ing ratio for each branch, of which a fraction is capable of being detected (N

) due to limitation of the detector, e.g. due to the detector geometry. This number is lower than the expected N

and by dividing them, ones gets the acceptance for the decay ( = N

/N

). Table 5 shows the result of the simu- lation for B

= 1.9 ∗ 10

.

Branch  = N

/N

η

→ π

π

η −

η

→ π

π

η → π

π

γγ 69.5(±1.6)%

η

→ π

π

η → π

π

π

π

π

51.8(±1.5)%

η

→ π

π

η → π

π

π

π

π

54.3(±1.7)%

η

→ π

π

η −

η

→ π

π

η → π

π

γγ 70.5(±1.7)%

η

→ π

π

η → π

π

π

π

π

53.0(±2.1)%

η

→ π

π

η → π

π

π

π

π

54.6(±2.4)%

Table 5: Accaptance for detecting a specific decay chain, based on the geometry of the KLOE detector.

One can see that the η

→ π

π

η → π

π

γγ and η

→ π

π

η → π

π

γγ branches has the highest acceptance due to the few particles in them. There- fore, those are the particles usually focused on in experiments. However, earlier experiment have achieved a much lower acceptance, around 2.5% [7].

Using this acceptance, B

= 1.9 ∗ 10

and Eq. (3), the expected number

6 Discussion and Conclusion

The goal of this thesis was to write the necessary code and simulate an experiment of the e

e

→ η

→ X event at DAΦNE and to, based on the simulation, discuss the possibility of detecting the decay of η

at DAΦNE.

By using ROOT, I created a simulation of the η

event and I limited myself to make a geometrical study of the event.

The results from my simulation (Table 3 and 4), shows that even with the low Branching Ratio for the decay, there will still be a reasonable number of η

particles being created and end products being detected. Therefore, an actual experiment at DAΦNE would be constructive, since there is quite a large chance for observing decay product of the η

mesons, if the high luminosity I have used in my simulation can be provided.

However, this depends on my simulation and, further, on my calculation of the acceptance being correct. The latest experiment [7] had an acceptance of 2.5%, much lower compared to the the ones in Table 5. As seen in the last part of Results (Section5), by using this acceptance the expected number of events being recorded became around 4, i.e. so small one would expect not to detect any η

events at all. Not detecting any η

events would at least constrain the Upper Limit closer to the theoretical value (Table 2).

My work was limited to just a geometrical study, which could explain the

difference between my calculated acceptances (Table 5) and the one used in

the experiment [7]. A continuation of this work would be to include more

of the constraints of the event and the detection to make a more realistic

simulation of the event. An other continuation would be to investigate the

large background (the e

e

collisions which do not create η

) and the impact

it has on the experiment.

References

[1] B.R. Martin and G. Shaw. Particle Physics (3rd ed.), Manchester Physics Series, John Wiley & Sons. (2008). ISBN 978-0-470-03294-7 [2] E. Abouzaid, et al. Measurement of the rare decay pi0 −− > e+e-.

Phys.Rev. D75 (2007) 012004 arXiv:hep-ex/0610072

[3] Marcin Berlowski. Search for η −− > e+e- decay with the WASA experiment. EPJ Web Conf. 37 (2012) 09007 arXiv:1301.6040 [nucl-ex]

[4] CELSIUS/WASA Collaboration (M. Berlowski et al.). Measurement of eta meson decays into lepton-antilepton pairs. Phys.Rev. D77 (2008) 032004 arXiv:0711.3531 [hep-ex]

[5] G. Agakishiev, et al. Inclusive dielectron spectra in p+p collisions at 3.5 GeV. Eur.Phys.J. A48 (2012) 64 arXiv:1112.3607 [nucl-ex]

[6] P. V. Vorobev et al. Upper limits of the electron widths of the C-even mesons η

, f

(975), f

(1270), f

(1300), a

(980), and a

(1320).

Sov.J.Nucl.Phys. 48 (1988) 273-276, Yad.Fiz. 48 (1988) 436-441

[7] CMD-3 Collaboration (E. Solodov, et al.). Search for e

e

→ η(958) reaction with CMD-3 at VEPP-2000. Unpublished as of July 2014 [8] A.E. Dorokhov. Rare decay π

→ e

e

as a test for Standard Model.

Phys.Part.Nucl.Lett. 7 (2010) 229-234 arXiv:0905.4577 [hep-ph]

[9] L. Bergstr¨ om. Rare Decay of a Pseudoscalar Meson Into a Lepton Pair

- A Way to Detect New Interactions?. Z.Phys.C 14, (1982) 129-134

[13] A. E. Dorokhov, M. A. Ivanov and S. G. Kovalenko. Complete structure dependent analysis of the decay P → l

l. Phys.Lett. B677 (2009) 145- 149 arXiv:0903.4249 [hep-ph]

[14] Franzini, P. and M. Moulson, The Physics of DAFNE and KLOE.

Ann.Rev.Nucl.Part.Sci. 56 (2006) 207-251 arXiv:hep-ex/0606033

[15] M. Boscolo, et al. Measurement of the luminosity at the DAFNE collider upgraded with the crab waist scheme. arXiv:0909.1913 [hep-ex]

[16] V.P. Druzhinin, S.I. Eidelman, S.I. Serednyakov and E.P. Solodov.

Hadron Production via e+e- Collisions with Initial State Radiation.

Rev.Mod.Phys. 83 (2011) 1545 arXiv:1105.4975 [hep-ex]

Appendix A : Code used in Simulation

The code used in my thesis, right now made to create the figures in Section 4.

Make note that each line in the code ends with a ”;”, exempt the ones starting with #include, and that texts following a ”//” are just comments.

// Some p e r f o r m a n c e i s s u e s was f o u n d . Recommendation i s t o b o x // i n w i t h ” / ∗ ” . . . ” ∗ / ” t h e d e c a y modes w h i c h a r e u n i n t e r e s t i n g // a t t h e moment t o r e d u c e t h e r u n t i m e .

#include ”TMath . h”

#include ”TF1 . h”

#include ”TH2 . h”

#include ”TRandom . h”

#include ” T L o r e n t z V e c t o r . h”

#include ” TGenPhaseSpace . h”

#include ”TCanvas . h”

#include ”TGraph . h”

#include ” T S t y l e . h”

#include ” T F i l e . h”

#include ” T P r o f i l e . h”

#include ” TTree . h”

#include ”TRandom3 . h”

#include ” T A t t F i l l . h”

// The d i f f e r e n t p a r t i c l e masses [ GeV ] double mw = 0 . 7 8 2 ;

double meta = 0 . 5 4 7 8 6 2 ; double mpi = 0 . 1 3 9 5 7 0 ; double metap = 0 . 9 5 7 7 8 ; double mpi0 = 0 . 1 3 4 9 7 7 ; i n t L = 1 e6 ; // L u m i n o s i t y

double M = ( 4 ∗ ( TMath : : Pi ( ) ) / ( metap ∗ metap ) ) ∗ 0 . 3 8 9 3 7 9 ∗ 1 e6 ; // P a r t

T L o r e n t z V e c t o r ∗ e t o t = new T L o r e n t z V e c t o r ( 0 , 0 , 0 , 0 ) ; e t o t −>SetPxPyPzE ( 0 . , 0 . , 0 . , metap ) ;

//−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−

// TH1Ds c r e a t e d f o r e v e r y t h i n g s a v e d .

// The t i t l e s i n d i c a t e s w i c h d e c a y and t h e [ ] b r a c k e t s // marks t h e p a r t i c l e s w h i c h i s d i s p l a y d i n t h e end .

// e+e− −> e t a p r i m −> [ p i +] p i − e t a ( e t a −> [ gamma ] gamma) TH1D ∗ h 1 1 1 p i = new TH1D( ”#p i ˆ{+} 111 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ;

TH1D ∗ h 1 1 1 p i u = new TH1D( ”#p i ˆ{+} 111 ( u ) ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 1 1 g = new TH1D( ”#gamma 111 ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ;

TH1D ∗ h 1 1 1 g u = new TH1D( ”#gamma 1 1 1 ( u ) ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 1 1 p i t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 111 ” , ”

” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 1 p i t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 1 1 1 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 1 g t h e t a = new TH1D( ” c o s (# t h e t a (#gamma) ) 111 ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 1 g t h e t a u = new TH1D( ” c o s (# t h e t a (#gamma) ) 1 1 1 ( u )

” , ” ” , 1 0 0 , − 1 , 1 ) ;

// e+e− −> e t a p r i m −> [ p i +] p i − e t a ( e t a −> p i+ p i − [ p i 0 ] ) TH1D ∗ h 1 1 2 p i = new TH1D( ”#p i ˆ{+} 112 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 1 2 p i u = new TH1D( ”#p i ˆ{+} 1 1 2 ( u ) ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 1 2 p i 0 = new TH1D( ”#p i ˆ{0} 112 ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 1 2 p i 0 u = new TH1D( ”#p i ˆ{0} 1 1 2 ( u ) ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 1 2 p i t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 112 ” , ”

” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 2 p i t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 1 1 2 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 2 p i 0 t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 112 ” ,

” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 2 p i 0 t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 1 1 2 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

// e+e− −> e t a p r i m −> [ p i +] p i − e t a ( e t a −> [ p i 0 ] p i 0 p i 0 ) TH1D ∗ h 1 1 3 p i = new TH1D( ”#p i ˆ{+} 113 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 1 3 p i u = new TH1D( ”#p i ˆ{+} 1 1 3 ( u ) ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 1 3 p i 0 = new TH1D( ”#p i ˆ{0} 113 ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 1 3 p i 0 u = new TH1D( ”#p i ˆ{0} 1 1 3 ( u ) ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 1 3 p i t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 113 ” , ”

” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 3 p i t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 1 1 3 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 3 p i 0 t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 113 ” ,

” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 1 3 p i 0 t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 1 1 3 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

// e+e− −>e t a p r i m −> [ p i 0 ] p i 0 e t a ( e t a −> [ gamma ] gamma) TH1D ∗ h 1 2 1 p i 0 = new TH1D( ”#p i ˆ{0} 121 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 2 1 p i 0 u = new TH1D( ”#p i ˆ{0} 1 2 1 ( u ) ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 2 1 g = new TH1D( ”#gamma 121 ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ;

TH1D ∗ h 1 2 1 g u = new TH1D( ”#gamma 1 2 1 ( u ) ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 2 1 p i 0 t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 121 ” ,

” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 1 p i 0 t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 1 2 1 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 1 g t h e t a = new TH1D( ” c o s (# t h e t a (#gamma) ) 121 ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 1 g t h e t a u = new TH1D( ” c o s (# t h e t a (#gamma) ) 1 2 1 ( u )

” , ” ” , 1 0 0 , − 1 , 1 ) ;

// e+e− −>e t a p r i m −> [ p i 0 ] p i 0 e t a ( e t a −> [ p i +] p i − p i 0 ) TH1D ∗ h 1 2 2 p i 0 = new TH1D( ”#p i ˆ{0} 122 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 2 2 p i 0 u = new TH1D( ”#p i ˆ{0} 1 2 2 ( u ) ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 2 2 p i = new TH1D( ”#p i ˆ{+} 122 ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 2 2 p i u = new TH1D( ”#p i ˆ{+} 1 2 2 ( u ) ” , ” ” , 1 0 0 , 0 . , 0 . 5 ) ; TH1D ∗ h 1 2 2 p i 0 t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 122 ” ,

” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 2 p i 0 t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ { 0 } ) ) 1 2 2 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 2 p i t h e t a = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 122 ” , ”

” , 1 0 0 , − 1 , 1 ) ;

TH1D ∗ h 1 2 2 p i t h e t a u = new TH1D( ” c o s (# t h e t a (# p i ˆ{+}) ) 1 2 2 ( u ) ” , ” ” , 1 0 0 , − 1 , 1 ) ;

// e+e− −>e t a p r i m −> [ p i 0 ] p i 0 e t a ( e t a −> p i 0 p i 0 p i 0 )

TH1D ∗ h 1 2 3 p i 0 1 = new TH1D( ”#p i ˆ{0} {1} 123 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ; TH1D ∗ h 1 2 3 p i 0 1 u = new TH1D( ”#p i ˆ{0} {1} 1 2 3 ( u ) ” , ” ”

, 1 0 0 , 0 , 0 . 5 ) ;

TH1D ∗ h 1 2 3 p i 0 2 = new TH1D( ”#p i ˆ{0} {2} 123 ” , ” ” , 1 0 0 , 0 , 0 . 5 ) ;

// C r e a t e s t h e p h a s e s p a c e f o r t h e d i f f e r e n t d e c a y s . // They a r e numbered i n o r d e r and how t h e y f o l l o w s , e . g . // d e c a y nr 121 ( p i 0 p i 0 e t a −> p i 0 p i 0 gamma gamma) s t a r t s // w i t h d e c a y nr 1 ( e+ e− −> e t a p r i m ) and f o l l o w s by d e c a y // nr 12 ( e t a p r i m −> p i 0 p i 0 e t a ) .

// Decay nr 1

TGenPhaseSpace ev1 ;

double dm1 [ 2 ] = { 0 . 0 , 0 . 0 } ; // Decay nr 11

TGenPhaseSpace ev11 ;

double dm11 [ 3 ] = {mpi , mpi , meta } ; // Decay nr 111

TGenPhaseSpace ev111 ;

double dm111 [ 2 ] = { 0 . 0 , 0 . 0 } ; // Decay nr 112

TGenPhaseSpace ev112 ;

double dm112 [ 3 ] = {mpi , mpi , mpi0 } ; // Decay nr 113

TGenPhaseSpace ev113 ;

double dm113 [ 3 ] = { mpi0 , mpi0 , mpi0 } ; // Decay nr 12

TGenPhaseSpace ev12 ;

double dm12 [ 3 ] = { mpi0 , mpi0 , meta } ; // Decay nr 121

TGenPhaseSpace ev121 ;

double dm121 [ 2 ] = { 0 . 0 , 0 . 0 } ; // Decay nr 122

TGenPhaseSpace ev122 ;

double dm122 [ 3 ] = {mpi , mpi , mpi0 } ; // Decay nr 123

TGenPhaseSpace ev123 ;

double dm123 [ 3 ] = { mpi0 , mpi0 , mpi0 } ;

// The B r a n c h i n g r a t i o s f o r t h e d i f f e r e n t decay , u s e d h e r e a s // l i m i t s f o r t h e l o o p s f u r t h e r down .

double d e c a y l i m i t 1 u p p e r = nev ∗ ( 0 . 4 2 9 + 0 . 2 2 2 ) ; double d e c a y l i m i t 1 1 u p p e r = nev ∗ 0 . 4 2 9 ;

double d e c a y l i m i t 1 1 1 u p p e r = nev ∗ 0 . 4 2 9 ∗ 0 . 3 9 4 1 ; double d e c a y l i m i t 1 1 2 l o w e r = d e c a y l i m i t 1 1 1 u p p e r ; double d e c a y l i m i t 1 1 2 u p p e r = nev

∗ 0 . 4 2 9 ∗ ( 0 . 3 9 4 1 + 0 . 2 2 9 2 ) ;

double d e c a y l i m i t 1 1 3 l o w e r = d e c a y l i m i t 1 1 2 u p p e r ; double d e c a y l i m i t 1 1 3 u p p e r = nev

∗ 0 . 4 2 9 ∗ ( 0 . 3 9 4 1 + 0 . 2 2 9 2 + 0 . 3 2 6 8 ) ;

double d e c a y l i m i t 1 2 l o w e r = d e c a y l i m i t 1 1 u p p e r ;

double d e c a y l i m i t 1 2 u p p e r = nev ∗ ( 0 . 4 2 9 + 0 . 2 2 2 ) ; double d e c a y l i m i t 1 2 1 u p p e r = nev ∗0.429+ nev

∗ 0 . 2 2 2 ∗ 0 . 3 9 4 1 ;

double d e c a y l i m i t 1 2 2 l o w e r = d e c a y l i m i t 1 2 1 u p p e r ; double d e c a y l i m i t 1 2 2 u p p e r = nev ∗0.429+ nev

∗ 0 . 2 2 2 ∗ ( 0 . 3 9 4 1 + 0 . 2 2 9 2 ) ;

double d e c a y l i m i t 1 2 3 l o w e r = d e c a y l i m i t 1 2 2 u p p e r ; double d e c a y l i m i t 1 2 3 u p p e r = nev ∗0.429+ nev

∗ 0 . 2 2 2 ∗ ( 0 . 3 9 4 1 + 0 . 2 2 9 2 + 0 . 3 2 6 8 ) ;

// The a n g l e l i m i t o f t h e d e t e c t o r KLOE, u s e d t o make // an a n g l e c u t among t h e p a r t i c l e s .

double t h e t a m i n = 2 2 ∗ (TMath : : Pi ( ) ) / 1 8 0 ; double thetamax = 1 5 8 ∗ (TMath : : Pi ( ) ) / 1 8 0 ;

//−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−

// The l o o p o f t h e d e c a y s . Every e t a p r i m e x p e c t e d a r e r e v i e w e d // and i s s i m u l a t e d t o d e t e r m i n e d i f i t becomes d e t e c t e d . // Again , b o x i n w i t h ” / ∗ ” . . . ” ∗ / ” t h e d e c a y modes w h i c h a r e // u n i n t e r e s t i n g a t t h e moment .

i n t n ; // Decay nr 1

f o r ( n = 1 ; n <= nev ; n++) {

i n t r = gRandom−>Uniform ( 1 , nev ) ; ev1 . SetDecay ( ∗ e t o t , 2 , dm1) ;

double w1 = ev1 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ g 1 = ev1 . GetDecay ( 0 ) ; T L o r e n t z V e c t o r ∗ g 2 = ev1 . GetDecay ( 1 ) ; e t a p r i m =(∗ g 1 ) +(∗ g 2 ) ;

//−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−

// Decay nr 11

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗ p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i p l u s 1 1 ) . E ( ) −(∗ p i p l u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i m i n u s 1 1 ) . Theta ( ) && ( ∗ p i m i n u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i m i n u s 1 1 ) . E ( ) −(∗ p i m i n u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ g 1 1 1 1 ) . Theta ( ) && ( ∗ g 1 1 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n i f ( t h e t a m i n <= ( ∗ g 1 1 1 2 ) . Theta ( ) && ( ∗ g 1 1 1 2 )

. Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n h 1 1 1 p i −> F i l l ( ( ∗ p i p l u s 1 1 ) . E ( ) , w1∗w11 ) ; h111 g−> F i l l ( ( ∗ g 1 1 1 1 ) . E ( ) , w1∗w11∗w111 ) ; }}}}}}

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗ p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 1 1 p i t h e t a −> F i l l ( c o s ( ( ∗ p i p l u s 1 1 ) . Theta ( ) ) , w1∗w11 ) ;

}

i f ( t h e t a m i n <= ( ∗ g 1 1 1 1 ) . Theta ( ) && ( ∗ g 1 1 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 1 1 g t h e t a −> F i l l ( c o s ( ( ∗ g 1 1 1 1 ) . Theta ( ) ) , w1∗w11∗w111 ) ;

}

h 1 1 1 p i u −> F i l l ( ( ∗ p i p l u s 1 1 ) . E ( ) , w1∗w11 ) ; h 1 1 1 g u −> F i l l ( ( ∗ g 1 1 1 1 ) . E ( ) , w1∗w11∗w111 ) ; h 1 1 1 p i t h e t a u −> F i l l ( c o s ( ( ∗ p i p l u s 1 1 ) . Theta ( ) ) ,

w1∗w11 ) ;

h 1 1 1 g t h e t a u −> F i l l ( c o s ( ( ∗ g 1 1 1 1 ) . Theta ( ) ) , w1∗

w11∗w111 ) ; }

// Decay nr 112

i f ( d e c a y l i m i t 1 1 2 l o w e r < r && r <=

d e c a y l i m i t 1 1 2 u p p e r ) ) {

ev112 . SetDecay ( ∗ e t a 1 1 , 3 , dm112 ) ; double w112 = ev112 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ p i p l u s 1 1 2 = ev112 . GetDecay ( 0 ) ; T L o r e n t z V e c t o r ∗ p i m i n u s 1 1 2 = ev112 . GetDecay ( 1 ) ; T L o r e n t z V e c t o r ∗ p i z e r o 1 1 2 = ev112 . GetDecay ( 2 ) ;

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗ p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i p l u s 1 1 ) . E ( ) −(∗ p i p l u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i m i n u s 1 1 ) . Theta ( ) && ( ∗ p i m i n u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i m i n u s 1 1 ) . E ( ) −(∗ p i m i n u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 2 ) . Theta ( ) && ( ∗ p i p l u s 1 1 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f (0.01 <= ( ( ∗ p i p l u s 1 1 2 ) . E ( ) −(∗ p i p l u s 1 1 2 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i m i n u s 1 1 2 ) . Theta ( ) && ( ∗ p i m i n u s 1 1 2 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i m i n u s 1 1 2 ) . E ( ) −(∗ p i m i n u s 1 1 2 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 1 2 ) . Theta ( ) && ( ∗ p i z e r o 1 1 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 1 2 ) . E ( ) −(∗ p i z e r o 1 1 2 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

h 1 1 2 p i −> F i l l ( ( ∗ p i p l u s 1 1 ) . E ( ) , w1∗w11 ) ; h 1 1 2 p i 0 −> F i l l ( ( ∗ p i z e r o 1 1 2 ) . E ( ) , w1∗w11∗w112

) ; }}}}}}}}}}

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗ p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 1 2 p i t h e t a −> F i l l ( c o s ( ( ∗ p i p l u s 1 1 ) . Theta ( )

h 1 1 2 p i 0 t h e t a u −> F i l l ( c o s ( ( ∗ p i z e r o 1 1 2 ) . Theta ( ) ) , w1∗w11∗w112 ) ;

}

// Decay nr 113

i f ( d e c a y l i m i t 1 1 3 l o w e r < r && r <=

d e c a y l i m i t 1 1 3 u p p e r ) {

ev113 . SetDecay ( ∗ e t a 1 1 , 3 , dm113 ) ; double w113 = ev113 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ p i z e r o 1 1 3 1 = ev113 . GetDecay ( 0 )

;

T L o r e n t z V e c t o r ∗ p i z e r o 1 1 3 2 = ev113 . GetDecay ( 1 )

;

T L o r e n t z V e c t o r ∗ p i z e r o 1 1 3 3 = ev113 . GetDecay ( 2 )

;

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗ p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i p l u s 1 1 ) . E ( ) −(∗ p i p l u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i m i n u s 1 1 ) . Theta ( ) && ( ∗ p i m i n u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f ( 0 . 0 1 <= ( ( ∗ p i m i n u s 1 1 ) . E ( ) −(∗ p i m i n u s 1 1 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) && ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 1 3 1 ) . E ( ) −(∗ p i z e r o 1 1 3 1 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 1 3 2 ) . Theta ( ) && ( ∗ p i z e r o 1 1 3 2 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 1 3 2 ) . E ( ) −(∗ p i z e r o 1 1 3 2 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 1 3 3 ) . Theta ( ) && ( ∗ p i z e r o 1 1 3 3 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 1 3 3 ) . E ( ) −(∗ p i z e r o 1 1 3 3 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

h 1 1 3 p i −> F i l l ( ( ∗ p i p l u s 1 1 ) . E ( ) , w1∗w11 ) ; h 1 1 3 p i 0 −> F i l l ( ( ∗ p i z e r o 1 1 3 1 ) . E ( ) , w1∗w11∗

w113 ) ; }}}}}}}}}}

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 1 ) . Theta ( ) && ( ∗

p i p l u s 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 1 3 p i t h e t a −> F i l l ( c o s ( ( ∗ p i p l u s 1 1 ) . Theta ( ) ) , w1∗w11 ) ;

}

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) && ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 1 3 p i 0 t h e t a −> F i l l ( c o s ( ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) ) , w1∗w11∗w113 ) ;

}

h 1 1 3 p i u −> F i l l ( ( ∗ p i p l u s 1 1 ) . E ( ) , w1∗w11 ) ;

h 1 1 3 p i 0 u −> F i l l ( ( ∗ p i z e r o 1 1 3 1 ) . E ( ) , w1∗w11∗w113 ) ;

h 1 1 3 p i t h e t a u −> F i l l ( c o s ( ( ∗ p i p l u s 1 1 ) . Theta ( ) ) , w1∗w11 ) ;

h 1 1 3 p i 0 t h e t a u −> F i l l ( c o s ( ( ∗ p i z e r o 1 1 3 1 ) . Theta ( ) ) , w1∗w11∗w113 ) ;

} }

//−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−

// Decay nr 12

i f ( d e c a y l i m i t 1 2 l o w e r < r && r <= d e c a y l i m i t 1 2 u p p e r ) {

ev12 . SetDecay ( e t a p r i m , 3 , dm12 ) ; double w12 = ev12 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ p i z e r o 1 2 1 = ev12 . GetDecay ( 0 ) ; T L o r e n t z V e c t o r ∗ p i z e r o 1 2 2 = ev12 . GetDecay ( 1 ) ; T L o r e n t z V e c t o r ∗ e t a 1 2 = ev12 . GetDecay ( 2 ) ; // Decay nr 121

i f ( r <= d e c a y l i m i t 1 2 1 u p p e r ) { ev121 . SetDecay ( ∗ e t a 1 2 , 2 , dm121 ) ; double w121 = ev121 . G e n e r a t e ( ) ;

( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ g 1 2 1 1 ) . Theta ( ) && ( ∗ g 1 2 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f ( t h e t a m i n <= ( ∗ g 1 2 1 2 ) . Theta ( ) && ( ∗ g 1 2 1 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

h 1 2 1 p i 0 −> F i l l ( ( ∗ p i z e r o 1 2 1 ) . E ( ) , w1∗w12 ) ; h121 g−> F i l l ( ( ∗ g 1 2 1 1 ) . E ( ) , w1∗w12∗w121 ) ; }}}}}}

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 2 1 p i 0 t h e t a −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 1 ) . Theta ( ) ) , w1∗w12 ) ;

}

i f ( t h e t a m i n <= ( ∗ g 1 2 1 1 ) . Theta ( ) && ( ∗ g 1 2 1 1 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

h 1 2 1 g t h e t a −> F i l l ( c o s ( ( ∗ g 1 2 1 1 ) . Theta ( ) ) , w1∗w12∗w121 ) ;

}

h 1 2 1 p i 0 u −> F i l l ( ( ∗ p i z e r o 1 2 1 ) . E ( ) , w1∗w12 ) ; h 1 2 1 g u −> F i l l ( ( ∗ g 1 2 1 1 ) . E ( ) , w1∗w12∗w121 ) ; h 1 2 1 p i 0 t h e t a u −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 1 ) . Theta

( ) ) , w1∗w12 ) ;

h 1 2 1 g t h e t a u −> F i l l ( c o s ( ( ∗ g 1 2 1 1 ) . Theta ( ) ) , w1∗

w12∗w121 ) ; }

// Decay nr 122

i f ( d e c a y l i m i t 1 2 2 l o w e r < r && r <=

d e c a y l i m i t 1 2 2 u p p e r ) {

ev122 . SetDecay ( ∗ e t a 1 2 , 3 , dm122 ) ; double w122 = ev122 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ p i p l u s 1 2 2 = ev122 . GetDecay ( 0 ) ; T L o r e n t z V e c t o r ∗ p i m i n u s 1 2 2 = ev122 . GetDecay ( 1 ) ; T L o r e n t z V e c t o r ∗ p i z e r o 1 2 2 = ev122 . GetDecay ( 2 ) ; i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 1 ) . Theta ( ) && ( ∗

p i z e r o 1 2 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 1 ) . E ( ) −(∗ p i z e r o 1 2 1 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 2 ) . Theta ( ) && ( ∗ p i z e r o 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle

f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 2 ) . E ( ) −(∗ p i z e r o 1 2 2 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 2 2 ) . Theta ( ) && ( ∗ p i p l u s 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f (0.01 <= ( ( ∗ p i p l u s 1 2 2 ) . E ( ) −(∗ p i p l u s 1 2 2 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i m i n u s 1 2 2 ) . Theta ( ) && ( ∗ p i m i n u s 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i m i n u s 1 2 2 ) . E ( ) −(∗ p i m i n u s 1 2 2 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 2 ) . Theta ( ) && ( ∗ p i z e r o 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 2 ) . E ( ) −(∗ p i z e r o 1 2 2 ) .M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

h 1 2 2 p i 0 −> F i l l ( ( ∗ p i z e r o 1 2 1 ) . E ( ) , w1∗w12 ) ; h 1 2 2 p i −> F i l l ( ( ∗ p i p l u s 1 2 2 ) . E ( ) , w1∗w12∗w122 )

; }}}}}}}}}}

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 2 2 p i 0 t h e t a −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 1 ) . Theta ( ) ) , w1∗w12 ) ;

}

i f ( t h e t a m i n <= ( ∗ p i p l u s 1 2 2 ) . Theta ( ) && ( ∗ p i p l u s 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle f o r

d e t e c t i o n

h 1 2 2 p i t h e t a −> F i l l ( c o s ( ( ∗ p i p l u s 1 2 2 ) . Theta ( ) ) , w1∗w12∗w122 ) ;

double w123 = ev123 . G e n e r a t e ( ) ;

T L o r e n t z V e c t o r ∗ p i z e r o 1 2 3 1 = ev123 . GetDecay ( 0 )

;

T L o r e n t z V e c t o r ∗ p i z e r o 1 2 3 2 = ev123 . GetDecay ( 1 )

;

T L o r e n t z V e c t o r ∗ p i z e r o 1 2 3 3 = ev123 . GetDecay ( 2 )

;

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 1 ) . E ( ) −(∗ p i z e r o 1 2 1 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 2 ) . Theta ( ) && ( ∗ p i z e r o 1 2 2 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 2 ) . E ( ) −(∗ p i z e r o 1 2 2 ) .M ( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 3 1 ) . E ( ) −(∗ p i z e r o 1 2 3 1 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 3 2 ) . Theta ( ) && ( ∗ p i z e r o 1 2 3 2 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 3 2 ) . E ( ) −(∗ p i z e r o 1 2 3 2 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 3 3 ) . Theta ( ) && ( ∗ p i z e r o 1 2 3 3 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

i f (0.01 <= ( ( ∗ p i z e r o 1 2 3 3 ) . E ( ) −(∗ p i z e r o 1 2 3 3 ) . M( ) ) ) { // S m a l l e s t d e t e c t a b l e e n e r g y

h 1 2 3 p i 0 1 −> F i l l ( ( ∗ p i z e r o 1 2 1 ) . E ( ) , w1∗w12 ) ; h 1 2 3 p i 0 2 −> F i l l ( ( ∗ p i z e r o 1 2 3 1 ) . E ( ) , w1∗w12∗

w123 ) ; }}}}}}}}}}

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 1 ) . Theta ( ) <= thetamax ) { // Angle f o r d e t e c t i o n

h 1 2 3 p i 0 1 t h e t a −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 1 ) . Theta ( ) ) , w1∗w12 ) ;

}

i f ( t h e t a m i n <= ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) && ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) <= thetamax ) { // Angle

f o r d e t e c t i o n

h 1 2 3 p i 0 2 t h e t a −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) ) , w1∗w12∗w123 ) ;

}

h 1 2 3 p i 0 1 u −> F i l l ( ( ∗ p i z e r o 1 2 1 ) . E ( ) , w1∗w12 ) ; h 1 2 3 p i 0 2 u −> F i l l ( ( ∗ p i z e r o 1 2 3 1 ) . E ( ) , w1∗w12∗

w123 ) ;

h 1 2 3 p i 0 1 t h e t a u −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 1 ) . Theta ( ) ) , w1∗w12 ) ;

h 1 2 3 p i 0 2 t h e t a u −> F i l l ( c o s ( ( ∗ p i z e r o 1 2 3 1 ) . Theta ( ) ) , w1∗w12∗w123 ) ;

} } }

//−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−

// C r e a t i n g h i s t o g r a m s g S t y l e −>S e t O p t S t a t ( 0 ) ;

// e+e− −>e t a p r i m −>[ p i +] p i −e t a ( e t a −>[gamma ] gamma) TCanvas ∗ c 1 1 = new TCanvas ( ” c 1 1 ” , ” c 1 1 ” ) ;

h 1 1 1 p i u −>S e t T i t l e ( ”#p i ˆ{+} Momentum and #gamma Energy ; E [ GeV ] ; ” ) ;

h 1 1 1 p i u −>Draw ( ) ;

h 1 1 1 p i −>S e t L i n e S t y l e ( 7 ) ; h 1 1 1 p i −>Draw ( ” same ” ) ;

h 1 1 1 g u −>S e t L i n e C o l o r ( kRed ) ; h 1 1 1 g u −>Draw ( ” same ” ) ;

h111 g−>S e t L i n e C o l o r ( kRed ) ; h111 g−>S e t L i n e S t y l e ( 7 ) ; h111 g−>Draw ( ” same ” ) ;

TCanvas ∗ c 1 2 = new TCanvas ( ” c 1 2 ” , ” c 1 2 ” ) ; c12−>D i v i d e ( 1 , 2 ) ;

c12−>cd ( 1 ) ;